Introduction

Obesity has become a global health problem, affecting approximately 13% of the world’s adult population, and over 340 million children and adolescents (WHO). This epidemic has profound implications not only for reproductive health but also for the well-being of subsequent generations. Previous studies have demonstrated that maternal obesity reduces the function of the hypothalamic-pituitary-ovarian (HPO) axis (X. Chen, Xiao, Cai, Huang, & Chen, 2022), deteriorates oocyte cytoplasmic quality and nuclear maturation, and disrupts genome methylation (Broughton & Moley, 2017). Reduced expression of Stella in oocytes induced by obesity results in global hypo-methylation in zygotes, which has an important contribution to the defective embryo development (Han et al., 2018). Furthermore, progeny of obese females have a higher risk of non-communicable diseases, such as obesity, diabetes, and cardiovascular diseases (Godfrey et al., 2017). Our previous study indicated that obesity altered the methylation status of Leptin, which might play a role in the metabolic disorders of female offspring. However, these methylation changes are not detected in the first- generation (F1) oocytes (Ge et al., 2014). The influence of maternal obesity on the genomic methylation of oocytes is still obscure. Thus, more studies are necessary to explore the role of DNA methylation in mediating the transgenerational transmission of metabolic syndrome induced by female obesity.

Obesity also perturbs glucose and lipid metabolism, which has negative effects on oocyte maturation and embryo development. Elevated levels of circulating free fatty acids in obese females contribute to oocyte lipotoxicity, while concomitantly diminishing mitochondrial function within oocytes (Broughton & Moley, 2017). The decreased melatonin levels are also reported in animals and humans by previous studies (Overberg et al., 2022; Virto et al., 2018). Melatonin has contributions to metabolism and DNA methylation in somatic cells (Davoodvandi, Nikfar, Reiter, & Asemi, 2022). These findings indicate that the disturbed metabolism induced by obesity is closely linked to the compromised oocyte quality and aberrant methylation patterns. This study aims to elucidate the effects of obesity-induced metabolic changes on oocyte genomic methylation and its inhereditary implications in a murine model.

Results

Obesity alters the genomic methylation in oocytes

The obese mouse model was induced via a high-fat diet (HFD) (Ge et al., 2014; Han et al., 2018), and mice fed with standard diet were used as a control (CD). The HFD group exhibited a significantly higher average body weight compared to CD group (n>86, Fig. S1A and B). Re-methylation in oocytes occurs in follicular development and is nearly complete at the germinal vesicle (GV) stage. The 5mC (5-methylcytosine) and 5hmC (5-hydroxymethylcytosine) levels in the GV oocytes of HFD mice were significantly higher than that in the CD group (n>30, Fig. 1A- C). To further explore the effects of maternal obesity on the oocyte methylation, we examined the genomic methylation of metaphase II (MII) oocytes using whole-genome bisulfite sequencing for small samples (WGBS, Novogene, Beijing, China). The information of reads count, mapping rate, conversion rate and sequencing depth was presented in Table S1. We found that global methylation in MII oocytes of HFD group was higher than that in the CD group (Fig. 1D). Methylated cytosine (C) can be classed into three types according to the context in the genome including CG, CHG and CHH (H=A, T or C). The methylated CG have more contributions to regulate gene expression. We found that the CG methylation level in MII oocytes of HFD was significantly higher than that in the CD group (Fig. 1E). Differentially methylated CG distributed across all chromosomes (Fig. S1C). To further analyze the distribution of methylation, each functional region of genes was equally divided into 20 bins, and then the average methylation levels in the functional regions were calculated, respectively. CGIs (CG islands) and CGI shores were predicted using cpgIslandExt and repeat sequences were predicted using RepeatMasker. Results showed that the hyper-methylation was distributed in the promoter, exon, upstream 2k, and downstream 2k regions (Fig. 1F and Fig. S1D) of genes in HFD oocytes. These findings suggest that maternal obesity results in hyper-genome methylation in oocytes.

Maternal obesity alters the DNA methylation of oocytes.

(A) Methylation levels of 5mC and 5hmC in oocytes (n>30). 5mC, 5- methylcytosine; 5hmC, 5-hydroxymethylcytosine; DAPI, chromatin.

(B, C) Relative fluorescence intensity of 5mC and 5hmC in GV oocytes.

(D) Genomic methylation level of MII oocytes examined by single-cell whole genome bisulfite sequencing. The control group (CD) has two replicates, and the obesity group (HFD) has three replicates.

(E) Average genomic CG methylation level in MII oocytes. CD, control group; HFD, obesity group; ** means p value < 0.01.

(F) CG methylation levels at different regions in MII oocytes. CGI, CpG island; utr5, 5’ untranslated region; utr3, 3’ untranslated region; repeat, repeat sequence.

(G) Total differentially methylated regions (DMRs) in oocytes of control and obesity groups. Hyper-DMRs, hypermethylated DMRs; hypo-DMRs, hypomethylated DMRs.

(H) Distribution of DMRs on chromosomes in MII oocytes. Outside-to-in: chromosomes, hyper-DMRs, TEs (transcription end regions), and gene, hypo-DMRs.

(I) KEGG pathway enrichment of genes with DMRs at the promoter regions, and the top 20 pathways are presented.

(J) Schedule of breeding. Female C57BL/6 mice fed with normal (CD) or high-fat diet (HFD) for 12 weeks were marked as F0. F1 was produced by F0 mated with normal males, respectively, and marked as CF1 and HF1; F2 was produced by female F1 mated with normal males and marked as CF2 and HF2, respectively.

Distribution of differentially methylated regions (DMRs)

We further analyzed the differentially methylated regions (DMRs) in oocytes, and identified 4340 DMRs between HFD and CD oocytes. These DMRs were defined by the following criteria: the number of CGs ≥ 4 and the absolute methylation difference ≥ 0.2. Among these, 2,013 were hyper- DMRs (46.38%), and 2327 were hypo-DMRs (53.62%) (Fig. 1G). These DMRs were distributed across all chromosomes (Fig. 1H). We then annotated DMRs into different genomic regions including promoter, exon, intron, CGI, CGI shore, repeat, TSS (transcription start site), TES (transcription end site), UTR3 (3’ end untranslated region), and UTR5 (5’ end untranslated region) regions (Fig. S2A), and the average methylation levels of DMRs in these regions were similar between HFD and CD oocytes (Fig. S2B). However, the hypo-DMRs were enriched in the UTR3, repeat, and intron regions compared with hyper-DMRs (Fig. S3).

