Skip to main content

Gene model-related m6A expression levels predict the risk of preeclampsia

Abstract

Background

This is the first study to explore the potential functions and expression patterns of RNA N6-methyladenosine (m6A) and potential related genes in preeclampsia.

Methods

We identified two m6A modification patterns through unsupervised cluster analysis and validated them by principal component analysis. We quantified the relative abundance of specific infiltrating immunocytes using single-sample gene set enrichment analysis (ssGSEA) and the Wilcoxon test. To screen hub genes related to m6A regulators, we performed weighted gene coexpression network analysis. Functional enrichment analysis was conducted for differential signalling pathways and cellular processes. Preeclampsia patients were grouped by consensus clustering based on differentially expressed hub genes and the relationship between different gene-mediated classifications and clinical features.

Results

Two m6A clusters in preeclampsia, cluster A and cluster B, were determined based on the expression of 17 m6A modification regulators; ssGSEA revealed seven significantly different immune cell subtypes between the two clusters. A total of 1393 DEGs and nine potential m6A-modified hub genes were screened. We divided the patients into two groups based on the expression of these nine genes. We found that almost all the patients in m6A cluster A were classified into hub gene cluster 1 and that a lower gestational age may be associated with more m6A-associated events.

Conclusions

This study revealed that hub gene-mediated classification is consistent with m6A modification clusters for predicting the clinical characteristics of patients with preeclampsia. Our results provide new insights into the molecular mechanisms of preeclampsia.

Peer Review reports

Background

Preeclampsia (PE) is one of the most serious hypertensive disorders of pregnancy, affecting 5–7% of all pregnant women and causing over 70,000 maternal deaths and 500,000 foetal deaths worldwide every year [1, 2]. PE presents as new-onset vascular dysfunction and may exacerbate stroke, kidney failure, pulmonary edema, liver rupture, and eclampsia during the third trimester. The exact aetiology of PE is not clear, although this is an interesting area of research. The widely held theory is that conditions such as impaired trophoblast invasion, incomplete spiral artery remodelling, endothelial dysfunction, oxidative stress and inappropriate immune responses may affect the expressed PE phenotype [3]. Nevertheless, its origin cannot be fully explained by our current understanding. In recent years, an increasing number of studies have focused on how epigenetic modification, including DNA methylation, noncoding RNA and, more marginally, histone posttranslational modification, impacts gene expression in the initiation and development of PE [4, 5]. Therefore, it is critical to explore the specific molecular mechanisms underlying the evolution of PE and identify new prognostic factors and potential therapeutic targets.

Canonical epigenetic modification research has focused on the DNA in chromatin, whereas with the deepening of research, epigenetic RNA modifications have recently begun to gain increasing attention. To date, over 160 kinds of RNA modifications have been confirmed, including N6-methyladenosine (m6A), 5-methylcytosine (m5C) and N1-methyladenosine (m1A), of which m6A modification is regarded as the most universal posttranscriptional modification in all eukaryotes, especially in determining the splicing, stability, export, decay and translation of mRNAs. Furthermore, similar to DNA and histone modification, m6A modification is dynamic and reversible in controlling gene expression; additionally, it is associated with three groups of enzymes, i.e., the “erasers” (demethylases) FTO and ALKBH5, the “readers” (binding proteins) YTHDFs and IGF2BPs, and the “writers” (RNA methyltransferases) METTL3, METTL14 and WTAP [6]. Emerging studies have shown that m6A modification and its regulators are involved in not only the pathogenesis and development of multiple diseases, including cancer [7], neurological diseases [8], infertility [9], obesity [10] and immunological diseases [11], but also the regulation of placental development and pathologies.

A previous study demonstrated that the level of global DNA methylation increased with gestational weeks, which means that gestational age may be responsible for interplacental variation and a programmed, reproducible change in DNA methylation [12]. A recent study indicated that compared to control placentas, low-birth-weight placentas were associated with significantly higher m6A expression in key genes related to lipid metabolism and angiogenesis by analysis of m6A modification levels in 44 placental samples from eight sows [10]. The study by Taniguchi revealed that m6A at the 5’-untranslated region (UTR) and near the stop codon in placental mRNA may play vital roles in foetal growth and PE [13]. Although the above discoveries provide new insights into the epigenetic mechanisms of PE, the general expression characteristics of m6A regulators in PE remain to be elucidated.

Herein, to systematically analyse the modification patterns of m6A regulators and related genes in PE, we carried out this study based on the public Gene Expression Omnibus (GEO) database. In addition, immune microenvironment infiltration characteristics were analysed. We found that m6A modification patterns were associated with some clinical features of PE and infiltration of multiple immune cells. We also revealed that m6A phenotype-related hub genes could well reflect the m6A level and be used to identify patients at high risk for PE.

