A newly defined risk signature, consisting of three m<sup>6</sup>A RNA methylation regulators, predicts the prognosis of ovarian cancer

Lili Fan; Ying Lin; Han Lei; Guang Shu; Liuer He; Zhipeng Yan; Hai Rihan; Gang Yin

doi:10.18632/aging.103811

Research Paper Volume 12, Issue 18 pp 18453—18475

A newly defined risk signature, consisting of three m⁶A RNA methylation regulators, predicts the prognosis of ovarian cancer

Lili Fan^{1,
,} , Ying Lin^{2,
,} , Han Lei^1, , Guang Shu^3, , Liuer He^3, , Zhipeng Yan^4, , Hai Rihan^3, , Gang Yin^1, ,

¹ Department of Pathology, Xiangya Hospital, School of Basic Medical Sciences, Central South University, Changsha, Hunan Province, China
² Department of Immunology, School of Basic Medical Sciences, Central South University, Changsha, Hunan Province, China
³ School of Basic Medical Sciences, Central South University, Changsha, Hunan Province, China
⁴ Hunan Cancer Hospital, the Affiliated Cancer Hospital of Xiangya School of Medicine, Central South University, Changsha, Hunan Province, China

* Equal contribution

Received: January 11, 2020 Accepted: July 20, 2020 Published: September 20, 2020

https://doi.org/10.18632/aging.103811
How to Cite

Copyright: © 2020 Fan et al. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Abstract

N6-methyladenosine (m⁶A) RNA methylation, involved in cancer initiation and progression, is dynamically regulated by the m⁶A RNA methylation regulators. However, the expression of m⁶A RNA methylation regulators in ovarian cancer and their correlation with prognosis remain elusive. Here, we demonstrated that the 18 central m⁶A RNA methylation regulators were expressed differently between ovarian cancer (OC) and normal tissues. By applying consensus clustering, all ovarian cancer patient cases can be divided into three subgroups (cluster1/2/3) based on overall expression levels of all 18 m⁶A RNA methylation regulators. We systematically analyzed the prognostic value of transcription levels of 18 m⁶A RNA methylation regulators in ovarian cancer and found that insulin-like growth factor 2 mRNA binding protein 1 (IGF2BP1), vir like m⁶A methyltransferase associated (VIRMA), and zinc finger CCCH-type containing 13 (ZC3H13) yield the highest scores for predicting the prognosis of ovarian cancer. Accordingly, we derived a risk signature consisting of transcription levels of these three selected m⁶A RNA methylation regulators as an independent prognostic marker for OC and validated our findings with data derived from a different ovarian cancer cohort. Moreover, by the Gene Set Enrichment Analysis (GSEA), we demonstrated that the three selected regulators were all correlated with pathways in cancer and WNT signaling pathways. In conclusion, m⁶A RNA methylation regulators are vital participants in ovarian cancer pathology; and IGF2BP1, VIRMA, and ZC3H13 mRNA levels are valuable factors for prognosis prediction and treatment strategy development.

Introduction

Ovarian cancer (OC) accounts for 2.5% of all malignancies among females, but for 5% of all cancer deaths due to its relatively high fatality rate, as about 80% of patients are diagnosed with the advanced disease [1]. At the time of diagnosis, most of the OC had metastasized to the uterus, bilateral appendage, omentum, and pelvic organs. With recent advances in surgery, chemotherapy, and novel immunotherapy, the overall survival of OC at every stage has been improved. However, there is still a lack of reliable prognostic indicators for OC.

Although the multiple layers of epigenetic regulation, such as modification of DNA and proteins, have been identified, the mechanism of RNA modification and its role in OC pathology remains unclear. Several common RNA modification forms have been documented including N6-methyladenosine (m⁶A), N1-methyladenosine (m¹A) and 5-methylcytosine (m⁵C) [2, 3]. As the most bountiful internal modification on mRNA, m⁶A occurs on adenosine and is enriched near stop codon and 3’ untranslated terminal region [4, 5]. Similar to DNA modification, the m⁶A RNA methylation is dynamically regulated by the m⁶A RNA methylation regulating elements, which are called “writers”, “erasers”, or “readers” depending on distinct functions. “Writers” introduce a methyl group on A of RRACH sequence (R = A or G, H = A, C or U) near the stop codon, 3’ untranslated region (UTR) and long internal exon through the methyltransferase complex (MTC) consisting of METTL3 as the core component and other related subunits including METTL14, WTAP, VIRMA, RBM15 and ZC3H13 [3, 6]. “Erasers” include FTO and ALKBH5, which mediate the demethylation reaction. These two layers of regulations make m⁶A a dynamic and reversible process. “Readers” are a group of RNA binding proteins that recognize the m⁶A methylation and perform corresponding functions. These proteins include YT521-B homology (YTH) domain-containing proteins, eukaryotic initiation factor 3 (eIF3), insulin-like growth factor 2 mRNA-binding proteins (IGF2BPs) and heterogeneous nuclear ribonucleoprotein (HNRNP) protein family members. Interestingly, recent studies have shown that the “writers” can also participate in post-methylation-regulation on target RNA and thereby partially function as “readers” [7, 8]. The m⁶A modification regulates a variety of critical steps in the RNA life cycle starting from transcription to degradation (such as transcription, splicing, exportation, translation, and degradation) and can influence cell process (such as cell cycle progression, cell proliferation, cell apoptosis, cell migration and invasion and cell differentiation) and physiological function (such as neural development, embryonic development and adipogenesis, etc.) through modulating the life cycle of target RNA [9–14].

Experimental evidence shows that m⁶A participates in cancer pathogenesis and development [15–17]. In OC, the expression of METTL3 is elevated and associated with poor patient survival. METTL3 enhances oncogene AXL expression, resulting in the promotion of epithelial to mesenchymal transition (EMT) [18, 19]. Moreover, IGF2BP1 stabilizes the mRNA of SRF to promote its expression, leading to enhanced expression of oncogene FOXK1 and PDLIM7 in tumor cells, and a more aggressive phenotype [20]. During the treatment of OC, elevated m⁶A level contributes to resistance to poly ADP-ribose polymerase inhibitors (PARPi) through upregulating the WNT signaling pathway via enhancing the stability of FZD10. Therefore, restraining the WNT signaling pathway in combination with PARPi represents a potential therapeutic strategy for OC [21]. Despite the accumulating data indicating vital roles of m⁶A in controlling physiological and pathological processes, our knowledge about the role of m⁶A in OC oncogenesis and prognosis is far from complete.

