Identification of m5C-related lncRNAs signature to predict prognosis and therapeutic responses in esophageal squamous cell carcinoma patients

Abstract

Esophageal squamous cell carcinoma (ESCC) has a dismal prognosis because of atypical early symptoms and heterogeneous therapeutic responses. 5-methylcytosine (m⁵C) modification plays an important role in the onset and development of many tumors and is widespread in long non-coding RNA (lncRNA) transcripts. However, the functions of m⁵C and lncRNAs in ESCC have not been completely elucidated. Herein, this study aimed to explore the role of m⁵C-related lncRNAs in ESCC. The RNA-seq transcriptome profiles and clinical information were downloaded from the TCGA-ESCC database. Pearson analysis was used to identify m⁵C-related lncRNAs. Then we established the m⁵C-related lncRNAs prognostic signature (m⁵C-LPS) using univariate Cox and least absolute shrinkage and selection operator (LASSO) regression analysis. Then, the prognostic value of m⁵C-LPS was evaluated internally and externally using the TCGA-ESCC and GSE53622 databases through multiple methods. We also detected the expression of these lncRNAs in ESCC cell lines and patient tissues. Fluorescence in situ hybridization (FISH) was used to detect the prognostic value of specific lncRNA. In addition, clinical parameters, immune status, genomic variants, oncogenic pathways, enrichment pathways, and therapeutic response features associated with m⁵C-LPS were explored using bioinformatics methods. We constructed and validated a prognostic signature based on 9 m⁵C-related lncRNAs (AC002091.2, AC009275.1, CAHM, LINC02057.1, AC0006329.1, AC037459.3, AC064807.1, ATP2B1-AS1, and UBAC2-AS1). The quantitative real-time polymerase chain reaction (qRT-PCR) revealed that most lncRNAs were upregulated in ESCC cell lines and patient tissues. And AC002091.2 was validated to have significant prognostic value in ESCC patients. A composite nomogram was generated to facilitate clinical practice by integrating this signature with the N stage. Besides, patients in the low-risk group were characterized by good clinical outcomes, favorable immune status, and low oncogenic alteration. Function enrichment analysis indicated that the risk score was associated with mRNA splicing, ncRNA processing, and DNA damage repair response. At the same time, we found significant differences in the responses to chemoradiotherapy between the two groups, proving the value of m⁵C-LPS in treatment decision-making in ESCC. This study established a novel prognostic signature based on 9 m⁵C-related lncRNAs, which is a promising biomarker for predicting clinical outcomes and therapeutic response in ESCC.

Figure 1

The expression pattern and interactive landscape of the m⁵C regulators in the TCGA-ESCC database. (A) Heatmap presenting the expression of 16 m⁵C regulators in normal esophageal (N) and ESCC (T) tissues from TCGA-ESCC database. p < 0.1, *p < 0.05, **p < 0.01, and ***p < 0.001. (B) Protein–protein interaction (PPI) network showing the interaction between m⁵C regulators. (C) Heatmap showing the Pearson correlation among 16 m⁵C regulators.

Figure 2

Identification and validation of the m⁵C-related lncRNAs prognostic signature (m⁵C-LPS) based on the cohort of TCGA-ESCC and GSE53622. (A,B) The minimum criterion of the LASSO regression algorithm was used to identify the most robust prognostic m⁵C-related lncRNAs. (C) Forest plot presenting the hazard ratio (HR) and 95% confidence interval (CI) of the 9 lncRNAs by the multivariate Cox regression. (D) The coefficients of the 9 lncRNAs contained in the m⁵C-LPS formula. (E) The distributions of the risk score, vital status, overall survival (OS), and expression levels of the 9 m⁵C-related lncRNAs in low- and high-risk groups in the cohort from TCGA-ESCC. (F) Kaplan–Meier (K–M) analysis demonstrated that patients with higher risk scores exhibited worse overall survival in the cohort from TCGA-ESCC. (G) K–M analysis demonstrated that patients with higher risk scores exhibited worse disease-free survival in the cohort from TCGA-ESCC. (H) The area under the curve (AUC) of the time-dependent ROC curves measures the predictive value of the risk score in the cohort from TCGA-ESCC. (I) The distributions of the risk score, vital status, overall survival (OS), and expression levels of the 9 m⁵C-related lncRNAs in low- and high-risk groups in the cohort from GSE53622. (J) K–M analysis demonstrated that patients with higher risk scores exhibited worse overall survival in the cohort from GSE53622. (K) AUC of the time-dependent ROC curves measuring the predictive value of the risk score in the cohort from GSE53622.

Figure 3

Co-expression status and expression level of m⁵C-related lncRNAs in TCGA-ESCC database and ESCC cell lines. (A) Sankey plot showing one-to-one matches between m⁵C genes, m⁵C-related lncRNAs, and their risk type. (B) Circle plot presenting the co-expression status of the 9 m⁵C-related lncRNAs with coefficients annotated. (C) The expression level of m⁵C-related lncRNAs in normal esophageal and ESCC tissues based on TCGA-ESCC database. (D) Representative Fluorescence in situ hybridization (FISH) images of AC002091.2 in ESCC tissues. Scale bars represent 50 μm. (E) K–M plot for overall survival grouped by AC002091.2 expression in 54 ESCC patients. *p < 0.05, **p < 0.01, ***p < 0.001.

Figure 4

The discrepancy in risk scores between different subgroups: T stage (A), M stage (B), stage (C), alcohol history (D), adjuvant postoperative pharmaceutical therapy (E), adjuvant postoperative radiotherapy (F), neoplasm status (G), disease free status (H), and overall survival status (I).

Figure 5

Verification of the independent prognostic value of m⁵C-LPS and construction of nomogram. Univariate (A) and multivariate (B) Cox regression analyses of the prognostic value of risk scores and other clinical parameters. (C) The ROC curves show the predictive value of the risk score and other clinical characteristics. (D) Nomogram composed of N stage and risk score was constructed to predict 1-, 2-, and 3-year survival rates. (E) Calibration plots were used to evaluate the nomogram for predicting 1-, 2-, and 3-year survival rates.

Figure 6

Investigation of immune status in different risk groups. (A) Heatmap revealing the immune and stromal cells infiltration in ESCC immune microenvironment. (B) Box plots showing the infiltration of the immune cells based on the Cibersort algorithm in different risk groups; *p < 0.05 and **p < 0.01. (C) Estimation of the coefficients for risk score with immune checkpoint genes. (D) Histogram showing the relationships between risk score and chemokines, receptors, MHC, immunoinhibitors, and immunostimulators.

Figure 7

Identification of therapeutic response features of different risk groups. (A) Heatmap revealing IC50 for chemotherapeutic agents and radiation-sensitivity index (RSI) for radiotherapy. Box plots showing the sensitivity of selected chemotherapeutic agents for patients in low- and high-risk groups based on CGP (B), CTRP (C), and GDSC (D) databases. (E) Boxplot showing the radiotherapy sensitivity for patients in low- and high-risk groups.