Methylation level at promoters strongly contributes to the regulation of gene expression. We then analyzed the enrichment of genes with DMRs at promoters in KEGG (Kyoto Encyclopedia of Genes and Genomes) pathways using KOBAS online. Results indicated that the genes with DMRs at promoters were significantly enriched in metabolic pathways including amino acid metabolism pathways, carbohydrate metabolism pathways, lipid metabolism pathways, and metabolism of cofactors and vitamins pathways (Fig. 1I, Table S2). A total of 35 genes with DMRs at promoters were included in metabolism pathways, and 19 of these genes were with hyper-DMRs and 16 of these genes were with hypo-DMRs (Table S3). These results suggest that the altered methylation in oocytes induced by maternal obesity may play a role in the metabolic disorders in offspring.

The disturbed methylation may be associated with the transgenerational inheritance of the metabolic disorders through females

Our recent study demonstrated that disturbed methylation in oocytes caused by uterine undernourishment can be partly transmitted to F2 oocytes via females, which may play a key role in the transgenerational inheritance of metabolic disorders (S. B. Tang et al., 2023). Here, we investigated the inheritance of altered methylation in HFD oocytes. F1 and F2 generations were produced as shown in the schedule in Fig. 1J: HF1 and CF1, female HFD and CD, respectively, were mated with control males; HF2 and CF2, female HF1 and CF1, respectively, were mated with control males. We examined the glucose and insulin tolerance (GTT and ITT), and found that the GTTs and ITTs of F0, F1, and F2 females were impaired (Fig. 2A-C). The inheritance of disrupted metabolism in females might be associated with the altered DNA methylation of oocytes. To address this question, we first examined the DNA methylation status of DMRs located at the promoters of Bhlha15 (also known as Mist1, basic helix-loop-helix, a transcription factor), Mgat1 (mannoside acetylglucosaminyltransferase 1), Taok3 (serine/threonine-protein kinase 3), Tkt (transketolase), Pik3cd (phosphatidylinositol-4, 5-bisphosphate 3- kinase catalytic subunit delta), and Pld1 (phospholipase D1) in the HFD and CD oocytes. We found that the methylation levels of hyper-DMRs including Bhlha15-DMR, Mgat1-DMR, and Taok3-DMR in HFD oocytes were significantly higher than in the CD group (Fig. 2D). The methylation level of hypo-DMRs, Pik3cd-DMR in HFD oocytes was significantly lower than that in the CD group (Fig. 2D). These results coincide with the genomic sequencing results. Whereas, the methylation level of Tkt-DMR (hypo-DMR) in HFD oocytes was higher than that in the CD group, and the methylation level of Pld1-DMR (hypo-DMR) was similar between the two groups, which contradicts the genomic sequencing results (Fig. S4A). These findings suggest that some regions are false positives in genomic sequencing. To exclude the effects of somatic cell contamination, we examined the methylation level of paternally imprinted gene H19, which was low in both HFD and CD oocytes (Fig. S4B). These findings indicate that the samples are not contaminated by somatic cells.

Transgenerational inheritance of metabolic disorders and altered DNA methylation.

(A-C) Glucose tolerance (GTT) and insulin tolerance ITT) were tested for female F0, F1, and F2, respectively. * p<0.05; ** p<0.01.

(D-F) DMR methylation at the promoter regions of Bhlha15, Mgat1, Taok3, Tkt, and Pid3cd in F0, F1, and F2 oocytes was respectively examined using bisulfite sequencing. At least 10 available clones from 80-100 oocytes were used to calculate the methylation level. White circle, unmethylated CG; black circle, methylated CG. * p<0.05; ** p<0.01.

(G) Inheritance of altered methylation in different generations was analyzed. * p<0.05; ** p<0.01.

We then examined the methylation of DMRs in F1 livers using bisulfite sequencing (BS). Ten livers from five litters were analyzed for each group. Results revealed that the methylation levels of Bhlha15-DMR and Mgat1- DMR were higher and the methylation level of Pik3cd-DMR was lower in HF1 livers than that in CF1 livers (Fig. S5A-C). The methylation level of Tkt-DMR in HF1 livers was lower than that in CF1 (Fig. S5D), although it was higher in HFD oocytes compared with CD. This contradiction may be associated with the uterine environment of obesity and the reprogramming during early embryo development. We further examined the expression of several genes with DMRs at promoters, including hyper-methylated genes, such as Bhlha15, Mgat1, Dgka, Pdpk1 and Taok3, and hypo-methylated genes, such as Igf1, Map3k8, Pld1, Tkt, Pik3cd and Sphk2. The expression levels of Bhlha15, Mgat1 and Pdpk1 were significantly lower and the expression levels of Map3k8, Tkt, Pik3cd, and Sphk2 were significantly higher in HF1 livers compared to that in CF1 (Fig. S5E and F). The expression trends of Bhlha15, Mgat1, Tkt, and Pik3cd were consistent with the methylation status at promoters in F1 livers. The expression of Dgka, Taok3, Igf1, and Pld1 was not affected in F1 livers (Fig. S4E and F). Mgat1 is associated with lipid metabolism and obesity (Jacobsson et al., 2012; Johansson et al., 2010), Tkt regulates glucose metabolism (Bartakova et al., 2016; Kang et al., 2018), and Pik3cd is involved in lipid metabolism and diabetes (Hood, Berger, Freedman, & Law, 2019; Wojcik et al., 2014). These results suggest that the altered methylation in HFD oocytes is partly transmitted to F1 livers via oocytes, and that the abnormal methylation may be a reason for the disturbed metabolism of HF1.

If the altered methylation in HFD oocytes was inherited by HF1 oocytes, it would be transmitted to F2 generation. Therefore, we examined the methylation of DMRs in F1 oocytes, and results revealed that Bhlha15- DMR, Mgat1-DMR, and Taok3-DMR were significantly hyper-methylated in HF1 oocytes compared with CF1 oocytes (Fig. 2E). The methylation level of Tkt-DMR was significantly lower (Fig. S6), and the methylation level of Pik3cd-DMR was slightly lower in the HF1 oocytes than that in the CF1 oocytes (Fig. 2E). These results indicate that at least a part of the altered methylation in HFD oocytes is transmitted to F1 oocytes via females.

