Three Novel Autophagy-Related lncRNAs as Prognostic Biomarkers for Lung Adenocarcinoma ()
1. Introduction
Lung cancer is currently the most common and lethal malignancy worldwide. The 2-year relative survival rate for non-small cell lung cancer (NSCLC) increased from 34% diagnostic rate between 2009 and 2010 to 42% between 2015 and 2016, including an absolute increase of 5% to 6% for each diagnostic stage [1]. Globocan 2020 data indicate the incidence and mortality of lung cancer are 37.0% worldwide and 39.8% in China [2]. Different histological subtypes have various origins, morphological appearance, epigenetic changes, clinical responses, and outcomes [3]. Therefore, investigating diagnostic and prognostic genes or markers for different subtypes could help implement precision therapies and improve patient survival.
Autophagy and apoptosis play important roles in lung cancer progression [4]. Eukaryotic cells respond to different stimuli and stress, such as starvation, hypoxia, drugs, and other harmful stimuli. To maintain cellular homeostasis, the degradation process of eukaryotic cells is regulated by autophagy-related genes (ARGs) [5]. Autophagy is a lysosomal degradation pathway physiologically involved in differentiation, survival, development and homeostasis, and pathologically in infection, cancer, neurodegeneration and aging as well as heart, liver, and kidney disease [6]. Notably, lncRNAs may act as competing endogenous RNAs (ceRNAs) which competitively bind with target miRNAs to regulate mRNAs [7]. In several studies, autophagy was shown to promote the tolerance of tumor cells to radiotherapy and chemotherapy, and lncRNAs play a crucial role in these processes [8] [9].
In the present study, the autophagy-related lncRNAs were mined from The Cancer Genome Atlas (TCGA) by searching for the lncRNA co-expressed genes of ARGs in LUAD and identified the autophagy lncRNAs associated with patient clinical data to investigate the molecular mechanisms of LUAD.
2. Materials and Methods
2.1. The Gene Expression and Clinical Data of LUAD Patients
The gene expression and clinical data of LUAD patients were downloaded from TCGA database (https://portal.gdc.cancer.gov/). In this study, among 551 LUADs, 497 tumor samples were compared with 54 adjacent tissues, and their gene expression and clinical data were collected. In a previous study, the M stage in the TNM staging system was suggested to be an unreliable independent prognostic factor [10], therefore, the missing and the M stage clinical data were excluded. In addition, after filtering, 497 cases of gene expression and 436 cases of clinical data were analyzed.
2.2. Autophagy-Related lncRNA Genes in LUAD
A list of 232 ARGs was collected from the Human Autophagy Database (HADb) (http://www.autophagy.lu/). After Spearman rank correlation test, 311 autophagy-related lncRNAs and their coefficient values were obtained using the following criteria: correlation > 0.7 and P < 0.001.
2.3. The Autophagy-Related lncRNA Prognostic Model
The autophagy-related lncRNAs were combined with survival data and then a univariate Cox’s proportional hazard analysis was performed as well as Kaplan-Meier analysis using the survival package in R. A P-value < 0.05 was considered statistically significant. To determine the optimum prognostic autophagy-related lncRNAs, the minimum Akaike information criterion (AIC) value was screened. Next, the OS curves of the optimum autophagy-related lncRNAs were drawn based on the gene expression levels. The LUAD patients were divided into high-risk and low-risk groups based on the median risk score calculated using the following formula: risk score = ΣCoef gene i × Gene i expression. Survival curves for each group were constructed based on Kaplan-Meier analysis.
2.4. Stability Evaluation of the Prognostic Model
To assess the stability of the constructed prognostic model, the clinical data were combined with the risk score of each LUAD patient. The clinical data included age, gender, clinical stage, T stage, and N stage. Then, univariate and multivariate prognostic analyses were performed to determine whether the three lncRNAs were independent prognostic factors using the survival package in R. The area under the ROC curve (AUC) of risk score and different clinical data were created to assess the stability for 1-year survival time using the survival ROC package in R.
2.5. Gene Set Enrichment Analysis (GSEA)
The Gene Set Enrichment Analysis (GSEA) was used to analyze the related pathways by choosing “c2.cp.kegg.v7.0.symbols.gmt” gene sets as reference gene sets. The above analysis was performed using GSEA 4.1.0 software.
3. Results
3.1. The Optimum Three Prognostic Autophagy-Related lncRNAs
The raw data contained 19,658 mRNAs and 14,142 lncRNAs. Spearman rank correlation test, univariate Cox’s proportional hazard analysis, and Kaplan-Meier analysis were used to screen 311 autophagy-related lncRNAs. Five prognostic autophagy-related lncRNAs were obtained: AC090948.1, AC011477.2, AC099850.3, TRG-AS1, and AC067852.3 (Table 1). Three optimum prognosis autophagy-related lncRNAs were acquired based on the minimum AIC value of 1551.55: AC011477.2, AC099850.3, and TRG-AS1 (Table 2).
3.2. Construction of Prognostic Model
The risk score of all LUAD patients was calculated based on the following formula: risk score = ΣCoef gene i × Gene i expression.
Risk score = −0.179635915 × AC011477.2 + 0.034555529 × AC099850.3 − 0.234455556 × TRG-AS1
Based on the median risk score (cutoff value = 1.02407), the LUAD patients were divided into high-risk and low-risk groups. A Sankey diagram of autophagy-related mRNA-lncRNA co-expression network was established (Table 3 and Figure 1).
To establish the optimal prognostic risk model, the survival curves of the three autophagy-related lncRNAs were drawn using survival package in R. The overall survival (OS) of patients with high TRG-AS1 and AC011477.2 expression levels was longer (Figure 2(a) and Figure 2(b)). Similarly, patients with a high AC099850.3 expression level correlated with poor OS (Figure 2(c)). Furthermore, all LUAD patients were divided into high-risk or low-risk groups based on the median risk score of 1.02407. Then, an OS curve based on Kaplan-Meier survival analysis was performed; LUAD patients in the low-risk score group had
![]()
Table 1. Five prognostic autophagy-related lncRNAs.
![]()
Table 2. Three selected autophagy-related lncRNAs.
![]()
Table 3. Network of mRNA and lncRNA.
![]()
Figure 1. The constructed autophagy-related mRNAs-lncRNAs co-expression network was visualized by Sankey diagram.
a better prognosis (Figure 2(d), P < 0.001). Based on the OS curve, the 5-year survival rate in the high-risk group was 26.51% (95% CI: 0.1842 - 0.382) versus 41.6% (95% CI: 0.307 - 0.563) in the low-risk group (P < 0.05).
To visualize the data, the risk curve (Figure 3(a)), survival status scatter plot (Figure 3(b)), and the heat map of the three lncRNAs’ expression (Figure 3(c)) were drawn.
3.3. Stability Validation of the Prognostic Model Composed of the Three Autophagy-Related lncRNAs
To evaluate the stability of the prognostic model composed of the three lncRNAs,
univariate and multivariate survival analyses of 436 patients’ clinical data and risk scores were performed using the survival package in R. Based on the forest plot (Figure 4(a) and Figure 4(b)), the hazard ratio (HR) of the risk score in the univariate Cox’s analysis was 2.041 (95% CI 1.537 - 2.711; P < 0.001) and 1.972 (95% CI 1.372 - 2.665; P < 0.001) in multivariate Cox’s analysis. Therefore, the three lncRNAs (AC011477.2, AC099850.3, TRG-AS1) could be used as independent prognostic factors for LUAD patients.
To further verify the stability of this model, ROC curves were created in the prognostic model. As shown in the ROC curves (Figure 5), AUC of the risk
![]()
(a)![]()
(b)![]()
(c)
Figure 3. Autophagy-related lncRNA risk score analysis of LUAD patients. (a) the risk curve based on the risk score of each sample; (b) the scatter plot based on the survival status of each sample; (c) Heatmap showed the expression profiles of the three autophagy-related lncRNAs in the low-risk and high-risk groups.
score was 0.700, indicating a higher diagnostic accuracy of risk score. Furthermore, the AUC of the clinical staging was 0.715, indicating that applying clinical staging for diagnosis was equally highly valuable.
In conclusion, the stability of this prognostic model was acceptable and could be applied for the prognosis of LUAD patients.
In addition, clinical data were divided into different groups based on correlation with risk score to perform differential analysis using t-test (Table 4). Women had lower risk than men, the risk of clinical stage I - II was less than of stage III - IV, the risk of T1 - T2 was less than of T3 - T4, and the risk of N0 was less than of N1 - N3 (P < 0.05).
![]()
(a)![]()
(b)
Figure 4. The univariate and multivariate Cox regression analysis of risk model score and clinical features regarding of prognostic value.
![]()
Figure 5. The ROC curves of risk score and the clinical data.
![]()
Table 4. The association between risk score and clinical characteristics.
3.4. GSEA
GSEA results showed different pathways were enriched in the high-risk and low-risk groups (Figure 6). KEGG reference gene sets were performed and results indicated enrichment in 36 pathways including pyrimidine metabolism, pentose phosphate pathway, citric acid cycle, and cell cycle in the high-risk group (P < 0.05) and enrichment in FC-EPSILON-RI signal pathway, intestinal immune network produced by IgA, and ABC transporters in the low-risk group (P < 0.05).
4. Discussion
Over the past few decades, researchers could investigate genomes more systematically due to application of high-throughput technologies [11]. Bioinformatics, interdisciplinary medicine, statistics, computer science, and applied mathematics, are disciplines that contributed to the development of the human genome project [12].
In several studies, autophagy has been shown to affect gene expression and even clinical treatment. Autophagy was shown to affect hepatocarcinogenesis and progression in several studies by regulating liver cancer immunity, oxidative stress, and cellular metabolism, including suppressing hepatocarcinogenesis at an early stage, promoting tumor growth at a progressive stage, and increasing chemoradiotherapy resistance of cancer cells [13] [14] [15]. Zhang et al. [16] found that lncRNA CASC2 might be useful for evaluating the expression of autophagy-related LC3-I, LC3-II, and p62. Several lncRNA-related models have been applied for the prognosis of patients with malignant tumors, such as lung, colon, liver, and breast cancers [17] [18] [19]. However, autophagy-related lncRNAs have not been extensively investigated.
![]()
Figure 6. KEGG analysis of the three autophagy-related lncRNAs by GSEA.
In the present study, a stable prognostic model of three autophagy-related lncRNAs (AC011477.2, AC099850.3, and TRG-AS1) in LUAD was established. AC099850.3 has been included in some hepatocellular carcinoma (HCC) prognosis studies, and the results showed that AC099850.3 may play an important role in the migration and proliferation of HCC cells. Wu et al. [20] found that AC099850.3 could promote the expression of cell cycle molecules BUB1, CDK1, PLK1, and TTK based on RT-PCR, and could inhibit CD155 and PD-L1 expression in the interference group based on Western blot analysis. The results indicated that AC099850.3 played a specific role in the immune response and warranted further research. Zhou et al. [21] showed that AC099850.3 was upregulated (logFC > 1) in tongue squamous cell carcinoma (TSCC).
He et al. [22] found that TRG-AS1 was highly expressed in TSCC patients and strongly associated with advanced TNM staging and poor OS. Sun et al. [23] demonstrated that TRG-AS1 was significantly overexpressed in HCC. Reportedly, TRG-AS1 silencing might inhibit the proliferation, migration, invasion, and epithelial mesenchymal transition (EMT). Xie et al. [24] found that TRG-AS1 regulated SUZ12 expression through the ceRNA of miR-877-5p, and TRG-AS1 might become a novel target for glioblastoma treatment. Notably, in the present study, TRG-AS1, a low-risk lncRNA for LUAD, had the opposite results in TSCC, liver cancer, and glioblastoma, which requires further research.
Currently, related studies in which lncRNA AC011477.2 has been investigated have not been performed. However, the present study can be used as a basis for future research on lncRNA AC011477.2 and clarify whether it has the potential to become a molecular diagnostic marker and tumor therapeutic target for LUAD.
GSEA enrichment analysis showed the high-risk group was enriched in the cell cycle, which could be associated with some of the above-mentioned effects of AC099850.3. However, further experimental proof is required. The results of the present study can be used for future prospective experiments to verify whether this prognostic model can be applicable for evaluating the prognosis of patients and identifying the specific mechanism(s) of these three lncRNAs in LUAD patients.
5. Conclusion
In the present study, an autophagy-related lncRNA prognostic model for LUAD consisting of AC011477.2, AC099850.3, and TRG-AS1 was established. The model showed acceptable stability for predicting the prognosis of LUAD patients and can be used for subsequent prospective experiments to improve clinical treatment.
Funding Information
Funding was provided by The National Natural Science Foundation of China (grant No. 81703001), the Natural Science Foundation of Hebei Province (grant No. H2020406050, H2021406021).