Epigenetic alternations of microRNAs and DNA methylation contribute to gestational diabetes mellitus

Abstract This study aimed to identify epigenetic alternations of microRNAs and DNA methylation for gestational diabetes mellitus (GDM) diagnosis and treatment using in silico approach. Data of mRNA and miRNA expression microarray (GSE103552 and GSE104297) and DNA methylation data set (GSE106099) were obtained from the GEO database. Differentially expressed genes (DEGs), differentially expressed miRNAs (DEMs) and differentially methylated genes (DMGs) were obtained by limma package. Functional and enrichment analyses were performed with the DAVID database. The protein‐protein interaction (PPI) network was constructed by STRING and visualized in Cytoscape. Simultaneously, a connectivity map (CMap) analysis was performed to screen potential therapeutic agents for GDM. In GDM, 184 low miRNA‐targeting up‐regulated genes and 234 high miRNA‐targeting down‐regulated genes as well as 364 hypomethylation–high‐expressed genes and 541 hypermethylation–low‐expressed genes were obtained. They were mainly enriched in terms of axon guidance, purine metabolism, focal adhesion and proteasome, respectively. In addition, 115 genes (67 up‐regulated and 48 down‐regulated) were regulated by both aberrant alternations of miRNAs and DNA methylation. Ten chemicals were identified as putative therapeutic agents for GDM and four hub genes (IGF1R, ATG7, DICER1 and RANBP2) were found in PPI and may be associated with GDM. Overall, this study identified a series of differentially expressed genes that are associated with epigenetic alternations of miRNA and DNA methylation in GDM. Ten chemicals and four hub genes may be further explored as potential drugs and targets for GDM diagnosis and treatment, respectively.

Functional and enrichment analyses were performed with the DAVID database. The protein-protein interaction (PPI) network was constructed by STRING and visualized in Cytoscape. Simultaneously, a connectivity map (CMap) analysis was performed to screen potential therapeutic agents for GDM. In GDM, 184 low miRNA-targeting upregulated genes and 234 high miRNA-targeting down-regulated genes as well as 364 hypomethylation-high-expressed genes and 541 hypermethylation-low-expressed genes were obtained. They were mainly enriched in terms of axon guidance, purine metabolism, focal adhesion and proteasome, respectively. In addition, 115 genes (67 up-regulated and 48 down-regulated) were regulated by both aberrant alternations of miRNAs and DNA methylation. Ten chemicals were identified as putative therapeutic agents for GDM and four hub genes (IGF1R, ATG7, DICER1 and RANBP2) were found in PPI and may be associated with GDM. Overall, this study identified a series of differentially expressed genes that are associated with epigenetic alternations of miRNA and DNA methylation in GDM. Ten chemicals and four hub genes