To confirm the transgenerational inheritance of the altered DNA methylation in HFD oocytes, we examined the methylation levels of DMRs in F2 livers. The methylation of Bhlha15-DMR and Mgat1-DMR was higher and the methylation level of Pik3cd-DMR was lower in HF2 livers than that in CF2 livers (Fig. S7A-C). The methylation of Tkt-DMR in HF2 livers was similar to that in CF2 group (Fig. S7D). In addition, the expression levels of Bhlha15, Pdpk1, and Mgat1 were significantly lower (Fig. S7E), and the expression levels of Pld1, Pik3cd, and Sphk2 were significantly higher in HF2 livers than that in CF2 (Fig. S7F). The expression of the other gene were similar in the livers between HF2 and CF2 (Fig. S7E and F). These results suggest that at least a part of the altered methylation in HF2 livers may be inherited from HFD oocytes, and that this alteration may be associated with the disrupted metabolism in F2 offspring.

We then examined the methylation of DMRs in F2 oocytes using BS, and found that the hyper-methylation of Bhlha15-DMR, Mgat1-DMR, and Taok3-DMR, and the hypo-methylation of Pik3cd-DMR and Tkt-DMR were maintained in HF2 oocytes (Fig. 2F).

To better understand the inheritance of the altered methylation, we analyzed the methylation level of DMRs among generations compared with that in HFD oocytes (Fig. 2G). For hyper-DMRs, the methylation of Bhlha15-DMR was maintained from HFD oocytes to HF2 oocytes. Compared to those in HFD oocytes, the methylation levels of Mgat1-DMR and Taok3-DMR were partly transmitted to HF2 oocytes. For hypo-DMRs, the methylation levels of Tkt-DMR and Pik3cd-DMR in HFD oocytes were inherited by HF2 (Fig. 2G). . We did not examine the methylation level of Pld1-DMR in F1 and F2 because it was similar in oocytes between CD and HFD (Fig. 2G). These results suggest that only a part of the altered methylation in HFD oocytes can be transmitted to F2 oocytes, and that disrupted methylation may be a reason for the inheritance of the metabolic disorders in F1 and F2.

Obesity alters the metabolomics of serum

Obesity alters the metabolism of glucose, fatty acids, and amino acids, which are essential for oogenesis. Thus, we suppose the altered metabolism may play a key role in the disturbed global methylation in HFD oocytes. We examined the metabolomics of serum using non-targeted approaches (BGI, Wuhan, China). We used LC-MS/MS to identify the variation in metabolites, including amino acids, carbohydrates, lipid, and phenols. The principal component analysis (PCA) showed that the PC1 and PC2, respectively, explained 44.6% and 9.96% of the total metabolite variation, respectively (Fig. 3A). The distribution of extracts between groups was distinguishable (Fig. 3A). We identified 538 differential metabolites based on the PLS-DA (robust orthogonal partial least squares-discriminant analysis) and T-test analysis with the following criteria: VIP (variable importance in projection) ≥ 1, fold change ≥ 1.2 or ≤ 0.83, and p- value<0.05, including 288 upregulated and 250 downregulated metabolites (Fig. 3B). The enrichment of differential metabolites was analyzed using KEGG, which revealed that differential metabolites were significantly enriched in tryptophan and vitamin B6 metabolism (Fig. 3C). The top 20 differential metabolites were presented in Fig. 3D. These results suggest that obesity disturbs the metabolomics of serum.

Maternal obesity alters metabolome of serum.

(A) Principal component analysis in CD and HFD mice.

(B) Differential metabolites in HFD serum compared with those in CD group. Red circles, upregulated metabolites; blue circles, downregulated metabolites.

(C) The enrichment of differential metabolites was analyzed using KEGG, and the top 10 enrichment terms are presented.

(D) Heat map of the top 20 differential metabolites in HFD serum.

(E-G) Comparison of the concentrations of pyridoxine, methionine, and tyrosine among the groups. * p<0.05; ** p<0.01; *** p<0.001.

(H-J) Concentrations of SAM, SAH and HCY in the livers were examined by ELISA. Ns, there was no statistical significance between groups.

(K) The concentration of SAM in oocytes was analyzed using ELISA. ** p<0.01.

(L) Relative concentration of melatonin in the serum. *** p<0.001.

(M) Genomic DNA methylation in oocytes was examined using immunofluorescence. CD, control group; HFD, obesity group; HFD+melatonin, obese mice were treated with exogenous melatonin for 14 days.

(N) Relative fluorescence intensity of 5mC was examined using Image J (CD, n=109; HFD, n=104; HFD+melatonin, n=96). * p<0.05; *** p<0.001.

Melatonin may play a key role in the genomic hyper-methylation of HFD oocytes

To investigate the associations between the metabolites and methylation in oocytes, we identified several differential metabolites in HFD mice, including pyridoxine (vitamin B6), L-methionine, melatonin, and L- tyrosine, which may be associated with the hyper-methylation of oocytes. The concentration of pyridoxine (vitamin B6), L-methionine, and L- tyrosine was significantly increased in HFD serum compared with CD group (Fig. 3D-G). In the methionine cycle, methionine adenosyltransferase catalyzes methionine and ATP into S-adenosyl methionine (SAM) (Shou, Pan, & Chin, 1969) which acts as a universal methyl donor (Lyon, Strippoli, Fang, & Cimmino, 2020). DNA methyltransferases (DNMTs) transfer the methyl group to cytosine and convert SAM to S-adenosyl homocysteine (SAH) which is further degraded to homocysteine (HCY) by SAH hydrolase (D. H. Chen, Wu, Hung, Hsieh, & Li, 2004). During these processes, vitamin B6 serves as a co-factor (Vaccaro & Naser, 2021). Excessive methionine and vitamin B6 intake has been reported to induce hyper-methylation (R. A. Waterland, 2006). These findings indicates that the hyper-methylation in HFD oocytes may be associated with the increased concentrations of methionine and pyridoxine. To confirm this hypothesis, we examined the concentrations of SAM, SAH, and HCY in livers and oocytes, which are crucial intermediate metabolites in the methionine cycle. Results showed that the concentrations of SAM, SAH, and HCY in HFD livers were similar to those in CD livers (Fig. 3H-J), and the concentration of SAM in HFD oocytes was lower than that in the CD group (Fig. 3K). These results suggest that the higher concentrations of methionine and pyridoxine in the serum of HFD mice may not be the main reason for the genomic hyper- methylation in oocytes.