The multiomics-based comprehensive analysis provides much more informative results in evaluating the expression and function of genes. In this study, we systematically analyzed the expression of 18 central m⁶A RNA methylation regulators in 379 OC with RNA sequencing data from The Cancer Genome Atlas (TCGA) datasets and in 88 normal samples with RNA sequencing data from the Genotype-Tissue Expression (GTEx) datasets. We aimed to evaluate the power of m⁶A RNA methylation regulators in predicting the prognosis of OC patients, and explore possible signaling pathways regulated by m⁶A RNA regulators in OC through comprehensive bioinformatical analysis (Figure 1). Our results indicated that the expressions of three m⁶A RNA methylation regulators, IGF2BP1, VIRMA, and ZC3H13, have strong power in predicting the prognosis of OC. We also built a risk signature gene set including these three selected m⁶A RNA methylation regulators and validated that the expression of this signature is highly correlated with the bad prognosis of OC patients using data derived from a different cohort.

Workflow chart of data generation and analysis. The study mainly incorporated two sections: comprehensive bioinformatics analysis in 18 m6A RNA methylation regulators (including Spearman correlation analysis, protein-protein interaction analysis, consensus cluster analysis, cluster survival analysis and so on) and in the three selected m6A RNA methylation regulators (including CCLE database methylation analysis, protein-protein interaction network analysis, univariate Cox, multivariate Cox and so on).

Figure 1. Workflow chart of data generation and analysis. The study mainly incorporated two sections: comprehensive bioinformatics analysis in 18 m⁶A RNA methylation regulators (including Spearman correlation analysis, protein-protein interaction analysis, consensus cluster analysis, cluster survival analysis and so on) and in the three selected m⁶A RNA methylation regulators (including CCLE database methylation analysis, protein-protein interaction network analysis, univariate Cox, multivariate Cox and so on).

Results

The different expression of 18 m⁶A RNA methylation regulators in normal ovarian and OC tissues

Given the critical functions of m⁶A RNA methylation regulators in tumorigenesis and development, we comprehensively explored the transcription of the 18 m⁶A RNA methylation regulators using the TCGA dataset. The RNA levels of m⁶A RNA methylation regulators were presented as heatmaps and box line diagrams (Figure 2A and 2B), which showed that the expression levels of m⁶A RNA methylation regulators in OC patients were significantly different from those of the normal controls. Based on the expression pattern, m⁶A RNA methylation regulators can be divided into two groups. One group (including IGF2BP2, IGF2BP1, IGF2BP3, ZC3H13, ALKBH5, RBM15, YTHDF3, YTHDF2, ELF3, and YTHDF1) is highly expressed in tumors, while the other (including WTAP, HNRNPC, METTL3, YTHDC3, YTHDC2, YTHDC1, VIRMA, METTL14, and FTO) is enriched in normal tissues.

Expression of m6A RNA methylation regulators and interaction among them. (A) The expression levels of 18 m6A RNA methylation regulators in normal controls (n = 88) and OC (n = 379) with agglomerative hierarchical clustering. (B) Box line diagram of 18 m6A RNA methylation regulators. (C) Spearman correlation analysis of the 18 m6A modification regulators. (D) The m6A modification-related interactions among the 18 m6A RNA methylation regulators. (E) Number of related nodes of m6A RNA methylation regulators (only show the number > 5).

Figure 2. Expression of m⁶A RNA methylation regulators and interaction among them. (A) The expression levels of 18 m⁶A RNA methylation regulators in normal controls (n = 88) and OC (n = 379) with agglomerative hierarchical clustering. (B) Box line diagram of 18 m⁶A RNA methylation regulators. (C) Spearman correlation analysis of the 18 m⁶A modification regulators. (D) The m⁶A modification-related interactions among the 18 m⁶A RNA methylation regulators. (E) Number of related nodes of m⁶A RNA methylation regulators (only show the number > 5).

To better understand the interactions among the 18 m⁶A RNA methylation regulators, we also inspected the correlation (Figure 2C) and interaction (Figure 2D and 2E) among these regulators. The protein-protein interaction network analysis indicated that WTAP, RBM15, METTL3, METTL14, and VIRMA have more connections with other regulators, while IGF2BP2, IGF2BP1, IGF2BP3 and ELF3 have less interaction with others (Figure 2D and 2E). Among all the “writer” genes, each gene, except for ZC3H13, has the same number of related nodes of m⁶A RNA methylation regulators (Figure 2D and 2E). The expression of ZC3H13 was also positively correlated with the “writer” gene RBM15 in OC, and negatively correlated with WTAP and METTL3, but no correlation with METTL14 and VIRMA (Figure 2C). Interestingly, the least interacting m⁶A RNA methylation regulators (ELF3, IGF2BP2, IGF2BP1, and IGF2BP3) are all “readers” (Figure 2C). We also noticed that FTO was predicted to interact with ALKBH5 (Figure 2D), and their expression levels were negatively correlated with each other in OC (Figure 2C). In addition, we further analyzed the experimentally determined interactions between these 18 m⁶A RNA methylation regulators and other proteins (Supplementary Figure 1A). It was clear that among the 18 m⁶A RNA methylation regulators, ZC3H13 had the most interactions with other proteins, mainly interacted with the cyclin-dependent kinases (CDKs) family, RNA polymerase II (POLR2) family and mediator (MED) complex family. Altogether, we concluded that most, but not all, of the 18 m⁶A RNA methylation regulators are closely associated with each other. Among them, ZC3H13 interacted mostly with CDKs, POLR2 and MED. Moreover, we also found that the frequencies of genetic changes (mutation or copy number change) of these 18 m⁶A RNA methylation regulators were relatively high (the maximum was 27%) in OC tissue compared with normal (Supplementary Figure 2), which might explain the altered expression of these genes in OC.

Consensus clustering of m⁶A RNA methylation regulators identified three clusters of OC

To further investigate the association between the expression profile of m⁶A RNA methylation regulators and the prognosis of these OC cases, we then focused on grouping all 379 OC cases according to the expression of all the 18 m⁶A RNA methylation regulators in an unbiased way through consensus clustering. With clustering stability increasing from k = 2 to 10 in the TCGA datasets, k = 3 seemed to be an acceptable selection based on the expression similarity of m⁶A RNA methylation regulators (Supplementary Figure 3A, 3B and Figure 3A). We then applied principal component analysis (PCA) to compare the transcriptional profile among cluster1, 2 and 3 groups. However, the results didn’t show a clear separation among them (Supplementary Figure 3C). Moreover, we observed that there was little difference in the overall survival (OS) rate among cluster1, 2 and 3 (Figure 3B). We further compared the clinicopathological characters of these three subgroups and found little difference among them (Figure 3C).

Divergent clinicopathological features and OS of OC in the cluster1/2/3 subgroups. (A) Consensus clustering matrix for k = 3. (B) Kaplan-Meier OS curves for 379 OC patients. (C) Heatmap and clinicopathologic characters of the three clusters (cluster1/2/3) defined by the m6A RNA methylation regulators’ consensus expression.

