Integrated miRNA-mRNA Expression Profile Analysis Reveals a Role of miR-146a-3p/TRAF6 in Plasma from Gestational Diabetes Mellitus Patients

BACKGROUND: By utilizing an integrative strategy, we screened noninvasive molecular markers for the early detection of GDM and constructed miRNA-mRNA regulatory networks associated with GDM. METHODS: A total of 3 microarray datasets (GSE98043, GSE19649 and GSE92772) of plasma samples comparing GDM pregnant women and healthy control pregnant women were downloaded from the GEO database. The GEO2R online platform was used to identify the differentially expressed genes (DEmRNAs) and the differentially expressed miRNAs (DEmiRNAs). The target genes of DEmiRNAs were identified using two independent and complementary types of information: computational target predictions and expression relationships. KEGG pathway annotation was performed for target genes and DEmiRNAs. An interaction network was constructed to identify hub genes of GDM. Another dataset (GSE92772) was used to externally verify the predictive ability of the hub genes. Gene set enrichment analysis (GSEA) was performed to explore the biological function of hub genes. RESULTS: A total of 264 DEmiRNAs and 1217 DEmRNAs were identified by comparing the microarray data of serum samples of GDM patients and healthy pregnant women. Hsa-miR-146a-3p ranked first because of its lowest P value and was selected for further analysis. A total of 47 target genes, including TRAF6, were shared between the computational target predictions and DEmRNAs and were identified as target genes of hsa-miR-146a-3p. Enrichment analysis indicated that GDM-related miRNAs were mainly enriched in the glypican pathway, proteoglycan syndecan-mediated signaling events, and syndecan-1-mediated signaling events. In addition, the glypican pathway was also one of the pathways regulated by hsa-miR-3 target genes of hsa-miR-146a-3p were identified using two independent and complementary types of information: computational target predictions and expression relationships. For computational target predictions, a total of 1366 mRNAs were predicted to be the target genes of hsa-miR-146a-3p in the mirDIP online platform. These 1366 target genes were further mapped to the 1217 DEmRNAs screened from the GSE19649 dataset. Finally, a total of 47 target genes were shared between the computational target predictions and DEmRNAs and were

146a-3p. The interaction network analysis of DEmRNAs identified TRAF6, CASP8, PTPN6, and CHD3 as hub genes involved in the pathophysiological process of GDM.
Next, these 4 hub genes were selected for independent external validation using the GSE19649 dataset. The expression of TRAF6, CASP8 and CHD3 in 8 pairs of GDM blood samples was confirmed to be higher than that in healthy pregnant women blood samples. However, the expression of PTPN6 in the blood samples of GDM patients was similar to that of healthy pregnant women blood samples. The AUC for predicting GDM was 0.813, 0.828, 0.813, and 0.703 for CHD3, PTPN6, and CASP8, respectively. GSEA showed that 9 hallmark gene sets of metabolism processes were enriched in TRAF6 function.
CONCLUSIONS: Three hub genes, TRAF6, CASP8, and CHD3, were identified and independently externally validated as potential GDM noninvasive serum markers.
Enrichment analysis indicated that GDM-related miRNAs were mainly enriched in the glypican pathway, proteoglycan syndecan-mediated signaling events, and syndecan-1-mediated signaling events. In addition, integrated miRNA-mRNA expression profile analysis showed that miR-146a-3p/TRAF6 might play a central role in the pathogenesis of gestational diabetes mellitus by involving the above pathways.