| INTRODUC TI ON
Gestational diabetes mellitus (GDM) is defined by 'the type of glucose intolerance that develops in the second and third trimester of pregnancy, resulting in hyperglycaemia of variable severity'. 1 GDM affects 2%-5% of pregnancies worldwide, with the significantly increased prevalence over the last decade. 2 Effects of GDM involve dysfunction of glucose metabolism and an increase of glucose accessibility by the foetus, which cause adverse pregnancy results or negatively impact the health of both mother and foetus. Moreover, GDM is related to adverse consequences not only during foetal development of pregnancy but also later in life. 3,4 Many literatures indicated genetic, epigenetic or environmental factors can contribute to GDM. 5 Moreover, increased insulin resistance stimulated by a genetic predisposition to impairment of pancreatic islet β-cell function may also contribute to GDM during pregnancy. 6 In the human umbilical vein and placental microvascular endothelia, GDM may also associate with foetoplacental vascular dysfunction characterized by increased NO synthesis and L-arginine transport activity. 7 While knowledge concerning the detailed processes governing the initiation and progression of GDM is still unknown and remains a key obstacle on the road of GDM treatment, the development of robust and accurate biomarkers will greatly facilitate the early detection and identification of biological features of GDM. Therefore, more potential biomarkers and some chemicals for GDM are urgent to be identified.
Currently, many studies suggested that epigenetic modification can play an important role in the pathogenesis of some multifactorial diseases, including GDM. 8,9 The epigenetic modification includes histone modification, DNA methylation and miRNA gene silencing.
All of them are closely linked with each other and influence protein synthesis patterns. 10 Disturbance of this complementary system may lead to dysfunction.
MicroRNAs are highly conserved small single-stranded non-coding RNA (18 ~ 25 nucleotides; nt) post-transcriptionally modulating gene expression in a variety of disease-related signalling processes and pathways. 11 They have emerged as promising diagnostic and therapeutic tools due to their association with GDM. For instance, miR-137 displays high expressions, whereas its target gene FNDC5 (fibronectin type III domain containing 5) is down expressed between women with GDM and normal in the placenta tissues, which can restrict the viability and migration of trophoblast cells. 12 Among the numerous aberrantly expressed miRNAs discovered in patients with GDM, miRNA-340 down-regulated PAIP1 (poly(A) binding protein-interacting protein 1), which involved in glucose and insulin regulation. 13 It has also been found that a low level of miR-21-3p in blood leucocytes of women may increase the risk of GDM. 14 In addition, associations of miR-21-3p with GDM were present only among women carrying male foetuses. 15 Furthermore, 32 different types of miRNA have been identified, differentially expressed in GDM women plasma compared to non-GDM women. Most of aberrantly expressed miRNAs are associated with glucose and insulin metabolism by affecting disrupting the MAPK (mitogen-activated protein kinase) signalling pathway or IRS genes. 16 Ding R, et al found miR-138-5p differentially expressed in GDM mouse placentas and significantly inhibited the proliferation and migration of trophoblasts (HTR-8/SVneo) by targeting the 3'-UTR of TBL1X (transducin beta-like 1X-linked) gene. 17 Tang X-W, et al reported miR-335-5p suppressed pancreatic islet β-cell secretion and enhanced insulin resistance by inhibiting VASH1 (vasohibin 1) in GDM mice, eventually activating the TGF-β pathway. 18 Aberrant DNA methylation, another common and best-studied epigenetic modification, makes a critical contribution in the regulation of genomic imprinting, genome stabilization gene expression and chromatin modification, which engage in placental development. 19 Numerous studies have revealed that aberrant DNA methylation has been reported to be engaged in adverse pregnancy outcomes, including GDM. 3 A study revealed that 56 pregnant women with GDM displayed a significant increase of global placental DNA methylation compared with 974 controls, independent of risk factors such as maternal age, BMI and recurrent miscarriages. 20 During the stages of embryogenesis, LEP (leptin) methylation profile was increased in placentas and umbilical cord blood in GDM. 5 Interestingly, DNA methylation of the LEP gene was evaluated in foetal tissues and associated with glycaemia after two hours of an oral glucose tolerance test (OGTT). 21 On the contrary, some genes with lower DNA methylation levels also take part in the GDM process. For instance, lower DNA methylation levels in the promoter of ADIPOQ on the foetal side of the placenta were correlated with higher maternal glucose levels and ADIPOQ is suspected to have insulin-sensitizing proprieties. In this study, data of mRNA expression profiling microarrays (GSE10 3552), miRNAs expression microarrays (GSE10 4297) and DNA methylation microarrays (GSE10 6099) were systematically analysed in order to identify core genes and pathways which contribute to GDM via epigenetic regulation.

| Microarray data
In this present study, data of mRNA and miRNA expression profiling microarray (GSE10 3552 and GSE10 4297) and gene methyla-