Figure 3. Divergent clinicopathological features and OS of OC in the cluster1/2/3 subgroups. (A) Consensus clustering matrix for k = 3. (B) Kaplan-Meier OS curves for 379 OC patients. (C) Heatmap and clinicopathologic characters of the three clusters (cluster1/2/3) defined by the m⁶A RNA methylation regulators’ consensus expression.

Development of a risk signature consisting of three m⁶A RNA methylation regulators

To better predict the clinical outcomes of OC with abnormal expression of m⁶A RNA methylation regulators, we engaged the least absolute shrinkage and selection operator (LASSO) Cox regression algorithm to the 18 regulators in the TCGA dataset, and obtained risk score through R packages LASSO regression analysis (Figure 4A, 4B). With the median risk score (median risk score = 5) as the cut-off point, we divided all patients into two groups and scrutinized notable differences in OS between the two groups (P < 0.05; Figure 4C). At the same time, we also examined the effect of different grades of OC patients on the prognosis. However, no statistical significance was found (Supplementary Figure 4). The receiver operating characteristic (ROC) analysis was used for testing if the survival prediction is sensitive and specific based on the risk score. The calculation of the area under the curve (AUC) values was carried out according to ROC curves (Figure 4D). And the ROC curve shows that our risk model yields supporting results as the AUC = 0.58. Hence, we compared the clinicopathological characters (including tumor pathological grade and age) of these two categories clustered by risk score (Supplementary Figure 5). We found that patients in the high-risk cluster and the low-risk cluster did not differ significantly in tumor pathological grades and age. Then, the distribution of risk score and survival status were also analyzed (Figure 4E, 4F). In Figure 4E, the risk score of each patient was arranged from low to high. Patients were divided into the low risk-group (blue dot) and the high-risk group (red dot) with the median risk score. From Figure 4F, we can see that the number of deaths of patients with high-risk scores is slightly larger than that of patients with low-risk scores.

Risk signature with 18 m6A RNA methylation regulators. (A) LASSO regression analysis of the 18 m6A RNA methylation regulators. (B) Tenfold cross-validation for tuning the parameter selection in the LASSO regression. The solid vertical lines indicate the partial likelihood deviance with standard error. The dotted vertical lines represent the optimal values of the tuning parameter (λ) by minimum criteria. (C) Kaplan-Meier OS curves for patients in the TCGA datasets designated to high- and low-risk groups depended on the risk score. (D) ROC curves demonstrated the predictive efficiency of the risk signature in OC of TCGA datasets. (E–F) Risk score and survival status for each patient.

Figure 4. Risk signature with 18 m⁶A RNA methylation regulators. (A) LASSO regression analysis of the 18 m⁶A RNA methylation regulators. (B) Tenfold cross-validation for tuning the parameter selection in the LASSO regression. The solid vertical lines indicate the partial likelihood deviance with standard error. The dotted vertical lines represent the optimal values of the tuning parameter (λ) by minimum criteria. (C) Kaplan-Meier OS curves for patients in the TCGA datasets designated to high- and low-risk groups depended on the risk score. (D) ROC curves demonstrated the predictive efficiency of the risk signature in OC of TCGA datasets. (E–F) Risk score and survival status for each patient.

We further explored the prognostic effect of individual m⁶A RNA methylation regulator in OC. We executed a univariate Cox regression analysis on the expression level of each m⁶A RNA methylation regulator in the TCGA dataset (Figure 5A). The results demonstrated that three out of 18 tested regulators were significantly correlated with OS (P < 0.1; Figure 5A–5C). These three genes (IGF2BP1, VIRMA, and ZC3H13) were all risky genes with HR > 1. Hence, we compared the clinicopathological characters (including tumor pathological grade and age) of the three regulators. It was clear that the expressions of the three selected regulators were all high in most high-risk patients (Figure 5D). To implement a quantitative method for superior outcome prediction, we established a nomogram that assimilated the three selected regulators associated with OC prognosis (Figure 5E). In this nomogram, a higher total point indicates a worse survival. Moreover, univariate and multivariate analyses for OS were executed to appraise whether clinicopathological characters (including age, stage and risk score) were independent prognostic factors of patient outcomes. Univariate analysis applying the Cox proportional hazards model for all variables demonstrated that risk score (P < 0.001, 95%CI HR 1.10-1.43) and age (P = 0.005, 95%CI HR 1.01-1.03) were all independent poor prognostic factors for OC patients (Figure 5F). Multivariate analysis applying the same variables as in the univariate analysis in the cohort supported that risk score (P = 0.002, 95%CI HR 1.08-1.41) and age (P = 0.014, 95%CI HR 1.00-1.03) were independent poor prognostic factors for OC patients (Figure 5G), but little meaningful association of OS was discovered with grade. Because different types of OC have different characteristics, and the OC samples in the TCGA database are all ovarian serous cystadenocarcinoma (OV), so we analyzed the expression of IGF2BP1, VIRMA and ZC3H13 in different types of OC in the Oncomine database (https://www.oncomine.org) (Supplementary Figure 6A–6C). The expression of IGF2BP1 in ovarian serous adenocarcinoma was higher than that in normal tissues; among ovarian clear cell adenocarcinoma, ovarian endometrioid adenocarcinoma, ovarian mucinous adenocarcinoma and ovarian serous adenocarcinoma, the expression of VIRMA was the highest in ovarian serous adenocarcinoma, and ZC3H13 was the highest in ovarian mucinous adenocarcinoma. In addition, we analyzed the expression profile of IGF2BP1, VIRMA and ZC3H13 in OC patients with different ages and tumor grades in the TCGA database on the UALCAN website (http://ualcan.path.uab.edu) [22] (Supplementary Figure 6D–6I). The expression of IGF2BP1 and ZC3H13 increased with the age of patients, while the expression of VIRMA in different age groups had no similar rule. Since there is only one sample in grade 1 and grade 4, we pay attention to the difference between grade 2 and grade 3. The expression of VIRMA in grade 3 is higher than that in grade 2, and there is statistical significance, while there is no significant difference in the expression of IGF2BP1 and ZC3H13 between grade 2 and grade 3. Moreover, the three selected m⁶A RNA methylation regulators also exist in other tumors, at least in many cell lines (Supplementary Figure 7). The mutation sites of the three selected m⁶A RNA methylation regulators in different cell lines have been shown in Supplementary Figure 7. This indicates that the m⁶A RNA methylation regulators may also play a role in other tumors too and have a certain correlation with prognosis, which needs further analysis. In particular, we also further explored the experimentally determined interactions between these three selected m⁶A RNA methylation regulators and the other proteins (Supplementary Figure 1B). Of interest, it shows that ZC3H13 interacts with the CDKs family, POLR2 family and MED complex family closely, which is consistent with the result in Supplementary Figure 1A, but IGF2BP1 had no interaction with any other proteins.

The selection of three m6A RNA methylation regulators, and their effect on OC prognosis and clinicopathological characteristics. (A) Cox univariate regression analyses were used to examine the associations between expression of 18 m6A RNA methylation regulators and prognosis. (B, C) The P-value of IGF2BP1, VIRMA and ZC3H13 D) Heatmap and clinicopathologic features of the three selected m6A RNA methylation regulators. (E) Nomogram for forecasting 1-year, 2-year and 3-year survival of clinically OC patients. The nomogram is used by adding up the points identified on the points scale for each variable. Based on the sum of these points projected on the bottom scales, it is used to predict the likelihood of individual patients surviving for 1-year, 2-year and 3-year. (F) Univariate analysis of the hazard ratios for risk score as independent prognostic elements to anticipate the OS. (G) Multivariate analysis of the hazard ratios for risk score as independent prognostic elements to predict the OS.