As presented in Fig. 3D, it is curious that the abundance of genistein, daidzein and dibutylpthalate were also altered in HFD serum compared with those in CD, which might have contributed to the altered DNA methylation in the oocytes. These metabolites might be from the diets or materials used to collect samples. To confirm these results, we examined the concentrations of genistein and dibutylpthalate using ELISA. The results revealed that the concentrations of these metabolites were similar between HFD and CD serum (Fig. S8), which suggest that these metabolites may have no effect on the altered methylation in oocytes.

As presented in Fig. 3D and L, the concentration of melatonin in the HFD serum was significantly lower than that in the CD group. Low concentration of melatonin has also been reported in obese rats (Virto et al., 2018) and humans (Overberg et al., 2022). Exogenous melatonin supplementation reduces body weight and improves lipid and glucose metabolism in both animals and humans (Guan, Wang, Cao, Dong, & Chen, 2021). Melatonin can also inhibit cancer by regulating DNA methylation status (Davoodvandi et al., 2022), improve DNA methylation reprogramming in the development of porcine cloned embryos (Qu et al., 2020), and affect DNA re-methylation in oocytes (Lan et al., 2018; Saeedabadi, Abazari-Kia, Rajabi, Parivar, & Salehi, 2018; Xiao et al., 2019). This indicates that the reduced melatonin may contribute to the hyper-methylation of HFD oocytes. Here, we found that if HFD mice were treated with exogenous melatonin for 14 days, the genomic hyper- methylation in HFD oocytes was significantly reduced (n≥96 Fig. 3M and N). These results suggest that the reduced melatonin concentration may be involved in regulating the hyper-methylation of HFD oocytes.

Melatonin regulates genomic methylation of oocytes by increasing the expression of DNMTs via the cAMP/PKA/CREB pathway

Melatonin receptors (MT1 and 2), which are coupled with inhibitor G- protein (Gi), have been identified in oocytes and granulosa cells (Jin et al., 2022; Wang et al., 2021). Activated Gi inhibits the activation of adenylyl cyclases (ADCYs), resulting in a decrease of cAMP (cyclic adenosine monophosphate), which regulates the activation of protein kinase A (PKA) and CREB (cAMP response element-binding protein) (Wongprayoon & Govitrapong, 2021). Elevated cAMP level increases the expression of DNMTs resulting in hyper-methylation in HL-1 cardiomyocytes (Fang et al., 2015). Hedrich et al. reported that CREMα induced hyper-methylation of CD8 cluster via increasing the expression of DNMT3a (Hedrich et al., 2014). We thus supposed that melatonin may regulate genomic methylation in oocytes via increasing the expression of DNMTs through the cAMP/PKA/CREB pathway (Fig. 4A). To confirm this hypothesis, female C57BL/6 mice fed with normal diet were treated with luzindole, an inhibitor of melatonin receptor, and the global methylation of 5mC and 5hmC was significantly increased in oocytes (n≥49 Fig. 4B-D). Luzindole did not affect the concentration of melatonin in the serum (Fig. S9A). However, excessive melatonin treatment significantly increased the concentration of melatonin in the serum (Fig. S9A) and decreased the methylation levels of 5mC and 5hmC in oocytes (Fig. 4B-D). These findings indicate that melatonin may regulate the re-methylation process in oocytes.

Melatonin regulates DNA methylation in oocytes.

(A) Schedule of the possible pathway by which melatonin regulates DNA methylation in oocytes. According to previous studies, we predicted that melatonin might regulate DNA methylation in oocytes via the cAMP/PKA/CREB pathway.

(B) Effects of melatonin and its inhibitor luzindole on oocyte methylation were examined using immunofluorescence.

(C, D) The relative fluorescence intensities of 5mC and 5hmC were analyzed using Image J (5mC: Control, n=81; Luzindole, n=83; Melatonin, n=86. 5hmC: Control, n=64; Luzindole, n=58; Melatonin, n=49). * p<0.05; *** p<0.001.

(E) The effects of melatonin and its inhibitor luzindole on the expression of adenylate cyclase (ADCY) in oocytes were examined by qPCR. * p<0.05; ** p<0.01.

(F) The concentration of cAMP in oocytes was examined by ELISA. * p<0.05; ** p<0.01.

To confirm whether melatonin regulates methylation in oocytes via the cAMP pathway, we examined the expression of ADCYs in oocytes using RT-PCR, and found that ADCY5, 6, and 9 were expressed in oocytes (Fig. S9B). The melatonin antagonist luzindole significantly increased the expression of ADCY6 and ADCY9 in oocytes, and melatonin reduced the expression of ADCY6 (Fig. 4E). However, the expression of ADCY5 was lower in the luzindole and melatonin groups compared with the control (Fig. 4E). In addition, the melatonin antagonist luzindole increased, while melatonin decreased the concentration of cAMP in oocytes, respectively (Fig. 4F). These results suggest that melatonin may regulate the synthesis of cAMP via ADCY6 in oocytes. To further confirm the role of cAMP in regulating the methylation of oocytes, we treated mice with SQ22536, an inhibitor of ADCYs, and found that this treatment significantly reduced the global methylation of 5mC and the concentration of cAMP in oocytes (n≥51 Fig. 5A-C). Whereas, the ADCYs activator forskolin significantly increased the cAMP concentration and global methylation of 5mC in oocytes (Fig. 5A-C). 8-Bromo-cAMP, a cAMP analogue, also increased the global methylation in oocytes (n≥41 Fig. 5D, E). cAMP functions by activating the downstream protein PKA. When we treated mice with H89 2HCL, a PKA antagonist, the global methylation of 5mC was significantly reduced in oocytes (n≥24 Fig. 5F, G). These results suggest that melatonin may mediate the methylation of oocytes via the cAMP/PKA pathway. cAMP activates PKA which further phosphorylates CREB to regulate gene expression. In oocytes, DNA re-methylation is regulated by DNMTs including DNMT3a, DNMT3l, and DNMT1. Therefore, we next investigated whether melatonin regulates DNMTs expression via the cAMP/PKA/CREB pathway in oocytes. We examined the expression of CREB1, CREM (cAMP responsive element modulator), CREB3l2 (cAMP responsive element binding protein 3 like 2), and ATF1 (activating transcription factor 1) in oocytes. Results showed that ADCYs activator forskolin treatment significantly increased the mRNA expression of CREB1 and CREM, and that the expression of CREB3l2 and ATF1 was slightly increased (Fig. 6A). In addition, the ADCYs inhibitor SQ22536 significantly reduced the expression of CREB1 and CREB3l2, although the expression of CREM and ATF1 was slightly decreased in oocytes (Fig. 6A). Furthermore, treatment with the ADCYs inhibitor SQ22536 significantly reduced the concentration of pCREB1, but it was increased by forskolin treatment in oocytes (n≥36 Fig. 6B, C). The pCREB1 level was also increased by 8-Bromo-cAMP and decreased by the PKA antagonist H89 2HCL in oocytes (n≥28 Fig. 6D-G). These results suggest that the expression and phosphorylation of CREB1 can be regulated by the cAMP/PKA pathway. Yang et al. demonstrated that CREB regulated DNMT3a expression in neurons of the dorsal root ganglion by binding to the promoter region (Yang et al., 2021). In the present study, the binding of pCREB1 with relative regions of DNMTs was examined using CUT & Tag assay. Each sample contained 500 GV oocytes, and two replicates were involved. The sequencing result revealed that five fragments including 10 pCREB1 binding motifs (predicted using online tool JASPAR, Table S4) were associated with DNMTs, including 3 fragments at intron 1 and distill intergenic regions of DNMT3a, 1 fragment at the promoter region of DNMT1, and 1 fragment at intron 13 of DNMT3l (Table S5). These results suggest that pCREB1 may have contributions to regulate the expression of DNMTs.