| BACKGROUND
In China, the incidence of gestational diabetes mellitus (GDM) has increased from approximately 5% to more than 16% since the implementation of diagnostic criteria for GDM in December 2011 [1]. Compared with healthy pregnant women, pregnant women with GDM have a higher incidence of pre-eclampsia and a cumulative incidence of type 2 diabetes within 5-16 years after delivery [2]. In addition, the risk of macrosomia, abortion, hypoglycemia, cardiac dysplasia, obesity, and metabolic syndrome in infants of mothers with GDM is also increased [3]. Due to the high prevalence of GDM and its adverse outcomes, it is necessary to explore noninvasive markers for the early detection of GDM.
In addition to common risk factors, some new environmental factors, such as passive smoking [4] and excessive fruit intake in mid-pregnancy[5], significantly influence the predisposition to GDM. Genomic studies have also revealed that GDM is a multigenic disease. At present, many genes, including protein-coding RNA and noncoding RNA, are considered to be involved in the pathogenesis of GDM. Serum CysC is significantly and independently associated with insulin resistance, which is helpful in identifying the risk of GDM in Chinese pregnant women [9]. In addition, the combined detection of plasma microRNA-16-5p, -17-5p and − 20a-5p is expected to provide early diagnosis for pregnant women with GDM[10]. However, these studies are only based on a single gene. To date, studies providing systematic insight into the composition and expression of the serum mRNA and miRNA transcriptome in GDM are still lacking. Therefore, in the present study, we utilized an integrative strategy to construct functional miRNA-mRNA regulatory networks by combining the reverse expression relationships between miRNAs and targets and computational predictions. Our study may provide clues to a better understanding of the mechanism of GDM and help in screening noninvasive molecular markers for the early detection of GDM.

| Screening of molecular markers
The GEO2R online platform was used to identify the differentially expressed genes (DEmRNAs) in the GSE19649 dataset and the differentially expressed miRNAs (DEmiRNAs) in the GSE98043 dataset between serum samples of GDM pregnant women and healthy control pregnant women. GEO2R is an interactive web tool that allows gene expression analysis of microarray datasets using the limma R package [12]. P < 0.05 and fold change > 1.5 or < 0.5 were set as the cutoff criteria.

| Functional enrichment analysis
The target genes of microRNAs were predicted by the mirDIP online platform [13], and the target genes with high scores (top 5%) were selected for further analysis.
The computationally predicted target genes and DEmRNAs were represented in a Venn diagram to identify real target genes. KEGG pathway annotation was performed for the real target genes and DEmiRNAs, respectively. All pathway analyses were performed using Funrich 3.1.3 software [14].

| External hub gene verification
Another dataset (GSE92772) was used as a testing set to verify the differential expression of hub genes between GDM women and healthy controls. Continuous variables were compared using paired t tests. A receiver operating characteristic (ROC) curve analysis was generated, and the area under the curve (AUC) was calculated to evaluate the ability of the hub gene to predict GDM.

| Gene set enrichment analysis (GSEA)
To explore the biological function of hub genes, GSEA was performed in patients with high and low hub gene expression levels in the GSE92772 dataset. The high and low expression groups were separated by the median gene expression.
Annotated gene sets c2.cp.kegg.v5.2.symbols. The GMT pathway database was used as the reference. P < 0.05 and |enrichment score (ES) | > 0.25 were set as the cutoff values.
A schematic overview of the proposed approach in this study is shown in Fig. 1.

