key: cord-0078507-wmqd2wfn
authors: Lu, Sha; Wang, Jiahao; Kakongoma, Nisile; Hua, Wen; Xu, Jiahui; Wang, Yunfei; He, Shutao; Gu, Hongcang; Shi, Jiantao; Hu, Wensheng
title: DNA methylation and expression profiles of placenta and umbilical cord blood reveal the characteristics of gestational diabetes mellitus patients and offspring
date: 2022-05-23
journal: Clin Epigenetics
DOI: 10.1186/s13148-022-01289-5
sha: 0535a53188573f3c3f14567c65a3ce4c2f7b4687
doc_id: 78507
cord_uid: wmqd2wfn

BACKGROUND: Gestational diabetes mellitus (GDM) is a common pregnancy-specific disease and is growing at an alarming rate worldwide, which can negatively affect the health of pregnant women and fetuses. However, most studies are limited to one tissue, placenta or umbilical cord blood, usually with one omics assay. It is thus difficult to systematically reveal the molecular mechanism of GDM and the key influencing factors on pregnant women and offspring. RESULTS: We recruited a group of 21 pregnant women with GDM and 20 controls without GDM. For each pregnant woman, reduced representation bisulfite sequencing and RNA-seq were performed using the placenta and paired neonatal umbilical cord blood specimens. Differentially methylated regions (DMRs) and differentially expressed genes (DEGs) were identified with body mass index as a covariate. Through the comparison of GDM and control samples, 2779 and 141 DMRs, 1442 and 488 DEGs were identified from placenta and umbilical cord blood, respectively. Functional enrichment analysis showed that the placenta methylation and expression profiles of GDM women mirrored the molecular characteristics of “type II diabetes” and “insulin resistance.” Methylation-altered genes in umbilical cord blood were associated with pathways “type II diabetes” and “cholesterol metabolism.” Remarkably, both DMRs and DEGs illustrated significant overlaps among placenta and umbilical cord blood samples. The overlapping DMRs were associated with “cholesterol metabolism.” The top-ranking pathways enriched in the shared DEGs include “growth hormone synthesis, secretion and action” and “type II diabetes mellitus.” CONCLUSIONS: Our research demonstrated the epigenetic and transcriptomic alternations of GDM women and offspring. Our findings emphasized the importance of epigenetic modifications in the communication between pregnant women with GDM and offspring, and provided a reference for the prevention, control, treatment, and intervention of perinatal deleterious events of GDM and neonatal complications. SUPPLEMENTARY INFORMATION: The online version contains supplementary material available at 10.1186/s13148-022-01289-5.

