m6A-Related lncRNAs Predict Overall Survival of Patients and Regulate the Tumor Immune Microenvironment in Osteosarcoma

Abstract

Background: m6A-related lncRNAs have demonstrated great potential tumor diagnostic and therapeutic targets. The goal of this work was to find m6A-regulated lncRNAs in osteosarcoma patients.

Method: The Cancer Genome Atlas (TCGA) database was used to retrieve RNA sequencing and medical information from osteosarcoma sufferers. The Pearson's correlation test was used to identify the m6A-related lncRNAs. A risk model was built using univariate and multivariable Cox regression analysis. Kaplan-Meier survival analysis and receiver functional requirements were used to assess the risk model's performance (ROC). By using the CIBERSORT method, the associations between the relative risks and different immune cell infiltration were investigated. Lastly, the bioactivities of high-risk and low-risk subgroups were investigated using Gene Set Enrichment Analysis (GSEA).

Result: A total of 531 m6A-related lncRNAs were obtained from TCGA. Seven lncRNAs have demonstrated prognostic values. A total of 88 OS patients were separated into cluster 1, cluster 2, and cluster 3. The overall survival rate of OS patients in cluster 3 was more favorable than that of those in cluster 1 and cluster 2. The average Stromal score was much higher in cluster 1 than in cluster 2 and cluster 3 (P < 0.05). The expression levels of lncRNAs used in the construction of the risk prediction model in the high-risk group were generally lower than those in the low-risk group. Analysis of patient survival indicated that the survival of the low-risk group was higher than that of the high-risk group (P < 0.0001) and the area under the curve (AUC) of the ROC curve was 0.719. Using the CIBERSORT algorithm, the results revealed that Macrophages M0, Macrophages M2, and T cells CD4 memory resting accounted for a large proportion of immune cell infiltration. By GSEA analysis, our results implied that the high-risk group was mainly involved in unfolded protein response, DNA repair signaling, and epithelial-mesenchymal transition signaling pathway and glycolysis pathway; meanwhile, the low-risk group was mainly involved in estrogen response early and KRAS signaling pathway.

Conclusion: Our investigation showed that m6A-related lncRNAs remained tightly connected to the immunological microenvironment of osteosarcoma tumors, potentially influencing carcinogenesis and development. The immune microenvironment and immune-related biochemical pathways can be changed by regulating the transcription of M6A modulators or lncRNAs. In addition, we looked for risk-related signaling of m6A-related lncRNAs in osteosarcomas and built and validated the risk prediction system. The findings of our current analysis will facilitate the assessment of outcomes and the development of immunotherapies for sufferers of osteosarcomas.

Figure 1

The correlation between 10 m6A-related genes and 7 lncRNA prognostic genes. ^∗P < 0.05; ^∗∗P < 0.01.

Figure 2

Consensus clustering of m6A-related prognostic lncRNAs. (a) TCGA osteosarcoma cohorts were grouped into three clusters according to the consensus clustering matrix (k = 3). (b) Overall survival analysis revealed a better overall survival of osteosarcoma patients in cluster 3 than those in cluster 1 and cluster 2. (c) The heatmap of the 3 clusters along with general information of patients.

Figure 3

(a) The comparison of the immune score between the low-risk and high-risk groups (p > 0.05). (b–e) Comparison of immune score: (b) Stromal score, (c) Immune score, (d) ESTIMATE Score, and (e) Tumor Purity.

Figure 4

(a) The risk score distribution; (b) survival time scatter diagram; (c) clinical and pathological characteristics and varied lncRNA expression patterns in high- and low groups are depicted in a heatmap; (d) the risks model's Kaplan–Meier survival line; (e) ROC curve analysis.

Figure 5

The violin plot of 22 tumor-infiltrating immune cell types in low- and high-risk groups. The infiltration of 22 immune cell types in osteosarcoma was analyzed by the CIBERSORT algorithm.

Figure 6

The correlations between immune infiltration and m6A-related lncRNAs. ^∗P < 0.05; ^∗∗P < 0.01.

Figure 7

Abnormally activated signaling pathways in the two subgroups after Gene Set Enrichment Analysis. (a–d) Performed in the high-risk group, including unfolded protein response (ES = 0.42, P=0.02, FDR = 0.110), DNA repair signaling pathway (ES = 0.41, P=0.0096, FDR = 0.070), epithelial-mesenchymal transition signaling pathway (ES = 0.55, P=0.032, FDR = 0.337), and glycolysis pathway (ES = 0.44, P=0.032, FDR = 0.087). (e–f) Performed in the low-risk group, including estrogen response early (ES = 0.40, P=0.038, FDR = 1.0) and KRAS signaling pathway (ES = 0.40, P=0.039, FDR = 1.0).