The prognostic value and immune landscaps of m6A/m5C-related lncRNAs signature in the low grade glioma

Abstract

Background: N6-methyladenosine (m6A) and 5-methylcytosine (m5C) are the main RNA methylation modifications involved in the oncogenesis of cancer. However, it remains obscure whether m6A/m5C-related long non-coding RNAs (lncRNAs) affect the development and progression of low grade gliomas (LGG).

Methods: We summarized 926 LGG tumor samples with RNA-seq data and clinical information from The Cancer Genome Atlas and Chinese Glioma Genome Atlas. 105 normal brain samples with RNA-seq data from the Genotype Tissue Expression project were collected for control. We obtained a molecular classification cluster from the expression pattern of sreened lncRNAs. The least absolute shrinkage and selection operator Cox regression was employed to construct a m6A/m5C-related lncRNAs prognostic signature of LGG. In vitro experiments were employed to validate the biological functions of lncRNAs in our risk model.

Results: The expression pattern of 14 sreened highly correlated lncRNAs could cluster samples into two groups, in which various clinicopathological features and the tumor immune microenvironment were significantly distinct. The survival time of cluster 1 was significantly reduced compared with cluster 2. This prognostic signature is based on 8 m6A/m5C-related lncRNAs (GDNF-AS1, HOXA-AS3, LINC00346, LINC00664, LINC00665, MIR155HG, NEAT1, RHPN1-AS1). Patients in the high-risk group harbored shorter survival times. Immunity microenvironment analysis showed B cells, CD4 + T cells, macrophages, and myeloid-derived DC cells were significantly increased in the high-risk group. Patients in high-risk group had the worse overall survival time regardless of followed TMZ therapy or radiotherapy. All observed results from the TCGA-LGG cohort could be validated in CGGA cohort. Afterwards, LINC00664 was found to promote cell viability, invasion and migration ability of glioma cells in vitro.

Conclusion: Our study elucidated a prognostic prediction model of LGG by 8 m6A/m5C methylated lncRNAs and a critical lncRNA regulation function involved in LGG progression. High-risk patients have shorter survival times and a pro-tumor immune microenvironment.

Keywords: Immune landscape; Long non-coding RNA; Low grade glioma; Prognostic signature; RNA methylation.

Fig. 1

Screening and clustering of m6A/m5C-related lncRNAs in Low Grade Gliomas. A The co-expression network of m6A/m5C regulators and related lncRNAs. B The heatmap of differential expression of 14 screened lncRNAs with prognostic value between 105 brain normal control tissues from the GTEx dataset and 509 LGG tumor tissues from the TCGA database. C Consensus matrix of the consensus clustering based on the 14 screened m6A/m5C-related lncRNAs for optimal k = 2 in 509 LGG tumor tissues from the TCGA cohort. D Principal component analysis (PCA) of cluster1 and cluster2. E Kaplan–Meier curve of overall survival time between cluster1 and cluster2 (p < 0.0001)

Fig. 2

Clinical characteristics, tumor immune landscapes and tumor mutational burden among clustering subgroup in the TCGA-LGG dataset. A Heatmap of the clinicopathological features between the two distinct clusters for the TCGA-LGG, the relative expression levels of respective lncRNAs ware plotted in heat map after normalization by using “scale” function. B Heatmap of the differential expression status of 29 immune checkpoint genes between cluster1 and cluster2. C, D Box plots showing differences in the infiltration of immune cells estimated using TIMER2 and CIBERSORT respectively. The Kruskal–Wallis test was used to determine the statistical significance of the difference between cluster1 and cluster2. E GSEA analysis displayed key pathways in cluster1 and cluster2. F Waterfall plots showing the distribution of the top 15 most frequent somatic mutation in the two clusters, the upper bar graph shows TMB

Fig. 3

Construction of m6A/m5C-related lncRNAs prognostic signature. A, B 8 m6A/m5C-related lncRNAs were constructed by using LASSO cox regression. C Bar plot shows the correlation coefficient of screened 8 m6A/m5C-related lncRNAs. (D, G) The distributions of risk scores, overall survival status and the expression of 8 m6A/m5C-related lncRNAs in the TCGA-LGG and CGGA-LGG. E, H Kaplan–Meier curve of OS time of high-risk group and low-risk group in the TCGA-LGG and CGGA-LGG cohort. F, I ROC curve of the risk score at 1-, 3- and 5-years’ follow-up in the TCGA-LGG and CGGA-LGG cohort. J Differences in the low- and high-risk group between various clinicalpathological features. The Wilcoxon test was used to determine the statistical significance

Fig. 4

Cox regression analysis and comprehensive nomogram for m6A/m5C-related lncRNAs prognostic signature. A–D Forest plots showing univariate and multivariate Cox regression analysis for OS in the TCGA-LGG patients and CGGA-LGG patients. E, F The ROC curve of the risk score and other clinical characteristics in the TCGA-LGG and CGGA-LGG cohort. G The comprehensive nomogram predicting the clinical outcomes of LGG patients with 1-, 3- and 5-year survival based on the TCGA-LGG cohort. H, I ROC curve and calibration plots for predicting the 1-, 3- and 5-year OS based on the TCGA-LGG cohort

Fig. 5

Tumor immune landscapes and tumor mutational burden between low- and high-risk subgroup in the TCGA-LGG dataset. A The sankey diagram of the relationship among the cluster, risk group and OS status in the TCGA-LGG patients. B Heatmap of the differential expression status of 29 immune checkpoint genes between low- and high-risk subgroup. C Dot plots show the correlation between ESTIMATE score, Stromal score, Immune score and risk score in the TCGA-LGG cohort. D, E Box plots showing differences in the infiltration of immune cells estimated using TIMER2 and CIBERSORT respectively. The Kruskal–Wallis test was used to determine the statistical significance of the difference between low- and high-risk subgroup. F Dot plots show the correlation between the infiltration levels of neutrophil, myeloid dendritic cell, macrophage M2 with risk score in the TCGA-LGG cohort. G Kaplan–Meier curve of OS time of low- and high-TMB patients in the TCGA-LGG cohort. H Kaplan–Meier curve of OS time of low- and high-TMB patients combined with risk score in the TCGA-LGG cohort. I The landscape of CNV in low- and high-risk groups

Fig. 6

Predictive value of risk scores in clinical treatment subgroups. A, B Kaplan–Meier curves of OS for low- and high-risk groups with or without TMZ treatment in the TCGA-LGG and CGGA-LGG cohort. C, D Kaplan–Meier curves of OS for low- and high-risk groups with or without radiotherapy in the TCGA-LGG and CGGA-LGG cohort

Fig. 7

Knockdown of LINC00664 inhibited cell viability, invation and migration in glioma cells. A, B LINC00664 expression levels in U87 and U251 cells transfected with siNC or two distinct specific LINC00664 siRNAs was detected by RT-PCR. Cell viability of transfected glioma cells was detected by CCK-8 assay at 0, 24, 48, 72 h after incubation. C Trans-well assay were used to detect the invasion ability of U87 and U251 cells. D wound healing assay were used to detect the migration ability of U87 and U251 cells. The migration area means the distance between the two edges at the 24 h relative to 0 h. ***p < 0.001, **p < 0.01