Dissecting the Role of N6-Methylandenosine-Related Long Non-coding RNAs Signature in Prognosis and Immune Microenvironment of Breast Cancer
Jinguo Zhang
Benjie Shan - Department of Medical Oncology, The First Affiliated Hospital of USTC, Division of Life Science and Medicine, University of Science and Technology of China, Hefei, China
Breast cancer (BC) represents a molecularly and clinically heterogeneous disease. Recent progress in immunotherapy has provided a glimmer of hope for several BC subtypes. The relationship between N6-methyladenosine (m6A) modification and long non-coding RNAs (LncRNAs) is still largely unexplored in BC. Here, with the intention to dissect the landscape of m6A-related lncRNAs and explore the immunotherapeutic value of the m6A-related lncRNA signature, we identified m6A-related lncRNAs by co-expression analysis from The Cancer Genome Atlas (TCGA) and stratified BC patients into different subgroups. Furthermore, we generated an m6A-related lncRNA prognostic signature. Four molecular subtypes were identified by consensus clustering. Cluster 3 preferentially had favorable prognosis, upregulated immune checkpoint expression, and high level of immune cell infiltration. Twenty-one m6A-related lncRNAs were applied to construct the m6A-related lncRNA model (m6A-LncRM). Survival analysis and receiver operating characteristic (ROC) curves further confirmed the prognostic value and prediction performance of m6A-LncRM. Finally, high- and low-risk BC subgroups displayed significantly different clinical features and immune cell infiltration status. Overall, our study systematically explored the prognostic value of the m6A-related LncRNAs and identified a high immunogenicity BC subtype. The proposed m6A-related LncRNA model might serve as a robust prognostic signature and attractive immunotherapeutic targets for BC treatment.
Introduction
The Global Cancer Statistics reported that female breast cancer (BC) surpassed lung cancer as the most diagnosed cancer, with an estimated 2.3 million new cases (11.7%) (Sung et al., 2021). BC is generally regarded as a heterogeneous disease in terms of its molecular features, histological composition, and clinical characteristics (McCart Reed et al., 2021; Sadeghalvad et al., 2021). With the evolution of high-throughput technologies, we tend to subtype the disease into clinically relevant molecular subtypes, including normal-like, luminal A and B, HER2-enriched, and basal-like or more intrinsic subtypes (Ochoa et al., 2020; Morgan et al., 2021). The five clinically molecular subtypes have drastically different treatment selection, prognosis, and tumor biology (Zubair et al., 2020). The purpose of BC molecular subtyping is to design personalized treatment strategies for patients based on the emerging evidence of BC classification and treatment responses (Pashayan et al., 2020). However, the subtyping of BC remains unexplored and challenging, which manifests as a heterogeneity of therapeutic responses and prognosis within the same clinical subtypes.
N6-Methyladenosine (m6A) is the most abundant internal epigenetic modification occurring in RNA molecules (Wang et al., 2020). The discovery of m6A RNA modification has added a new layer of regulatory mechanism controlling gene expression (Nombela et al., 2021). Increasing evidence implicates that m6A regulators (“writers,” “erasers,” and “readers”) exert an essential role in multiple types of cancer by regulation of the “epi-transcriptome” in cancers (Huang et al., 2020; Liu et al., 2020). Remarkably, both m6A “writers” and “erasers” regulators are abnormally overexpressed and perform an oncogenic role in BC (Wei et al., 2020). The long non-coding RNAs (lncRNAs), a large class of conserved endogenous RNAs, are characterized as transcripts longer than 200 nucleotides and no protein-coding potential (Statello et al., 2021). m6A modification has been observed not only in messenger RNA but also in lncRNAs, affecting the fate of the modified RNA molecules (Lence et al., 2019). The antibody-based approach m6A-seq (MeRIP-seq) analysis has revealed that lncRNAs were also involved in such modification (He et al., 2020). Currently, m6A modification of long non-coding RNAs is still largely unexplored in BC.
N6-methyladenosine modification on RNA was involved in various stages of the RNA life cycle, including RNA transcription, processing, splicing, degradation, and RNA translation (Wei et al., 2020). m6A regulators act as either oncogenes or tumor suppressors to mediate the development and progression of BC. Aberrant expression of METTL3 and METTL14 has been reported to be associated with BC cell proliferation, migration, and metastasis (Wu et al., 2019; Shi et al., 2020). Demethylase FTO was highly expressed in BC, and upregulation of FTO enhanced the aggressiveness properties of BC, especially for HER2-overexpressing BC (Tan et al., 2015). Another m6A demethylase, ALKBH5, was recently reported to play a vital role in tumor formation and self-renewal of BC stem cells (Zhang et al., 2016a,b). The oncogenic role of “readers” proteins such as YTHDF1, YTHDF1, HNRNPC, HNRNPA2B1, and IGF2BP2 has been previously investigated in BC (Zheng et al., 2021). Regarding the association of m6A regulators and clinical features of BC, m6A regulators have been reported to significantly correlate with the clinicopathological characteristics, survival outcomes, and antitumor immune response in BC (He et al., 2021). Presently, increased attention was being paid to lncRNA signature models for its prognostic and predictive potential in cancer treatment (Li et al., 2019; Xu et al., 2021). Thus, dissecting the role of m6A-related lncRNA signature in BC may help researchers identify significant potential biomarkers for clinical applications. The emergence of immune checkpoint inhibitors (ICIs) has contributed to a revolutionary shift for traditional cancer treatment. Recent preclinical and clinical trials have supported that ICIs are also promising approaches in treating BC, particularly for triple-negative breast cancer (TNBC). PD-L1 checkpoint inhibitor atezolizumab in combination with nanoparticle albumin-bound (nab)-paclitaxel was approved in PD-L1-positive advanced or metastatic TNBC (Schmid et al., 2018). In the KEYNOTE-355 study, progression-free survival benefit was observed in pembrolizumab combined with chemotherapy compared with chemotherapy alone (Cortes et al., 2020). Depending on the results of the KEYNOTE-355 study, pembrolizumab combined with chemotherapy is currently FDA-approved for unresectable locally advanced or metastatic PD-L1-positive TNBC. In the present study, we aimed to depict the landscape of m6A-related lncRNAs and explore the immunotherapeutic value of the m6A-related lncRNA signature in the TCGA breast cancer (TCGA-BRCA) cohort. Furthermore, an m6A-related lncRNA prognostic model was constructed using least absolute shrinkage and selection operator (LASSO) regression. The role of immune infiltration and microenvironment heterogeneity in different BC subtypes and high-risk and low-risk subgroups was also investigated.
Materials and Methods
Data Processing of the Cancer Genome Atlas-Breast Cancer Dataset
A total of 1,109 BC patients and 113 normal paired tissues from the TCGA-BRCA (The Cancer Genome Atlas-Breast Cancer) program were enrolled in our study. The public RNA sequencing and clinical information of patients were obtained from TCGA1. The expression of 24 m6A-related genes was extracted from the TCGA-BRCA, including expression data on writers (METTL3, METTL14, METTL16, WTAP, VIRMA or KIA1499, RBM15, RBM15B, and ZC3H13), readers (YTHDC1, YTHDC2, YTHDF1, YTHDF2, YTHDF3, HNRNPC, HNRNPA2B1, IGF2BP1, IGF2BP2, IGF2BP3, FMR1, LRPPRC, and RBMX), and erasers (FTO and ALKBH5). The lncRNAs in the TCGA dataset were identified based on the annotation of Genome Reference Consortium Human Build 38 (GRCh38), and 14,086 lncRNAs were evaluated in the TCGA-BRCA transcriptome matrix. In this study, lncRNAs we identified were composed of eight types of transcripts (lincRNA, sense overlapping, retained intron, antisense, processed transcript, sense intronic, 3′ overlapping ncRNA, and non-coding and macro lncRNA). In this study, BC samples without complete survival data were excluded from all analyses. The clinical features of TCGA-BRCA samples are presented in Supplementary Table 1.
Identification of N6-Methyladenosine-Related Long Non-coding RNAs
The Spearman correlation analysis was applied to screen m6A-related lncRNAs with the criteria of correlation coefficient > 0.4 and p < 0.001 in the TCGA-BRCA dataset. The network of m6A-related genes and m6A-related lncRNAs was visualized by Cytoscape software 3.5.1. To explore the prognostic value of m6A-related lncRNAs, univariate Cox regression analysis was conducted to select the m6A-related lncRNAs that were associated with patients’ overall survival (OS). The HR value and Cox p value were calculated. The Wilcoxon signed-rank test was used to compare the expression of m6A-related lncRNAs in BRCA tissues versus normal tissues.
Breast Cancer Subtype Defined by N6-Methyladenosine-Related Long Non-coding RNAs
To systematically assess the roles and functions of m6A-related lncRNAs, cluster analysis was performed to subgroup BC patients into different groups using the ConsensusClusterPlus package of R software (Wilkerson and Hayes, 2010). Based on the expression level of prognostic-related lncRNAs, the consistent clustering algorithm was applied to determine the clustering number of samples. Then survival analysis was conducted to explore the prognostic value of different clusters. Comparison of clinical traits and different clusters was evaluated by the Chi-square test. The expression of checkpoint members including PDCD1 (PD-1), CD274 (PD-L1), CTLA4, LAG3, ICOS, and IDO1 was estimated in different BRCA clusters using the “limma” package.
Construction of a Risk Model Based on N6-Methyladenosine-Related Long Non-coding RNAs
Based on univariate Cox regression, LASSO regression analysis was analyzed to explore the association between m6A-prognosis lncRNAs and BC risk. The TCGA BRCA patients were randomly divided into training cohort and validation cohort. Twenty-one m6A-related lncRNAs were selected to construct the best risk score model with the “glmnet” package. We calculated the riskScore of each sample with the formula: RiskScore = ExpressionLncRNA1 × CoefficientLncRNA1 + ExpressionLncRNA2 × CoefficientmLncRNA2 + … ExpressionLncRNAn × CoefficientLncRNAn. The samples were divided into high-risk and low-risk groups with the median value of the riskScore. Subsequently, Kaplan–Meier curves and receiver operating characteristic (ROC) curves of 1 year using “timeROC” package were adopted to evaluate the predictive accuracy of the risk model. In addition, univariate Cox regression analysis and multivariate Cox regression analysis were performed to assess the independence of the risk model by regarding risk score as a single characteristic factor. To evaluate the applicability of the model, we combined the clinical traits with riskScore to plot the survival curves.
Evaluation of Immune Infiltration
To specifically analyze the immune infiltration, the CIBERSORT analytical tool was adopted to identify the difference of 22 types of immune cells in different BRCA clusters or high- and low-risk groups using the Wilcoxon rank-sum test (Newman et al., 2015). A comparison of the fraction of 22 immune cells in BRCA clusters was produced with “limma” packages in R. We performed ESTIMATE (Estimation of Stromal and Immune cells in MAlignant Tumors using Expression data) analysis to evaluate the level of immune cell infiltration, stromal level, and estimate score for each BRCA sample (Yoshihara et al., 2013). The immune, stromal, and estimate scores in different BRCA clusters were compared.
Gene Set Enrichment Analysis
In our analysis, GSEA was used to determine the biological pathway exchanges in different BRCA clusters. The transcriptome expression data of TCGA-BRCA were analyzed. Gene sets with nominal p value < 0.05 were considered significant enrichments. Kyoto Encyclopedia of Genes and Genomes (KEGG) and Gene Ontology (GO) gene sets were performed with GSEA 4.1.0 tool.
Statistical Analysis
All analyses were performed using R software (version 4.0.2). The Wilcoxon signed-rank test was adopted for two-group comparison, and the Kruskal–Wallis test was applied for comparison of more than two groups. With the best cutoff value for subgroup stratification, Kaplan–Meier curves were plotted to compare the OS among various subgroups and the log-rank p value was calculated. Univariate and multivariate Cox regression analyses were conducted to evaluate the independent prognostic value of the risk model based on m6A-related lncRNAs. p-value < 0.05 was considered statistically significant.
Results
Selection of N6-Methyladenosine-Related Long Non-coding RNAs in the Cancer Genome Atlas-Breast Cancer Cohort
We downloaded the RNA-Seq data of the TCGA-BRCA dataset and identified 14,086 lncRNAs using GRCh38 annotation. The expression of 24 m6A-related regulators was extracted, and then Spearman correlation analysis was performed. Using a threshold of correlation coefficient > 0.4 and p < 0.001, 491 lncRNAs were significantly correlated with m6A-related regulators and selected for further analysis. To directly visualize the correlation between lncRNAs and m6A-related regulators, a gene network was plotted by Cytoscape. The co-expression network of m6A-related regulators and lncRNAs showed that 14 m6A-related regulators (RBM15, RBMX, YTHDC1, YTHDC2, YTHDF1, YTHDF3, METTL3, METTL16, METTL14, IGF2BP2, IGF2BP3, HNRNRA2B1, VIRMA, and FTO) were significantly related to some lncRNAs, of which RBM15 was the hub regulator (Figure 1A). To further investigate the clinical significance of m6A-related lncRNAs, univariate Cox regression was performed. The forest plot displayed that 51 m6A-related lncRNAs were significantly associated with the OS of BC patients (Figure 1B). A heatmap plot compared the expression of 51 m6A-related lncRNAs for TCGA BC tissues versus normal tissues, and we found that all prognostic lncRNAs were abnormally expressed in BC (Figure 1C).
Figure 1. Identification of prognostic m6A-related lncRNAs in BC patients. (A) Construction of the m6A-regulator–lncRNA co-expression network using Cytoscape. Within networks, m6A-related genes are displayed in red nodes and m6A-related lncRNAs in green nodes. (B) A forest map showed 51 prognostic m6A-related lncRNAs identified by univariate Cox proportional hazard regression. (C) A heatmap displayed differentially expressed 51 m6A-related lncRNAs in the TCGA-BRCA dataset using the “limma” method. *p < 0.05; **p < 0.01; and ***p < 0.001.
Identification of N6-Methyladenosine-Related Long Non-coding RNAs and Immune Subtypes of Breast Cancer
Molecular subtyping classification provides a fundamental basis for finding the optimal treatment for a particular patient. Therefore, we constructed consensus clustering analysis based on the expression profiles of 51 m6A-related prognostic lncRNAs. Intrigued by cumulative distribution, we chose k = 4 to subgroup BC patients where the sample cluster remained stable and robust (Figure 2A and Supplementary Figure 1). Kaplan–Meier analysis of OS revealed that Cluster 3 had the best prognosis, whereas Cluster 4 had the worst survival outcome (Figure 2B). We also explored the association of clinical characteristics with different BC clusters. The m6A-related lncRNA subgroup was only correlated with age at diagnosis (p < 0.05) (Supplementary Figure 2). Recently, immune therapies have evolved into a promising strategy for solid tumors. Therefore, we further estimated the expression of immune checkpoints in different BC clusters. As shown in Figure 2C, Cluster 3 patients displayed a significantly higher expression of PD-L1 compared to Cluster 1, Cluster 2, and Cluster 4. As expected, the same expression patterns of other immune checkpoint members including CTLA-4, PD-1, LAG3, and IDO1 in Cluster 3 were also observed (Figures 2D–H). Then, the correlation between PD-L1 and 51 m6A-related prognostic lncRNAs was analyzed. We observed a significant correlation between PD-L1 and 16 m6A-related lncRNAs (AL359076.1, ECE1-AS1, AL365330.1, MIATNB, AC005104.1, U73166.1, AP000919.3, AC018926.2, AC004846.2, AL592301.1, AL138789.1, AL513218.1, AC084876.1, LRRC8C-DT, AL513190.1, and AC012442.2) (Figure 2I). Moreover, the correlation of other immune checkpoint members and 51 m6A-related prognostic lncRNAs is presented in Supplementary Figure 3. The above findings indicated that BC patients with cluster 3 showed higher immunogenicity and might respond to the ICI therapy.
Figure 2. Consensus clustering of m6A-related lncRNAs and immune checkpoint member expression in different BC subtypes. (A) Consensus clustering matrix for k = 4. (B) Kaplan–Meier curves of OS for four clusters in TCGA-BRCA. (C–H) The expression of immune checkpoints PD-L1, PD-1, CTLA-4, ICOS, LAG3, and IDO1 in different BC subtypes. (I) Association of PD-L1 and 51 prognostic m6A-related lncRNAs was determined using Spearman correlation analysis. ns, no significance; *p < 0.05; **p < 0.01; and ***p < 0.001.
Immune Landscape of N6-Methyladenosine-Related Long Non-coding RNAs Subtypes of Breast Cancer
To thoroughly explore immune infiltration status in different BC clusters, CIBERSORT algorithms were used to assess the composition of the immune microenvironment. BC patients in cluster 3 tended to have a higher immune score and estimate score and a lower stroma score (Figures 3A–C). We further calculated individual immune infiltration in different BC clusters using the ESTIMATE program. Relative to other BC clusters, higher levels of activated CD4 memory T cells and follicular helper T cells were found in cluster 3 BC samples (Figures 3D,E). Cluster 3 BC samples also showed a higher infiltration of CD8 T cells compared to cluster 1 and cluster 4 (Figure 3F). High infiltration of M1 macrophages and low infiltration of M0 macrophages and M2 macrophages were observed in cluster 3 BC samples (Figures 3G–I). In addition, cluster 3 also presented the lowest level of resting mast cells and the highest level of dendritic cells among different BC clusters (Figures 3J,K). The above results suggested that cluster 3 might contain immune-hot tumors.
Figure 3. Immune characteristics among four m6A-related lncRNA BC clusters. (A–C) Different expression of immune score (A), stromal score (B), and ESTIMATE score (C) in four BC clusters. (D–K) Differences in the levels of infiltration of CD4 memory T cells, follicular helper T cells, CD8 T cells, M1 macrophages, M0 macrophages, M2 macrophages, resting mast cells, and dendritic cells in four BC clusters.
Gene Set Enrichment Analysis
To annotate the potential biological process and pathways between cluster 3 and other clusters, GSEA analysis was computed using gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) signatures based on individual sample expression profiles. Our results showed that various immune response processes were significantly enriched in cluster 3 samples (Figure 4). GO signature revealed that MHC protein binding, natural cell-mediated immunity, T cell receptor signaling pathway, T cell activation, T help immune response, and so on were positively enriched in cluster 3. KEGG signature further revealed that cluster 3 was positively correlated with antigen processing and presentation and natural killer cell-mediated cytotoxicity.
Figure 4. Gene set enrichment analysis (GSEA) indicating that various immune response processes were enriched in the BC cluster 3 patients compared to other BC subtypes. The above results demonstrated that biological processes such as lymphocyte activation, innate immune response, natural killer cell-mediated immunity, MHC protein binding, T cell activation, antigen processing and presentation, and so on were positively enrichment in BC cluster 3.
Generation of an N6-Methyladenosine-Related Long Non-coding RNAs Model for Prognosis
To construct an m6A-related lncRNA model for predicting the OS of BC patients, LASSO regression was performed based on the 51 m6A-related prognostic lncRNAs in the TCGA BRCA dataset. Twenty-one lncRNAs were identified as powerful prognostic lncRNAs via the LASSO analysis (Figures 5A,B). We then randomly divided BC patients into two sets: the training cohort (N = 536) and testing cohort (N = 533). The risk score for each group was calculated using the coefficient of 21 lncRNAs (Figure 5C). The detailed information of the 21 m6A-related lncRNAs for constructing the prognostic signature was presented in Table 1. According to the median risk score, each cohort was defined into high-risk group and low-risk group. Kaplan–Meier curves displayed that BC patients with higher risk scores had shorter OS than those with lower risk scores in the training cohort and testing cohort (Figures 5D,E). Moreover, ROC curves further demonstrated that m6A-LncRM had a promising prediction performance in the training cohort and testing cohort (Figures 5F,G). The ROC curves for 3 or 5 years were presented in Supplementary Figure 4. The risk score distribution and OS status for the training cohort was shown in Figure 5H. The scatter plots displayed the survival status of the training cohort and indicated that the death number of BC patients was increased with increasing risk score. Similar results were also observed in the testing cohort (Figure 5I). The expression heatmap of the 21 lncRNAs is shown in Supplementary Figure 5.
Figure 5. Risk model from m6A-related lncRNAs. (A) Elucidation for LASSO coefficient profiles of prognostic m6A-related lncRNAs. (B) The least absolute shrinkage was performed and selection operator (LASSO) regression model for OS. (C) Bar plot displayed the coefficients of 21 selected m6A-related lncRNAs. (D,E) Kaplan–Meier curves for OS of TCGA-BRCA patients based on the risk stratification in the training cohort (D) and test cohort (E). (F,G) ROC analysis for OS prediction of TCGA-BRCA patients in the training cohort (F) and test cohort (G). (H,I) Distribution of risk score and OS status in the TCGA training cohort (H) and TCGA test cohort (I).
Table 1. The detailed information of the 21 m6A-related lncRNAs used to construct the prognostic signature.
Prognostic Value of the N6-Methyladenosine-Related Long Non-coding RNAs Model Signature
To evaluate whether the effect of m6A-LncRM on survival outcome was an independent risk factor for BC patients, we conducted univariate and multivariate Cox analyses for variables including age, tumor size, tumor stage, node status, and m6A-LncRM risk score. Univariable Cox regression analysis showed that m6A-LncRM was significantly correlated with OS of patients with BC [training cohort: HR = 1.626, 95% CI (1.343, 1.968), p < 0.001, Figure 6A; test cohort: HR = 2.967, 95% CI (1.413, 6.228), p = 0.004, Figure 6B]. Multivariate Cox analysis further indicated that m6A-LncRM remained an independent prognostic factor [training cohort: HR = 1.591, 95% CI (1.251, 2.022), p < 0.001, Figure 6C test cohort: HR = 2.575, 95% CI (1.263, 5.251), p = 0.009, Figure 6D]. In addition, clinical stratification analyses were performed after being adjusted by clinical factors including age, tumor size, node status, and tumor stage, and the results indicated that the high-risk group had worse OS than the low-risk group across all clinically stratified subgroups (Figures 6E–L).
Figure 6. Evaluation of the prognostic role of m6ALncM signature. (A,B) The univariate Cox forest map of risk model score and clinical features in the TCGA training cohort (A) and TCGA test cohort (B). (C,D) The multivariate Cox forest plot of risk model score and clinical characteristics in the training cohort (C) and test cohort (D). (E–L) The m6A-LncRM signature retained its prognostic value in multiple subgroups of BC patients. Survival analysis in low- and high-risk groups adjusted by clinical factors including age, tumor size, node status, and tumor stage.
Relationship of the N6-Methyladenosine-Related Long Non-coding RNAs Model Signature and Clinical Features
We further investigated whether clinical features were associated with the m6A-LncRM risk score. Overall, the heatmap displayed that BC clusters, immuneScore, patients’ age, and node status between high- and low-risk groups were significantly different (Figure 7A). More specifically, the risk scores of cluster 2 and cluster 3 patients were significantly decreased relative to those of cluster 1 and cluster 4 (Figure 7B). BC patients with high immuneScore harbored significantly low risk scores (Figure 7C). Compared with stage I–II BC patients, patients with stage III–IV had higher risk scores (Figure 7D). BC patients with N1–N3 metastasis had higher risk scores versus the N0 group (Figure 7E). Moreover, metastatic disease also tended to have a higher risk score compared with that for BC patients without metastasis (Figure 7F). However, patients’ age and tumor size had no impact on risk scores (Figures 7G,H).
Figure 7. Clinical evaluation by the m6A-LncRM signature. (A) A strip heatmap along with the scatter diagram showed that clinical stage, BC clusters, and ImmuneScore were significantly associated with the riskScore. (B–H) BC patients with different clinicopathological features including m6A-LncRNA clusters, ImmuneScore, tumor stage, node status, tumor size, patient age, and metastasis status had different levels of risk scores.
Correlations Between N6-Methyladenosine-Related Long Non-coding RNAs Model Signature and Immune Cell Infiltration
We then investigated the differences of immune cell infiltration between high-risk and low-risk groups. Firstly, the expression of immune checkpoints was estimated. As shown in Figure 8A, low-risk BC patients displayed a significantly higher expression of PD-L1 compared to high-risk BC patients. The expression trend was nearly identical for other immune checkpoint members including PD-1, CTLA4, ICOS, and LAG3 (Figures 8B–E). CIBERSORT analysis showed that the risk score exhibited a negative relation to the infiltration of three types of immune cells: CD8 T cells (r = −0.18, p < 0.001), naive B cells (r = −0.23, p < 0.001), and plasma cells (r = −0.09, p = 0.008) (Figures 8F–H). The infiltration level of M0 macrophages and M2 macrophages was positively correlated with the risk score (Figures 8I,J).
Figure 8. Relationship of immune cell infiltration and the m6A-LncRM signature. (A–E) The expression of immune checkpoint members including PD-L1, PD-1, CTLA-4, ICOS, and LAG3 in low-risk and high-risk BC groups. (F–J) m6ALncM signature negatively related to CD8 T cells, naive B cells, and plasma cells, while positive relationships were observed between m6ALncM and M0 macrophages and M2 macrophages.
Discussion
In the current study, we aimed to identify the prognostic value of m6A-related lncRNAs. The expression of 24 m6A-related genes and 14,086 lncRNAs was extracted from the TCGA-BRCA cohort. Spearman correlation analysis and univariate Cox regression selected 51 m6A-related lncRNAs for further analysis. Consensus clustering analysis revealed a BC subtype (Cluster 3) with best prognosis, and a high expression of immune checkpoints was observed in Cluster 3. CIBERSORT and ESTIMATE analyses further revealed that a high level of immune cell infiltration was found in Cluster 3. GSEA analysis showed that cluster 3 samples were positively enriched in various immune response processes. Furthermore, 21 m6A-related lncRNAs were selected to establish m6A-LncRM for predicting the OS of BC patients. Survival analysis and ROC curves confirmed a promising prognostic value and prediction performance of m6ALncRM. Multivariate Cox regression analysis showed that m6A-LncRM was able to serve as an independent risk factor. Based on the median risk score, the high-risk group had worse survival outcome than the low-risk group across all clinically stratified subgroups. Finally, low- and high-risk BC subgroups displayed significantly different clinical features and immune cell infiltration status. In the present study, we identified a high immunogenicity BC subtype and developed a novel signature m6A-LncRM based on m6A-related lncRNA profiles of BC, which might serve as a predictive biomarker. Furthermore, we also developed a scoring system termed as m6A-LncRM risk score to evaluate individual m6A-related lncRNA patterns, which might be useful in identifying BC patients with high immunogenicity. Our study aimed to investigate the role of m6A in BC from a systems biology point of view, but not focusing on the mechanism of m6A modification in lncRNAs.
Existing classification methods of BC subtypes are restricted to mRNA expression analysis (Zhang et al., 2020). The PAM50 intrinsic subtypes of BC, based on the gene expression pattern, represents one of the most widely used protein-coding gene panels in clinical practice (Parker et al., 2009). Protein-coding transcripts merely represent about 2% of the human genome, while the rest are non-coding RNAs without protein-coding potential (Zhang et al., 2013). Recent bioinformatic research has identified an immune-related lncRNA signature model to predict the risk of 5-year recurrence-free survival in BC (Lai et al., 2020). Other studies have depicted the lncRNA landscape using transcriptome microarrays in BC and identified deregulated lncRNA expression patterns across different molecular subtypes (Cedro-Tanda et al., 2020). However, the overall perspective of m6A regulators involved in the dysregulation of lncRNAs in BC is not fully understood. As the most extensively studied RNA modification, m6A has been recognized as an epigenetic regulator affecting RNA splicing, stability, and translation. Currently, the field of lncRNA m6A modifications in BC is still in its relative infancy; thus, investigation of m6A-related lncRNAs from a “big data” perspective is warranted. LncRNAs could also recruit m6A regulators to exert their function. For instance, LINC00942 directly interacted with METTL14 protein by recognizing specific sequences, supporting its posttranscriptional m6A methylation modification of downstream targets in BC (Sun et al., 2020). Hypoxia-induced lncRNA KB-1980E6.3 was reported to recruit “readers” protein IGF2BP1 that enhanced c-Myc mRNA stability in BC (Zhu et al., 2021). The m6A modifications in lncRNAs were also widely observed. Recent findings indicated that METTL3-mediated m6A modification induced the high expression of LINC00958 by facilitating its RNA transcript stability in BC (Rong et al., 2021). The m6A modification of lncRNA Pvt1 was reported to participate in skin tissue homeostasis and wound repair (Lee et al., 2021). MeRIP-seq analysis detected multiple m6A sites of MALAT1, and its m6A residues could recruit YTHDC1 to nuclear speckles (Wang et al., 2021). In this study, we identified m6A-related lncRNAs using the correlation analysis and found that RBM15 had a complex network with numerous lncRNAs. The TCGA database and immunohistochemistry staining analysis indicated that RBM15 expression was upregulated in BC specimens (Liu et al., 2019). Up until now, no studies have yet analyzed its function and m6A-dependent mechanism in BC.
Since infiltration of various immune cells in the tumor microenvironment was observed in BC, BC was no longer considered as an immunological quiescent tumor type (Gil Del Alcazar et al., 2020). The advent of immunotherapy has attracted great attention of researchers and clinicians, which encompasses vaccines, oncolytic viruses, chimeric antigen receptor-modified T cells (CAR-T), relevant nanotechnology, and immune checkpoint blockade (Adams et al., 2019). Cancer vaccines have been under investigation in numerous ongoing clinical trials of BC such as dendritic cells, HER2-based vaccines, and peptide-targeting neoantigen vaccines (Hamilton et al., 2012; Fennemann et al., 2019). However, CAR-T cell therapy for BC is still in the early phase of clinical trials (Jin et al., 2020; Ferraro et al., 2021). Recent advances in nanotechnology allowed nanoparticles to co-deliver immunomodulatory agents that exhibited promising efficacy in BC (Li et al., 2021). Remarkably, the PD-1/PD-L1 monoclonal antibody was gradually coming to light and became a promising strategy for advanced TNBC patients (Schmid et al., 2020). Previously, BC was not believed to be an immunogenic tumor. Recent research progress has confirmed that some TNBC and HER2-positive subtypes were immunogenic (de Melo Gagliato et al., 2020). Numerous studies have reported that PD-L1 expression in breast cancer or inflammatory cells (Ghebeh et al., 2006; Joneja et al., 2017) and PD-L1 expression of cancer cells have long been suggested as a reliable indicator for anti-PD1/PD-L1 drugs (Reck et al., 2016). In our analysis, we identified a BC subtype that has a high expression abundance of immune checkpoint members, suggesting that this BC cluster might be suitable for immunotherapy approaches. To select most appropriate patients for ICI therapies, biomarkers that predict patients’ response to immunotherapy in BC are warranted. Infiltrating lymphocytes (TILs) have emerged as potential biomarkers of immunotherapy response. Increasing evidence revealed that the existence of TIL in the tumor microenvironment is strongly associated with better response rates for ICI therapies in “hot” or “inflamed” tumors (Rooney et al., 2015; Wein et al., 2017; Molinero et al., 2019). Apart from this, a high level of TIL infiltration also predicted a favorable response to neoadjuvant chemotherapy in BC (Wimberly et al., 2015). In this study, a higher level of TILs was found in cluster 3 BC samples relative to other BC clusters, which indicated that cluster 3 patients might benefit from ICI therapies. In concordance with immune infiltration, survival analysis also displayed a favorable outcome in cluster 3 patients. Finally, the GSEA results further suggested that immune-related biological processes were significantly enriched in cluster 3.
To improve the accuracy and efficacy on predicting outcome, we constructed a 21-m6A-related lncRNA model and found that this 21-lncRNA signature would independently predict OS in BC patients. Previous studies have attempted to identify an lncRNA-related panel for predicting the risk of tumor recurrence, neoadjuvant treatment response, and survival outcome (Guo et al., 2016; Wang et al., 2017; Li et al., 2018). Among the 21 candidate lncRNAs, only OTUD6B-AS1 and COL4A2-AS1 had been reported as prognostic factors in BC (Yao et al., 2019; Ma et al., 2020); the rest of lncRNAs were identified as prognostic signatures in BC for the first time. In the proposed model, we classified patients into low- and high-risk groups; we then reassessed the survival outcome, univariate and multivariate analyses, tumor immune infiltration, and expression of immune checkpoints. The low-risk group was associated with favorable survival, high level of TIL infiltration, high ImmuneScore, and high expression of immune checkpoint members, which indicated that this modeling algorithm functioned well.
Although m6A modifications of LncRNAs are still largely unexplored in BC, our analysis may provide the overall profile of m6A-related lncRNAs in BC. To our knowledge, this is the first report of the m6A-related lncRNA signature in BC. Undeniably, several drawbacks of our study should be noticed that needs to be further explored. Further in vitro and in vivo experiments are required to confirm the relationship of m6A-related regulators and lncRNA expressions. To explore the suitability and applicability of our model, more independent BC cohorts are needed to verify the prognostic value of the m6A-related lncRNAs model. Since the analyzed data were derived from databases, another limitation is that it is retrospective in nature. Additionally, the mechanism of m6A-related lncRNAs needs to be further investigated to improve the immunotherapy efficacy of BC.
Conclusion
This study systematically explored the prognostic value of the m6A-related lncRNAs in BC. A high immunogenicity BC subtype was identified that could be a potential candidate for antitumor immunotherapy. The 21-lncRNA model is a robust biomarker that can independently predict OS of TCGA-BRCA datasets and immunotherapy efficacy.
Data Availability Statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.
Author Contributions
JZ performed most of the analyses and wrote the manuscript. XH and JC designed the study. BS, LL, JD, QS, and QZ performed some of the analyses and edited the manuscript. All authors contributed to the article and approved the submitted version.
Funding
This work was supported by grants from the Key Science and Technology Program of Anhui province, China (1401042007) and Fundamental Research Funds for the Central Universities (WK9110000067).
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2021.711859/full#supplementary-material
Footnotes
References
Adams, S., Gatti-Mays, M. E., Kalinsky, K., Korde, L. A., Sharon, E., Amiri-Kordestani, L., et al. (2019). Current landscape of immunotherapy in breast cancer: a review. JAMA Oncol. 5, 1205–1214. doi: 10.1001/jamaoncol.2018.7147
Cedro-Tanda, A., Ríos-Romero, M., Romero-Córdoba, S., Cisneros-Villanueva, M., Rebollar-Vega, R. G., Alfaro-Ruiz, L. A., et al. (2020). A lncRNA landscape in breast cancer reveals a potential role for AC009283.1 in proliferation and apoptosis in HER2-enriched subtype. Sci. Rep. 10:13146. doi: 10.1038/s41598-020-69905-z
Cortes, J., Cescon, D. W., Rugo, H. S., Nowecki, Z., Im, S. A., Yusof, M. M., et al. (2020). Pembrolizumab plus chemotherapy versus placebo plus chemotherapy for previously untreated locally recurrent inoperable or metastatic triple-negative breast cancer (KEYNOTE-355): a randomised, placebo-controlled, double-blind, phase 3 clinical trial. Lancet 396, 1817–1828. doi: 10.1016/s0140-6736(20)32531-9
de Melo Gagliato, D., Buzaid, A. C., Perez-Garcia, J., and Cortes, J. (2020). Immunotherapy in breast cancer: current practice and clinical challenges. Biodrugs 34, 611–623. doi: 10.1007/s40259-020-00436-9
Fennemann, F. L., de Vries, I. J. M., Figdor, C. G., and Verdoes, M. (2019). Attacking tumors from all sides: personalized multiplex vaccines to tackle intratumor heterogeneity. Front. Immunol. 10:824. doi: 10.3389/fimmu.2019.00824
Ferraro, E., Drago, J. Z., and Modi, S. (2021). Implementing antibody-drug conjugates (ADCs) in HER2-positive breast cancer: state of the art and future directions. Breast Cancer Res. 23:84. doi: 10.1186/s13058-021-01459-y
Ghebeh, H., Mohammed, S., Al-Omair, A., Qattan, A., Lehe, C., Al-Qudaihi, G., et al. (2006). The B7-H1 (PD-L1) T lymphocyte-inhibitory molecule is expressed in breast cancer patients with infiltrating ductal carcinoma: correlation with important high-risk prognostic factors. Neoplasia 8, 190–198. doi: 10.1593/neo.05733
Gil Del Alcazar, C. R., Alečković, M., and Polyak, K. (2020). Immune escape during breast tumor progression. Cancer Immunol. Res. 8, 422–427. doi: 10.1158/2326-6066.cir-19-0786
Guo, W., Wang, Q., Zhan, Y., Chen, X., Yu, Q., Zhang, J., et al. (2016). Transcriptome sequencing uncovers a three-long noncoding RNA signature in predicting breast cancer survival. Sci. Rep. 6:27931. doi: 10.1038/srep27931
Hamilton, E., Blackwell, K., Hobeika, A. C., Clay, T. M., Broadwater, G., Ren, X. R., et al. (2012). Phase 1 clinical trial of HER2-specific immunotherapy with concomitant HER2 kinase inhibition [corrected]. J. Transl. Med. 10:28. doi: 10.1186/1479-5876-10-28
He, R. Z., Jiang, J., and Luo, D. X. (2020). The functions of N6-methyladenosine modification in lncRNAs. Genes Dis. 7, 598–605. doi: 10.1016/j.gendis.2020.03.005
He, X., Tan, L., Ni, J., and Shen, G. (2021). Expression pattern of m(6)A regulators is significantly correlated with malignancy and antitumor immune response of breast cancer. Cancer Gene Ther. 28, 188–196. doi: 10.1038/s41417-020-00208-1
Huang, H., Weng, H., and Chen, J. (2020). m(6)A modification in coding and non-coding RNAs: roles and therapeutic implications in cancer. Cancer Cell 37, 270–288. doi: 10.1016/j.ccell.2020.02.004
Jin, J., Cheng, J., Huang, M., Luo, H., and Zhou, J. (2020). Fueling chimeric antigen receptor T cells with cytokines. Am. J. Cancer Res. 10, 4038–4055.
Joneja, U., Vranic, S., Swensen, J., Feldman, R., Chen, W., Kimbrough, J., et al. (2017). Comprehensive profiling of metaplastic breast carcinomas reveals frequent overexpression of programmed death-ligand 1. J. Clin. Pathol. 70, 255–259. doi: 10.1136/jclinpath-2016-203874
Lai, J., Chen, B., Zhang, G., Li, X., Mok, H., and Liao, N. (2020). Molecular characterization of breast cancer: a potential novel immune-related lncRNAs signature. J. Transl. Med. 18:416. doi: 10.1186/s12967-020-02578-4
Lee, J., Wu, Y., Harada, B. T., Li, Y., Zhao, J., He, C., et al. (2021). N(6) -methyladenosine modification of lncRNA Pvt1 governs epidermal stemness. EMBO J. 40:e106276. doi: 10.15252/embj.2020106276
Lence, T., Paolantoni, C., Worpenberg, L., and Roignant, J. Y. (2019). Mechanistic insights into m(6)A RNA enzymes. Biochim. Biophys. Acta Gene Regul. Mech. 1862, 222–229. doi: 10.1016/j.bbagrm.2018.10.014
Li, H., Gao, C., Liu, L., Zhuang, J., Yang, J., Liu, C., et al. (2019). 7-lncRNA assessment model for monitoring and prognosis of breast cancer patients: based on Cox regression and co-expression analysis. Front. Oncol. 9:1348. doi: 10.3389/fonc.2019.01348
Li, J., Wang, W., Xia, P., Wan, L., Zhang, L., Yu, L., et al. (2018). Identification of a five-lncRNA signature for predicting the risk of tumor recurrence in patients with breast cancer. Int. J. Cancer 143, 2150–2160. doi: 10.1002/ijc.31573
Li, Y., Miao, W., He, D., Wang, S., Lou, J., Jiang, Y., et al. (2021). Recent progress on immunotherapy for breast cancer: tumor microenvironment, nanotechnology and more. Front. Bioeng. Biotechnol. 9:680315. doi: 10.3389/fbioe.2021.680315
Liu, L., Liu, X., Dong, Z., Li, J., Yu, Y., Chen, X., et al. (2019). N6-methyladenosine-related genomic targets are altered in breast cancer tissue and associated with poor survival. J. Cancer 10, 5447–5459. doi: 10.7150/jca.35053
Liu, L., Wang, Y., Wu, J., Liu, J., Qin, Z., and Fan, H. (2020). N(6)-Methyladenosine: a potential breakthrough for human cancer. Mol. Ther. Nucleic Acids 19, 804–813. doi: 10.1016/j.omtn.2019.12.013
Ma, W., Zhao, F., Yu, X., Guan, S., Suo, H., Tao, Z., et al. (2020). Immune-related lncRNAs as predictors of survival in breast cancer: a prognostic signature. J. Transl. Med. 18:442. doi: 10.1186/s12967-020-02522-6
McCart Reed, A. E., Kalinowski, L., Simpson, P. T., and Lakhani, S. R. (2021). Invasive lobular carcinoma of the breast: the increasing importance of this special subtype. Breast Cancer Res. 23:6. doi: 10.1186/s13058-020-01384-6
Molinero, L., Li, Y., Chang, C. W., Maund, S., Berg, M., Harrison, J., et al. (2019). Tumor immune microenvironment and genomic evolution in a patient with metastatic triple negative breast cancer and a complete response to atezolizumab. J. Immunother. Cancer 7:274. doi: 10.1186/s40425-019-0740-8
Morgan, A. J., Giannoudis, A., and Palmieri, C. (2021). The genomic landscape of breast cancer brain metastases: a systematic review. Lancet Oncol. 22, e7–e17. doi: 10.1016/s1470-2045(20)30556-8
Newman, A. M., Liu, C. L., Green, M. R., Gentles, A. J., Feng, W., Xu, Y., et al. (2015). Robust enumeration of cell subsets from tissue expression profiles. Nat. Methods 12, 453–457. doi: 10.1038/nmeth.3337
Nombela, P., Miguel-López, B., and Blanco, S. (2021). The role of m(6)A, m(5)C and Ψ RNA modifications in cancer: novel therapeutic opportunities. Mol. Cancer 20:18. doi: 10.1186/s12943-020-01263-w
Ochoa, S., de Anda-Jáuregui, G., and Hernández-Lemus, E. (2020). Multi-omic regulation of the PAM50 gene signature in breast cancer molecular subtypes. Front. Oncol. 10:845. doi: 10.3389/fonc.2020.00845
Parker, J. S., Mullins, M., Cheang, M. C., Leung, S., Voduc, D., Vickery, T., et al. (2009). Supervised risk predictor of breast cancer based on intrinsic subtypes. J. Clin. Oncol. 27, 1160–1167. doi: 10.1200/jco.2008.18.1370
Pashayan, N., Antoniou, A. C., Ivanus, U., Esserman, L. J., Easton, D. F., French, D., et al. (2020). Personalized early detection and prevention of breast cancer: ENVISION consensus statement. Nat. Rev. Clin. Oncol. 17, 687–705. doi: 10.1038/s41571-020-0388-9
Reck, M., Rodríguez-Abreu, D., Robinson, A. G., Hui, R., Csöszi, T., Fülöp, A., et al. (2016). Pembrolizumab versus chemotherapy for PD-L1-positive non-small-cell lung cancer. N. Engl. J. Med. 375, 1823–1833. doi: 10.1056/NEJMoa1606774
Rong, D., Dong, Q., Qu, H., Deng, X., Gao, F., Li, Q., et al. (2021). m(6)A-induced LINC00958 promotes breast cancer tumorigenesis via the miR-378a-3p/YY1 axis. Cell Death Discov. 7:27. doi: 10.1038/s41420-020-00382-z
Rooney, M. S., Shukla, S. A., Wu, C. J., Getz, G., and Hacohen, N. (2015). Molecular and genetic properties of tumors associated with local immune cytolytic activity. Cell 160, 48–61. doi: 10.1016/j.cell.2014.12.033
Sadeghalvad, M., Mohammadi-Motlagh, H. R., and Rezaei, N. (2021). Immune microenvironment in different molecular subtypes of ductal breast carcinoma. Breast Cancer Res. Treat. 185, 261–279. doi: 10.1007/s10549-020-05954-2
Schmid, P., Adams, S., Rugo, H. S., Schneeweiss, A., Barrios, C. H., Iwata, H., et al. (2018). Atezolizumab and nab-paclitaxel in advanced triple-negative breast cancer. N. Engl. J. Med. 379, 2108–2121. doi: 10.1056/NEJMoa1809615
Schmid, P., Rugo, H. S., Adams, S., Schneeweiss, A., Barrios, C. H., Iwata, H., et al. (2020). Atezolizumab plus nab-paclitaxel as first-line treatment for unresectable, locally advanced or metastatic triple-negative breast cancer (IMpassion130): updated efficacy results from a randomised, double-blind, placebo-controlled, phase 3 trial. Lancet Oncol. 21, 44–59. doi: 10.1016/s1470-2045(19)30689-8
Shi, Y., Zheng, C., Jin, Y., Bao, B., Wang, D., Hou, K., et al. (2020). Reduced expression of METTL3 promotes metastasis of triple-negative breast cancer by m6A methylation-mediated COL3A1 up-regulation. Front. Oncol. 10:1126. doi: 10.3389/fonc.2020.01126
Statello, L., Guo, C. J., Chen, L. L., and Huarte, M. (2021). Gene regulation by long non-coding RNAs and its biological functions. Nat. Rev. Mol. Cell Biol. 22, 96–118. doi: 10.1038/s41580-020-00315-9
Sun, T., Wu, Z., Wang, X., Wang, Y., Hu, X., Qin, W., et al. (2020). LNC942 promoting METTL14-mediated m(6)A methylation in breast cancer cell proliferation and progression. Oncogene 39, 5358–5372. doi: 10.1038/s41388-020-1338-9
Sung, H., Ferlay, J., Siegel, R. L., Laversanne, M., Soerjomataram, I., Jemal, A., et al. (2021). Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J. Clin. 71, 209–249. doi: 10.3322/caac.21660
Tan, A., Dang, Y., Chen, G., and Mo, Z. (2015). Overexpression of the fat mass and obesity associated gene (FTO) in breast cancer and its clinical implications. Int. J. Clin. Exp. Pathol. 8, 13405–13410.
Wang, G., Chen, X., Liang, Y., Wang, W., and Shen, K. (2017). A long noncoding RNA signature that predicts pathological complete remission rate sensitively in neoadjuvant treatment of breast cancer. Transl. Oncol. 10, 988–997. doi: 10.1016/j.tranon.2017.09.005
Wang, T., Kong, S., Tao, M., and Ju, S. (2020). The potential role of RNA N6-methyladenosine in cancer progression. Mol. Cancer 19:88. doi: 10.1186/s12943-020-01204-7
Wang, X., Liu, C., Zhang, S., Yan, H., Zhang, L., Jiang, A., et al. (2021). N(6)-methyladenosine modification of MALAT1 promotes metastasis via reshaping nuclear speckles. Dev. Cell 56, 702–715.e8. doi: 10.1016/j.devcel.2021.01.015
Wei, M., Bai, J. W., Niu, L., Zhang, Y. Q., Chen, H. Y., and Zhang, G. J. (2020). The complex roles and therapeutic implications of m(6)A modifications in breast cancer. Front. Cell Dev. Biol. 8:615071. doi: 10.3389/fcell.2020.615071
Wein, L., Savas, P., Luen, S. J., Virassamy, B., Salgado, R., and Loi, S. (2017). Clinical validity and utility of tumor-infiltrating lymphocytes in routine clinical practice for breast cancer patients: current and future directions. Front. Oncol. 7:156. doi: 10.3389/fonc.2017.00156
Wilkerson, M. D., and Hayes, D. N. (2010). ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics 26, 1572–1573. doi: 10.1093/bioinformatics/btq170
Wimberly, H., Brown, J. R., Schalper, K., Haack, H., Silver, M. R., Nixon, C., et al. (2015). PD-L1 expression correlates with tumor-infiltrating lymphocytes and response to neoadjuvant chemotherapy in breast cancer. Cancer Immunol. Res. 3, 326–332. doi: 10.1158/2326-6066.cir-14-0133
Wu, L., Wu, D., Ning, J., Liu, W., and Zhang, D. (2019). Changes of N6-methyladenosine modulators promote breast cancer progression. BMC Cancer 19:326. doi: 10.1186/s12885-019-5538-z
Xu, F., Huang, X., Li, Y., Chen, Y., and Lin, L. (2021). m(6)A-related lncRNAs are potential biomarkers for predicting prognoses and immune responses in patients with LUAD. Mol. Ther. Nucleic Acids 24, 780–791. doi: 10.1016/j.omtn.2021.04.003
Yao, Y., Zhang, T., Qi, L., Zhou, C., Wei, J., Feng, F., et al. (2019). Integrated analysis of co-expression and ceRNA network identifies five lncRNAs as prognostic markers for breast cancer. J. Cell. Mol. Med. 23, 8410–8419. doi: 10.1111/jcmm.14721
Yoshihara, K., Shahmoradgoli, M., Martínez, E., Vegesna, R., Kim, H., Torres-Garcia, W., et al. (2013). Inferring tumour purity and stromal and immune cell admixture from expression data. Nat. Commun. 4:2612. doi: 10.1038/ncomms3612
Zhang, A., Zhou, N., Huang, J., Liu, Q., Fukuda, K., Ma, D., et al. (2013). The human long non-coding RNA-RoR is a p53 repressor in response to DNA damage. Cell Res. 23, 340–350. doi: 10.1038/cr.2012.164
Zhang, C., Samanta, D., Lu, H., Bullen, J. W., Zhang, H., Chen, I., et al. (2016a). Hypoxia induces the breast cancer stem cell phenotype by HIF-dependent and ALKBH5-mediated m6A-demethylation of NANOG mRNA. Proc. Natl. Acad. Sci. U.S.A. 113, E2047–E2056. doi: 10.1073/pnas.1602883113
Zhang, C., Zhi, W. I., Lu, H., Samanta, D., Chen, I., Gabrielson, E., et al. (2016b). Hypoxia-inducible factors regulate pluripotency factor expression by ZNF217- and ALKBH5-mediated modulation of RNA methylation in breast cancer cells. Oncotarget 7, 64527–64542. doi: 10.18632/oncotarget.11743
Zhang, J., Wang, L., Xu, X., Li, X., Guan, W., Meng, T., et al. (2020). Transcriptome-based network analysis unveils eight immune-related genes as molecular signatures in the immunomodulatory subtype of triple-negative breast cancer. Front. Oncol. 10:1787. doi: 10.3389/fonc.2020.01787
Zheng, F., Du, F., Zhao, J., Wang, X., Si, Y., Jin, P., et al. (2021). The emerging role of RNA N6-methyladenosine methylation in breast cancer. Biomark. Res. 9:39. doi: 10.1186/s40364-021-00295-8
Zhu, P., He, F., Hou, Y., Tu, G., Li, Q., Jin, T., et al. (2021). A novel hypoxic long noncoding RNA KB-1980E6.3 maintains breast cancer stem cell stemness via interacting with IGF2BP1 to facilitate c-Myc mRNA stability. Oncogene 40, 1609–1627. doi: 10.1038/s41388-020-01638-9
Keywords: breast cancer, m6A-related LncRNAs, immune infiltration, gene signature, prognosis
Citation: Zhang J, Shan B, Lin L, Dong J, Sun Q, Zhou Q, Chen J and Han X (2021) Dissecting the Role of N6-Methylandenosine-Related Long Non-coding RNAs Signature in Prognosis and Immune Microenvironment of Breast Cancer. Front. Cell Dev. Biol. 9:711859. doi: 10.3389/fcell.2021.711859
Received: 19 May 2021; Accepted: 10 September 2021;
Published: 06 October 2021.
Edited by:
Jawed Akhtar Siddiqui, University of Nebraska Medical Center, United StatesReviewed by:
Kristina Andrijauskaite, The University of Texas Health Science Center at San Antonio, United StatesSumeet Nayak, KSQ Therapeutics, United States
Zili Zhen, Fudan University, China
Copyright © 2021 Zhang, Shan, Lin, Dong, Sun, Zhou, Chen and Han. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Jian Chen, jchen@fsyy.ustc.edu.cn; orcid.org/0000-0003-0759-7201; Xinghua Han, hxhmail@ustc.edu.cn; orcid.org/0000-0002-0632-6955