Identification and validation of a novel prognostic signature with DNA methylation-driven genes for survival prediction of breast cancer

doi:10.21203/rs.3.rs-2070655/v1

Download PDF

Research Article

Identification and validation of a novel prognostic signature with DNA methylation-driven genes for survival prediction of breast cancer

https://doi.org/10.21203/rs.3.rs-2070655/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Background:

Research has been accruing to demonstrate that DNA methylation plays a crucial role in the diagnosis of breast cancer mainly through regulating mRNA expression. Our study aims to construct a risk signature based on the methylation-driven genes (MDGs) to predict patients’ prognoses and identify tumors’ underlying molecular mechanisms.

Methods:

The data included in this study were downloaded from TCGA and GEO databases. Subsequently, univariate Cox regression and LASSO Cox regression analyses were constructed to identify prognostic MDGs and construct a risk signature. We have also used the ROC curve and Kaplan-Meier analysis to assess the predictive performance of the signature. Multivariate Cox regression analysis was used to identify the independent prognostic factor, and a nomogram was built to facilitate the use of the signature in clinical. Finally, GSVA, TISIDB, CIBERSORT, and drug-sensitive analyses were used to explore the potential mechanisms, and an eRNA network was constructed to identify potential regulators of the risk signature.

Results:

A total of 288 MDGs were identified in breast cancer, and 19 prognosis-related MDGs were included in the risk signature to predict patients’ overall survival with satisfactory performance. We identified that the 19-gene risk signature is an independent prognostic factor and could stratify patients into low- and high-risk groups with different prognoses. Furthermore, patients under different risk situations have diverse proportions of infiltrating immune cells, frequently mutated genes, and sensitive drugs. Nomogram integrating risk signature and clinicopathological factors achieve excellent predictive ability. A ceRNA network consisting of 9 lncRNA, 38 miRNA, and 10 mRNA was constructed based on the MDGs identified in the risk signature.

Conclusions:

We have successfully constructed an MDG-based prognostic risk signature in breast cancer and established a corresponding nomogram model. Additionally, we uncovered the underlying molecular mechanisms and therapeutic targets in tumors with different risks.

Methylation-driven genes (MDGs)

Breast cancer

Prognosis

Risk signature

Biomarker

Breast cancer is the most common type of cancer worldwide, responsible for the death of over half a million people every year. In developing countries, the incidence and mortality rate of breast cancer cases keep increasing, seriously threatening women’s health[1]. However, in recent years, early detection and efficient systemic therapies could successfully improve the prognosis of some patients in developed countries. Nevertheless, despite these advancements in clinical diagnoses, how to estimate cancer risks and how to avoid overtreatment and undertreatment remain unclear. Therefore, it is essential to identify effective markers in staging primary breast cancer and guiding appropriate treatments.

DNA methylation is a kind of epigenetic modification that is widely investigated in cancer. Proper DNA methylation is crucial to maintaining genome stability, and aberrant DNA methylation is a factor that leads to dysregulation of gene expression, contributing to cancer development[2]. A number of studies have delineated the pattern and function of DNA methylation in breast cancer[3, 4]. Martin Widschwendter et.al proposed that the methylation of tumor-associated genes can contribute to six new capabilities of breast carcinogenesis, such as p16^Ink4A in disturbing the “limitless replicative potential of cells”, ERα and ERβ in regulating “self-sufficiency in growth signals”[5]. Aside from single methylated genes investigated with MSP, genome-wide DNA methylation profiles identified groups of hypo/hyper-methylated genes that could assist breast cancer’s prognosis[6]. Meanwhile, predictive signatures with a cluster of methylated genes/ sites were also proposed in breast cancer[7–11].

In general, hypermethylation is closely correlated to gene silencing, and always indicates gene activation. Currently, many studies have shown that hypo/hypermethylation genes accompanied with correspondingly irregular mRNA expression are common in cancer, and could promote cancer development through mediating cancer-related pathways[12]. Hence, these methylation genes are identified as methylation-driven genes (MDGs) and are considered potential biomarkers in assessing cancer risk and guiding therapies[13, 14]. In search of effective and sensitive predictive factors for cancer, the integration of several molecular biomarkers into a single model is widely used for prognostic evaluation. Since DNA methylation can simultaneously influence multiple genes rather than a single gene, predictive signatures consisting of a cluster of MDGs are more valuable in staging breast cancer.

To identify effective risk signature of breast cancer based on MDGs, we extracted data from the TCGA database and constructed a 19-gene risk signature that could be used to stratify patients into low-risk and high-risk groups with diverse prognoses. Furthermore, we verified the signature with both datasets, and achieved satisfactory predictive performance. Subsequently, we explored the molecular mechanism and investigated sensitive drugs in patients with different risks, aiming to provide different insights on precise treatment and effective mechanisms for breast cancer patients. Additionally, we outlined a nomogram with the signature and clinicopathological characteristics that could help clinicians and patients to make accurate decisions. At last, we built a ceRNA network based on the hub MDGs to identify potential regulators about the signature.

1. Data acquisition and procession