Role of cAMP in DNA methylation in oocytes.

(A) Female mice were respectively treated with the ADCY inhibitor SQ22536 or activator forskolin. Oocyte methylation was examined using immunofluorescence.

(B) The relative intensity of fluorescence in oocytes was analyzed using Image J (Control, n=107; SQ22536, n=51; Forskolin, n=57). ** p<0.01; *** p<0.001.

(C) cAMP concentration in oocytes was examined using ELISA. * p<0.05; ** p<0.01.

(D) Female mice were treated with the cAMP analogue 8-Bromo-cAMP, and oocyte methylation was examined using immunofluorescence. The relative fluorescence intensity of 5mC was analyzed using Image J (E) (Control, n=41; 8-Bromo-cAMP, n=42). ** p<0.01.

(F, G) Female mice were treated with the PKA (protein kinase A) antagonist H 89 2HCL, and then oocyte methylation was examined using immunofluorescence. The relative fluorescence intensity of 5mC was analyzed using Image J (G) (Control, n=24; H 89 2HCl, n=25). * p<0.05.

Effects of cAMP on CREB1.

(A) The mRNA expression of cAMP-response element binding (CREB) proteins in oocytes was examined by qPCR. * p<0.05.

(B and C) Phosphorylated CREB1 (pCREB1) in oocytes was examined using immunofluorescence, and the relative fluorescence intensity of pCREB1 was examined by Image J (C) (Control, n=36; SQ22536, n=48; Forskolin, n=41). * p<0.05; ** p<0.01; *** p<0.001.

(D and E) After treatment with the cAMP analogue 8-Bromo-cAMP, pCREB1 in oocytes was examined using immunofluorescence. The relative fluorescence intensity was analyzed using Image J (E) (Control, n=28; 8-Bromo-cAMP, n=28). *** p<0.001.

Next, we investigated the expression of DNMT1, DNMT3a, and DNMT3l in oocytes. The melatonin inhibitor luzindole (slightly) and the ADCYs activator forskolin (significantly) increased the expression of DNMT1 and DNMT3a, respectively (Fig. 7A, B). Melatonin and the ADCYs inhibitor SQ22536 significantly and slightly reduced the expression of DNMT1 and DNMT3a in oocytes, respectively (Fig. 7A, B). The protein level of DNMT3a in GV oocytes was also significantly increased by 8-Bromo- cAMP and decreased by the PKA antagonist H 89 2HCL, respectively (n≥48 Fig. 7C-F). Although DNMT1 is well known as a maintenance methyltransferase, it also contributes to de novo methylation in oocytes (Li et al., 2018). Therefore, we examined the localization of DNMT1 in oocytes, and found that 8-Bromo-cAMP treatment significantly increased the localization of DNMT1 in the nucleus of oocytes, but it was reduced by the PKA antagonist H89 2HCL (n≥22 Fig. 7G-J). When the activation of DNMTs was inhibited by 5-azacytidine, the methylation level in GV oocytes was significantly decreased (n≥43 Fig. S10). These results suggest that melatonin may influence the genomic methylation of oocytes via regulating the expression of DNMT1 and DNMT3a mediated via the cAMP/PKA/CREB pathway.

Role of the melatonin/cAMP/PKA pathway in the expression of DNMTs.

(A) The expression levels of DNMT1, DNMT3a and DNMT3l in oocytes were examined using qPCR after the treatment with SQ22536 and forskolin. * p<0.05.

(B) The relative expressions of DNMT1, DNMT3a and DNMT3l in oocytes were examined using qPCR after the treatment with luzindole and melatonin. * p<0.05.

(C) After 8-Bromo-cAMP treatment, the relative expression of DNMT3a in oocytes was examined using immunofluorescence and calculated by Image J (D) (Control, n=54; 8-Bromo-cAMP, n=70). ** p<0.01.

(E and F) Treatment with the PKA antagonist H 89 2HCL treatment significantly reduced the level of DNMT3a in oocytes examined using immunofluorescence (Control, n=62; H 89 2HCl, n=48). ** p<0.01.

(G and H) DNMT1 localization in the oocyte nucleus was examined using immunofluorescence after 8-Bromo-cAMP treatment (Control, n=30; 8- Bromo-cAMP, n=31). *** p<0.001.

(I and J) The localization of DNMT1 in oocyte nucleus was reduced by the treatment with the PKA antagonist H 89 2HCL (Control, n=22; H 89 2HCl, n=28). ** p<0.01.

Increased DNMTs mediate hyper-methylation of HFD oocytes via cAMP/PKA/CREB pathway