| RESULTS
By comparing the microarray data of serum samples of GDM patients and those of healthy pregnant women, 264 DEmiRNAs were identified in the GSE98043 dataset, of which 137 genes were upregulated and 127 genes were downregulated. In addition, a total of 1217 DEmRNAs were identified in the GSE19649 dataset. The top 10 DEmiRNAs and DEmRNAs are listed in Tables 1-2. Hsa-miR-146a-3p ranked first because of its lowest P value and was selected for further analysis. The target genes of hsa-miR-146a-3p were identified using two independent and complementary types of information: computational target predictions and expression relationships. For computational target predictions, a total of 1366 mRNAs were predicted to be the target genes of hsa-miR-146a-3p in the mirDIP online platform. These 1366 target genes were further mapped to the 1217 DEmRNAs screened from the GSE19649 dataset. Finally, a total of 47 target genes were shared between the computational target predictions and DEmRNAs and were identified as real target genes of hsa-miR-146a-3p (Fig. 2).
Enrichment analysis of the signaling pathways of DEmiRNAs indicated that GDMrelated miRNAs were mainly enriched in the glypican pathway (P < 0.001), proteoglycan syndecan-mediated signaling events (P < 0.001), and syndecan-1mediated signaling events (P < 0.001) (Fig. 3, Table 3). In addition, analysis of the signaling pathways of the real target genes of hsa-miR-146a-3p found that the glypican pathway was also one of the pathways regulated by hsa-miR-146a-3p (P < 0.001). In addition, the DEmRNAs TRAF6, MEF2A, TNFAIP3, SKIL, FLT1 and IRS2 were also targeted by hsa-miR-146a-3p in the present study and were enriched in the glypican pathway. To explore the hub DEmRNAs involved in GDM, we further performed an interaction network analysis of DEmRNAs screened from the GSE19649 dataset (Fig. 4). The results indicated that TRAF6, with a total of 16 connecting nodes, was the top hub gene involved in the pathophysiological process of GDM. Other hub genes included CASP8 (number of connecting nodes: 6), PTPN6 (number of connecting nodes: 6), and CHD3 (number of connecting nodes: 6). Next, these 4 hub genes were selected for independent external validation using the GSE19649 dataset. The expression of TRAF6, CASP8 and CHD3 in 8 pairs of GDM blood samples was confirmed to be higher than that in healthy pregnant women blood samples. However, the expression of PTPN6 in the blood samples of GDM patients is similar to that of healthy pregnant women blood samples (Fig. 5). The AUC for predicting GDM was 0.813, 0.828, 0.813, and 0.703 for CHD3, PTPN6, TRAF6, and CASP8, respectively ( Fig. 6, Table 4). Furthermore, TRAF6 was involved in the following three core pathways of GDM: glypican pathway, proteoglycan syndecan-mediated signaling events, and syndecan-1-mediated signaling events. To determine the potential mechanism of TRAF6 in GDM, GSEA was conducted (Table 5). A total of 9 hallmark gene sets of metabolism processes were enriched, including cysteine and methionine metabolism, pyruvate metabolism, sphingolipid metabolism, beta alanine metabolism, propanoate metabolism, galactose metabolism, starch and sucrose metabolism, butanoate metabolism, and arginine and proline metabolism processes. noninvasive markers of GDM in early and mid-term pregnancy based on miRNAomics. The expression of several plasma miRNAs including miRNA-29a and miRNA-16 was found to be associated with GDM, and both were also differentially expressed in our present study (miRNA-29a, P = 0.003; miRNA-16, P = 0.008) [15,17].
The results of the present study also provide insights into the molecular mechanisms that underlie GDM. miRNAs are involved in many biological processes through posttranscriptional regulation [18]. They regulate gene expression and transcription by complementary binding to the 3' end of target genes [18]. However, to date, little is known about the regulatory process between miRNAs and mRNAs in GDM. Thus, we performed an integrative analysis to construct functional miRNAgene regulatory networks by combining the expression relationships between miRNAs and targets and computational predictions. First, enrichment analysis of the involved signaling pathways was carried out based on DEmiRNAs. The metabolic pathways glypican pathway, proteoglycan syndecan-mediated signaling events, and syndecan-1-mediated signaling events, were the core pathways regulated by miRNAs in GDM. Glypicans belong to the family of membrane-bound heparan sulfate proteoglycans that are linked to the cell surface and are involved in the regulation of growth factor activity [19]. Proteoglycan families are involved in a wide range of diseases, including obesity and diabetes [20]. Impaired glucose metabolism in GDM may alter the expression of proteoglycans, which may impair the biological functions of the placenta [21]. A previous study demonstrated that the expression of a member of the proteoglycan family, endocan, was increased in the human placenta from obese women with GDM and in response to proinflammatory stimuli [20]. In addition, the metabolic effect of high glucose and high osmotic pressure of GDM patients may contribute to the increased proteoglycan perlecan expression in diabetic placentas [22]. Leelalertlauw et al. [23] found that serum glypican 4 levels were increased with increasing degrees of obesity. Furthermore, Ussar et al. [24] demonstrated that circulating glypican-4 levels correlate with BMI and insulin sensitivity in humans. Moreover, glypican-4 interacts with the insulin receptor, enhances insulin receptor signaling, and increases adipocyte differentiation. These results suggest that abnormal expression of proteoglycan family genes in GDM is helpful for understanding the pathogenesis of GDM.
In the present study, hsa-miR-146a-3p, ranked first in the altered DEmiRNAs of GDM because of its lowest P value, was downregulated with a log (fold change) of -12.76 in the blood samples of GDM patients. Changes in miR-146a in diabetic patients and animals were reported in a previous study. Liu et al. [25] demonstrated that treatment of diabetic mice with miR-146a mimics robustly reduced diabetic peripheral neuropathy. Then, the target genes of hsa-miR-146a-3p were identified based on the parallel DEmRNAs in the present study. The results indicated that TRAF6 was not only identified as a target of hsa-miR-146a-3p but also a hub gene involved in the pathophysiological process of GDM. TRAF6 has been demonstrated to be a target gene of miR-146a at the experimental level in previous immune studies [26]. Regarding diabetes, Kamali et al found that miR-146a antagomir significantly increased TRAF6 mRNA levels in human umbilical vein endothelial cells under hyperglycemic conditions. In addition, previous research revealed that culture of aortic endothelial cells in high glucose medium significantly increased TRAF6 expression in a time-dependent manner. In addition, TRAF6 mediated high glucoseinduced endothelial dysfunction via NF-KB and AP-1-dependent signaling [27]. To date, the role of miR-146a and its target gene TRAF6 in the pathogenesis of GDM has never been explored. In the present study, GSA revealed that several metabolic processes were enriched in the mechanism of TRAF6 in GDM. Moreover, TRAF6 was involved in all three core pathways regarding glypican metabolism in GDM based on pathway analysis of DEmiRNAs. Furthermore, we confirmed the altered expression of TRAF6 in 8 pairs of GDM blood samples to be higher than that in healthy pregnant blood samples in another independent patient group. Thus, serum TRAF6 may serve as a novel biomarker for GDM.
Other core genes, such as CASP8 and CHD3, were confirmed to be higher in 8 pairs of GDM blood samples than in healthy pregnant women blood samples. At present, there is no research on the relationship between CASP8 or CHD3 and GDM. Wu et al. [27] found that maternal type 2 diabetes mellitus causes heart defects in the developing embryo manifested by excessive apoptosis in heart cells, with increasing apoptosis markers of cleaved caspase 3 and 8. Whether CASP8 can be used as a marker of fetal heart defects in pregnant women with GDM needs to be further studied.
There are limitations to this study. The serum GDM markers screened from the GSE19649 and GSE98043 datasets were independently externally validated using the GSE19649 dataset. The sample size was relatively small. However, the clinical information, such as age and BMI, of GDM and healthy pregnant women were matched in the GSE19649 dataset, which improved the test effectiveness. In addition, other GDM noninvasive markers explored in this study still need to be further verified by expanding the sample size in the future.
In summary, a series of potential GDM noninvasive serum markers were identified by integrated miRNA-mRNA expression profile analysis. Enrichment analysis indicated that GDM-related miRNAs were mainly enriched in the glypican pathway, proteoglycan syndecan-mediated signaling events, and syndecan-1-mediated signaling events. In addition, miR-146a-3p/TRAF6 might play a central role in the pathogenesis of GDM through the above pathways.

5.
Huang W-Q, Lu Y, Xu M, Huang J, Su Y-X, Zhang C-X: Excessive fruit consumption during the second trimester is associated with increased likelihood of gestational diabetes mellitus: a prospective study.
Scientific reports 2017, 7:43620.  Figure 1 Schematic overview of the proposed approach in this study.

Figure 1
Schematic overview of the proposed approach in this study.

Figure 2
Analysis of the real target genes of hsa-miR-146a-3p. The computationally predicted target g Figure 2 Analysis of the real target genes of hsa-miR-146a-3p. The computationally predicted target g Figure 3 Enrichment analysis of signaling pathways involved in DEmiRNAs.

Figure 3
Enrichment analysis of signaling pathways involved in DEmiRNAs.

Figure 4
Interaction network analysis of DEmRNAs.

Figure 4
Interaction network analysis of DEmRNAs. The receiver operating characteristic (ROC) curve analysis validated that all four hub genes w Figure 6 The receiver operating characteristic (ROC) curve analysis validated that all four hub genes w