The DNA methylation data, mRNA expression data and corresponding clinicopathological features of TCGA-BRCA were retrieved from UCSC Xena database (http://xena.ucsc.edu/). Methylation data was processed with the R package “ChAMP”, genes with “Tumor AVG < 0.3 & Control AVG > 0.2 & adj.P < 0.05” and “Tumor AVG > 0.3 & Control AVG < 0.2 & adj.P < 0.05” were defined as differential methylated genes (DMGs). Differential expression of mRNA in 1217 patients were analyzed with the R package “limma”, and genes with “|Log₂FC| > 1 & adj.P < 0.05” were defined as differential mRNA expressed genes (DEGs).

For external validation, 327 and 252 patients with integrally mRNA expression data and survival data were separately obtained from the GEO database (GSE20685, GSE21653). Detailed information of the GSE20685 and GSE21653 datasets were both obtained from the GPL570 platform. Data were processed with R package “limma”, and mRNAs with “|Log₂FC| > 1 & adj.P < 0.05” were defined as DEGs.

2. Construction and validation of the risk signature and nomogram

Overlapping genes with opposite methylation and mRNA expression were defined as methylation-driven genes (MDGs). 288 MDGs were brought into univariate COX regression analysis with the R package “survival” and the LASSO regression analysis with the R package “glmnet”. Through selection, MDGs closely associated with patient’s prognosis (P < 0.05) were identified as candidate genes. Furthermore, 19 hub MDGs obtained from LASSO analysis were used to construct the risk signature. The riskscore of each patient that based on the signature was calculated using the formula of “Risk score = \({\sum }_{i=1}^{n}Expi*i\)”. n represents the number of hub genes in the risk signature, βi was the regression coefficient and Expi was the relative expression value of each hub gene.

Based on the risk signature, the riskscore of each patient was calculated, and the median riskscore was act as a cutoff value to stratify patients into low-risk and high-risk groups. The association between risk signature and overall survival was further analyzed according to the risk groups. Survival curves were analyzed with Kaplan-Meier and calculated with log-rank. Time-dependent receiver operating characteristic (ROC) curve was constructed and the area under the curve (AUC) was calculated to assess the predictive accuracy of the risk signature. The nomogram for predicting the overall survival of breast cancer was constructed with the riskscore and clinicopathological characteristics (age, gender, T stage, N stage, M stage, TNM stage). Calibration curve and C-index were used to assess the performance of the prognostic nomogram.

3. Bioinformatic analyses

3.1 GO and KEGG analysis

GO and KEGG analyses were performed by the R package “ClusterProfiler (Version3.14.3)”, Q value < 0.05 were considered statistically significant for GO and KEGG analysis.

3.2 PPI analysis

Protein-protein interaction (PPI) among MDGs was analyzed with Cytoscape software (Version3.9)

3.3 GSVA analysis

GSVA was performed to evaluate differential pathways between low-risk and high-risk groups with the R package of “GSVA”, Q value < 0.05 was considered as statistically significant.

3.4 TISIDB analysis

TISIDB is a website facilitating in comprehensive investigation of the interaction between tumors and immune system. Thus, TISIDB was utilized in our study to explore the relationship between signature gene sets and the immune system, related data were downloaded from the TISIDB website (http://cis.hku.hk/TISIDB/index.php).

3.5 Investigation of tumor infiltrating immune cells

The CIBERSORT was used to analyze the RNA-seq data of breast cancer in different groups, thereby to investigate the relative proportion of 22 infiltrating immune cells. Relative data used in CIBERSORT were downloaded from CIBERSORT web portal (http://cibersort.stanford.edu/). Furthermore, the correlations among each immune cell were analyzed by Spearman correlation. The differences of immune infiltrating cells between low-risk and high-risk groups were analyzed by Wilcox test, and P < 0.05 was considered as statistically significant. The Spearman correlation analysis was used to analyze the correlations between riskscore and infiltrating immune cells.

3.6 Investigation of tumor mutation and neoantigen

Tumor mutation, neoantigens, and tumor mutational burden (TMB) data were all obtained from the TCGA data portal. The analysis of top mutated genes in low-risk and high-risk groups using Wilcox test.

3.7 Exploration Of Potentially Effective Therapeutic Drugs

Based on the data in “The Genomics of Drug Sensitivity in Cancer (GDSC)” database, predictive analysis with the R package of “pRRophetic” was performed to investigate the common administrating drugs’ sensitivity of each tumor. The differences in the IC₅₀ between low-risk and high-risk groups were compared with Wilcox test.

3.8 Construction Of Cerna Network

The interaction network of mRNA-miRNA was analyzed with FunRich (Version 3.1.3), and the interaction network of miRNA-lncRNA was analyzed with the DIANA-LncBase database. The miRNA-lncRNA interactions were further screened with a threshold of MIC > 0.4. Finally, a ceRNA network was established by integrating mRNA-miRNA interaction and lncRNA-miRNA interactions.

4. Statistical analyses

Spearman correlation analysis was used to analyze the correlation between DNA methylation and mRNA expression as well as each subtype of infiltrating immune cells. Multivariate Cox regression analysis was used to identify the hazard ratio (HR) of each prognostic factor. All the statistical analyses were performed by R software (Version3.6.2), P < 0.05 was considered statistically significant.

1. Identification and function analysis of Methylation-driven genes in breast cancer

To identify an effective predictive model based on DNA methylation-driven genes (MDGs) that could evaluate the risk of breast cancer, we combined different analyses, including hub genes’ identification, risk signature construction and validation, mechanisms exploration, predictive nomogram’s construction, and ceRNA network’s construction (Fig. 1A).

Furthermore, to identify potentially valuable methylated genes that can be used in evaluating breast cancer’s prognosis, we performed the genomic analysis of differential methylated genes (DMGs) and differential mRNA expressed genes (DEGs) based on the data extracted from the TCGA database. We have identified 2788 methylated genes exhibiting significantly 5-methylation (5mC) change between breast cancer tissues and normal tissues. Among these, 2505 differential methylated genes (DMGs) exhibited an evident 5mC increase (Tumor AVG > 0.3 & Control AVG < 0.2 & adj.P < 0.05), and 283 DMGs had a considerable 5mC decrease (Tumor AVG < 0.3 & Control AVG > 0.2 & adj.P < 0.05) (Fig. 1B).

Moreover, 1722 significantly differential mRNA expressed genes (DEGs) were identified between breast cancer and normal tissues (|Log₂ FC| > 1 & adj.P < 0.05), including 615 up-regulated DEGs and 1107 down-regulated DEGs in cancer tissues (Fig. 1C). Through the above genomic analysis, we observed 288 methylation-driven genes (MDG) with opposite 5mC and mRNA change in breast cancers, which can be used as candidate genes to establish an effective risk signature to evaluate patient’s prognosis (Fig. 1D). Subsequently, results of the GO and KEGG analyses showed that the candidate 288 genes were significantly enriched in some cancer-related pathways, such as extracellular structure organization, collagen − containing extracellular matrix, amide binding pathways in GO analysis (Fig. 1E), focal adhesion, and JAK − STAT signaling pathways in KEGG analysis (Fig. 1F). The protein-protein interaction (PPI) network of the 288 candidate DMGs is illustrated in detail in Figure S1.

2. Establishment of a risk signature with 19 methylation-driven genes for predicting the overall survival of breast cancer

To identify eligible genes in establishing risk signature, we performed univariate Cox regression analysis in combination with LASSO regression analysis based on the 288 candidate MDGs (Fig. 2A−2B). In univariate Cox regression analysis, 33 genes associated with patient’s prognosis (P < 0.05) were selected and enrolled in the following LASSO regression analysis. Subsequently, we randomly divided patients into training (n = 855) and testing (n = 214) sets (Table 1), and selected the best model for predicting overall survival of breast cancer patients through LASSO analysis in training set. Finally, 19 out of 33 MDGs were included in the risk signature, and the corresponding coefficients were computed as shown in Fig. 2C. Based on the coefficient value and the expression of each candidate MDGs, we constructed a risk signature as follows: Risk Score = (−0.12824∗MDFIC) + (−0.10812∗ STAT5A) + (−0.09103∗ TNFRSF18) + (−0.07636∗ NDRG2) + (−0.0627∗EGR2) + (−0.05786∗ NRG1) + (−0.05401∗KRT19) + (−0.05311∗N4BP2L1) + (−0.05093 ∗BATF) + (−0.04251∗CXCL14) + (−0.02905∗CCND2) + (−0.01227∗C2orf40) + (−0.00543∗SFRP1) + (−0.00011∗PPP1R15A) + (0.008829∗MAL2) + (0.016824∗NPR3) + (0.10699∗C2CD2) + (0.114751∗EMP1) + (0.164026∗NT5E).

To evaluate the efficiency of the risk signature, we classified patients into low-risk and high-risk groups according to their median score. Results of the Kaplan-Meier analysis showed that in the training set, patients assessed with low risks had better overall survival rates than those with high risks, and similar survival benefits were observed in the individuals assigned to the testing set (Fig. 2D). Additionally, in both sets, the AUC of the survival assessment model ranged from 0.7 to 0.78 (Fig. 2E), indicating that the risk signature is valuable in predicting the prognosis of breast cancer. Alongside the model assessment, we performed a correlation analysis and observed that the expression of the 19 MDGs was mostly associated with each other (P < 0.05) (Fig. 2F). In addition, we identified a considerable negative correlation between DNA methylation and mRNA expression in NDRG2, C2CD2, SFRP1, C2orf40, and STAT5A (|Cor| > 0.5) (Fig. 2G−2K).

3. Investigation of the potential molecular regulators in patients with different risks

To investigate molecular pathways interfering patients’ prognoses, we performed gene set variation analysis (GSVA) between low-risk and high-risk groups. GSVA revealed that the G2M checkpoint, mitotic spindle, and TNFα signaling pathways were significantly enriched in the high-risk group, and fatty acid metabolism, DNA repair, angiogenesis, and P53 pathways were enriched in the low-risk group (Fig. 3A).

Aiming to investigate whether the immune system regulators are associated with our risk signature, we performed a correlation analysis based on data extracted from the TISIDB database. In comparison to the low-risk group, chemokine and immune receptor-associated genes were lowly expressed in the high-risk group (Fig. 3B). The relationship between risk score and immune checkpoint molecules was further evaluated by using TISIDB data, and the results indicated that PDCD1, TGFBR1, TGFB1, and CD244 are closely correlated with patient’s risk score (P < 0.05) (Fig. 3C). To further explore the effect of infiltrating immune cells in breast cancer tissues, we analyzed the 22 immune cell subpopulations in each sample with CIBERSORT algorithm (Figure S2) and performed a spearman correlation analysis considering the quantity of each infiltrating immune cell (Fig. 3D). The results showed that tumors in the high-risk group have a great number of macrophages M0 and macrophages M2 cells, while tumors in the low-risk group have more B cell naïve, T cells CD8, and T cells follicular helper cells (Fig. 3E). The outcomes of further correlation analysis revealed that the risk score is positively correlated with the level of macrophages M0 and macrophages M2 cells, and negatively correlated with T cells CD8 and dendritic cells resting (Fig. 3F). These results suggest that immune-cell characteristics are associated with our risk signature.

We also assessed the mutation scenario between patients with different risks. As shown in Fig. 3G, 49% and 20% of patients in the high-risk group had TP53 and TTN mutation, respectively. In the low-risk group, these rates decreased to 27% and 16% of patients. However, 39% and 20% of patients in the low-risk group had PIK3CA and CDH1 mutation compared to 34% and 10% of individuals in the high-risk group. Subsequently, we calculated the tumor mutational burden (TMB) and the level of neoantigen in each patient. In comparison to the low-risk group, both indicators were significantly higher in patients in the high-risk set (Fig. 3H).

4. Exploring the potential therapeutic agents in clinical treatment based on the risk model

To explore potential therapeutic strategies, we attempted to analyze the efficiency of some standard chemotherapeutics in the low-risk and high-risk groups, aiming to identify potential therapeutic agents. In our analysis, we observed that a high-risk score was associated with a higher half inhibitory centration (IC₅₀) for docetaxel (P < 0.05) (Fig. 4A) and paclitaxel (P < 0.05) (Fig. 4B), and a lower IC₅₀ for cisplatin (P < 0.05) (Fig. 4C), which indicated that patients in the low-risk group may be sensitive to docetaxel and paclitaxel-based chemotherapy and those in the high-risk group may be more sensitive to cisplatin. No significant differences of doxorubicin were observed between patients with different risks (Fig. 4D). Similar to the common chemotherapeutics, high-risk scores were also associated with part of the epidermal growth factor receptor and tyrosine kinase inhibitors, such as gefitinib (P < 0.0001) (Fig. 4E) and lapatinib (P < 0.05) (Fig. 4F). No significant differences were observed in imatinib (Fig. 4G). The efficacy of therapeutic drugs such as gemcitabine, cytarabine, metformin, methotrexate, and bexarotene also exhibit significant differences between low-risk and high-risk groups (Fig. 4H−4L). Taken together, the predictive model may act as a potential indicator for therapeutic decision making.

5. Validation of the prognostic 19-gene risk signature in external cohorts

To ensure the usefulness and effectiveness of the 19-gene risk signature in different populations, we performed an external validation assessment with the data of GSE20658 and GSE21653 in the GEO database. Based on the risk signature, patients were stratified into low-risk and high-risk groups and the survival rates of those in the high-risk group were significantly lower than those in the low-risk group (Fig. 5A−5B). The AUC values in the GSE20658 dataset for predicting the overall survival after 3-, 5- and 7-years were 0.64, 0.65, and 0.67, respectively (Fig. 5C). Similar to the results of the GSE20658 dataset, the AUC values in the GSE21653 dataset were separately 0.67, 0.65, and 0.66 (Fig. 5D).

6. Construction of the prognostic nomogram with the 19-gene risk signature for survival analysis of breast cancer

Based on the risk signature and clinicopathological characteristics of patients, we constructed a prognostic nomogram to provide clinicians with an effective quantitative method to predict patients’ prognoses (Fig. 6A). Subsequently, we validated the nomogram with a calibration plot and concordance index analysis (C-index). The calibration plot displayed excellent performance, and the C-index of the nomogram was 0.778, indicating that the prediction efficiency of the prognostic nomogram was satisfactory (Fig. 6B). To evaluate the independent effect of the 19-gene risk signature on breast cancer survival, we performed a multivariate COX regression analysis with the predictive risk score and clinicopathological parameters, including age, gender, TNM stage, T stage, N stage, and M stage. As shown in Fig. 6C, age and the risk score of the signature were independent factors for predicting overall survival (Age: HR = 1.04, 95%CI = 1.01−1.1; Risk score: HR = 4.58, 95%CI = 1.78−11.7). We have also observed that the predictive risk score exhibited differences among patients with different fustats and T stages (P < 0.05) (Fig. 6D).

7. Construction of the lncRNA-miRNA-mRNA ceRNA network

Considering long non-coding RNAs (lncRNAs) can act as competitive endogenous RNAs (ceRNAs) to bind to microRNAs (miRNAs), and consequently interfere with the expression of target MDGs. We performed reverse prediction with FunRich analysis and identified 123 pairs of mRNA-miRNA based on the 19 genes involved in the risk signature. Subsequently, according to the miRNA, 22057 pairs of mRNA-miRNA-lncRNA were identified in the DIANA-LncBase database. To further distinguish crucial lncRNAs that may regulate target MDGs, Maximal Information Coefficient (MIC) algorithm was used to identify the robustness of the pairwise relationship between lncRNA and miRNA-mRNA (MIC > 0.4). Finally, a ceRNA network with 73 pairs of lncRNA-miRNA-mRNA was constructed, which comprised of lncRNA, 38 miRNAs, and 10 mRNAs, and may be a potential regulatory mechanism involved in the prediction of the 19-gene risk signature (Fig. 7A)

Accumulating evidence has demonstrated that aberrant DNA methylation plays a crucial role in promoting breast cancer development through regulating gene expression, which is named methylation-driven genes (MDGs). In this study, we identified 288 MDGs in breast cancer and constructed a risk signature based on 19 prognosis-associated MDGs with a moderate predictive ability of patients’ overall survival. Furthermore, we explored the possible molecular features of immune regulation, mutation patterns, and drug sensitivity in patients with different risks. We have also established a nomogram model that combined the risk signature and clinicopathological characteristics, suitable to serve as a predictive biomarker in clinical applications. At last, we built a ceRNA network with 73 pairs of lncRNA-miRNA-mRNA to uncover the underlying molecular in regulating the proposed risk signature.

In recent years, with the development of genomics in cancer, predictive signatures with a cluster of methylated genes and methylation sites of breast cancer have been proposed by many studies. Chuntao Tao et al. constructed a novel prognostic biomarker for breast cancer with 7 DNA methylation sites[7]. In another similar study conducted by Guoyang et al.’s, 4 dmCpG sites were identified and included in the OS-associated promoter CpG sites signature for predicting breast cancer’s prognosis[8]. Xiongdong Zhong et al. screened out 10 prognosis-associated methylation-driven genes found in breast cancer and constructed a predictive nomogram model with 3 out of 10 genes and several clinical features[9]. Yinqi Gao et al. established a prognostic signature in the TCGA-TNBC cohort with 5 differentially methylated sites (DMSs)[10]. Ming Zhang et al. reported a predictive model built with 7 DMSs, which could be used to distinguish breast cancer from normal tissues[11]. Compared to previous predictive models constructed with the abundance of DNA methylation, our predictive signature was based on prognosis-associated methylation-driven genes’ mRNA expression, and we verified the efficiency of the model with multiple validations. Additionally, we noticed that previous studies only proposed predictive models with several DMSs, so our study further performed mechanism exploration and potential drugs’ investigation. Moreover, we constructed a predictive nomogram and mRNA-miRNA-lncRNA network that have important implications in clinical application and identifying novel biomarkers for further research. Our study not only established a new method in breast cancer’s classification, but also shed light on new perspectives for the further molecular investigation and potential treatments for breast cancer.

In our study, all the methylation-driven genes involved in the risk signature were significantly associated with breast cancer prognosis. Amid these genes, 12 out of 19 were previously reported to be associated with breast cancer, a result that could be confirmed via molecular investigation. In previous studies, KRT19[15], CCND2[16], SFRP1[17], C2ORF40[18], and NDRG2[19] could act as tumor suppressors in breast cancer, while TNFRSF18[20], BATF[21], CXCL14[22], and MAL2 [23]could promote the progression of breast cancer. In addition, studies have formerly identified that NRG1[24] and STAT5A[25] play a dual role as both candidate oncogene and tumor suppressor, but the function of NT5E is still controversial[26, 27]. We observed that some of the target genes might promote cancer development through altering their methylation, such as NT5E, NRG1, NDRG2, CCND2[16, 28], and C2orf40. In line with previous reports, a significant correlation between methylation level and mRNA expression of NDRG2, CCND2, and C2orf40 was also observed in our study. Aside from tumor cell regulation, some of the genes, including CXCL14, MAL2[29], NT5E, and TNFRSF18 [20, 30]were demonstrated to regulate breast cancer through tumor immunity and be potential targets of immunotherapy.

The tumor’s immune microenvironment plays a key role in tumor development, clinical therapeutic efficacy, and overall survival of diagnosed patients. Therefore, we speculated that tumor immunity might be the potential mechanism by which the risk signature can be associated with patients’ prognoses. Through immune analysis with CIBERSORT algorithm and data extracted from the TISIDB database, different subpopulations of immune cells were observed in both the low- and high-risk groups. The level of macrophage M0/M2 was positively correlated with the risk score and significantly enriched in the high-risk group, while in the low-risk group, the proportion of CD8⁺ T and B naïve cells was higher. Previous studies have demonstrated that tumor-associated macrophages M0 and M2 are generally associated with poor prognosis[31], while CD8⁺ T and B cells are correlated with better breast cancer prognosis[32, 33]. The survival advantages/ disadvantages of the low-risk/ high-risk patients indicated that the subpopulation of infiltrating immune cells might be the potential factors that affect prognosis and may develop as effective immunotherapeutic targets in breast cancer classified with our risk signature.

Aside from tumor immunity, tumor mutation is also pivotal to cancer development. In our research, we identified that TP53 mutation was more frequently observed in patients from the high-risk group, which has a higher proportion of macrophage M0 and M2. Differently, in patients from the low-risk group, the proportion of tumors with PIK3CA mutation is higher. Our conclusions regarding the relationship between tumor mutation and immune cells are consistent with Lin Li et.al’s results, who reported that TP53-mutated cancers have higher level of antitumor immune signatures[34]. Moreover, previous reports have uncovered that TP53 mutation is commonly observed in triple-negative breast cancer (TNBC)[35], and PIK3CA mutation is frequently detected in Luminal and Her2 breast cancer[36, 37]. These conclusions, combined with the results of our analysis, reinforce the idea that more patients with high risks are prone to develop molecular features of TNBC, and patients with low risk are more likely to have Luminal features.

Apart from the molecular investigation, we also expect to screen out candidate therapeutic targets and suitable drugs for patients with different risks. We have observed higher TMB and neoantigen load were observed in patients from the high-risk group, and many of them are also more sensitive to cisplatin. Considering previous studies have suggested that high TMB and neoantigen load are positively correlated with immunotherapy’s sensitivity[38], so we concluded that patients with high risks included in our study may benefit from immunotherapy and chemotherapy based on the use of cisplatin. In contrast, patients in the low-risk group were sensitive to common taxane chemo drugs (paclitaxel and docetaxel) other than anthracycline (doxorubicin). Moreover, low-risk patients are also more sensitive to some antineoplastic drugs and anti-HER2 drugs, which have been proven to affect breast cancer, but have not been widely used in clinical treatment as much as cytarabine, gemcitabine, gefitinib, lapatinib, and methotrexate.

It is widely known that lncRNAs may act as ceRNAs to regulate mRNAs expression, thereby, we performed a ceRNA network analysis to uncover the potential regulators of our risk signature. In our analysis, we identified 38 miRNAs and 9 lincRNAs associated with the target DMGs, and some of the lincRNAs were reported to play a role in cancer development early on. Studies have previously detected that HOXA-AS2 [39]and LINC00707 [40]could promote tumor invasion and metastasis. Other related studies also clarified that the high expression of MAGI2-AS3[41, 42], ADAMTS9-AS2[43], and LINC01140[44] could inhibit the proliferation of breast cancer cells. Until now, no studies have identified the relationship between SOX9-AS1, LINC01354, LINC01198, and breast cancer, but their oncogenic roles in other tumors (e.g., liver cancer, lung cancer, and colorectal cancer) have been well demonstrated.

Our study constructed a risk signature with 19 methylation-driven genes, and the signature is predictable in multiple validations. However, there are still several limitations to our study. Firstly, we only performed bioinformatic analysis without experimental validation. Secondly, further prospective studies are expected to evaluate the value of the proposed risk signature and nomogram in clinical practice. Furthermore, experimental exploration is still needed to elucidate the underlying mechanism of the risk signature and identify the effective therapeutic targets.

MDGs

methylation-driven genes

DMGs

differential methylated genes

DEGs

differential mRNA expressed genes

DMSs

differentially methylated sites

ROC

receiver operating characteristic

AUC

area under the curve

PPI

Protein-protein interaction

TMB

tumor mutational burden

GDSC

The Genomics of Drug Sensitivity in Cancer

hazard ratio

5mC

5-methylation

IC₅₀

half inhibitory centration

C-index

concordance index

lncRNAs

long non-coding RNAs

ceRNAs

endogenous RNAs

miRNAs

microRNAs

MIC

Maximal Information Coefficient

TNBC

triple-negative breast cancer

TPR

true positive rate

FPR

false positive rate

LRisk

low-risk group

HRisk

high-risk group

female

male

Non-M

non-metastasic breast cancer

breast cancer with metastasis

Ethics approval and consent to participate

Our study was approved by the Ethics Committee of China Medical University

Consent for publication

Not applicable

Availability of data and materials

The database supporting the conclusions of this article is available in the TCGA database and GEO database.

Competing interests

The authors declare that they have no competing interests

Funding

This work was supported by grants from National Natural Science Foundation of China under Grant 82073282

Authors' contributions

SLW and XMY downloaded and collected the data. SLW contributed to data analysis and drafting the manuscript. SLW, XMY, XYM and FJ contributed to the study design and data interpretation. All the authors contributed to the article and approved the manuscript.

Acknowledgements

Not applicable

Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, Bray F: Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin 2021.
Klutstein M, Nejman D, Greenfield R, Cedar H: DNA Methylation in Cancer and Aging. Cancer Res 2016, 76(12):3446–3450.
Gyorffy B, Bottai G, Fleischer T, Munkacsy G, Budczies J, Paladini L, Borresen-Dale AL, Kristensen VN, Santarpia L: Aberrant DNA methylation impacts gene expression and prognosis in breast cancer subtypes. Int J Cancer 2016, 138(1):87–97.
Szyf M: DNA methylation signatures for breast cancer classification and prognosis. Genome Med 2012, 4(3):26.
Widschwendter M, Jones PA: DNA methylation and breast carcinogenesis. Oncogene 2002, 21(35):5462–5482.
Ennour-Idrissi K, Dragic D, Durocher F, Diorio C: Epigenome-wide DNA methylation and risk of breast cancer: a systematic review. BMC Cancer 2020, 20(1):1048.
Tao C, Luo R, Song J, Zhang W, Ran L: A seven-DNA methylation signature as a novel prognostic biomarker in breast cancer. J Cell Biochem 2020, 121(3):2385–2393.
Guo Y, Mao X, Qiao Z, Chen B, Jin F: A Novel Promoter CpG-Based Signature for Long-Term Survival Prediction of Breast Cancer Patients. Front Oncol 2020, 10:579692.
Zhong X, Zhong G: Prognostic biomarker identification and tumor classification in breast cancer patients by methylation and transcriptome analysis. FEBS Open Bio 2021.
Gao Y, Wang X, Li S, Zhang Z, Li X, Lin F: Identification of a DNA Methylation-Based Prognostic Signature for Patients with Triple-Negative Breast Cancer. Med Sci Monit 2021, 27:e930025.
Zhang M, Wang Y, Wang Y, Jiang L, Li X, Gao H, Wei M, Zhao L: Integrative Analysis of DNA Methylation and Gene Expression to Determine Specific Diagnostic Biomarkers and Prognostic Biomarkers of Breast Cancer. Front Cell Dev Biol 2020, 8:529386.
McMahon KW, Karunasena E, Ahuja N: The Roles of DNA Methylation in the Stages of Cancer. Cancer J 2017, 23(5):257–261.
Fu B, Du C, Wu Z, Li M, Zhao Y, Liu X, Wu H, Wei M: Analysis of DNA methylation-driven genes for predicting the prognosis of patients with colorectal cancer. Aging (Albany NY) 2020, 12(22):22814–22839.
Lu T, Chen D, Wang Y, Sun X, Li S, Miao S, Wo Y, Dong Y, Leng X, Du W et al: Identification of DNA methylation-driven genes in esophageal squamous cell carcinoma: a study based on The Cancer Genome Atlas. Cancer Cell Int 2019, 19:52.
Saha SK, Kim K, Yang GM, Choi HY, Cho SG: Cytokeratin 19 (KRT19) has a Role in the Reprogramming of Cancer Stem Cell-Like Cells to Less Aggressive and More Drug-Sensitive Cells. Int J Mol Sci 2018, 19(5).
Hung CS, Wang SC, Yen YT, Lee TH, Wen WC, Lin RK: Hypermethylation of CCND2 in Lung and Breast Cancer Is a Potential Biomarker and Drug Target. Int J Mol Sci 2018, 19(10).
Clemenceau A, Diorio C, Durocher F: Role of Secreted Frizzled-Related Protein 1 in Early Mammary Gland Tumorigenesis and Its Regulation in Breast Microenvironment. Cells 2020, 9(1).
Lu J, Wen M, Huang Y, He X, Wang Y, Wu Q, Li Z, Castellanos-Martin A, Abad M, Cruz-Hernandez JJ et al: C2ORF40 suppresses breast cancer cell proliferation and invasion through modulating expression of M phase cell cycle genes. Epigenetics 2013, 8(6):571–583.
Kloten V, Schlensog M, Eschenbruch J, Gasthaus J, Tiedemann J, Mijnes J, Heide T, Braunschweig T, Knuchel R, Dahl E: Abundant NDRG2 Expression Is Associated with Aggressiveness and Unfavorable Patients' Outcome in Basal-Like Breast Cancer. PLoS One 2016, 11(7):e0159073.
Zhu MMT, Burugu S, Gao D, Yu J, Kos Z, Leung S, Horst BA, Nielsen TO: Evaluation of glucocorticoid-induced TNF receptor (GITR) expression in breast cancer and across multiple tumor types. Mod Pathol 2020, 33(9):1753–1763.
Zhang Z, Lin M, Wang J, Yang F, Yang P, Liu Y, Chen Z, Zheng Y: Calycosin inhibits breast cancer cell migration and invasion by suppressing EMT via BATF/TGF-beta1. Aging (Albany NY) 2021, 13(12):16009–16023.
Westrich JA, Vermeer DW, Colbert PL, Spanos WC, Pyeon D: The multifarious roles of the chemokine CXCL14 in cancer progression and immune responses. Mol Carcinog 2020, 59(7):794–806.
Bhandari A, Shen Y, Sindan N, Xia E, Gautam B, Lv S, Zhang X: MAL2 promotes proliferation, migration, and invasion through regulating epithelial-mesenchymal transition in breast cancer cell lines. Biochem Biophys Res Commun 2018, 504(2):434–439.
Chua YL, Ito Y, Pole JC, Newman S, Chin SF, Stein RC, Ellis IO, Caldas C, O'Hare MJ, Murrell A et al: The NRG1 gene is frequently silenced by methylation in breast cancers and is a strong candidate for the 8p tumour suppressor gene. Oncogene 2009, 28(46):4041–4052.
Koptyra M, Gupta S, Talati P, Nevalainen MT: Signal transducer and activator of transcription 5a/b: biomarker and therapeutic target in prostate and breast cancer. Int J Biochem Cell Biol 2011, 43(10):1417–1421.
Jeong YJ, Oh HK, Choi HR, Park SH: Methylation of the NT5E Gene Is Associated with Poor Prognostic Factors in Breast Cancer. Diagnostics (Basel) 2020, 10(11).
Lo Nigro C, Monteverde M, Lee S, Lattanzio L, Vivenza D, Comino A, Syed N, McHugh A, Wang H, Proby C et al: NT5E CpG island methylation is a favourable breast cancer biomarker. Br J Cancer 2012, 107(1):75–83.
Zhu W, Qin W, Hewett JE, Sauter ER: Quantitative evaluation of DNA hypermethylation in malignant and benign breast tissue and fluids. Int J Cancer 2010, 126(2):474–482.
Fang Y, Wang L, Wan C, Sun Y, Van der Jeught K, Zhou Z, Dong T, So KM, Yu T, Li Y et al: MAL2 drives immune evasion in breast cancer by suppressing tumor antigen presentation. J Clin Invest 2021, 131(1).
Zhou Y, Heitmann JS, Clar KL, Kropp KN, Hinterleitner M, Engler T, Koch A, Hartkopf AD, Zender L, Salih HR et al: Platelet-expressed immune checkpoint regulator GITRL in breast cancer. Cancer Immunol Immunother 2021, 70(9):2483–2496.
Qiu SQ, Waaijer SJH, Zwager MC, de Vries EGE, van der Vegt B, Schroder CP: Tumor-associated macrophages in breast cancer: Innocent bystander or important player? Cancer Treat Rev 2018, 70:178–189.
Fu C, Jiang A: Dendritic Cells and CD8 T Cell Immunity in Tumor Microenvironment. Front Immunol 2018, 9:3059.
Tsou P, Katayama H, Ostrin EJ, Hanash SM: The Emerging Role of B Cells in Tumor Immunity. Cancer Res 2016, 76(19):5597–5601.
Li L, Li M, Wang X: Cancer type-dependent correlations between TP53 mutations and antitumor immunity. DNA Repair (Amst) 2020, 88:102785.
Duffy MJ, Synnott NC, Crown J: Mutant p53 in breast cancer: potential as a therapeutic target and biomarker. Breast Cancer Res Treat 2018, 170(2):213–219.
Martinez-Saez O, Chic N, Pascual T, Adamo B, Vidal M, Gonzalez-Farre B, Sanfeliu E, Schettini F, Conte B, Braso-Maristany F et al: Frequency and spectrum of PIK3CA somatic mutations in breast cancer. Breast Cancer Res 2020, 22(1):45.
Miricescu D, Totan A, Stanescu S, II, Badoiu SC, Stefani C, Greabu M: PI3K/AKT/mTOR Signaling Pathway in Breast Cancer: From Molecular Landscape to Clinical Aspects. Int J Mol Sci 2020, 22(1).
Wang P, Chen Y, Wang C: Beyond Tumor Mutation Burden: Tumor Neoantigen Burden as a Biomarker for Immunotherapy and Other Types of Therapy. Front Oncol 2021, 11:672677.
Wang J, Su Z, Lu S, Fu W, Liu Z, Jiang X, Tai S: LncRNA HOXA-AS2 and its molecular mechanisms in human cancer. Clin Chim Acta 2018, 485:229–233.
Yuan RX, Bao D, Zhang Y: Linc00707 promotes cell proliferation, invasion, and migration via the miR-30c/CTHRC1 regulatory loop in breast cancer. Eur Rev Med Pharmacol Sci 2020, 24(9):4863–4872.
Du S, Hu W, Zhao Y, Zhou H, Wen W, Xu M, Zhao P, Liu K: Long non-coding RNA MAGI2-AS3 inhibits breast cancer cell migration and invasion via sponging microRNA-374a. Cancer Biomark 2019, 24(3):269–277.
Yang Y, Yang H, Xu M, Zhang H, Sun M, Mu P, Dong T, Du S, Liu K: Long non-coding RNA (lncRNA) MAGI2-AS3 inhibits breast cancer cell growth by targeting the Fas/FasL signalling pathway. Hum Cell 2018, 31(3):232–241.
Deva Magendhra Rao AK, Patel K, Korivi Jyothiraj S, Meenakumari B, Sundersingh S, Sridevi V, Rajkumar T, Pandey A, Chatterjee A, Gowda H et al: Identification of lncRNAs associated with early-stage breast cancer and their prognostic implications. Mol Oncol 2019, 13(6):1342–1355.
Li D, Li L, Cao Y, Chen X: Downregulation of LINC01140 is associated with adverse features of breast cancer. Oncol Lett 2020, 19(2):1157–1164.

Table 1 is available in the Supplementary Files section.

No competing interests reported.

Table1.jpg
Table 1. Characteristics of the patients in TCGA training and TCGA testing groups.
FigureS1.jpg
Figure S1. protein-protein interaction network (PPI) of the 288 candidate MDGs.
FigureS2.jpg
Figure S2. The landscape of immune infiltration in breast cancer.

Download PDF

Version 1

posted

You are reading this latest preprint version

Identification and validation of a novel prognostic signature with DNA methylation-driven genes for survival prediction of breast cancer

Status:

Version 1

Abstract

Background:

Methods:

Results:

Conclusions:

Figures

Introduction

Material And Methods

1. Data acquisition and procession

2. Construction and validation of the risk signature and nomogram

3. Bioinformatic analyses

3.1 GO and KEGG analysis

3.2 PPI analysis

3.3 GSVA analysis

3.4 TISIDB analysis

3.5 Investigation of tumor infiltrating immune cells

3.6 Investigation of tumor mutation and neoantigen

3.7 Exploration Of Potentially Effective Therapeutic Drugs

3.8 Construction Of Cerna Network

4. Statistical analyses

Results

1. Identification and function analysis of Methylation-driven genes in breast cancer

2. Establishment of a risk signature with 19 methylation-driven genes for predicting the overall survival of breast cancer

3. Investigation of the potential molecular regulators in patients with different risks

4. Exploring the potential therapeutic agents in clinical treatment based on the risk model

5. Validation of the prognostic 19-gene risk signature in external cohorts

6. Construction of the prognostic nomogram with the 19-gene risk signature for survival analysis of breast cancer

7. Construction of the lncRNA-miRNA-mRNA ceRNA network

Discussion

List Of Abbreviations

Declarations

References

Tables

Additional Declarations

Supplementary Files

Status:

Version 1