| Data process
Raw gene expression profiles data were pre-processed by R and Bioconductor packages. After background correction, logarithm transformation and normalization were conducted, differentially expressed genes (DEGs), differentially expressed miRNAs (DEMs) and differentially methylation probes (DMPs) were screened by using the limma package in R (version 3.6.0). 23 After converting from probe level to gene level using the IlluminaHumanv4.db (version 1.26.0), DEGs, DEMs and DMPs were screened with P < 0.05 and | t |> 2 as the cut-off criteria. DMPs located in the gene region were assigned to the corresponding genes, which were defined as differentially methylation genes (DMGs). To find out epigenetic alternations in GDM, jvenn online software (http://jvenn.toulo use.inra.fr/app/ examp le.html) was adopted to identify overlapping genes from the DEG in GES103552 and DMGs in GSE10 6099. Subsequently, DEGs and potential targets of DEMs were overlapped to obtained low miRNA-targeting up-regulated genes and high miRNA-targeting down-regulated genes. Besides, hypomethylation-high-expression genes and hypermethylation-low-expression genes were summarized and obtained through overlapping aberrant methylated and expressed genes.

| Prediction of miRNA Target genes and Construction of miRNA-mRNA Network
The target genes of DEMs were predicted by using miRWalk online database (http://zmf.umm.uni-heide lberg.de/apps/zmf/mirwa lk/index.html), which includes five different databases to predict miRNA target genes (miRanda, miRDB, Targetscan, RNA22 and miR-Walk). Predicted genes which were fitted at least 3 databases were considered as the target gene of DEMs. After aligning DEMs and DEGs, we used the Cytoscape tool (v 3.7.1) to visualize the entire miRNA-mRNA network.

| Protein-protein Interaction (PPI) Network Construction and Module Analysis
The Search Tool for the Retrieval of Interacting Genes (STRING) online tool was used to construct a PPI network of hypomethylationhigh expression genes and hypermethylation-low-expression genes, respectively. After PPI was visualized (combined score >0.7), modules within the PPI network was obtained by the Molecular Complex Detection (MCODE) in Cytoscape software. MCODE score >3 and the number of nodes >4. The functional enrichment analysis of the genes in individual modules was achieved by DAVID with a significance threshold of P<0.05. Hub genes were screened with connection degree >10.

| Drug discovery in CMap
CMap (Connectivity Map) database (https://www.broad insti tute.org) is an open database that can be used to identify connections among small molecules which sharing a mechanism of chemicals, physiological processes and action, and then predict potential drugs in silicon. 24 CMap analysis is used to predict potential small molecular compounds which can induce or reverse the altered expression of DEGs in cell lines. The link between the chemicals and query genes was measured via a connectivity score ranged from −1 to 1 and P < 0.05.

| Molecular docking analysis
Molecular docking can more intuitively show and predict the interaction between compounds and target proteins encoded by miRNAassociated aberrant methylation DEGs via Autodock (version 4.2.6).
Protein crystal structures were downloaded from PDB (Protein Data Bank, https://www.wwpdb.org) and chemical structures were obtained from zinc15 online database (http://zinc.docki ng.org). First, the protein crystal structures were imported into Atutodock tools.
Following the removal of irrelevant water molecules and ions, the addition of polar hydrogen atoms and assignment AD4 type, the proteins were prepared for docking. Compounds in the mol2 format were then imported into the software and protein-ligand docking was run using a genetic algorithm with an optimized genetic algorithm. 25

| The microarray data information and identification of aberrantly methylated-differentially expressed genes in GDM
The characteristics of the studies based on the GEO dataset are presented in Table S1. In GSE10 3552, a total of 5212 DEGs were identi-  Red indicates that the expression of genes is relatively up-regulated or the level of methylation is hypermethylated, blue indicated that the expression of genes is relatively down-regulated or the level of methylation is hypomethylated As for genes methylation microarray, 7387 hypomethylated CpG sites located within 4553 genes and 10010 hypermethylated CpG sites located within 5607 genes were found in GSE10 6099. All differentially methylated CpG sites from each autosomal chromosome are shown in the circus plot ( Figure 1A). The distribution of differentially methylated CpG sites in six different genomic subregions is shown in Figure 1B. Additionally, differentially methylation genes (DMGs) are evenly distributed on autosomes ( Figure 1C).  Figure 1D and 1E, which suggested that they can be distinguished between GDM and normal.

| Down-regulated genes and highexpression miRNA
Enrichment analysis of high miRNA-targeting down-regulated genes suggested that 45 GO terms were identified with the thresholds of Pvalue <0.05, which were mainly associated with cell division and protein binding (Table S3, Figure 3C). The most enriched KEGG analysis terms were Purine metabolism, Metabolic pathways, Biosynthesis of antibiotics, Steroid biosynthesis and Mismatch repair ( Figure 3E).
From the miRNA-mRNA network, a total of 14 genes were regulated by two miRNAs and TNFRSF9, ST8SIA4, IPMK and RSBN1L were regulated by three miRNAs ( Figure 3D). No overlapping has-miR-191-3p significant target genes were predicted in miRWalk.

| High-Expression and Hypomethylation Genes
Total 364 hypomethylation-high-expression genes were enriched in 116 GO enrichment terms with the thresholds of P-value <0.05, such as signal transduction and cell migration (Table S4, Figure 4A).

KEGG pathway analysis indicated enrichment of Pathways in cancer,
Focal adhesion and ECM-receptor interaction ( Figure 4A). In total, 532 nodes and 1153 edges are shown in the PPI network ( Figure 4B). The top ten genes ranked by degree were identified as hub genes, included ABL1, COL4A2, FBN1, LAMB2, LTBP1, POLR2A, RHOT2, SDC2, TGOLN2 and VHL (Table S5). In these 10 hub genes, POLR2A is presented with the highest degree (degree = 31). Moreover, in order to study the important modules found in the PPI network, the top two significant modules were selected by MCODE plug-in with 8.00 and 6.27 scores ( Figure 4E), respectively, and the functional annotation of the genes involved in the modules was analysed using PANTHER classification system (http://geneo ntolo gy.org). GO analysis showed that module 1 and module 2 were mainly related to TGFbeta in the extracellular matrix and protein polyubiquitination, respectively (data was not shown). Furthermore, Reactome pathway analysis enriched module 1 and module 2 genes in pathways of post-translational protein phosphorylation and neddylation.

| Low-expression and Hypermethylation Genes
For hypermethylation-low-expression genes, 117 GO terms were recognized with the thresholds of P-value <0.05 in the DAVID database (Table S4, Figure 4C). KEGG pathway analysis identified enriched pathways of Proteasome, Porphyrin and chlorophyll metabolism, Spliceosome and DNA replication. In total, 359 nodes and 248 edges were shown in the PPI network ( Figure 4D).
Moreover, a total of 359 nodes were analysed using the plug-in MCODE. The top 2 significant modules with 20.00 and 13.35 scores were selected ( Figure 4F). Enrichment analysis of Modules 1 and 2 indicated that hypermethylation targeting down-regulated genes participate in the termination of RNA polymerase II transcription and proteasomal ubiquitin-independent protein catabolic process, respectively. Reactome pathway enrichment analysis revealed that modules 1 and 2 were significantly enriched in pathways, including mRNA Splicing and Regulation of ornithine decarboxylase.

| DEGs associated with Both DNA Methylation and Aberrant miRNA
Interestingly, several DEGs were regulated by both aberrant alternations of DNA methylation and miRNA, which might demonstrate more vital and valuable function underlying GDM. Sixty-seven genes such as ABCE1, ANKRD46, ANP32E and ATG7 were up-regulated under the modulation of both hypomethylation and decreased miRNA ( Figure 5A). Simultaneously, forty-eight genes, including SPATS2, SERPINE1, TACC1, ADD1 and NEK6 were down-regulated under modulation of both hypermethylation and increased miRNA ( Figure 5B). The DNA methylation site and its relation to CpG island, as well as the specific regulatory miRNA and binding site, are summarized in Table S6. Moreover, the results of functional enrichment analysis of up-or down-regulated genes in GDM are shown in Figure 5C (Table 1) 29 In the TA B L E 1 Ten chemicals were predicted as putative therapeutic agents for GDM  Since there are no effective drugs for gestational diabetes mellitus, the online database was used to aid the prediction of some drugs.
At present, CMap is a tool of practical value for exploring new drugs and reusing existing drugs, and its effectiveness has been confirmed by many studies. 37 From the CMap database, 10 chemicals, includ- Ranbp2 is a vital, large, mosaic, pleiotropic nucleoporin and localized at the cytoplasmic peripheral side of the nuclear pore complex, which commands proteostasis of selective substrates and the nuclear-cytoplasmic trafficking in a cell-type dependent manner. 43,44 In cancer cells, IGF-1R first binds to the dynactin subunit p150 (glued), which transports the receptor to the nuclear pore complex where it co-localizes with importin-β and then binds to Ranbp2. 45 However, Ranbp2 has never been linked directly to GDM, with only some preliminary explorations in the embryo. 46 ATG7, as an important autophagy-related protein, is involved in activating ubiquitin-like protein (UBL), which is essential for the formation of autophagosome in the recognized pathway. 47 Although the role of autophagy on GDM was controversial, autophagy was activated in GDM placentas 48 and the birth weight of foetuses was significantly decreased in labyrinth layer-specific ATG7 knockout mouse models. 49 There are some deficiencies in our survey. Due to data availability, this study did not analyse the association of clinical data such as clinical parameters and prognosis with epigenetic changes.
In addition, the effects of abnormal methylation and expression In addition, the effects of abnormal methylation and expression of miRNAs on gene expression were not validated in experiments.
Thus, further evaluations in clinical trials are required to validate these genes.

| CON CLUS ION
This study indicated a series of aberrantly methylated-differ-

ACK N OWLED G EM ENT
We are grateful to all the members for their generous participation.

CO N FLI C T O F I NTE R E S T
The authors declare that they have no conflict of interests.

AUTH O R S' CO NTR I B UTI O N
ZWQ and SYP were responsible for the statistical analysis and wrote the manuscript. LJW, ZZF, LM, FXP, ZQX and XJH contributed to review and revision of the manuscript. LSQ, ZWJ, ZMH and DJ were responsible for the study design. All authors read and approved the final manuscript.

E TH I C S A PPROVA L A N D CO N S E NT TO PA RTI CI PATE
Not applicable.

CO N S E NT FO R PU B LI C ATI O N
Not applicable.

DATA AVA I L A B I L I T Y S TAT E M E N T
The authors declare that the data supporting the findings of this study are available within the article.