Methods

Data acquisition

Gene expression data and the clinical data were obtained from the GEO database (https://www.ncbi.nlm.nih.gov/geo/). GSE75010 (n = 80) and GSE60438 (n = 60) were gathered in the present study for further analysis. To normalize the expression data from different datasets, the R package “sva” was used to remove batch effects [14]. A total of 21 regulators related to m6A modification were included with reference to previous studies in the follow-up analyses. These 21 m6A regulators included eight writers (METTL3, RBM15, METTL14, RBM15B, CBLL1, WTAP, KIAA1429, ZC3H13), eleven readers (YTHDF1, YTHDF2, YTHDF3, YTHDC1, YTHDC2, IGF2BP1, FMR1, HNRNPA2B1, ELAVL1, HNRNPC, LRPPRC) and two erasers (ALKBH5, FTO). Among them, 17 m6A methylase regulators were supported by expression profile data.

Sample classification

To determine distinct m6A modification patterns in PE, we conducted unsupervised consensus clustering analysis based on the 17 m6A regulators by using the R package Consensus Cluster Plus, and the consensus clustering algorithm ran 1,000 repetitions to guarantee the stability of classification. Principal component analysis (PCA) was conducted to further verify the expression patterns of the 17 m6A regulators in different modification patterns.

Correlation between m6A regulators and immune characteristics

To quantify the relative abundance of specific infiltrating immunocytes in each sample, a single-sample gene-set enrichment analysis (ssGSEA) algorithm was used. The Wilcoxon test was used to compare the enrichment scores among the groups.

Identification and functional analysis of m6A-mediated genes

To identify m6A regulator-mediated genes, the empirical Bayesian approach of the ‘limma’ R package was applied to determine differentially expressed genes (DEGs) between two distinct m6A subtypes. The false discovery rate (FDR) was applied for multiple testing correction of raw P values through the Benjamini–Hochberg method. A |log2| fold change (FC)|> 1 and an FDR < 0.01 were set as the thresholds for identifying DEGs.

The R software package weighted gene coexpression network analysis (WGCNA) was used to identify the gene modules that were significantly correlated with different m6A modifications. The analysis setting included the following main steps: (1) construct the similarity matrix; (2) select the weighting coefficient, β, and transform the similarity matrix into an adjacency matrix; (3) transform the adjacency matrix into a topological overlap matrix (TOM); and (4) perform hierarchical clustering to identify modules. Finally, we selected the minimum module eigengenes (MEs) with > 25 genes and a cutting height of 0.25. Relevant modules of high importance for clinical traits were identified (Additional file 1: Table S1).

In addition, the correlations between modules and m6A modification patterns were identified by mean absolute values of gene significance (GS) through Pearson correlation. The larger the mean absolute value was, the more strongly the module was related to the m6A modification patterns.

To further study the function of these m6A-related module genes, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis and Gene Ontology (GO) annotation analysis were performed using the R package ClusterProfiler through GO (http://www.geneontology.org) and KEGG database (www.kegg.jp/kegg/kegg1.html). P < 0.05 was chosen as the cut-off value.

Identification of the m6A-related hub genes for classifying PE patients

We defined the genes with a GS over 0.5 and a module membership (MM) over 0.75 in the green gene module as hub genes. To classify PE RNA sequence results into different gene-mediated classifications, we used k-means unsupervised clustering analyses based on different expression patterns of nine hub genes. Through methods of fast t-distributed stochastic neighbour embedding (t-SNE) and dimension reduction, PE samples were clearly classified into different subtypes.

Result

Characteristics of datasets and patients

A total of 140 PE RNA sequence results were obtained from two different database (GSE75010 and GSE60438). The median age of the patients was 32 years. Sixty-two patients had a gestational age of more than 32 weeks.

We reviewed the literature and collated a list of 21 m6A RNA methylation regulators. There were 17 m6A regulators expressed in both datasets. These 17 m6A regulators included METTL3, RBM15, RBM15B, WTAP, KIAA1429, CBLL1, ALKBH5, FTO, YTHDC1, YTHDC2, YTHDF1, YTHDF2, YTHDF3, IGF2BP1, FMR1, LRPPRC and ELAVL1. Their chromosomal locations are shown in Fig. 1a.

Fig. 1
figure 1

Bioinformatics analysis of the expression and genetic variation of m6A regulators in PE. a Overview of m6A gene locus and gene information. b Consensus clustering identified two subgroups according to the expression of m6A regulators. c Heatmap of the expression of 19 m6A regulators in two distinct m6A clusters. PE: preeclampsia

m6A molecular subtypes

According to the expression similarity of m6A RNA methylation regulators, k = 2 seemed to be an adequate selection with clustering stability increasing from k = 2 to 6 in the cohort (Additional file 2: Fig. S1, Additional file 3: Fig. S2, Additional file 4: Fig. S3). Based on unsupervised clustering, we identified two distinct modification patterns, including 65 cases in m6A cluster A and 75 cases in m6A cluster B (Fig. 1b). Most m6A RNA methylation regulator expression showed clear distinctions and significant differences in the two cluster subgroups (Fig. 1c). PE samples could be completely distinguished into two clusters based on PCA. Furthermore, the clinical features between the two subtypes were compared. Due to the limited clinical data, maternal age with a cut-off value of 32 was the only significantly different feature (P = 0.011) (Fig. 2a, Additional file 5: Fig. S4).

Fig. 2
figure 2

m6A methylation modification patterns are correlated with the clinical characteristics and immune cell infiltration. a, b clinicopathological features between the two clusters were compared. c Heatmap of immune cell infiltration in two distinct m6A clusters

Immune status evaluation among m6A subtypes

After ssGSEA, 28 kinds of immune cell subtypes were identified, and seven of them were significantly different between the two clusters according to the Wilcoxon test (Fig. 2b). In m6A cluster A, monocytes, activated dendritic cells, plasmacytoid dendritic cells, type 17 T helper cells, type 1 T helper cells, mast cells and neutrophils were upregulated.

Heterogeneity of other biological processes among subtypes

We identified 1393 DEGs with a fold change |logFC|> 0.15 and an FDR < 0.05. After WGCNA, we further analysed 135 samples (five were outliers) and eventually identified 9 modules. From the heatmap of module-trait correlations, the green module, including 87 genes, was most highly related to the m6A cluster (Additional file 6: Fig. S5, Fig. 3a–d). To further explore module biological function, GO enrichment analyses were conducted, and the results are shown in Fig. 4a–c. The biological functions were mainly enriched in covalent chromatin modification, endomembrane system organization and protein localization to the Golgi apparatus. The results of KEGG analysis are shown in Fig. 4d.

Fig. 3
figure 3

Biological characteristics of key m6A module. a Dendrogram of all DEGs obtained by average linkage hierarchical clustering. b Heatmap of the Correlation between module eigengenes and m6A modification patterns. Each cell contains the correlation coefficient and P-value. c Bar plots of mean GS across modules. We selected the ME green module for subsequent analysis. d Scatter plot of module eigengenes related to m6A modification patterns in the green module. DEGs: differentially expressed genes; GS: gene significance; ME: minimum module eigengene

Fig. 4
figure 4

ac Genes from green module in the biological process, molecular function, and cellular component. d Enrichment plot conducted via KEGG analysis. KEGG: Kyoto Encyclopedia of Genes and Genomes

Hub gene screening and molecular subtypes

We further selected nine hub genes from the green module by setting MM > 0.7 and GS > 0.5. These hub genes included DHX30, VPS4A, RNF26, ACTN4, GTF3C1, MORC2, MYO1C, CLPTM1 and GATAD2B. We found that these nine hub genes were significantly different among the m6A clusters (Fig. 5a; P < 0.001). Furthermore, we explored the relationship between the expression of the nine hub genes and gestational age and found that the expression levels of RNF26, MORC2 and MYO1C were significantly associated with a low gestational age (Fig. 5b; P < 0.05).

Fig. 5
figure 5

a Expression of nine hub genes in two distinct m6A clusters. b Expression of nine hub genes at different gestational age of patients with preeclampsia. The asterisks represented the statistical P-value (ns: P > 0.05; *P < 0.05; **P < 0.01; ***P < 0.001; ****P < 0.0001)

According to the results obtained, we used the K-means unsupervised clustering method to classify the PE samples (Additional file 6: Figure S5). As shown in Fig. 6a, the PE samples could be well divided into two clusters based on m6A hub gene expression. We found that the expression levels of the nine genes were significantly different between the two clusters, with lower gene expression in gene cluster 1 and higher gene expression in gene cluster 2 (P < 0.05, Fig. 6b, c), which could be regarded as signatures for different PE classification. For a better understanding of the association of m6A clusters, gestational age and hub gene clusters, we conducted a correlation analysis and found that almost all patients in m6A cluster A were classified into hub gene cluster 1. The expression of all nine m6A-related hub genes was significantly higher in gene cluster 1 than in cluster 2 (Fig. 6d). However, no significant difference in gestational age was observed between m6A cluster A and cluster B.

Fig. 6
figure 6

Different PE subgroups identified by unsupervised cluster analysis based on m6A hub gene expression. a Consensus clustering identified two gene clusters according to the expression of m6A regulators. b, c Heat map and violin diagram of levels of nine hub genes in cluster 1 and cluster 2; d Sankey diagram showing the relationship of gestational age, m6A clusters, and hubgene clusters. PE: preeclampsia

There are limited clinical data in GSE60438, while GSE75010 has more details. The common clinical features of these three sets include gestational weeks, infant sex and maternal age. We further analysed the differences in clinical features between the two PE groups. We found that when using 32 gestational weeks as a cut-off, cluster 1 had more patients with at least 32 gestational weeks (P = 0.086). There were no differences in infant sex or maternal age. Because GSE75010 has more details, we analysed the clinical characteristics of the 1-min Apgar score, the 5-min Apgar score, diagnosis of HELLP, diagnosis of intrauterine growth restriction (IUGR), maternal BMI, previous history of hypertensive pregnancy, miscarriage and nulliparity. We found that patients in cluster 1 were more likely to present with IUGR (P = 0.026).

Discussion

PE is a unique hypertensive disorder in human pregnancy, and placental trophoblast dysfunction plays significant roles in the pathogenesis of this pregnancy disorder. The exact pathophysiology that causes PE remains unclear; however, genetic, immunological, endocrine, and environmental factors have been implicated in its pathogenesis [15]. The placenta plays an essential role in the development of PE. While the mechanism of abnormal placentation is controversial, animal models have demonstrated that uteroplacental ischaemia drives the hypertensive, multiorgan failure response observed in maternal preeclamptic syndrome [16]. Previous studies found increased m6A expression and increased m6A RNA methylation in placental trophoblasts in PE patients, suggesting that aberrant m6A RNA methylation may contribute to trophoblast dysfunction in this pregnancy disorder. Further determining whether increased m6A RNA methylation is associated with altered clinical features and identifying different subtypes according to m6A methylation levels will be helpful for understanding this disease.

In this study, we first analysed the expression of 21 m6A factors and identified the correlations between them. According to the expression levels of m6A, 140 samples were well divided into two groups. We then investigated the clinical characteristics of such groups and found that when using 32 as the cut-off for maternal age, cluster A was enriched in the high maternal age group (P = 0.011). Cluster A was also related to a high rate of IUGR (P = 0.055). The number of reports regarding the profiling of m6A modification and its potential role in the placenta of PE is small. Several studies have investigated the clinical features in PE patients and healthy controls and found no difference in maternal age in the whole group [17, 18]. Currently, no study has analysed the correlation between maternal age and m6A levels in PE patients, and our results need to be further explored.

Given the indispensable function of RNA m6A modification in various bioprocesses, it is reasonable to speculate that deregulation of m6A modification may also be associated with PE. Recent studies showed that METTL3 and METTL14 were upregulated in PE and verified that the significant increase in HSPA1A mRNA and protein expression was regulated by m6A modification, suggesting that m6A plays a key role in the regulation of gene expression [18]. Increased METTL3 expression and m6A RNA methylation were also observed by Taniguchi et al., and they suggested that aberrant m6A modification may contribute to trophoblast dysfunction in PE [17]. Furthermore, m6A at the 5’-UTR and nearby stop codon in placental mRNA may play important roles in foetal growth and PE through MeRIP-Seq conducted on human placentas obtained from mothers of infants of various birth weights [13]. For example, in the present study, infants in the high m6A expression group were more likely to have IUGR, indicating that m6A may affect foetal growth.

To identify the differences between the two m6A clusters, we analysed the infiltrating immune cells and found that several types of immune cells were correlated with m6A levels. PE is associated with chronic immune activation that leads to an increased production of inflammatory cytokines by proinflammatory T cells and a decrease in regulatory and anti-inflammatory cytokines, which further promotes an inflammatory state during PE [19, 20]. Macrophages, natural killer (NK) cells, dendritic cells (DCs), T cells, and T regulatory cells (Tregs) are present in the decidua and are required for the normal invasion of trophoblasts during placentation [21]. The Th1 and Th17 subclasses are responsible for promoting inflammation during PE, while Tregs and Th2 cells are decreased and are therefore unable to properly control the inflammation associated with increased inflammatory T-cell populations [22]. It has been proven that m6A methylation of mRNA controls T-cell homeostasis [23]. In the present study, Th17 and Th1 cells were enriched in m6A cluster A, which indicated that low levels of m6A promoted the activation of Th17 and Th1 cells.

Oxidative stress is regarded as an important immune event during PE and contributes to the development of inflammation in the vasculature, which includes the production of inflammatory IL-17 [24]. Another role of IL-17 that contributes to hypertension during pregnancy is stimulating B-cell production of agonistic autoantibody to the Ang II, type 1 receptor [25, 26]. Moreover, increased inflammatory cytokines, such as TNF-α, IL-6 and IL-17, are observed in the circulation of RUPP rats, while regulatory cytokines, IL-10 and IL-4, are decreased [27]. IL-10 and Tregs also improve the pathophysiology of PE through Th1 regulation and by exhibiting anti-inflammatory effects in general. m6A modification has been demonstrated to be a major posttranscriptional regulator of immune responses in cells [23]. The immune environment in PE is complex, and m6A regulation may contribute to the development of PE. Further evaluation of the correlation between m6A and infiltrating immune cells in PE is important.

We identified 1393 DEGs between the two m6A clusters, and then we performed WGCNA. According to the results, 87 hub genes were selected for further analysis. GO and KEGG analyses revealed that these DEGs were significantly enriched in different pathways and functions, including vasopressin-regulated water reabsorption, the AMPK signalling pathway, covalent chromatin modification, positive regulation of intracellular transport, and organelle localization. AMPK activation is required for placental differentiation and vasodilation of uterine artery blood flow [28]. Lack of AMPK induces malplacentation, which results in angiogenic imbalance. The increase in serum AMPK in severely preeclamptic women suggests a compensatory mechanism for the angiogenic imbalance [28]. AMPK activators ameliorate preeclamptic symptoms, which indicates that AMPK is a potential therapeutic target for PE. As described before, Treg and Th17 cells respond to the development of PE, and AMPK activation restores the normal balance between Treg and Th17 cells and cures such an imbalance [29]. Future studies are recommended to verify whether AMPK signalling may provide a prospective therapeutic target for the prevention and treatment of PE.

α-Actinins (ACTNs) are cytoskeletal proteins that maintain cytoskeleton integrity and control cell movement. Several studies have indicated that ACTN4 may participate in endothelial cell regulation. Moreover, recent studies have shown that ACTN4 is involved in cell apoptosis, which plays an important role in the progression of PE [30, 31]. Zhao et al. observed decreased expression of ACTN4 in severe PE endotheliocytes and demonstrated that dysregulated ACTN4 expression may be associated with PE due to its effects on endothelial cell apoptosis via the p38-MAPK/p53 apoptosis pathway [32]. GATAD2B expression was significantly decreased in pregnant mouse and human myometrium during labour. Chen et al. suggested that GATAD2B serves as an important mediator of P4–PR suppression of proinflammatory and contractile genes during pregnancy [32]. Other genes identified in this study may also be related to PE progression and have not been well investigated in PE. Microrchidia family CW-type zinc finger 1 (MORC1) is a highly conserved nuclear protein that is increasingly recognized as an epigenetic regulator [33]. MORC2 is a newly identified chromatin remodelling enzyme with an emerging role in the DNA damage response, but the underlying mechanism remains largely unknown [34]. RNF26 belongs to the RING domain family of proteins, which plays a key role in organizing the endosomal pathway for efficient cargo transport by mediating the ubiquitination of SQSTM1 [34]. Moreover, RNF26 has been reported to limit the type I interferon response by promoting the autophagic degradation of IRF3 [35]. All three genes were differentially expressed between the two gene groups and different maternal age groups in the current study. The correlation between PE and the nine hub genes, especially the potential methylation, needs to be further revealed.

Previous studies have illustrated gestational age as one of the important factors that influence methylation, and some studies have reported that women undergoing spontaneous preterm delivery may possess differences in methylation compared with those with an additional week of gestation for a total of at least 37 weeks [12, 36]. Furthermore, the augmented m6A levels at the 5’-UTR in mRNAs of small-for-date placental samples were dominant compared to the reduction in m6A levels, whereas a decrease in m6A in the vicinity of stop codons was common in heavy-for-date placental samples [13]. There is a lack of studies investigating whether there are some differences in the expression level of m6A at different PE gestational ages. According to our results, gestational age did not affect the m6A level, but when compared in different gene groups, the patients with high gene expression were enriched in gestational weeks less than 32 weeks. This indicated that gestational age is an important confounder when studying variation in placental DNA methylation. Fewer gestational weeks may be associated with more m6A-associated events.

There are some limitations in this study. First, due to the limited clinical information and the lack of foetal outcome data in GSE60438, the differences between two m6A clusters as well as gene clusters could not be well identified. Furthermore, studies focusing on m6A methylation in PE are limited, and some genes revealed in this study are newly recognized and have not been well researched in PE. It is difficult to explain the mechanism of such genes identified in this study.

Conclusions

This is the first study to evaluate the m6A level and potential related genes in PE patients. Our results revealed that a high m6A level is associated with lower maternal age and IUGR occurrence. Several immune cells may contribute to different m6A clusters and clinical features. Furthermore, nine hub genes that mediated classification could well reflect the m6A level. The high hub gene expression group is correlated with a high IUGR occurrence rate and gestational weeks. Our results indicated that the m6A level plays an essential role in PE development and outcome for both the mother and foetus. Currently, there is no gene risk model for PE patients. The nine hub gene model established in this study reflected that different m6A levels can be used to well identify patients at high risk for PE.

Availability of data and materials

Microarray datasets (GSE75010 and GSE60438) for this study are openly available in Gene Expression Omnibus database at https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE75010 and https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE60438, respectively.

Abbreviations

PE:

Preeclampsia

m6A:

N6-methyladenosine

m5C:

5-Methylcytosine

m1A:

N1-methyladenosine

UTR:

Untranslated region

GEO:

Gene expression omnibus

PCA:

Principal component analysis

ssGSEA:

Single-sample gene-set enrichment analysis

DEGs:

Differentially expressed genes

FDR:

False discovery rate

FC:

Fold change

WGCNA:

Weighted gene coexpression network analysis

TOM:

Topological overlap matrix

MEs:

Minimum module eigengenes

GS:

Gene significance

KEGG:

Kyoto encyclopedia of genes and genomes

GO:

Gene ontology

MM:

Module membership

t-SNE:

T-distributed stochastic neighbour embedding

IUGR:

Intrauterine growth restriction

NK:

Natural killer

DCs:

Dendritic cells

Tregs:

T regulatory cells

ACTNs:

α-Actinins

MORC1:

Microrchidia family CW-type zinc finger 1

References

  1. Hogan MC, Foreman KJ, Naghavi M, Ahn SY, Wang M, Makela SM, et al. Maternal mortality for 181 countries, 1980–2008: a systematic analysis of progress towards Millennium development goal 5. Lancet. 2010;375(9726):1609–23.

    Article  Google Scholar 

  2. Wanderer JP, Leffert LR, Mhyre JM, Kuklina EV, Callaghan WM, Bateman BT. Epidemiology of obstetric-related ICU admissions in Maryland: 1999–2008*. Crit Care Med. 2013;41(8):1844–52.

    Article  Google Scholar 

  3. Phipps EA, Thadhani R, Benzing T, Karumanchi SA. Pre-eclampsia: pathogenesis, novel diagnostics and therapies. Nat Rev Nephrol. 2019;15(5):275–89.

    Article  Google Scholar 

  4. Bianco-Miotto T, Mayne BT, Buckberry S, Breen J, Rodriguez Lopez CM, Roberts CT. Recent progress towards understanding the role of DNA methylation in human placental development. Reproduction. 2016;152(1):R23-30.

    CAS  Article  Google Scholar 

  5. Apicella C, Ruano CSM, Méhats C, Miralles F, Vaiman D. The role of epigenetics in placental development and the etiology of preeclampsia. Int J Mol Sci. 2019;20(11):2837.

    CAS  Article  Google Scholar 

  6. Zaccara S, Ries RJ, Jaffrey SR. Reading, writing and erasing mRNA methylation. Nat Rev Mol Cell Biol. 2019;20(10):608–24.

    CAS  Article  Google Scholar 

  7. Sun T, Wu R, Ming L. The role of m6A RNA methylation in cancer. Biomed Pharmacother. 2019;112:108613.

    CAS  Article  Google Scholar 

  8. Han M, Liu Z, Xu Y, Liu X, Wang D, Li F, et al. Abnormality of m6A mRNA methylation is involved in Alzheimer’s disease. Front Neurosci. 2020;14:98.

    Article  Google Scholar 

  9. Vanova I, Much C, Di Giacomo M, Azzi C, Morgan M, Moreira PN, et al. The RNA m6A reader YTHDF2 is essential for the post-transcriptional regulation of the maternal transcriptome and oocyte competence. Mol Cell. 2017;67(6):1059–67.

    Article  Google Scholar 

  10. Song T, Lu J, Deng Z, Xu T, Yang Y, Wei H, et al. Maternal obesity aggravates the abnormality of porcine placenta by increasing N6-methyladenosine. Int J Obes (Lond). 2018;42(10):1812–20.

    CAS  Article  Google Scholar 

  11. Zhang X, Zhang S, Yan X, Shan Y, Liu L, Zhou J, et al. m6A regulator-mediated RNA methylation modification patterns are involved in immune microenvironment regulation of periodontitis. J Cell Mol Med. 2021;25(7):3634–45.

    CAS  Article  Google Scholar 

  12. Novakovic B, Saffery R. The ever growing complexity of placental epigenetics–role in adverse pregnancy outcomes and fetal programming. Placenta. 2012;33(12):959–70.

    CAS  Article  Google Scholar 

  13. Taniguchi K, Kawai T, Kitawaki J, Tomikawa J, Nakabayashi K, Okamura K, et al. Epitranscriptomic profiling in human placenta: N6-methyladenosine modification at the 5ʹ-untranslated region is related to fetal growth and preeclampsia. FASEB J. 2020;34(1):494–512.

    CAS  Article  Google Scholar 

  14. Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. 2012;28(6):882–3.

    CAS  Article  Google Scholar 

  15. ACOG Practice Bulletin No. 202 summary: gestational hypertension and preeclampsia. Obstet Gynecol. 2019;133(1):1.

    Google Scholar 

  16. Palei AC, Spradley FT, Warrington JP, George EM, Granger JP. Pathophysiology of hypertension in pre-eclampsia: a lesson in integrative physiology. Acta Physiol (Oxf). 2013;208(3):224–33.

    CAS  Article  Google Scholar 

  17. Gu Y, Chu X, Morgan JA, Lewis DF, Wang Y. Upregulation of METTL3 expression and m6A RNA methylation in placental trophoblasts in preeclampsia. Placenta. 2021;103:43–9.

    CAS  Article  Google Scholar 

  18. Wang J, Gao F, Zhao X, Cai Y, Jin H. Integrated analysis of the transcriptome-wide m6A methylome in preeclampsia and healthy control placentas. PeerJ. 2020;8:e9880.

    Article  Google Scholar 

  19. Bennett WA, Lagoo-Deenadayalan S, Stopple JA, Barber WH, Hale E, Brackin MN, et al. Cytokine expression by first-trimester human chorionic villi. Am J Reprod Immunol. 1998;40(5):309–18.

    CAS  Article  Google Scholar 

  20. Formby B. Immunologic response in pregnancy. Its role in endocrine disorders of pregnancy and influence on the course of maternal autoimmune diseases. Endocrinol Metab Clin North Am. 1995;24(1):187–205.

    CAS  Article  Google Scholar 

  21. Nagamatsu T, Schust DJ. The immunomodulatory roles of macrophages at the maternal-fetal interface. Reprod Sci. 2010;17(3):209–18.

    CAS  Article  Google Scholar 

  22. Raghupathy R. Cytokines as key players in the pathophysiology of preeclampsia. Med Princ Pract. 2013;22(Suppl 1):8–19.

    Article  Google Scholar 

  23. Liu C, Yang Z, Li R, Wu Y, Chi M, Gao S, Sun X, Meng X, Wang B. Potential roles of N6-methyladenosine (m6A) in immune cells. J Transl Med. 2021;19(1):251.

    CAS  Article  Google Scholar 

  24. Yiyenoğlu ÖB, Uğur MG, Özcan HÇ, Can G, Öztürk E, Balat Ö, et al. Assessment of oxidative stress markers in recurrent pregnancy loss: a prospective study. Arch Gynecol Obstet. 2014;289(6):1337–40.

    Article  Google Scholar 

  25. Lockwood CJ, Yen CF, Basar M, Kayisli UA, Martel M, Buhimschi I, et al. Preeclampsia-related inflammatory cytokines regulate interleukin-6 expression in human decidual cells. Am J Pathol. 2008;172(6):1571–9.

    CAS  Article  Google Scholar 

  26. LaMarca B, Parrish MR, Wallace K. Agonistic autoantibodies to the angiotensin II type I receptor cause pathophysiologic characteristics of preeclampsia. Gend Med. 2012;9(3):139–46.

    Article  Google Scholar 

  27. Gadonski G, LaMarca BB, Sullivan E, Bennett W, Chandler D, Granger JP. Hypertension produced by reductions in uterine perfusion in the pregnant rat: role of interleukin 6. Hypertension. 2006;48(4):711–6.

    CAS  Article  Google Scholar 

  28. Kumagai A, Itakura A, Koya D, Kanasaki K. AMP-activated protein (AMPK) in pathophysiology of pregnancy complications. Int J Mol Sci. 2018;19(10):3076.

    Article  Google Scholar 

  29. Gualdoni GA, Mayer KA, Göschl L, Boucheron N, Ellmeier W, Zlabinger GJ. The AMP analog AICAR modulates the Treg/Th17 axis through enhancement of fatty acid oxidation. FASEB J. 2016;30(11):3800–9.

    CAS  Article  Google Scholar 

  30. Liu QY, Lei JX, LeBlanc J, Sodja C, Ly D, Charlebois C, et al. Regulation of DNaseY activity by actinin-alpha4 during apoptosis. Cell Death Differ. 2004;11(6):645–54.

    CAS  Article  Google Scholar 

  31. Lomert E, Turoverova L, Kriger D, Aksenov ND, Nikotina AD, Petukhov A, et al. Co-expression of RelA/p65 and ACTN4 induces apoptosis in non-small lung carcinoma cells. Cell Cycle. 2018;17(5):616–26.

    CAS  PubMed  PubMed Central  Google Scholar 

  32. Zhao J, Peng W, Ran Y, Ge H, Zhang C, Zou H, et al. Dysregulated expression of ACTN4 contributes to endothelial cell injury via the activation of the p38-MAPK/p53 apoptosis pathway in preeclampsia. J Physiol Biochem. 2019;75(4):475–87.

    CAS  Article  Google Scholar 

  33. Li DQ, Nair SS, Kumar R. The MORC family: new epigenetic regulators of transcription and DNA damage response. Epigenetics. 2013;8(7):685–93.

    CAS  Article  Google Scholar 

  34. Jongsma ML, Berlin I, Wijdeven RH, Janssen L, Janssen GM, Garstka MA, et al. An ER-associated pathway defines endosomal architecture for controlled cargo transport. Cell. 2016;166(1):152–66.

    CAS  Article  Google Scholar 

  35. Qin Y, Zhou MT, Hu MM, Hu YH, Zhang J, Guo L, et al. RNF26 temporally regulates virus-triggered type I interferon induction by two distinct mechanisms. PLoS Pathog. 2014;10(9):e1004358.

    Article  Google Scholar 

  36. Wilson SL, Leavey K, Cox BJ, Robinson WP. Mining DNA methylation alterations towards a classification of placental pathologies. Hum Mol Genet. 2018;27(1):135–46.

    CAS  Article  Google Scholar 

Download references

Acknowledgements

We would like to acknowledge AJE for its linguistic assistance during the preparation of this manuscript.

Funding

This work was supported by The construction Fund of Medical Key Disciplines of Hangzhou and Medical and Health Technology Project of Hangzhou (Z20220026).

Author information

Authors and Affiliations

Authors

Contributions

HW: Conceptualization, Funding acquisition, Project administration, Supervision, Writing—review and editing. ZY: Conceptualization, Data curation, Formal analysis, Investigation, Writing—review and editing. LY: Conceptualization, Data curation, Formal analysis, Methodology, Writing—original draft, Writing—review and editing. CC: Data curation, Methodology, Software, Writing—original draft, Writing—review and editing. DM: Data curation, Formal analysis, Visualization, Writing—review and editing. WY: Data curation, Formal analysis, Writing—original draft, Writing—review and editing. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Ying Zhu or Wei Hu.

Ethics declarations

Ethics approval and consent to participate

No ethical approval or consent was required as the data in the article was taken from a database and did not involve human/animal tissue.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

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

Supplementary Information

Additional file 1: Table S1.

The number of genes in each module.

Additional file 2: Fig. S1.

Heatmap of the matrix of co-occurrence proportions for preeclampsia samples.

Additional file 3: Fig. S2.

Consensus clustering cumulative distribution function for k = 2–6.

Additional file 4: Fig. S3.

Evaluation of the pairwise correlations among 17 m6A regulators’ expression.

Additional file 5: Fig. S4.

Difference in maternal age between two m6A clusters.

Additional file 6: Fig. S5.

Analysis of the scale-free ft index and the mean connectivity for various soft-thresholding powers.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Li, Y., Chen, C., Diao, M. et al. Gene model-related m6A expression levels predict the risk of preeclampsia. BMC Med Genomics 15, 103 (2022). https://doi.org/10.1186/s12920-022-01254-4

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12920-022-01254-4

Keywords

  • Preeclampsia
  • N6-methyladenosine
  • Hub gene
  • Prognosis
  • Bioinformatic gene analysis
  • Biomarkers