To explore how melatonin regulates the hyper-methylation of HFD oocytes, we examined the expression of DNMTs. Results showed that maternal obesity in mice significantly increased the expression of DNMT1, DNMT3a, and DNMT3l in oocytes (Fig. 8A). Thus, we aimed to confirm whether the cAMP/PKA/CREB pathway mediated the increased expression of DNMTs in HFD oocytes. We examined the concentration of cAMP in oocytes, and found that maternal obesity significantly increased the concentration of cAMP in oocytes compared with that in CD group (Fig. 8B). The mRNA expression of CREB1, but not CREM, in the HFD oocytes was significantly increased compared with that in CD group (Fig. 8C), and the pCREB1 level in HFD oocytes was also significantly increased (n≥49 Fig. 8D, E). However, the increased level of pCREB1 was reduced by exogenous melatonin treatment (Fig. 8D, E). When obese females were treated with the PKA antagonist H89 2HCL, both 5mC and pCREB1 levels were significantly reduced in oocytes (n≥17 Fig. 8F-I). These results suggest that reduced melatonin in obese mice may increase the expression of DNMTs via the cAMP/PKA/CREB pathway.

Melatonin regulates DNMTs expression via cAMP/PKA/CREB pathway in HFD oocytes.

(A) The relative expression of DNMT1, DNMT3a, and DNMT3l in HFD oocytes was examined using qPCR. * p<0.05; ** p<0.01.

(B) The concentration of cAMP in HFD oocytes was examined using ELISA. ** p<0.01.

(C) The relative expressions CREB1 and CREM in HFD oocytes were tested using qPCR. * p<0.05; ** p<0.01.

(D and E) The level of pCREB1 in oocytes was examined using immunofluorescence, and the relative fluorescence intensity was calculated by Image J (E) (CD, n=69; HFD, n=49; HFD+melatonin, n=61). HFD, oocytes from obese mice; CD, oocytes from control mice; HFD + melatonin, oocytes from obese mice treated with exogenous melatonin. * p<0.05; *** p<0.001.

(F and G) Treatment with the PKA antagonist H89 2HCL reduced the methylation level of HFD oocytes (CD, n=48; HFD, n=31; HFD+H 89 2HCl, n=27). ** p<0.01; *** p<0.001.

(H and I) The level of pCREB1 in HFD oocytes was also decreased by the treatment with the PKA antagonist H89 2HCL (CD, n=17; HFD, n=17; HFD+H 89 2HCl, n=22). * p<0.05; ** p<0.01; ns, no statistical significance between groups. (J and K) Treatment with the PKA antagonist H89 2HCL reduced the localization of DNMT1 in HFD oocytes (CD, n=24; HFD, n=29; HFD+H 89 2HCl, n=25). ** p<0.01; *** p<0.001; ns, no statistical significance between groups.

The increased expression of DNMT3a (Fig. 8A) may have contributed to the genomic hyper-methylation of HFD oocytes because the primary function of DNMT3a is de novo DNA methylation. The de novo methylation function of DNMT1 is usually prevented by Stella (also known as Dppa3 or PGC7) in oocytes (Bostick et al., 2007). Nevertheless, maternal obesity significantly reduces the expression of Stella in oocytes (Han et al., 2018), which indicates that DNMT1 may also contribute to the hyper-methylation of the HFD oocytes. We found that the localization of DNMT1 in the nuclei of HFD oocyte was significantly increased, which could be reduced by the PKA antagonist H89 2HCL treatment (n≥24 Fig. 8J, K). These results suggest that the decreased melatonin induced by maternal obesity increases the expression of DNMTs via the cAMP/PKA/CREB pathway, which results in the hyper-methylation in HFD oocytes.

Discussion

Maternal obesity has negative effects on oocyte quality and offspring health, but the underlying mechanisms are not well understod. In the present study, we found that maternal obesity induced hyper-methylation in oocytes, and that the abnormal methylation, at least in part, is transmitted to F2 oocytes in females, which may be associated with the occurrence and inheritance of metabolic disorders. Maternal obesity induced metabolic changes in the serum. The decreased melatonin in serum may be involved in regulating the hyper-methylation of HFD oocytes by increasing the expression of DNMTs, which is mediated by the cAMP/PKA/CREB pathway.

Transgenerational epigenetic inheritance is common in plants, but related investigations in mammals are hindered by epigenetic reprogramming events during gametogenesis and early embryo development (Schmitz & Ecker, 2012; Xavier, Roman, Aitken, & Nixon, 2019). In mammals, isogenic agouti viable yellow (Avy) and axin-fused (AxinFu) mice, whose phenotypes are regulated by the DNA methylation level of IAP (intra-cisternal A particle long terminal repeat) respectively located at the upstream and in intron 6, are solid evidence confirming the transgenerational inheritance of epigenetic modifications (Morgan, Sutherland, Martin, & Whitelaw, 1999; Rakyan, Blewitt, Druker, Preis, & Whitelaw, 2002). Nevertheless, epigenetic modifications can be affected by environmental factors such as metabolic diseases and diet. The intergenerational inheritance of phenotypes and epigenetic changes induced by maternal environmental factors have been confirmed by previous studies (B. Chen et al., 2022; Ge et al., 2014), but there are still many debates about the transgenerational inheritance of epigenetic changes induced by environment. Rats from stressed mothers are more likely to be stressed, and that can be transmitted across generations, but this transgenerational inheritance is not mediated by gametes (Weaver et al., 2004). Females cannot mediate the transgenerational inheritance of hyper- methylation induced by diet in Avy mice (Robert A. Waterland, Travisano, & Tahiliani, 2007), but another study reported that if the Avy allele is from the father, the hyper-methylation induced by diet during pregnancy can be retained in germ cells (Cropley, Suter, Beckman, & Martin, 2006). Anway et al. reported that the abnormal spermatogenesis induced by the exposure of vinclozolin during pregnancy can be transmitted across four generations via sperm (Anway, Cupp, Uzumcu, & Skinner, 2005). We previously reported that maternal obesity disturbed DNA methylation status of imprinted genes in oocytes (Ge et al., 2014), but the transgenerational inheritance was not observed. We also demonstrated that disturbed methylation in oocytes induced by undernourishment in utero could be inherited, at least partly, by F2 oocytes via females (S.-B. Tang et al., 2023). Recently, Takahashi et al. edited DNA methylation of promoter-associated CGIs, and found that the edited DNA methylation, associated with disrupted metabolism, was stably inherited by multiple generations (Takahashi et al., 2023). In the present study, we found that maternal obesity induced genomic hyper-methylation in oocytes, and that a part of the abnormal methylation was transmitted to F2 via female gametes. Moreover, the transmission of metabolic disorders has also been observed across two generations. These results suggest that the transgenerational inheritance of abnormal methylation induced by maternal obesity, at least in part, can be mediated by oocytes, which may be a reason for the inheritance of metabolic disorders.