Gestational diabetes mellitus (GDM) is a common pregnancy-specific disease, characterized by glucose intolerance. GDM usually occurs in the middle or late pregnancy [1] [2] [3] . Increased insulin resistance and impaired insulin secretion are the main causes of GDM [4] . Due to different lifestyles, ethnicities, diagnostic methods, the median (interquartile range) prevalence (%) of GDM varies significantly from the lowest in Europe (6.1, 1. to the highest in the Middle East and North Africa (15.2, 8.8-20) during 2005-2018 [5] . Independent studies show that factors affecting the occurrence of GDM include maternal weight, gestational weight gain, diet and exercise, ethnicity, maternal age, and family history of diabetes [6] [7] [8] [9] [10] [11] . Women with GDM are prone to pregnancy-related complications or perinatal adverse events, such as hypertension, preeclampsia, preterm birth, shoulder dystocia, perinatal morbidity, and death [12] . Fetuses exposed to GDM have an increased risk of macrosomia, fetal distress, congenital malformations, and jaundice [13] . Meanwhile, the offspring of GDM women has an increased risk of metabolic diseases and cardiovascular diseases, such as obesity and type II diabetes. GDM is also related to immune function loss and cognitive impairment in the offspring [14] . The molecular mechanism of insulin resistance and insulin secretion defect in pregnant women is still under intensive investigation.

From a fertilized egg to the complete fetus, human genome remains relatively stable while the epigenome undergoes regulated reprogramming including dramatic DNA methylation changes and systemic histone modification alternations [15, 16] . Moreover, the epigenome of fetus is sensitive to the pregnant women's intrauterine environment, especially high-glucose levels, which eventually leads to phenotypic changes [17] . How the exposure of fetal to maternal GDM contributes to increasing offspring risk remains to be explored. Results from animal studies suggest that the exposure to intrauterine high-glucose environment may cause disorders of hypothalamic neuropeptide neurons, impair blood glucose homeostasis and lead to abnormal renal development and islet β-cell dysfunction [18, 19] . Epidemiological studies have also confirmed that exposure of intrauterine high-glucose environment may cause the transfer of lipid matrix molecules to the fetus, resulting in disorders associated with adipocyte dysfunction and fatty acid accumulation [20] . Using bisulfite pyrosequencing assays, Houde and colleagues identified three candidate genes, LRP1B, BRD2 and CACNA1D, whose DNA methylation alterations in placenta and umbilical cord blood were related to energy metabolism and maternal glucose levels [21] . Through methylated DNA immunoprecipitation (MeDIP) assay and next-generation sequencing (NGS), Rong and colleagues revealed genome-wide changes in the placenta of GDM patients and transcriptional variations of four GDM-related genes [22] . However, the genes and pathways involved in GDM have not been systematically characterized and thus remain elusive. Therefore, genome-scale DNA methylation analysis and transcriptome studying of the placenta and paired umbilical cord blood from women with or without GDM can potentially uncover the underlying molecular mechanism.

Reduced representation bisulfite sequencing (RRBS) is an accurate and cost-effective method that can provide genome-wide mapping 5-methylcytosines, particularly enriching for coverage in promoter and enhancer regions that are crucial for transcriptional regulation [23, 24] . RNA-seq allows to profile thousands of gene expressions precisely and simultaneously [25] . Here, we employed the RRBS and RNA-seq methods to systematically analyze the genome-wide methylation and expression profiles of placenta and paired umbilical cord blood from a group of patients with GDM and control pregnant women without GDM. The aim of this study is to reveal critical epigenetic changes associated with GDM in pregnant women and to identify differentially expressed genes (DEGs) which may potentially contribute to GDM.

Forty-one pregnant women were recruited into this study, including 21 pregnant women with GDM and 20 without GDM as control. Those pregnant women delivered 41 babies, including 4 babies with macrosomia, 2 of which were from the GDM group. Anthropometric and metabolic parameters of pregnant women and newborns are shown in Table 1 . We observed that both the maternal weight and body mass index (BMI) were significantly higher in the GDM group compared with those in the control group (Wilcoxon rank-sum test, P value < 0.05), while other characteristics, such as maternal age, gestational weight gain, maternal height, gestational weeks, and fetal birth weight, did not show significant differences between two groups. Given the significant differences of maternal BMI between GDM and control groups, BMI will be used as a covariate in the subsequent search for differentially methylated regions (DMRs) and DEGs.

DMRs are genomic regions with differential methylation status across multiple CpG sites, which are considered to play an important role in gene imprinting and regulating the expression of nearby genes that were called differentially methylated genes (DMGs). In placenta samples, after excluding the effect of maternal BMI, a total of 2779 DMRs were identified, including 1446 hyper-methylated DMRs (hyper-DMRs) associated with 1293 hyper-methylated DMGs (hyper-DMGs) and 1333 hypo-methylated DMRs (hypo-DMRs) associated with 1197 hypo-methylated DMGs (hypo-DMGs) (Fig. 1A , Additional file 2: Table S2 ). Among them, 209 hyper-DMRs and 189 hypo-DMRs were located in promoter regions, respectively. The mean methylation of hyper-DMR and hypo-DMR in GDM and control groups is shown (Fig. 1B) . As expected, compared with the control group, the distribution of mean methylation levels of all hyper-DMRs or hypo-DMRs in GDM group showed significant difference (Wilcoxon rank-sum test, P value < 0.05). Additionally, hyper-DMRs and hypo-DMRs widely spread across different chromosomes (Additional file 1: Fig. S1C ), indicating that the placenta has extensive and significant methylation changes in a genome-wide manner.

To further characterize DMRs identified, we performed functional enrichment analysis using Gene Ontology (GO) terms and identified 401 and 288 significantly enriched terms for hyper-DMGs and hypo-DMGs (FDR < 0.05), respectively (Additional file 2: Table S3, S5) . For hyper-DMGs, the enriched pathways mainly included neural and developmental pathways, endocrine-, and immune-related pathways (Additional file 1: Fig. S1D ). For hypo-DMGs, the enriched GO terms were similar to that in hyper-DMRs, but also included other ones, such as "cognition, " "long-term memory, " and epigenetic modification-related pathways (Additional file 1: Fig. S1E ). Meanwhile, enrichment analysis with Kyoto Encyclopedia of Genes and Genomes (KEGG) database identified 31 and 3 pathways for hyper-DMGs and hypo-DMGs, respectively ( Fig. 1C , D, Additional file 2: Table S4 , S6, FDR < 0.05). Of the pathways associated with hyper-DMGs, "oxytocin signaling" (29 genes, FDR = 7.14e-05), "calcium signaling" (28 genes, FDR = 5.47e-06), and "insulin signaling" (19 genes, FDR = 8.63e-03) ranked in the top of the list (Fig. 1C) . Interestingly, the "type II diabetes mellitus" pathway was significantly enriched in hypo-DMGs (FDR = 0.033), suggesting that GDM may have the epigenetic characteristics of type II diabetes mellitus (Fig. 1D) .

To further explore the pathogenesis of GDM, we next performed RNA-seq analysis using placenta samples from the GDM and control groups. After adjustments for multiple testing, we obtain 113 up-regulated and 9 down-regulated genes (FDR < 0.05). Given this number of DEGs is so small, it is challenging to interpret these genes in terms of biological processes and signaling pathways. However, this is not totally uncommon and has been observed in several previous type II diabetes mellitus-related studies [26, 27] . For this reason, researchers have developed the Gene Set Enrichment Analysis (GSEA) method, which was designed to detect modest but coordinate changes in the expression of groups of functionally related genes [28] . This also suggests that top-ranking genes are still good candidates for subsequent analysis, even though individual one of which is not statistically significant. For this reason, we used a P value-based threshold (|Fold changes|> 2, P value < 0.05) to select DEGs. With this threshold, 1001 genes were upregulated, and 441 genes were down-regulated ( Fig. 2A , Additional file 2: Table S7 , S8). The up-regulated genes were annotated to multiple endocrine and metabolicrelated pathways, including "endocrine resistance, " "cortisol synthesis, and secretion, " "dilated cardiomyopathy, " and "insulin secretion" (Fig. 2B , Additional file 2: Table S9 and S10). For down-regulated genes, 341 and 32 terms were enriched in GO and KEGG functional enrichment analysis, respectively (Additional file 2: Table S11 and S12). Examples of enriched KEGG terms include "antigen processing and presentation" (FDR = 8.89e-12), and "type I diabetes mellitus" (FDR = 2.81e-10) (Additional file 1: Fig. S2A ). GO enrichment analysis showed that pathways related to immune regulation, cell proliferation, and glycolipid metabolism were most enriched (Additional file 1: Fig. S2B ). We next conducted GSEA using expression profile of placenta to detect modest but coordinate changes in the expression of groups of functionally related genes. Of note, 28/40 and 31/32 of KEGG positive and negative enriched pathways were validated by GSEA results, respectively (FDR < 0.25), which show that DGE-based results are largely consistent with that from GSEA in the pathway-level (Additional file 2: Table S13, S14). Interestingly, we discovered two significant insulin-related pathways, "insulin secretion" (FDR = 0.0898, Fig. 2C ) and "insulin resistance" (FDR = 0.111, Additional file 1: Fig. S2D ). Other positively enriched pathways include the "notch signaling pathway, " "autophagy, " "growth hormone synthesis, secretion and action, " "prolactin signaling pathway, " and "aldosterone synthesis and secretion pathways" ( Table 2 ). Most notably, pathway "type I diabetes mellitus" ranked first in the negative enriched pathways (FDR < 0.001, Additional file 1: Fig. S2E ). Other negatively enriched pathways, such as "one carbon pool by folate" (FDR = 3.65e-3), "primary bile acid biosynthesis" (FDR = 6.96e-4), and "folate biosynthesis pathway" (FDR = 0.012), are listed in Table 3 .

Next, to investigate whether methylation will affect the expression of genes that may attribute to GDM, we combined methylation data and transcriptome data for an integrated analysis. The analysis indicated that 155 genes showed both methylation changes and expression alternations between the GDM and the control group (P value < 2.2e-16, odds ratio = 4.69, Fisher's Exact Test, Fig. 2D ). Though not statistically significant after adjustment for multiple comparisons, the top-ranking pathways enriched in the overlapping genes include "Glutamatergic synapse" (ADCY5/ADCY7/GRIK4/PLA2G4A/ PLCB1/SHANK2/SHANK3/SLC1A6, P value = 5.32e-06), "cortisol synthesis and secretion" (ADCY5/ADCY7/ LDLR/PLCB1, P value = 2.29e-03), "insulin secretion" (ADCY5/ADCY7/PLCB1/RYR2, P value = 6.28e-3), "insulin signaling pathway" (HK2/PRKCZ/SH2B2/SHC2/ SOCS3, P value = 6.34e-3), "type II diabetes" (HK2/ PRKCZ/SOCS3, P value = 7.10e-3) and "bile secretion" (ADCY5/ADCY7/LDLR, P value = 0.040) (Fig. 2E ). We then selected 10 genes of interest and showed their mean methylation difference, expression change, and associated pathways in Table 4 . Besides, from the 155 overlapped gene, we screened out four gene, SFN, PLCB1, TFAP2A and DAB2IP, whose gene expression level and gene mean methylation level show significant correlation and functional for glycolipid metabolism (Fig. 2F ). Meanwhile, using the GSEA method, placenta CpG methylation generally had a negative regulatory effect on gene expression (P value < 0.001, NES = -0.255, Additional file 1: Fig. S2F ). In summary, we found that some genes involved in glycolipid metabolism, insulin secretion, and insulin resistance had been epigenetically modified, which might contribute directly or indirectly to GDM.

To explore whether the exposure to increased glucose level could cause epigenetic changes in fetus, we next conducted RRBS and RNA-seq analyses using the umbilical cord blood paired with placenta described above. Due to the existence of placental barrier, the changes of fetal epigenome were subtler comparing with that in placenta. Genome-wide methylation profiling by RRBS identified far less numbers of DMRs in contrast to that from the paired placenta. Specifically, after excluding the effect of maternal BMI, we detected 75 hyper-DMRs and 66 hypo-DMRs which annotated to 74 hyper-DMGs and 65 hypo-DMGs, respectively (Fig. 3A , B; Additional file 1: Fig. S3A ). Next, KEGG pathway enrichment analysis was conducted to identify pathways of differentially methylated genes. Though not statistically significant after adjustment for multiple comparisons, some enlightening pathways were identified in hypo-DMRs (Fig. 3B , red bar), including "cAMP signaling pathway" (ATP2B3/DRD5), "type II diabetes mellitus" (PRKCD), "cholesterol metabolism" (OSBPL5), "aldosterone synthesis and secretion" (ATP2B3), "insulin resistance" (PRKCD). On the other hand, pathways "endocrine resistance" (JAG2/RB1), "PI3K-Akt signaling pathway" (BCL2L11/LPAR2/TCL1A), and "bile secretion" (SLC10A2) were found in hyper-methylated genes (Fig. 3B , blue bar). RNA-seq analyses of umbilical cord blood samples revealed 271 up-regulated genes and 217 down-regulated genes (|Fold change|> 2, P value < 0.05, Additional file 2: Table S15, S16) by comparing GDM patients to the controls (Fig. 3C) . KEGG enrichment analysis showed that the up-and down-regulated genes were associated with 21 and 6 pathways, respectively (P value < 0.05, Fig. 3D ). Interestingly, two pathways, "human papillomavirus infection" and "ovarian steroidogenesis" were overlapped for up-regulated genes between placenta and umbilical cord blood, showing that there is a certain relationship of expression profiles between pregnant women and offspring.

Due to minor changes in gene expression levels, GSEA was then performed as a cutoff-free method for pathway enrichment analysis (Additional file 2: Table S17, Fig. 3 Alterations for umbilical cord blood were related to insulin secretion and resistance. A Genomic annotation of DMRs in umbilical cord blood. Blue bar represents hyper-DMRs and red bar represents hypo-DMRs. B KEGG pathway enrichment analysis of hypo-DMGs (red bar) and hyper-DMGs (blue bar) of umbilical cord blood. C Volcano plot of umbilical cord blood expression profile. D KEGG pathway enrichment analysis of up-DEGs (red bar) and hyper-DEGs (blue bar) in umbilical cord blood S18). In details, "collecting duct acid secretion" (FDR = 0.235), "FC gamma R-mediated phagocytosis" (FDR = 0.244), and "staphylococcus aureus infection" (FDR = 0.249) were significantly activated. By contrast, "systemic lupus erythematosus" (FDR = 0.0816), "one carbon pool by folate" (FDR = 0.0844), "fatty acid biosynthesis" (FDR = 0.203), and "primary bile acid biosynthesis" (FDR = 0.239) were significantly repressed. It is worth noting that pathway "one carbon pool by folate" were also enriched in the placenta samples. These results preliminarily show that there is a strong correlation between the current genomic changes and the offspring, and these changes are related to glucose metabolism and insulin, suggesting the possibility of metabolic diseases in the offspring.

In umbilical cord blood, we found one overlapping gene ZNF423 (DNA methylation delta = − 0.4, log 2 (Fold Change) = − 2.15), whose gene methylation and expression were both altered. Expression of ZNF423 has been reported directly to be proportional to the size of adipocytes [29] . The gene product of ZNF423 promotes the transformation of non-adipocytes into adipocytes and inhibits subcutaneous adipogenesis, leading to adipocyte hypertrophy and inflammation, and finally developing to obesity and insulin resistance [30] . The significant down-regulation of ZNF423 in umbilical cord blood may be related to the higher possibility of suffering obesity in offspring exposed to GDM.

After discovering altered pathways shared by umbilical cord blood and placenta samples, we further explored the correlation between placenta and umbilical cord blood through an integrated analysis. The analysis revealed that 36 DMRs co-existed in the placentas and paired umbilical cord bloods (P value < 2.2e-16, odds ratio = 9.04, Fisher's exact test, Fig. 4A ). The top-ranking pathways enriched in Fig. 4 Significant correlation between the changes of placenta and umbilical cord blood, and involves glycolipid metabolism. A The Venn Diagram of DMGs in placenta and umbilical cord blood. Significance was tested by Fisher's Exact Test. B The Venn Diagram of DEGs in placenta and umbilical cord blood. Significance was tested by Fisher's Exact Test. C Mean methylation level of OSBPL5 DMR for placenta and umbilical cord blood in different group, log 2 TPM of each group is annotated on the right of each panel. Yellow bar represents GDM samples and green bar represents control samples the shared DMRs include "glycosphingolipid biosynthesis-globo and isoglobo series" (A4GALT, P value = 0.012), "N-glycan biosynthesis" (ALG10, P value = 0.042) and "cholesterol metabolism" (OSBPL5, P value = 0.042).

We then investigated altered gene expressions in placentas and paired umbilical cord blood samples by comparing GDM patients to the controls. The results demonstrated that 23 genes were overlapped between placenta and paired umbilical cord blood (P value = 2.3e-03, odds ratio = 2.04, Fisher's exact test, Fig. 4B ). The top-ranking pathways enriched in the shared genes include "ECM-receptor interaction" (FRAS1/SPP1, P value = 4.92e-03), "growth hormone synthesis, secretion and action" (ADCY5/SOCS3, P value = 8.84e-03) and "type II diabetes mellitus" (SOCS3, P value = 0.055).

The gene OSBPL5 encode a protein which belongs to the oxysterol-binding protein family, and proteins in this family are critical to maintenance of cholesterol balance in human body [31] . Comparison of the methylation levels of OSBPL5 across multiple specimens is illustrated in Fig. 4C . The DMR in placenta was annotated into one of the exons in the OSBPL5 gene, including four CpG sites (DNA methylation delta = 0.35; Fig. 4C left panel) . Likewise, the DMR in umbilical cord blood located in an intron, including six CpG sites (DNA methylation delta = − 0.38; Fig. 4C right panel) . Profiling of transcriptome showed that expression of OSBPL5 (log 2 TPM) in the placenta of GDM patient was 2.4-fold of that in the control group (GDM 4.034 vs. control 2.767), but the difference was small in umbilical cord blood (GDM 0.757 vs. control 0.863). These results suggest that the methylation and expression profile of umbilical cord blood are closely related to that in placenta, and the genes in glucose metabolism-related pathway are more likely to be regulated by DNA methylation.

To the best of our knowledge, this is the first study that use epigenomic and transcriptomic assays to characterize placenta and the paired umbilical cord blood samples from GDM patients and controls. Our data demonstrated that placenta undergoes extensive methylation changes in genomic regions that are related to glucose metabolismrelated pathways. Remarkably, methylation-altered genes in umbilical cord blood were associated with pathways insulin resistance and insulin secretion. And we also found that DMGs and DEGs were significantly overlapped between placenta and umbilical cord blood, indicating that the GDM conditions could affect fetus.

An early study found that insulin resistance in GDM patients led to a compensatory increase of the synthesis and secretion of insulin, which can promote the absorption and metabolism of blood glucose in islet β-cells [32] .

However, when insulin resistance increases to a certain level, over-secreted insulin fails to maintain blood glucose levels in a normal range [33] and GDM patients could display hyperglycemia and hyperinsulinemia simultaneously [34] . GSEA result of placental RNA-seq data showed that several genes associated with the insulin secretion and the insulin resistance pathways were significantly up-regulated. More importantly, genes with both DNA methylation alternations and expression changes were identified, including PRKCZ, ADCY5, CACNA1C, PLCB1. The DNA methylation alternations in those genes might contribute to expression changes of these gene, which further aggravates insulin resistance and insulin secretion.

At the beginning of pregnancy, the placenta releases hormones into the blood of pregnant women, raising the blood glucose concentration to ensure adequate nutrition for the fetus [35, 36] . Our data showed that many genes related with hormone production and secretion were highly expressed in placenta of GDM patients, which involved in "cortisol synthesis and secretion, " "growth hormone synthesis and secretion" and "aldosterone synthesis and secretion pathway" (Table 2 ). However, it is reported that excessive cortisol levels in pregnant women's blood can promote anxiety [37] , increase the risk of obesity in offspring [38] and reduce children's cognitive ability [39] . Our investigation discovered decreased methylation levels in multiple genes of the cortisol synthesis and secretion pathway, suggesting that cortisol levels could be altered in GDM patients through placenta. This hypothesis is further supported by the GSEA result of RNA-seq data, which listed pathway "cortisol synthesis and secret" as one of the most positively enriched pathways (Additional file 1: Fig. S2C ). The genes of which both methylation statuses and expression levels varied in this pathway included CACNA1C, LDLR, and ADCY9. The CACNA1C gene encodes an alpha-1 subunit of a voltage-dependent calcium channel and is related to the exocytosis of many hormones. Xu et al. [40] reported that miR-153 worked as a negative regulator for the expression of CACNA1C, which further increased the secretion of insulin. An independent GWAS study revealed that CACNA1C gene is associated with diabetic cataract [41] . Besides, the co-binding of low-density lipoprotein receptor (LDLR) and insulin receptor (IR) reduces the ability of LDLR to clean up LDLR in blood, and insulin can destroy this co-binding [42] .

Moreover, compared with control pregnancy, GDM patients had significantly lower serum total bile acids at 24 weeks of gestation (Additional file 1: Fig. S4A ) and at 40 weeks of gestation (Additional file 1: Fig. S4B ) (Wilcoxon rank-sum test, P value < 0.05). Previous studies have shown that bile acids can activate glucose induced insulin secretion [43] and have a strong correlation with insulin resistance and nonalcoholic fatty liver [44] . Based on the technique of liquid chromatography, Li et al. reported that lower levels of serum bile acids in early pregnancy were independently associated with an increased risk of GDM in Chinese pregnant women [45] . At the same time, our GSEA analysis also showed that several genes associated with the primary bile acid biosynthesis pathway were significantly down-regulated in placental tissues (Additional file 1: Fig. S4C ). Expression of CH25H, CYP7B1 and HSD3B7 were significantly down-regulated in core enrichment genes, which played an important role in catalytic synthesis of bile acids. CYP7B1 is also involved in the synthesis of dehydroepiandrosterone and pregnenolone [46] . Both dehydroepiandrosterone and pregnenolone are related to hippocampus-associated memory and learning, which may reveal the relationship between GDM patients and the risk of memory decline. We believe that this is the first report to connect bile acid metabolism with GDM using placental high-throughput sequencing data. Notably, bile acid has been reported for the treatment of patients with obesity and diabetes [47] , suggesting that bile acids can potentially be used in pregnant women to prevent and treat GDM in the future.

In placental tissues, 155 genes were differentially methylated and differentially expressed (partially displayed in Table 4 ), and some of these genes were involved in key metabolic pathway. From those genes, we manifested that the expression of three critical genes in the type II diabetes pathway was altered, including HK2 (downregulated), PRKCZ (up-regulated) and SOCS3 (downregulated). HK2 is the key gene of glycolysis, and its down-regulation has been shown directly leading to the reduction of glucose metabolism [48] . SOCS3 can inhibit insulin secretion, and down-regulation of this gene promote insulin secretion, resulting in insulin resistance [49] . The activation of PRKCZ depends on insulin, and over-expression of this gene affects the expression of the insulin-like growth factor 1 receptor (IGF1R) gene [50] . Additionally, previous research show that genes illustrated in Fig. 2F play key roles in insulin secretion, lipid droplet formation and are associated with type II diabetes [51] [52] [53] . Taken together, those results suggested that the pathogenesis of GDM may be similar to type II diabetes, characterized with excessive insulin secretion and insulin resistance, and DNA methylation plays an important role in the process.

Previous studies have indicated that the methylation levels of metabolism-related genes or genes involved in pathways could be changed after exposing to persistent high-glucose levels in pregnant women with GDM, which is termed "metabolic programming" [54] [55] [56] . In present results, several pathways showed descending methylation levels in placenta of GDM patients, including the "aldosterone synthesis and secretion pathway, " "insulin resistance" and "type II diabetes mellitus pathway". Besides, for up-regulated genes in placenta of GDM patients, enriched pathways included "EGFR tyrosine kinase inhibitor resistance, " "insulin secretion and "growth hormone synthesis, secretion and action". Our results indicated that the glucose and lipid metabolism of offspring were probably affected by intrauterine exposure to GDM, both in epigenetic and expression levels. Furthermore, there was a significant relevance between the changes of placenta and umbilical cord blood, regardless differentially methylated genes (Fig. 4A ) or differentially expressed genes (Fig. 4B) . Most notably, genes in the one carbon pool by folate pathway were significantly down-regulated in placentas and umbilical cord bloods, indicating that both pregnant women and fetuses with GDM may have the characteristics of folate deficiency. However, a meta-analysis regarding to the correlation between folic acid and the pathogenesis of GDM provides controversial results [57] . Larger cohort studies are necessary to explore the correlations. The above results show that there are apparent associations between the methylome and the transcriptome in placenta and umbilical cord blood, which share the common characteristics of insulin resistance and type II diabetes.

Still, there are few limitations in this study. Firstly, due to the complexity of placenta dissection structure, it is inevitable to obtain a variety of tissues and cell types during sampling. Using single-cell sequencing technology or spatial transcriptome technology will make a more in-depth study on the pathogenesis and mother-child interaction of GDM comparing with traditional bulk RNA-seq. In addition, in differentially gene expression analysis, we did not correct P value by FDR method, because very few genes (placenta 122, umbilical cord blood 0) had significant differences with FDR cutoff of 5%. Therefore, we used GSEA method to support our point of view, which results are largely consistent with DGE-based results. Finally, as confounding factor, BMI weights much in present research, and we did not find a suitable tool to exclude the influence of BMI during calling DMRs. Instead, we divided samples into two groups with equal number of samples, BMI-high and BMI-low, and used DSS to identify BMI-related DMRs, which were removed from subsequent analysis.

In summary, our study systematically profiled the methylome and transcriptome of placenta and umbilical cord blood of pregnant women with or without GDM, and investigated the relationship between methylation and gene expression. Our results provide detailed molecular changes among GDM patients compared with the controls, which pave the way for the further investigation of GDM pathogenicity.

Pregnant women were recruited at Hangzhou Women's Hospital (Hangzhou Maternity and Child Health Care Hospital, China) from January 1, 2020 to May 31, 2020. All participants signed a written consent forms and this study was approved by the ethics committee of Hangzhou Women's Hospital. GDM is diagnosed according to the International Association of Diabetes and Pregnancy Study Groups (IADPSG) guidelines. Specifically, at least one of the following conditions must be met: 1. Fasting blood glucose ≥ 5.1 mmol/L; 2. Blood glucose ≥ 10.0 mmol/L one hour after the oral glucose tolerance test (OGTT); 3. Blood glucose ≥ 8.5 mmol/L two hours post the OGTT [58] . Pregnant women who met the following conditions were enrolled in this study: (1) age 24-39; (2) singleton pregnancy; (3) full-term birth (37-40 weeks of gestation); (4) detailed information for the OGTT; (5) records showing that the perinatal examination, delivery and infant physical examination were conducted in our hospital. Pregnant women with normal blood glucose and no history of adverse diseases were recruited in the control group. We excluded the pregnant women with the following conditions from our study: (1) having chronic diabetes complicated with pregnancy or preeclampsia. (2) suffering from liver or kidney dysfunctions, or other chronic diseases which required longterm drug treatment; (3) showing a mental disorder or a serious infection; (4) smoking, drinking and drug abuse during pregnancy; (5) fetal chromosome abnormalities; (6) lacking registration information in the Women's and Children's Insurance Systems. The difference of physiological parameters of pregnant women and newborns between GDM group and control group was performed by Wilcoxon rank-sum test.

Within 15 min after delivery, 5 ml of neonatal umbilical cord blood was collected using two PAXGenRNA sampling tubes (2.5 ml/tube), and collected samples were stored at -80 ℃ [59] . The placental tissues were collected from both the maternal and fetal side. Specifically, using a clean scalpel, one full-thick placental tissue was incised 5 cm from the perimeter of the placenta (a size of about 1 × 1 × 2 cm, including maternal and fetal side). Then, the tissues were rinsed with phosphate buffer (PBS) or normal saline until no visible blood and stored at -80 ℃ in a cryopreservation tube [60] . DNA or RNA was extracted after homogenizing the tissues together.

Genomic DNA of placenta and umbilical cord blood specimens was extracted by using the TIANamp genomic DNA kit (Cat. No. DP304-02, TIANGEN Biotech CO. LTD., Beijing, China) per the manufacture's recommendations. 50-300 ng of purified genomic DNA was utilized to generate RRBS libraries as previously described [61] . Indexed RRBS libraries were pooled accordingly and sequenced in a NovaSeq 6000 sequencer (Illumina, USA) with 100-base paired-end reads. Similarly, total RNAs were isolated using a miRNeasy Mini kit (QIAGEN, Germany) according to the manufacture's recommendations. We used the Ribo-off rRNA depletion kit (Cat. No. N406-01, Vazyme Biotech Co. Ltd, Nanjing, China) to remove ribosomal RNA from total RNAs following the manufacturer's instructions. To generate RNA-seq libraries, 10-100 ng rRNA-depleted RNAs were utilized by using the VAHTS total RNA-seq (H/M/R) library prep kit for Illumina (Cat. No. NR603-01, Vazyme Biotech Co. Ltd, Nanjing, China) per the manufacturer's instructions. We sequenced the RNA-seq libraries in the NovaSeq 6000 sequencer with 100-base paired-end reads.

The sequencing data volume of a single sample is around 8 gigabases (Gb) and the average sequencing depth of the target region of the enzyme digestion fragment is > 30 × , the coverage proportion of CpG region is > 40%, and the total data volume is 100 Gb. Similarly, transcriptome profile for each placenta and umbilical cord blood sample was also conducted. VAHTSTM Total RNA-seq (H/M/R) Library Prep Kit for Illumina ® based on Ribo-zero method was used to constructed transcriptome RNA library and high-throughput sequencing was performed subsequently. 10 Gb volume of transcriptome data for each sample was necessary.

For RRBS data, FastQC (version 0.11.9) [62] software was used to evaluate the quality of the original RRBS sequencing reads, and TrimGalore (version 0.6.5) [63] was further applied to trim adapters and low-quality bases. The trimmed reads were then aligned to the hg19 genome using BSMAP [64] software. Methylation calling of BAM file obtained in the previous step was conducted by MethylDackel [65] . R package DSS [66] was used to discover de novo DMRs. In this research, DMRs were defined as regions whose mean methylation level across CpG sites increased or decreased > 0.1, width > 50 bp, CpG number > 3, and P value < 0.05. Specifically, statistical tests for differential methylation at each CpG site were performed by function "DMLtest" with default parameters, which output was then processed using function "callDMR" (parameters: delta = 0.1, p.threshold = 0.05) to call DMRs-GDM by comparing GDM group with control group. Meanwhile, to eliminate the influence of confounding factors, we divided samples into two groups with equal number of samples, BMIhigh and BMI-low, and used DSS to identify BMI-related DMRs, which were removed from subsequent analysis. Finally, DMRs were annotated using R package ChIPseeker [67] with default parameters to find differentially methylated genes (DMGs) and genomic distribution of DMRs.

As for RNA-seq data, quality control and base trimming were similarly conducted by using FastQC and TrimGalore, respectively. STAR (2.6.1d) [68] was used to map the quality-controlled reads to the hg19 genome. Then, the transcripts were quantified by using software kallisto [69] , and the expression abundance of transcripts was generated accordingly. RNASeQC is used to count the quality of genome alignment. Subsequently, using the output of kallisto, gene expression levels were estimated from the abundance of transcription level by function "tximport" of R package tximport [70] . Gene differential expression analysis between GDM group and control group was performed by function DESeq of R package DESeq2 [71] using default parameters. Briefly, dispersion was estimated by fitting a dispersion-mean relation via a robust gamma-family GLM (General linear model); library size was estimated by the standard median ratio method introduced in DESeq; finally, the GLM coefficients were tested for significance using Wald test. Meanwhile, BMI was added as a covariate in the difference analysis to exclude the influence of confounding factors.

The unique mapping ratio of RRBS and RNA-seq data were illustrated in Additional file 1: Fig. S1A , B, respectively. All quality control table of RRBS and RNA-seq data were listed in Additional file 2: Table S1 .1-S1.4.

The genomic annotation of DMRs was conducted by R package ChIPseeker [67] and consequently were grouped to 6 classes: promoter, exon, intron, intergenic, UTR (untranslated region) and downstream regions. To investigate the enriched functional pathways of DMGs and DEGs, functional enrichment analysis of GO (Gene Ontology) biological process and KEGG pathway annotation was performed by the R package clusterProfiler [72] (parameters: pvalueCutoff = 0.05, qvalueCutoff = 0.05, pAdjustMethod = "fdr").

Gene Set Enrichment Analysis (GSEA) was used to detect modest but coordinate changes in the expression of groups of functionally related genes [73] . Java command line tool gsea (http:// www. gsea-msigdb. org/ gsea/ msigdb/ downl oad_ file. jsp? fileP ath=/ resou rces/ softw are/ gsea2-2. 2.4. jar) was used to find the key pathways in gene expression profiles of placenta and umbilical cord blood (parameters: xtools.gsea.GseaPreranked -scor-ing_scheme weighted -collapse true -mode Max_probe -norm None -nperm 1000 -include_only_symbols true make_sets true -plot_top_x 20 -rnd_seed timestamp -set_max 3000 -set_min 15 -zip_report false -gui false). Specifically, all genes in gene expression profile were firstly sorted by Z-scores, which were converted from P values and associated log2 fold change using function qnorm. Then, enrichment significance test of each KEGG pathway was performed by GSEA, and P values of all pathways were adjusted for multiple testing. According to official advice, pathways with FDR < 0.25 were considered statistically significant.

To explore the effect of DNA methylation on gene expression under GDM environment, we compared differentially methylated genes and differentially expressed genes in a tissue of interest, including placenta and umbilical cord blood. A gene was counted as an overlap if it was differentially expressed and near DMR. The statistical significance of overlapping was tested by Fisher's exact test with all genes in gene expression profile as a background. Similarly, we compared DMGs in placenta and umbilical cord blood, as well as DEGs in placenta and umbilical cord blood. Spearman's rank correlation test was used to assess the significance of correlation between mean methylation level of DMGs and mean gene expression level of DEGs across samples. The biological function of overlapped genes was then evaluated by GO/ KEGG functional enrichment analysis. 

The online version contains supplementary material available at https:// doi. org/ 10. 1186/ s13148-022-01289-5.

Additional file1: Fig. S1 Placenta shows genome-wide methylation alteration associated with glucose metabolism in GDM patients. Unique mapping ratio of placenta and blood samples for A RRBS data and B RNA-seq data. C Genomic distribution of DMRs in placenta. GO pathway enrichment result of placenta for D hyper-DMGs and E hypo-DMGs. The number in brackets represents the number of enriched genes. Fig. S2 Methylation contributes to expression change of genes associated with insulin signaling pathway. A KEGG pathway enrichment analysis of down-DMGs in placenta. B "Cortisol synthesis and secretion" C "Adherens junction, " D "Insulin resistance" and E "Type I diabetes mellitus" pathway enrichment result of placenta expression profile GSEA analysis. F GSEA analysis of placenta DMGs to placenta expression profile. Fig. S3 Alterations for umbilical cord blood were related to insulin secretion and resistance A Genomic distribution of DMRs in umbilical cord blood. B Differences of DMRs mean methylation levels between GDM and control umbilical cord blood samples. Yellow box represents GDM samples and green box represents control samples. *** P value < 0.001, Wilcoxon rank-sum test. Table S9 GO pathway enrichment of UP-DEGs of placenta. Table S10 KEGG pathway enrichment of UP-DEGs of placenta. Table S11 GO pathway enrichment of DN-DEGs of placenta. Table S12 KEGG pathway enrichment of DN-DEGs of placenta. Table S13 Positive enriched pathways of GSEA for placenta gene expression profile. Table S14 Negative enriched pathways of GSEA for placenta gene expression profile. Table S15 Up-regulated differentially expressed genes of umbilical cord blood. Table S16 Down-regulated differentially expressed genes of umbilical cord blood. Table S17 Positive enriched pathways of GSEA for umbilical cord blood gene expression profile. Table S18 Negative enriched pathways of GSEA for umbilical cord blood gene expression profile.

A American Diabetes. 2. Classification and diagnosis of diabetes: standards of medical care in diabetes-2020

Gestational diabetes mellitus

Gestational diabetes: diagnosis classification, and clinical care

The pathophysiology of gestational diabetes mellitus

Gestational diabetes mellitus

The association of general obesity, central obesity and visceral body fat with the risk of gestational diabetes mellitus: evidence from a systematic review and meta-analysis

The association between weight gain at different stages of pregnancy and risk of gestational diabetes mellitus

A lifestyle intervention of weightgain restriction: diet and exercise in obese women with gestational diabetes mellitus

The prevalence of gestational diabetes mellitus recurrence-effect of ethnicity and parity: a metaanalysis

Maternal age and the risk of gestational diabetes mellitus: a systematic review and meta-analysis of over 120 million participants

Family history of diabetes and the risk of gestational diabetes mellitus in Iran: a systematic review and meta-analysis

Gestational diabetes mellitus: mechanisms treatment, and complications

Gestational diabetes mellitus epigenetically affects genes predominantly involved in metabolic diseases

Gestational diabetes mellitus

DNA methylation: roles in mammalian development

The mammalian epigenome

Developmental plasticity and developmental origins of non-communicable disease: theoretical considerations and epigenetic mechanisms

Programming' of orexigenic and anorexigenic hypothalamic neurons in offspring of treated and untreated diabetic mother rats

Altered nephrogenesis due to maternal diabetes is associated with increased expression of IGF-II/mannose-6-phosphate receptor in the fetal kidney

Maternal serum triglyceride at 24-32 weeks' gestation and newborn weight in nondiabetic women with positive diabetic screens

BRD2 and CACNA1D: new candidate genes in fetal metabolic programming of newborns exposed to maternal hyperglycemia

DNA methylation profiles in placenta and its association with gestational diabetes mellitus

Genome-scale DNA methylation mapping of clinical samples at single-nucleotide resolution

Preparation of reduced representation bisulfite sequencing libraries for genome-scale DNA methylation profiling

RNA-seq: a revolutionary tool for transcriptomics

Gene expression profile in skeletal muscle of type 2 diabetes and the effect of insulin treatment

Microarray profiling of skeletal muscle tissues from equally obese, non-diabetic insulin-sensitive and insulinresistant pima Indians

PGC-1alpha-responsive genes involved in oxidative phosphorylation are coordinately downregulated in human diabetes

Epigenetic modifications of the Zfp/ ZNF423 gene control murine adipogenic commitment and are dysregulated in human hypertrophic obesity

Adipogenesis and metabolic health

The OSBP-related protein family in humans

Diabetes mellitus: complications and therapeutics

Mechanisms of insulin action and insulin resistance

Hyperglycemia-and hyperinsulinemia-induced insulin resistance causes alterations in cellular bioenergetics and activation of inflammatory signaling in lymphatic muscle

Elevated maternal cortisol leads to relative maternal hyperglycemia and increased stillbirth in ovine pregnancy

Intermediate metabolism in normal pregnancy and in gestational diabetes

Salivary cortisol response to infant distress in pregnant women with depressive symptoms

Maternal cortisol during pregnancy and infant adiposity: a prospective investigation

Elevated maternal cortisol levels during pregnancy are associated with reduced childhood IQ

The Ia-2beta intronic miRNA, miR-153, is a negative regulator of insulin and dopamine secretion through its effect on the Cacna1c gene in mice

A genome-wide association study provides new evidence that CACNA1C gene is associated with diabetic cataract

The association between insulin and low-density lipoprotein receptors

Bile acids acutely stimulate insulin secretion of mouse beta-cells via farnesoid X receptor activation and K(ATP) channel inhibition

Bile acids and insulin resistance: implications for treating nonalcoholic fatty liver disease

Bile acid metabolites in early pregnancy and risk of gestational diabetes in Chinese women: a nested case-control study

CYP7B1: one cytochrome P450, two human genetic diseases, and multiple physiological functions

Bile acid control of metabolism and inflammation in obesity, type 2 diabetes dyslipidemia, and nonalcoholic fatty liver disease

HK2/hexokinase-II integrates glycolysis and autophagy to confer cellular protection

SOCS3 inhibits insulin signaling in porcine primary adipocytes

Atypical protein kinase C zeta: potential player in cell survival and cell migration of ovarian cancer

Keratin 17 Is induced in oral cancer and facilitates tumor growth

Phospholipase C-beta1 potentiates glucose-stimulated insulin secretion

A network-based bioinformatics approach to identify molecular biomarkers for type 2 diabetes that are linked to the progression of neurological diseases

Metabolic programming of MEST DNA methylation by intrauterine exposure to gestational diabetes mellitus

Epigenetic dysregulation of the IGF system in placenta of newborns exposed to maternal impaired glucose tolerance

Gestational diabetes alters offspring DNA methylation profiles in human and rat: identification of key pathways involved in endocrine system disorders insulin signaling, diabetes signaling, and ILK signaling

Vitamin B12 and folate levels during pregnancy and risk of gestational diabetes mellitus: a systematic review and meta-analysis

Frequency of gestational diabetes mellitus at collaborating centers based on IADPSG consensus panel-recommended criteria the hyperglycemia and adverse pregnancy outcome (HAPO) study

Gestational diabetes mellitus

Epigenetics in cancer

Global delay in nascent strand DNA methylation

BSMAP: whole genome bisulfite sequence MAPping program

• fast, convenient online submission • thorough peer review by experienced researchers in your field • rapid publication on acceptance • support for research data, including large and complex data types • gold Open Access which fosters wider collaboration and increased citations maximum visibility for your research: over 100M website views per year submit your research ? Choose BMC

A Bayesian hierarchical model to detect differentially methylated loci from single nucleotide resolution sequencing data

ChIPseeker: an R/Bioconductor package for ChIP peak annotation, comparison and visualization

STAR: ultrafast universal RNA-seq aligner

Near-optimal probabilistic RNA-seq quantification

Differential analyses for RNA-seq: transcript-level estimates improve gene-level inferences

Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2

Clusterprofiler: an R package for comparing biological themes among gene clusters

Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations

We thank all the participants of this research study and the staff helping with this study. We also thank the referees for making important suggestions for analysis of our data.

Author contributions WH contributed to study design; SL, NK, WH, JX, and YW contributed to data collection; JW, SH, HG, and JS involved in data analysis; WH and JS involved in study supervision; JW, SL, and HG contributed to manuscript writing. All authors have approved the submitted version of the manuscript. All authors agree to be personally accountable for the manuscript's content. All authors read and approved the final manuscript. 

The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

This study was performed in accordance with the Declaration of Helsinki and approved by the ethics committee of Hangzhou Women's Hospital. All participants provided written informed consent. Ethics Number: No. 201800901.

Not applicable.

The authors declare that they have no competing interests.