Figure 5. The selection of three m⁶A RNA methylation regulators, and their effect on OC prognosis and clinicopathological characteristics. (A) Cox univariate regression analyses were used to examine the associations between expression of 18 m⁶A RNA methylation regulators and prognosis. (B, C) The P-value of IGF2BP1, VIRMA and ZC3H13 < 0.1. (D) Heatmap and clinicopathologic features of the three selected m⁶A RNA methylation regulators. (E) Nomogram for forecasting 1-year, 2-year and 3-year survival of clinically OC patients. The nomogram is used by adding up the points identified on the points scale for each variable. Based on the sum of these points projected on the bottom scales, it is used to predict the likelihood of individual patients surviving for 1-year, 2-year and 3-year. (F) Univariate analysis of the hazard ratios for risk score as independent prognostic elements to anticipate the OS. (G) Multivariate analysis of the hazard ratios for risk score as independent prognostic elements to predict the OS.

To further understand the biology of the three genes, we analyzed the protein expression of IGF2BP1, VIRMA and ZC3H13 in OC in the Clinical Proteomic Tumor Analysis Consortium (CPTAC) samples [18] (Supplementary Figure 8). The results suggest that the mRNA and protein levels of ZC3H13 are both high in cancer tissues compared with the normal control (Figure 2A and Supplementary Figure 8C). However, the protein expression of IGF2BP1 and VIRMA in the normal control group and tumor group was not consistent with the mRNA expression levels (Figure 2A and Supplementary Figure 8A–8B). This suggests that the mRNA may serve as a reservation. Only under certain circumstances, such as hypoxia or immune stimulation, can the protein be translated. A similar mode of action could be observed in the production of some cytokines. Besides, the mRNA may play a role in regulating the expression of other proteins by generating microRNAs. The specific details need further experiments and analysis.

Validation of the risk signature using data collected from a different OC cohort, and the exploration of signaling pathways that they involve

To figure out the prognostic importance of each gene of the signature composed of IGF2BP1, VIRMA, and ZC3H13, the OS of patients with a high expression level of any gene of the signature was compared to that of patients with low expression. We noticed that OC patients with high VIRMA expression had a shorter median OS than those with low expression (P < 0.05, Figure 6B). Unexpectedly, for IGF2BP1 and ZC3H13, the OS rate of patients didn’t associate with their expression levels (Figure 6A and 6C). However, through univariate analysis, the association between the risk signature genes developed in this study can be validated using data derived from a different cohort downloaded from the Gene Expression Omnibus (GEO) database (Figure 6D–6F). This result further confirmed the efficiency of the risk signature to predict prognosis developed in this study. Using the Gene Set Enrichment Analysis (GSEA), we found that IGF2BP1, VIRMA, and ZC3H13 are all associated with pathways in cancer and WNT signaling pathway (Figure 6G–6I, Supplementary Table 1).

The OS analysis, GSEA analysis, and validation in the GEO database of the three selected m6A RNA methylation regulators. (A–C) OS survival curve of OC patients based on the three selected m6A RNA methylation regulators levels. (D–F) The validation of the three selected m6A RNA methylation regulators using the GEO database through univariate analysis. (G–I) Enrichment of genes in the Kyoto Encyclopedia of Genes and Genomes (KEGG) different pathways by GSEA.

Figure 6. The OS analysis, GSEA analysis, and validation in the GEO database of the three selected m⁶A RNA methylation regulators. (A–C) OS survival curve of OC patients based on the three selected m⁶A RNA methylation regulators levels. (D–F) The validation of the three selected m⁶A RNA methylation regulators using the GEO database through univariate analysis. (G–I) Enrichment of genes in the Kyoto Encyclopedia of Genes and Genomes (KEGG) different pathways by GSEA.

Discussion

Previous studies have pointed out that up- or down-regulation of specific RNA m⁶A methylation regulators are associated with the oncogenesis of many different tumors. Additionally, the same m⁶A methylation regulators might have different functions in different tumors. At present, the main work has been devoted to the study of the mechanism of m⁶A in promoting tumorigenesis [23–28].

OC is one of the deadliest gynecological malignancies. Most patients have stage III~IV at the moment of diagnosis, and the prognosis is poor [29, 30]. Classical epigenetics, restricted to DNA or protein modification, plays critical roles in OC initiation, malignant progression, and prognosis [31–33]. In our study, we found that the expression of another layer of epigenetic regulators, namely m⁶A RNA methylation, is also firmly associated with the malignancy and prognosis of OC.

Here, we firstly analyzed the expression of m⁶A RNA methylation regulators in OC and normal tissues and the relationship between their expression. Then we applied the consensus clustering to divide all OC samples into three clusters and analyzed the expression of m⁶A RNA methylation regulators and different clinicopathological variables according to the clustering. However, the results of PCA did not show a clear distinction among cluster1, 2 and 3, and there was almost no significant difference in terms of prognosis and clinical case characteristics among the three groups. The possible reasons could be that the sample size is not big enough to reflect the difference, or the clustering algorithm is not sensitive enough for these data.

Next, we explored the prognostic value of each m⁶A RNA methylation regulators and developed a risk signature applying three chose m⁶A RNA methylation regulators, IGF2BP1, VIRMA and ZC3H13, which are selected by the Cox univariate analysis and LASSO Cox regression analysis. Based on this signature, we established a nomogram that assimilated the three selected regulators associated with OC prognosis and used univariate analysis and multivariate analysis to assess the prognostic value of the three m⁶A RNA methylation regulators. At last, the OS analysis, GSEA analysis were applied to data collected from a different OC cohort to validate the prognostic value of the three selected m⁶A RNA methylation regulators. This work provided a different biomarker other than the tumor stage for predicting the prognosis of OC.