During the methylation process, the methyl group is donated by SAM which is generated from homocysteine, 5-methyltetrahydrofolate (5mTHF), and methionine. 5mTHF is an intermediate of one-carbon metabolism (Mentch & Locasale, 2016). One-carbon units such as folate and vitamin B12 are crucial for the establishment of methylation (Mentch & Locasale, 2016). Disturbed glucose and lipid metabolism also has a negative influence on DNA methylation (Keller et al., 2014). These indicate that abnormal metabolism induced by maternal obesity (King, 2006) may play a key role in the genomic hyper-methylation in oocytes. In the present study, we found that the metabolomics of serum from HFD mice was distinguishable from that of serum from CD mice. Although the concentrations of vitamin B6 and methionine were higher in HFD serum than that in CD, this may not be an important reason for the genomic hyper- methylation in oocytes because the concentrations of SAM, SAH, and HCY in the livers and oocytes was similar between CD and HFD mice. In humans, obesity reduces the melatonin level in circulation (Overberg et al., 2022; Virto et al., 2018). In the present study, we also found that maternal obesity induced by high-fat diet reduced the concentration of melatonin in the serum. Melatonin not only can decrease body weight, but also regulates DNA methylation of somatic cells and germ cells (Davoodvandi et al., 2022; Lan et al., 2018). Nevertheless, the molecular mechanism by which melatonin regulates DNA methylation in oocytes is still unclear. Melatonin has two receptors, MT1 and MT2, both of which have been identified in oocytes (Wang et al., 2021). Melatonin receptors interact with the inhibitor G-protein and can regulate gene expression via the cAMP/PKA/CREB pathway (Wongprayoon & Govitrapong, 2021). CREM and cAMP mediate DNA methylation in somatic cells by regulating the expression of DNMT3a (Fang et al., 2015; Hedrich et al., 2014). In neurons, CREB interacts with the promoter of DNMT3a to regulate its expression and DNA methylation (Yang et al., 2021). In the present study, we found that melatonin, which is mediated by the cAMP/PKA/CREB pathway, regulates methylation in oocytes by increasing the expression of DNMT1 and DNMT3a. Similar results were also observed in HFD oocytes. Hyper-methylation of HFD oocytes can be reduced by exogenous melatonin and PKA inhibitors. These suggest that decreased melatonin levels are involved in regulating the genomic hyper-methylation of HFD oocytes by increasing the expression of DNMTs, which is mediated by the cAMP/PKA/CREB pathway.

During follicular development, re-methylation in oocytes is catalyzed mainly by DNMT3a and DNMT3l (Kibe et al., 2021). DNMT1 is usually responsible for maintaining DNA methylation, but DNMT1 also contributes to CG methylation in oocytes (Shirane et al., 2013). During normal oocyte development, DNMT1 is mainly prevented from the nuclei by Stella. When the Stella level is knocked out, Uhrf1 (Ubiquitin-like containing PHD Ring Finger 1) moves to the nucleus from the cytoplasm and recruits DNMT1 to chromatin, resulting in hyper-methylation (Li et al., 2018). Maternal obesity significantly deceases the expression of Stella in oocytes (Han et al., 2018). These suggest that maternal obesity may induce hyper-methylation in oocytes. In the present study, we found that maternal obesity increased the genomic DNA methylation of GV and MII oocytes, the expression of DNMT1 and its localization at chromatin in GV oocytes. These suggest that reduced Stella in HFD oocytes may recruit more DNMT1 into chromatin resulting in hyper-methylation. Nevertheless, Han LS et al. reported that maternal obesity has no significant influence on whole genome methylation of GV oocytes and Hou YJ et al. reported that maternal obesity reduces whole genome methylation of NSN (no Hoechst- positive rim surrounding the nucleolus) GV oocytes. This contradiction may be associated with the methods used to examine the methylation level and the sample size.

In summary, we found that maternal obesity induced genomic hyper- methylation in oocytes and that at least a part of the altered methylation can be transmitted to F2 oocytes, which may be a reason for the inheritance of metabolic disorders. Furthermore, reduced melatonin in HFD mice was involved in regulating the genomic hyper-methylation of oocytes by increasing the expression of DNMTs, and this process was mediated by the cAMP/PKA/CREB pathway. However, there are some limitations for the present study: there is not enough evidence to confirm the role of altered DNA methylation in metabolic disorders in the offspring of obese mothers The molecular mechanisms by which DNA methylation escapes reprogramming in oogenesis have not been elucidated. There are may be other mechanisms involved in regulating genomic hyper-methylation in HFD oocytes. Therefore, more studies are needed in the short future.

Materials and methods Mice

C57BL/6 mice were purchased from Jinan Pengyue Company (Jinan, China). Mice were housed in the controlled room with 12 h light and 12 h dark cycle, and at 23-25℃. The Animal Ethics Committee of Qingdao Agricultural University supported all procedures (QAU201900326).

For the obesity model, female C57BL/6 mice at the age of 4 weeks were randomly divided into two groups fed with high-fat diet (HFD, Research Diets, D12492, USA) and normal diet (CD) for 12 weeks, respectively. We examined the body weight every week. The formulation of the diets is presented in Table S6.

Offspring were produced according to the schedule in Fig.1J. For F1 offspring, female HFD and CD mice were mated with normal adult male C57BL/6 mice, respectively, and the offspring were marked as HF1 and CF1. To avoid the effects of males on methylation, the same males were used to produce F2: female HF1 and CF1 mating with normal males, marked as HF2 and CF2.

Immunofluorescence

Briefly, oocytes were fixed with 4% PFA (paraformaldehyde), permeabilized with 0.5% TritonX-100, and blocked with 1% BSA (bovine serum albumin). After that, the oocytes were incubated with primary antibodies overnight at 4℃. The secondary antibodies were stained for 1 h at room temperature. For 5mC and 5hmC staining, the oocytes were treated with 4 N HCl for 10 min after permeabilization, and then transferred into 100mM Tris-HCl for 10 min. After washed 3 times using PBS/BSA with 0.05% Tween 20, the oocytes were blocked with PBS/BSA for 1 h, and incubated with primary antibodies overnight at 4℃. Secondary antibodies were stained for 1 h at room temperature. Then, oocytes were transferred into DAPI with mounting and sealed. Fluorescence signal was examined using a laser scanning confocal microscope (Leica SP5, Germany). The relative fluorescence intensity was examined using Image J.

Antibodies