For the three selected regulators, IGF2BP1, VIRMA, and ZC3H13, the GSEA study result indicated that they were also all correlated with pathways in cancer and WNT signaling pathways. Corresponding to the results in the GSEA, Felicite K. Noubissi et al. [34] indicated that IGF2BP1 plays a central role in carcinoma development. They also mentioned that IGF2BP1 was a direct target of the WNT/β-catenin signaling. Consistently, another analysis using TCGA data indicated that the expression of IGF2BP1 was negatively associated with survival in prostate cancer [35]. Compared with these regulators, there are relatively few studies on the signal pathway mediated by VIRMA and ZC3H13. More relevant studies are needed in the future to reveal the signaling pathways that these genes are involved in and their physiological and pathological mechanisms both in vitro and in vivo.

Another interesting finding is, when we analyzed the experimentally determined protein-protein interactions between these three selected m⁶A RNA methylation regulators and the other proteins (Supplementary Figure 1B), we found ZC3H13 interacted with CDKs family, POLR2 family and MED complex family closely, which is consistent with the result in Supplementary Figure 1A. In mammalian cells, cyclin-dependent kinases (CDKs) control critical cell cycle checkpoints and RNA polymerase II-dependent transcriptional events in response to extracellular and intracellular signals leading to proliferation [36, 37]. Mediator (MED) is a large multiprotein complex conserved in all eukaryotes that plays an essential role in transcriptional regulation. The mediator comprises 30 subunits in humans that form three main modules and a separable four-subunit kinase module [38]. The mediator complex interacts with DNA-binding gene-specific transcription factors to modulate transcription by RNA polymerase II [39]. This suggested that ZC3H13 was closely involved in the transcription process, which was consistent with the known function of the m⁶A RNA methylation regulators.

As for the inconsistency between IGF2BP1 and VIRMA protein levels and mRNA levels (Figure 2A and Supplementary Figure 8A, 8B), this is a phenomenon that has been observed before in studying the expression of other genes. Yanqing Liu et al. used Pearson’s correlation analysis of scatter plots to reveal an inconsistent relationship between HuR protein levels and mRNA levels in colorectal cancer tissues, and pointed out that the inconsistency between HuR protein and mRNA levels indicated that some post-transcriptional gene regulation mechanisms were involved in the control of HuR expression [40]. Pernilla Israelsson et al. evaluated the relative mRNA expression and the corresponding protein expression of the cytokines IL-6, IL-8, TNF-α and TNF-β/LTA at seven consecutive time points in the kinetic experiments with ovarian cancer cell lines OVCAR-3 and SKOV-3 and compared with T cell line Jurkat, which served as a control [41]. They also observed inconsistency between protein and mRNA expression. In addition, Xinyu Ren et al. found that PD1/PDL1 mRNA and protein expressions were inconsistent too in triple-negative breast cancer [42]. In short, the inconsistency of protein and mRNA expression levels observed in this study is not a rare phenomenon, the specific details of which require more experiments and analysis.

In conclusion, our studies comprehensively manifested the expression and prognostic value of m⁶A RNA methylation regulators in OC (Figure 7). The three selected m⁶A RNA methylation regulators, which were OC prognosis-associated factors, were also enriched in the biological processes and signaling pathways that drive the malignant progression of OC. In brief, our study provides novel markers for evaluating OC prognosis and furnishes significant proof for future research on the role of RNA m⁶A methylation in OC.

Outline for the expression changes, mechanism, and potential functions of m6A RNA methylation regulators in OC.

Figure 7. Outline for the expression changes, mechanism, and potential functions of m⁶A RNA methylation regulators in OC.

Materials and Methods

Datasets