Primary antibodies used in the present study included anti-5mC antibody (Abcam, ab73938), anti-5hmC antibody (Abcam, ab214728), anti-pCREB antibody (Cell Signaling Technology, 9198S), anti-DNMT3a antibody (Active motif, 61478), and anti-DNMT1 antibody (Active motif, 39204).

Whole-genome bisulfite sequencing and analysis

Metaphase II (MII) oocytes were collected from the oviduct. For each sample, 100 oocytes from at least 10 mice were pooled together and transferred to lysis buffer. Genomic DNA was fragmented, and the end was repaired. Then, fragmentations were ligated with adapters. Bisulfite treatment was performed using EZ DNA Methylation-Direct (Zymo Research). Lambda DNA was used as a control. After that, the sequencing library was established and sequenced using Illumina HiSeq/NovaSeq (Novogene, China). The raw data quality was evaluated using FastQC, and low quality data and adapters were trimmed using fastp. Clean data were compared to the reference genome mm10. Methylated C site calling was performed using Bismark. Differentially methylated regions (DMRs) were identified using DSS-single. The enrichment of genes in the KEGG pathway was carried out using the online tool KOBAS.

qPCR

Total RNA was extracted from oocytes or tissues using RNAprep Pure Micro Kit (Tiangen, DP420) or RNA Easy Fast Kit (Tiangen, DP451). cDNA was synthesized using Hifair III 1st Strand cDNA Synthesis Kit (Yeasen, China). cDNA was used as templates to examined the relative expression of genes. Housekeeping genes Ppia and Gapdh were used as references. Relative expression was calculated as 2-⊿⊿Ct.

Bisulfite sequencing (BS)

Each sample included 5 oocytes and at least 20 samples were used for each DMR. The samples were treated as described in a previous study (Ge et al., 2014). Briefly, samples were treated with lysis buffer and 0.3 M NaOH, respectively. After that, samples were embedded in 2% low melting point agarose (Sigma), which was treated with fresh bisulfite solution (2.5 M sodium metabisulfite, Merck; 125 mM hydroquinone, Sigma; pH 5) for 4 h. Treated DNA was used as templates to amplify the target fragment using nest-PCR. PCR products were cloned to T-vector and sequenced. Methylation status was analyzed using BiqAnalyzer, which can remove the low conversion rate (<97%) and possible repeat sequences. At least 10 available clones was used for each DMR.

Glucose and insulin tolerance

Glucose and insulin tolerance (GTT and ITT) were examined as previously reported (S. B. Tang et al., 2023). Briefly, mice were treated with glucose at 2 g/kg body weight or insulin (Actrapid®, Novo Nordisk) at 10 IU/kg body weight after 16 h or 4 h fasting, respectively. After that, blood glucose was measured by tail blood at 0, 30, 60, 90, and 120 min, respectively.

ELISA

Concentrations of cAMP, SAM, SAH, and HCY were examined using ELISA kits (Jinma Biotechnology Co. Ltd, Shanghai, China) according to the handbook. A standard curve was produced using four-parameters logistics.

Non-target metabolomics in serum

Metabolites in the serum were examined using LC-MS/MS (BGI, Wuhan, China). Raw date were treated with Compound Discoverer 3.1 (Thermo Fisher Scientific, USA). After that, preprocessing of the exported data was performed using metaX. Metabolites were identified according to the databases of BMDB (BGI), mzCloud and ChemSpider (HMDB, KEGG, LipidMaps). Identified metabolites were annotated according to KEGG and HMDB. Differential metabolites were scanned using PCA and PLS- DA combined with fold changes and Student’s t test.

Chemicals

Inhibitors used in the present study included luzindole (Sigma), SQ22536 (Selleck), forskolin (Selleck), H89 2HCL (Selleck), and azacitidine (Selleck). 8-Bromo-cAMP was purchased from Selleck. Melatonin (Sigma) was injected by tail vein for 14 days. The other chemicals were administrated via intraperitoneal injection for 14 days. The control group was injected with the appropriate solutions.

CUT & Tag and sequencing

For each sample, 500 oocytes were pooled together for this experiment, and two replicates were performed. Library was established using Hyperactive® Universal CUT&Tag Assay Kit (Vazyme, China) according to the manufacturer’s instructions. Briefly, oocytes were washed with washing buffer and transferred into tubes with ConA Beads Pro for 10 min at room temperature. The tubes were placed on a magnetic rack for 2 min, and then discarded the supernatant. After that, 50 μl of precooled antibody buffer and primary antibody were added and incubated overnight at 4 ℃.

Then, tubes were put on a magnetic rack for 2 min, discarded supernatant, and added 50 μl Dig-wash buffer with secondary antibody (1:100) into tubes incubated at room temperature for 60 min. Samples were washed three times, and then incubated with pA/G-Tnp pro for 1 h at room temperature. After washing, 50 μl of TTBL was added to the samples and incubated at 37℃ for 1 h. Then, 2 μl 10% SDS and DNA Spike-in were added to the samples, which were incubated at 55 ℃ for 10 min. Supernatant was transferred into a new tube after being putted on a magnetic rack for 2 min. 25 μl DNA Extract Beads Pro was added into supernatant, incubated 20 min at room temperature, and then washed two times with B&W buffer. After that, the DNA Extract Beads Pro was re- suspended in 15 μl ddH2O, and amplified at 60℃ for 20 cycles. Library quality was examined using Qubit, AATI, and QPCR, and then sequenced using NovaSeq 6000 (Novogene, China). Adaptor and low-quality reads were removed from the raw data, and clean data were used for further analysis. Reads were mapped to the mouse reference genome mm39 using Bowtie2. Peak calling was performed using MACS2 at q-value <0.05. The peak was subsequently annotated into related gene regions.

Statistical analysis

Average data are presented as mean ± SE (standard error), and the statistical difference was calculated using two-tail independent-samples t test. Methylation level is presented as a percentage, and the statistical difference was calculated using the chi-square test. If p value < 0.05, the statistical difference was considered to be significant.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (31872312), National R&D Program of China (2022YFC2703500), the Breeding Plan of Shandong Provincial Qingchuang Research Team (Innovation Team of Farm Animal Cloning 012–1622001), and the Doctor Foundation of Qingdao Agricultural University (6631116008).

Conflicts of interest statement

There are no conflicts of interest to declare.

Data availability

The raw data of sequencing are submitted to the database BGI Sub with No. CRA011654.