In March 2019, we obtained the RNA-seq transcriptome data of 379 OC patients and the corresponding clinicopathological information of 587 OC patients from the TCGA database (http://cancergenome.nih.gov/) and obtained the RNA-seq transcriptome data of 88 normal human ovarian tissues from GTEx database (https://www.gtexportal.org/home/datasets). For the RNA-seq data, TCGA samples (n = 379) were normalized by fragment per kilobase of exon model per Million (FPKM, namely Fragment Per Kilobase Million, which is defined in this way [43]).

Selection and differential expression analysis of m⁶A RNA methylation regulators

We collected a list of 18 m⁶A RNA methylation regulators from published literature [44, 45]. Next, we systematically contrasted the expression of these m⁶A RNA methylation regulators in ovarian with different clinicopathological characters. All data were processed using the R software (version 3.4.0). The “limma” package was used for identifying DEGs between the OC samples and matched non-cancerous samples. The screening conditions for the differential genes were: P < 0.05(“*”), P < 0.01(“**”), P < 0.001(“***”). Heat maps of differential genes were drawn using the R-package, “pheatmap”.

Bioinformatic analysis

To evaluate the prognostic value of m⁶A RNA methylation regulators, we executed univariate Cox regression analyses of their expression in the TCGA dataset, from which we selected three regulators virtually associated with survival (P < 0.1), which we chose for further functional research and development of a potential risk signature with the LASSO Cox regression algorithm [46, 47]. Finally, three regulators and their coefficients were decided by the minimum criteria, choosing the best penalty parameter λ associated with the TGGA datasets. The risk score for the signature was counted applying the formula:

$R i s k s c o r e = \sum_{i = 1}^{n} C o d f i * x_{i}^{[48]},$

where Coefi is the coefficient, and x_i is the z-score transformed relative expression value of each selected regulator. This formula was applied to count a risk score for each patient in TGGA datasets. The high-risk subtype (samples with the risk score higher than 5) and the low-risk subtype (samples with risk score lower than 5) were defined in OC cases based on the risk score of its tumor samples. The ROC analysis was used for testing if the survival prediction is sensitive and specific based on the risk score.

PPI network construction of m⁶A RNA methylation regulators

Protein-protein interactions (PPI) analysis was conducted to reveal the molecular mechanisms of a list of 18 m⁶A RNA methylation regulators in ovarian cancer. We utilized the Search Tool for the Retrieval of Interacting Genes (STRING) protein database 11.0 (http://string-db.org/) to construct the PPI networks. An interaction score > 0.4 was regarded as the cut-off criterion.

Gene set enrichment analysis (GSEA) of the three selected m⁶A RNA methylation regulators

Gene Set Enrichment Analysis (GSEA) of prognosis-related MeDEGs was performed using GSEA 3.0 software with gene set c2 (cp.kegg.v.6.2.symbols.gmt). High throughput RNA expression of 379 ovarian cancer genes from TCGA was utilized as the dataset. Each sample was defined as either “H” or “L”, depending on whether it was greater than the median mRNA expression value of prognosis-related DEGs or not. The number and type of permutations were set at “1000” and “phenotype,” respectively. An enrichment score >0.4 and P < 0.05 were regarded as statistically significant.

Tumor subgroup gene expression and survival analyses

UALCAN (http://ualcan.path.uab.edu/index.html) is a portal for facilitating tumor subgroup gene expression and survival analyses for analyzing cancer OMICS data. Using TCGA transcriptome and clinical patient data, compare across different tumor subgroups as defined by the patient’s age and tumor grade through the expression level of the gene. Finally, the primary tumor samples were categorized using OC patient clinical data, and boxplots were generated of the expression level of each gene across various subgroups.

Survival Analysis of IGF2BP1, VIRMA and ZC3H13

Download the datum matrix containing all ovarian patients’ prognosis information from the GEO database. The overall survival analysis was conducted using the only patient with survival data and gene expression data from RNA-seq. For each gene, a tab-separated input file was created with columns for TCGA sample id, Time (days_to_death), Status (Alive or Dead), and Expression level. Samples were categorized into two groups according to the level of gene expression: (1) High expression (with values above the average) and (2) Low expression (with values below the average). The data was processed by R packages “survival”. Kaplan-Meier curves were drawn to demonstrate the relationship between the patient’s overall survival and gene expression levels of m⁶A RNA methylation regulators. The relationship was tested by the log-rank test.

Validation of the IGF2BP1, VIRMA and ZC3H13

To determine the robustness of this model, we used the same coefficients from the training set to validate in the validation sets including GSE30161 (n = 58), GSE9891 (n = 285), GSE63885 (n = 101), GSE19829 (n = 70), GSE18520 (n = 63), GSE26193 (n = 107), GSE27651 (n = 49), GSE14764 (n = 80), GSE3149 (n = 153), GSE23554 (n = 28), and TCGA (n = 427) dataset. We used multivariate Cox proportional analysis to determine a panel of prognostic genes. The calculation of the patient’s risk score in the training set was performed according to the formulate obtained from the multivariate Cox proportional model. The forest plots were used to display the multivariable Cox results, including all the above variables. The “forestplot” R packages were used to draw forest plots.

Independent prognostic analysis

The univariate and multivariable Cox regression analyses were utilized to access the prognostic value of the risk score generated from the multivariate model. The demographics and clinical information, including age and grade, were used for model correction.

Statistical analysis

One-way ANOVA was applied to contrast the expression levels of m⁶A RNA methylation regulators in ovaries with normal patients group (GTEx datasets) and tumor patients group (TCGA datasets), and t-tests were applied to contrast the expression levels in OC patients for grade, age and survival status.

Patients were grouped into three clusters by consensus expression of m⁶A RNA methylation regulators or were separated into low-risk and high-risk groups applying the median risk score (came from the risk signature) as the cut-off value. Chi-square tests were applied to contrast the distribution of patients’ age, survival status and the grade between the two risk groups.

To contrast the risk scores of the signature for ovaries with different clinicopathologic, a one-way ANOVA or t-test was conducted to contrast the risk scores in patients divided by clinical or molecular-pathological features. Univariate and multivariate Cox regression analyses were conducted to evaluate the prognostic value of the risk score and various clinical and molecular-pathological features. The prediction efficiency of the risk signature was tested with the ROC curve.

The Kaplan-Meier method with a two-sided log-rank test was applied to contrast the OS of the patients in the cluster 1/2/3 groups or in the low- and high-risk groups. All statistical analyses were executed utilizing R v3.4.1 (https://www.r-project.org/) and Prism 8 (GraphPad Software Inc., La Jolla, CA).

Supplementary Materials

Supplementary Figures

Supplementary Table 1

Abbreviations

m⁶A: N6-methyladenosine; OC: ovarian cancer; IGF2BP1: insulin-like growth factor 2 mRNA binding protein 1; VIRMA: vir like m6A methyltransferase associated; ZC3H13: zinc finger CCCH-type containing 13; GSEA: Gene Set Enrichment Analysis; m¹A: N1-methyladenosine; m⁵C: 5-methylcytosine; UTR: untranslated region; MTC: methyltransferase complex; YTH: YT521-B homology; eIF3: eukaryotic initiation factor 3; IGF2BPs: Insulin-like growth factor 2 mRNA-binding proteins; HNRNP: heterogeneous nuclear ribonucleoprotein; EMT: epithelial to mesenchymal transition; PARPi: poly ADP-ribose polymerase inhibitors; TCGA: The Cancer Genome Atlas; GTEx: Genotype-Tissue Expression; CDKs: cyclin dependent kinases; POLR2: RNA polymerase II; MED: mediator; PCA: principal component analysis; OS: overall survival; CDF: cumulative distribution function; LASSO: the least absolute shrinkage and selection operator; ROC: receiver operating characteristic; AUC: area under the curve; CPTAC: Clinical Proteomic Tumor Analysis Consortium; OV: ovarian serous cystadenocarcinoma; CCLE: Cancer Cell Line Encyclopedia; KEGG: Kyoto Encyclopedia of Genes and Genomes; GEO: Gene Expression Omnibus; PPI: protein-protein interactions; STRING: Search Tool for the Retrieval of Interacting Genes.

Author Contributions

Lili Fan conceived and designed the study, also crafted figures and tables, and was responsible for the critical reading of the manuscript. Also, she contributed to the data collection and analysis. Ying Lin, Han Lei and Liuer He was responsible for the writing and helped with the figures and tables. Han Lei, Zhipeng Yan, Liuer He, Guang Shu, Hai Rihan and Gang Yin supervised and contributed to the critical reading of the manuscript.

Acknowledgments

The authors gratefully acknowledge contributions from the TCGA network, the GTEx network and the CCLE Network.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Funding

This study was supported by a grant from the National Key R&D Program of China, Stem Cell and Translation Research (No. 2016YFA0102000).

References

1. Torre LA, Trabert B, DeSantis CE, Miller KD, Samimi G, Runowicz CD, Gaudet MM, Jemal A, Siegel RL. Ovarian cancer statistics, 2018. CA Cancer J Clin. 2018; 68:284–96. https://doi.org/10.3322/caac.21456 [PubMed]
2. Frye M, Harada BT, Behm M, He C. RNA modifications modulate gene expression during development. Science. 2018; 361:1346–49. https://doi.org/10.1126/science.aau1646 [PubMed]
3. Roundtree IA, Evans ME, Pan T, He C. Dynamic RNA modifications in gene expression regulation. Cell. 2017; 169:1187–200. https://doi.org/10.1016/j.cell.2017.05.045 [PubMed]
4. Meyer KD, Saletore Y, Zumbo P, Elemento O, Mason CE, Jaffrey SR. Comprehensive analysis of mRNA methylation reveals enrichment in 3’ UTRs and near stop codons. Cell. 2012; 149:1635–46. https://doi.org/10.1016/j.cell.2012.05.003 [PubMed]
5. Hu Y, Wang S, Liu J, Huang Y, Gong C, Liu J, Xiao Y, Yang S. New sights in cancer: component and function of N6-methyladenosine modification. Biomed Pharmacother. 2020; 122:109694. https://doi.org/10.1016/j.biopha.2019.109694 [PubMed]
6. He L, Li H, Wu A, Peng Y, Shu G, Yin G. Functions of N6-methyladenosine and its role in cancer. Mol Cancer. 2019; 18:176. https://doi.org/10.1186/s12943-019-1109-9 [PubMed]
7. Choe J, Lin S, Zhang W, Liu Q, Wang L, Ramirez-Moya J, Du P, Kim W, Tang S, Sliz P, Santisteban P, George RE, Richards WG, et al. mRNA circularization by METTL3-eIF3h enhances translation and promotes oncogenesis. Nature. 2018; 561:556–60. https://doi.org/10.1038/s41586-018-0538-8 [PubMed]
8. Pendleton KE, Chen B, Liu K, Hunter OV, Xie Y, Tu BP, Conrad NK. The U6 snRNA m⁶A methyltransferase METTL16 regulates SAM synthetase intron retention. Cell. 2017; 169:824–35.e14. https://doi.org/10.1016/j.cell.2017.05.003 [PubMed]
9. Zhang Z, Wang M, Xie D, Huang Z, Zhang L, Yang Y, Ma D, Li W, Zhou Q, Yang YG, Wang XJ. METTL3-mediated N⁶-methyladenosine mRNA modification enhances long-term memory consolidation. Cell Res. 2018; 28:1050–61. https://doi.org/10.1038/s41422-018-0092-9 [PubMed]
10. Meng TG, Lu X, Guo L, Hou GM, Ma XS, Li QN, Huang L, Fan LH, Zhao ZH, Ou XH, OuYang YC, Schatten H, Li L, et al. Mettl14 is required for mouse postimplantation development by facilitating epiblast maturation. FASEB J. 2019; 33:1179–87. https://doi.org/10.1096/fj.201800719R [PubMed]
11. Wang T, Kong S, Tao M, Ju S. The potential role of RNA N6-methyladenosine in cancer progression. Mol Cancer. 2020; 19:88. https://doi.org/10.1186/s12943-020-01204-7 [PubMed]
12. Liu S, Zhuo L, Wang J, Zhang Q, Li Q, Li G, Yan L, Jin T, Pan T, Sui X, Lv Q, Xie T. METTL3 plays multiple functions in biological processes. Am J Cancer Res. 2020; 10:1631–46. [PubMed]
13. Du K, Zhang L, Lee T, Sun T. m⁶A RNA methylation controls neural development and is involved in human diseases. Mol Neurobiol. 2019; 56:1596–606. https://doi.org/10.1007/s12035-018-1138-1 [PubMed]
14. Zhong H, Tang HF, Kai Y. N6-methyladenine RNA modification (m6A): an emerging regulator of metabolic diseases. Curr Drug Targets. 2020. [Epub ahead of print]. https://doi.org/10.2174/1389450121666200210125247 [PubMed]
15. Cheng M, Sheng L, Gao Q, Xiong Q, Zhang H, Wu M, Liang Y, Zhu F, Zhang Y, Zhang X, Yuan Q, Li Y. The m⁶A methyltransferase METTL3 promotes bladder cancer progression via AFF4/NF-κB/MYC signaling network. Oncogene. 2019; 38:3667–80. https://doi.org/10.1038/s41388-019-0683-z [PubMed]
16. Yang F, Jin H, Que B, Chao Y, Zhang H, Ying X, Zhou Z, Yuan Z, Su J, Wu B, Zhang W, Qi D, Chen D, et al. Dynamic m⁶A mRNA methylation reveals the role of METTL3-m⁶A-CDCP1 signaling axis in chemical carcinogenesis. Oncogene. 2019; 38:4755–72. https://doi.org/10.1038/s41388-019-0755-0 [PubMed]
17. Zhou Z, Lv J, Yu H, Han J, Yang X, Feng D, Wu Q, Yuan B, Lu Q, Yang H. Mechanism of RNA modification N6-methyladenosine in human cancer. Mol Cancer. 2020; 19:104. https://doi.org/10.1186/s12943-020-01216-3 [PubMed]
18. Antony J, Tan TZ, Kelly Z, Low J, Choolani M, Recchi C, Gabra H, Thiery JP, Huang RY. The GAS6-AXL signaling network is a mesenchymal (mes) molecular subtype-specific therapeutic target for ovarian cancer. Sci Signal. 2016; 9:ra97. https://doi.org/10.1126/scisignal.aaf8175 [PubMed]
19. Hua W, Zhao Y, Jin X, Yu D, He J, Xie D, Duan P. METTL3 promotes ovarian carcinoma growth and invasion through the regulation of AXL translation and epithelial to mesenchymal transition. Gynecol Oncol. 2018; 151:356–65. https://doi.org/10.1016/j.ygyno.2018.09.015 [PubMed]
20. Müller S, Glaß M, Singh AK, Haase J, Bley N, Fuchs T, Lederer M, Dahl A, Huang H, Chen J, Posern G, Hüttelmaier S. IGF2BP1 promotes SRF-dependent transcription in cancer in a m⁶A- and miRNA-dependent manner. Nucleic Acids Res. 2019; 47:375–90. https://doi.org/10.1093/nar/gky1012 [PubMed]
21. Fukumoto T, Zhu H, Nacarelli T, Karakashev S, Fatkhutdinov N, Wu S, Liu P, Kossenkov AV, Showe LC, Jean S, Zhang L, Zhang R. N⁶-methylation of adenosine of FZD10 mRNA contributes to PARP inhibitor resistance. Cancer Res. 2019; 79:2812–20. https://doi.org/10.1158/0008-5472.CAN-18-3592 [PubMed]
22. Chandrashekar DS, Bashel B, Balasubramanya SA, Creighton CJ, Ponce-Rodriguez I, Chakravarthi BV, Varambally S. UALCAN: a portal for facilitating tumor subgroup gene expression and survival analyses. Neoplasia. 2017; 19:649–58. https://doi.org/10.1016/j.neo.2017.05.002 [PubMed]
23. Yang X, Wang J. Precision therapy for acute myeloid leukemia. J Hematol Oncol. 2018; 11:3. https://doi.org/10.1186/s13045-017-0543-7 [PubMed]
24. Cui Q, Shi H, Ye P, Li L, Qu Q, Sun G, Sun G, Lu Z, Huang Y, Yang CG, Riggs AD, He C, Shi Y. m⁶A RNA methylation regulates the self-renewal and tumorigenesis of glioblastoma stem cells. Cell Rep. 2017; 18:2622–34. https://doi.org/10.1016/j.celrep.2017.02.059 [PubMed]
25. Visvanathan A, Patil V, Arora A, Hegde AS, Arivazhagan A, Santosh V, Somasundaram K. Essential role of METTL3-mediated m⁶A modification in glioma stem-like cells maintenance and radioresistance. Oncogene. 2018; 37:522–33. https://doi.org/10.1038/onc.2017.351 [PubMed]
26. Du M, Zhang Y, Mao Y, Mou J, Zhao J, Xue Q, Wang D, Huang J, Gao S, Gao Y. MiR-33a suppresses proliferation of NSCLC cells via targeting METTL3 mRNA. Biochem Biophys Res Commun. 2017; 482:582–89. https://doi.org/10.1016/j.bbrc.2016.11.077 [PubMed]
27. Ma JZ, Yang F, Zhou CC, Liu F, Yuan JH, Wang F, Wang TT, Xu QG, Zhou WP, Sun SH. METTL14 suppresses the metastatic potential of hepatocellular carcinoma by modulating N⁶-methyladenosine-dependent primary MicroRNA processing. Hepatology. 2017; 65:529–43. https://doi.org/10.1002/hep.28885 [PubMed]
28. Zhang C, Samanta D, Lu H, Bullen JW, Zhang H, Chen I, He X, Semenza GL. Hypoxia induces the breast cancer stem cell phenotype by HIF-dependent and ALKBH5-mediated m⁶A-demethylation of NANOG mRNA. Proc Natl Acad Sci USA. 2016; 113:E2047–56. https://doi.org/10.1073/pnas.1602883113 [PubMed]
29. Nash Z, Menon U. Ovarian cancer screening: current status and future directions. Best Pract Res Clin Obstet Gynaecol. 2020; 65:32–45. https://doi.org/10.1016/j.bpobgyn.2020.02.010 [PubMed]
30. Bonifácio VD. Ovarian cancer biomarkers: moving forward in early detection. Adv Exp Med Biol. 2020; 1219:355–63. https://doi.org/10.1007/978-3-030-34025-4_18 [PubMed]
31. Natanzon Y, Goode EL, Cunningham JM. Epigenetics in ovarian cancer. Semin Cancer Biol. 2018; 51:160–69. https://doi.org/10.1016/j.semcancer.2017.08.003 [PubMed]
32. Singh A, Gupta S, Sachan M. Epigenetic biomarkers in the management of ovarian cancer: current prospectives. Front Cell Dev Biol. 2019; 7:182. https://doi.org/10.3389/fcell.2019.00182 [PubMed]
33. Moufarrij S, Dandapani M, Arthofer E, Gomez S, Srivastava A, Lopez-Acevedo M, Villagra A, Chiappinelli KB. Epigenetic therapy for ovarian cancer: promise and progress. Clin Epigenetics. 2019; 11:7. https://doi.org/10.1186/s13148-018-0602-0 [PubMed]
34. Noubissi FK, Yedjou CG, Spiegelman VS, Tchounwou PB. Cross-talk between Wnt and hh signaling pathways in the pathology of basal cell carcinoma. Int J Environ Res Public Health. 2018; 15:1442. https://doi.org/10.3390/ijerph15071442 [PubMed]
35. Huang X, Zhang H, Guo X, Zhu Z, Cai H, Kong X. Insulin-like growth factor 2 mRNA-binding protein 1 (IGF2BP1) in cancer. J Hematol Oncol. 2018; 11:88. https://doi.org/10.1186/s13045-018-0628-y [PubMed]
36. Parua PK, Fisher RP. Dissecting the pol II transcription cycle and derailing cancer with CDK inhibitors. Nat Chem Biol. 2020; 16:716–24. https://doi.org/10.1038/s41589-020-0563-4 [PubMed]
37. Sánchez-Martínez C, Lallena MJ, Sanfeliciano SG, de Dios A. Cyclin dependent kinase (CDK) inhibitors as anticancer drugs: recent advances (2015-2019). Bioorg Med Chem Lett. 2019; 29:126637. https://doi.org/10.1016/j.bmcl.2019.126637 [PubMed]
38. Verger A, Monté D, Villeret V. Twenty years of mediator complex structural studies. Biochem Soc Trans. 2019; 47:399–410. https://doi.org/10.1042/BST20180608 [PubMed]
39. Soutourina J. Transcription regulation by the mediator complex. Nat Rev Mol Cell Biol. 2018; 19:262–74. https://doi.org/10.1038/nrm.2017.115 [PubMed]
40. Liu Y, Chen X, Cheng R, Yang F, Yu M, Wang C, Cui S, Hong Y, Liang H, Liu M, Zhao C, Ding M, Sun W, et al. The jun/miR-22/HuR regulatory axis contributes to tumourigenesis in colorectal cancer. Mol Cancer. 2018; 17:11. https://doi.org/10.1186/s12943-017-0751-3 [PubMed]
41. Israelsson P, Dehlin E, Nagaev I, Lundin E, Ottander U, Mincheva-Nilsson L. Cytokine mRNA and protein expression by cell cultures of epithelial ovarian cancer-methodological considerations on the choice of analytical method for cytokine analyses. Am J Reprod Immunol. 2020; 84:e13249. https://doi.org/10.1111/aji.13249 [PubMed]
42. Ren X, Wu H, Lu J, Zhang Y, Luo Y, Xu Q, Shen S, Liang Z. PD1 protein expression in tumor infiltrated lymphocytes rather than PDL1 in tumor cells predicts survival in triple-negative breast cancer. Cancer Biol Ther. 2018; 19:373–80. https://doi.org/10.1080/15384047.2018.1423919 [PubMed]
43. Leek JT. Svaseq: removing batch effects and other unwanted noise from sequencing data. Nucleic Acids Res. 2014; 42:e161. https://doi.org/10.1093/nar/gku864 [PubMed]
44. Deng X, Su R, Weng H, Huang H, Li Z, Chen J. RNA N⁶-methyladenosine modification in cancers: current status and perspectives. Cell Res. 2018; 28:507–17. https://doi.org/10.1038/s41422-018-0034-6 [PubMed]
45. Yang Y, Hsu PJ, Chen YS, Yang YG. Dynamic transcriptomic m⁶A decoration: writers, erasers, readers and functions in RNA metabolism. Cell Res. 2018; 28:616–24. https://doi.org/10.1038/s41422-018-0040-8 [PubMed]
46. Sauerbrei W, Royston P, Binder H. Selection of important variables and determination of functional form for continuous predictors in multivariable model building. Stat Med. 2007; 26:5512–28. https://doi.org/10.1002/sim.3148 [PubMed]
47. Hu X, Martinez-Ledesma E, Zheng S, Kim H, Barthel F, Jiang T, Hess KR, Verhaak RG. Multigene signature for predicting prognosis of patients with 1p19q co-deletion diffuse glioma. Neuro Oncol. 2017; 19:786–95. https://doi.org/10.1093/neuonc/now285 [PubMed]
48. Zhou Z, Huang R, Chai R, Zhou X, Hu Z, Wang W, Chen B, Deng L, Liu Y, Wu F. Identification of an energy metabolism-related signature associated with clinical prognosis in diffuse glioma. Aging (Albany NY). 2018; 10:3185–209. https://doi.org/10.18632/aging.101625 [PubMed]

Research Paper Volume 12, Issue 18 pp 18453—18475