Abstract

Objective. To construct a predictive signature based on autophagy-associated lncRNAs for predicting prognosis in lung adenocarcinoma (LUAD). Materials and Methods. Differentially expressed autophagy genes (DEAGs) and differentially expressed lncRNAs (DElncRNAs) were screened between normal and LUAD samples at thresholds of and value < 0.05. Univariate Cox regression analysis was conducted to identify overall survival- (OS-) associated DElncRNAs. The total cohort was randomly divided into a training group () and a validation group () at a ratio of 1 : 1. Multivariate Cox regression analysis was used to build prognostic models in the training group that were further validated by the area under curve (AUC) values of the receiver operating characteristic (ROC) curves in both the validation and total cohorts. Results. A total of 30 DEAGs and 2997 DElncRNAs were identified between 497 LUAD tissues and 54 normal tissues; however, only 1183 DElncRNAs were related to the 30 DEAGs. A signature consisting of 13 DElncRNAs was built to predict OS in lung adenocarcinoma, and the survival analysis indicated a significant OS advantage of the low-risk group over the high-risk group in the training group, with a 5-year OS AUC of 0.854. In the validation group, survival analysis also indicated a significantly favorable OS for the low-risk group over the high-risk group, with a 5-year OS AUC of 0.737. Univariate and multivariate Cox regression analyses indicated that only positive surgical margin (vs negative surgical margin) and high-risk group (vs low-risk group) based on the predictive signature were independent risk factors predictive of overall mortality in LUAD. Conclusions. This study investigated the association between autophagy-associated lncRNAs and prognosis in LUAD and built a robust predictive signature of 13 lncRNAs to predict OS.

1. Introduction

Lung cancer remains a significant public health problem threatening life, with 142,670 estimated deaths in the United States in 2019 and over 1.6 million deaths worldwide annually [1, 2]. Lung cancer generally consists of small cell lung cancer (SCLC) and non-small-cell lung cancer (NSCLC), with lung adenocarcinoma (LUAD) accounting for almost 50% of NSCLC cases [35]. Although various therapeutic approaches have been introduced for LUAD, there were still no obvious improvements in ameliorating unfavorable prognoses, especially in patients with metastatic disease. Metastases of LUAD to the nervous system, bone, liver, adrenal gland and even urethra tend to indicate poor therapeutic outcomes, and only some selected cases may benefit from systematic therapy [68]. The TNM staging system provides a relatively reliable predictive model for prognosis and remains the most frequently applied predictor of survival [9]. However, a comprehensive investigation of the underlying molecular mechanisms and cellular pathways may be effective potential diagnostic tools and therapeutic targets for LUAD. Whole-exome sequencing and immune profiling analyses of LUAD indicated that molecular and immune phenotypes were associated with survival and response to adjuvant therapy in the clinical outcomes and personalized immune-based therapy of LUAD [10].

Autophagy, a highly evolutionarily conserved catabolic process, recycles and degrades cellular components via lysosomes to provide material for biomolecule synthesis [11, 12]. Malfunctions in autophagy are involved in a wide range of diseases, including cancer, neurodegeneration, and autoimmune diseases [1316]. Autophagy is a double-edged sword with survival-supporting effects or cell death promotion in cancer cells, and it affects cancer cell responses to cytotoxic drugs [14]. Increasing evidence indicates that the interplay of autophagy and apoptosis is crucial in the pathophysiology of LUAD [17]. Long noncoding RNAs (lncRNAs), characterized by their noncoding function and their greater than 200 nucleotide length, are involved in carcinogenesis, cancer progression, and metastasis and can serve as robust diagnostic and predictive biomarkers in a variety of cancers [1822].

Considering the significance of autophagy and lncRNAs in cancer biology, this study is aimed at investigating autophagy genes and autophagy-associated lncRNAs in LUAD from the TCGA (The Cancer Genome Atlas) database and building an effective signature based on autophagy-associated lncRNAs to predict prognosis in LUAD.

2. Materials and Methods

2.1. Data Collection

We retrieved the FPKM (fragments per kilobase of transcript per million fragments mapped) (level 3) sequencing profiles of mRNAs and lncRNAs from the TCGA data portal (https://tcga-data.nci.nih.gov/tcga/) and clinical information from the cBio Cancer Genomics Portal (http://cbioportal.org) in August 2019. The autophagy genes were collected from the Human Autophagy Database (HADb; http://www.autophagy.lu/project.html).

2.2. Identification of Differentially Expressed RNAs

The differentially expressed autophagy genes (DEAGs) and differentially expressed lncRNAs (DElncRNAs) were screened between LUAD and normal tissues by the “limma” package in R, with thresholds set as and value < 0.05. Heatmaps of the DEAGs were plotted by the “pheatmap” package.

2.3. Functional Enrichment Analysis of the DEAGs

Functional enrichment analysis of the DEAGs was conducted using DAVID, including biological functions, cellular components, and molecular functions, and the Kyoto Encyclopedia of Genes and Genomes (KEGG, http://www.kegg.jp/) database was searched for significant pathways ( and enrichment ). Autophagy-related DElncRNAs (ARDElncRNAs) were identified by evaluating the expression correlations between the DElncRNAs and DEAGs using correlation in the total LUAD cohort.

2.4. Prognosis-Associated ARDElncRNAs

Univariate Cox regression analysis was conducted to explore overall survival- (OS-) associated ARDElncRNAs with the “survival” package in R. The interaction network between DEAGs and OS-associated ARDElncRNAs was constructed and visualized by Cytoscape v3.7.0 software.

2.5. Construction and Validation of the ARDElncRNA-Based Signature

The total TCGA cohort was randomly divided into a training group () and validation group () in a ratio of 1 : 1 by the “caret” package in R. Multivariate Cox regression analysis was performed to build a prognostic model in the training group using the “survival” package in R. The predictive value of the ARDElncRNA-based signature was evaluated by the area under the curve (AUC) values of the receiver-operator characteristic (ROC) curves in the training, validation, and total groups using the “survivalROC” package.

2.6. Independent Risk Factors for Overall Mortality in LUAD

The clinical variables included the age at diagnosis, sex, T status, N status, M status, AJCC TNM stage, surgical margin resection status, and risk score of the ARDElncRNA-based signature. The patients were divided by age into the following groups: ≤18-60, 60-80, >80 years old, and unknown. Sex was classified as male or female. T status was categorized as T1, T2, T3/4, and unknown. N status was negative (N0), positive (including N1, N2, and N3), and unknown. M status included negative, positive, and unknown. Stage was categorized into four types: I, II, III/IV, and unknown. Surgical margin resection status was divided into negative (R0), positive (R1/2), and unknown. The risk scores based on the ARDElncRNA-based signature were classified into low- and high-risk score groups. Univariate and multivariate Cox regression analyses were performed to explore independent predictors of OS in LUAD.

3. Results

3.1. Functional Enrichment Analysis of the DEAGs in LUAD

A total of 232 autophagy genes were collected from the HADb database, and we identified 30 DEAGs between 497 tumor tissues and 54 normal tissues from TCGA at the threshold of and value < 0.05. As demonstrated in the heatmap (Figure 1(a)), volcano plot (Figure 1(b)), and boxplot (Figure 1(c)), 18 autophagy genes were upregulated, while 12 were downregulated in LUAD compared with normal tissues. Gene functional enrichment analysis indicated that autophagy-associated mechanisms were frequently implicated. “Autophagy”, “process utilizing autophagic mechanism”, and “macroautophagy” were among the top 20 biological processes; “autophagosome” and “autophagosome membrane” were the enriched cellular components; and “autophagy-animal” was the significant KEGG pathway (Figures 2(a), 2(b), and 2(d)). Moreover, we found that these DEAGs were also enriched in carcinogenesis and immunology, such as “HIF-1 signaling pathway”, “bladder cancer”, “platinum drug resistance”, “EGFR tyrosine kinase inhibitor resistance”, “PD-L1 expression and PD-1 checkpoint pathway in cancer”, “T cell receptor signaling pathway”, and “Th17 cell differentiation” in KEGG analysis (Figure 2(d)).

3.2. DElncRNAs in LUAD

A total of 2997 DElncRNAs were screened between normal samples and lung adenocarcinoma samples at the threshold of and value < 0.05, among which 2346 DElncRNAs were upregulated while 651 DElncRNAs were downregulated. Based on the expression profiles of 497 lung adenocarcinoma tissues, 1183 ARDElncRNAs were identified with correlation and value < 0.05.

3.3. Prognosis-Associated ARDElncRNAs in LUAD

After the exclusion of LUAD patients with unavailable follow-up information or an OS duration less than 1 month, a total of 457 LUAD patients were included for univariate Cox regression analysis. Only 78 of 1183 ARDElncRNAs were found to be associated with the OS of LUAD patients, and expression correlations between DEAGs and OS-associated lncRNAs are summarized in Table 1. An interaction network was built on the 78 OS-associated ARDElncRNAs and 16 corresponding associated DEAGs, including ATIC, BIRC5, DAPK2, EIF4EBP1, ERBB2, ERO1A, FOS, GAPDH, GRID1, IFNG, MAP1LC3C, NLRC4, PARP1, PTK6, TMEM74, and VMP1, as shown in Figure 3.

3.4. OS-Associated ARDElncRNA-Based Signature for OS in LUAD

The 457 lung adenocarcinoma patients were randomly divided into a training group () and a validation group (). The 78 identified OS-associated ARDElncRNAs were first analyzed by multivariate Cox regression in the training group to build a prognostic model. A signature consisting of 13 ARDElncRNAs was built to predict the OS of LUAD, including AC010999.2, AC084117.1, AC016747.2, LINC01338, AL356608.1, AL162632.3, SEMA6A-AS2, AC012213.2, LINC01691, AC034102.8, TMPO-AS1, AL121829.2, and AC069542.1. The detailed signature was as follows: The training group was further divided into a low-risk group and a high-risk group by the median risk score of the 13-ARDElncRNA signature. The Kaplan-Meier survival analysis indicated a significant OS advantage of the low-risk group over the high-risk group, as shown in Figure 4(a) (), and the sensitivity and specificity of this model in predicting OS was favorable with a 5-year OS AUC of 0.854 (Figure 4(b)).

3.5. Validation of the ARDElncRNA-Based Signature for OS in LUAD

The predictive value of the ARDElncRNA-based signature was further evaluated in the validation and total groups. The survival curve analysis indicated significantly a favorable OS for the low-risk group compared with the high-risk group (Figure 4(c), ), with a 5-year OS AUC of 0.737 (Figure 4(d)) in the validation group. Moreover, in the total LUAD group combining the training and validation groups, there was still a significant difference between the low-risk and high-risk groups (Figure 4(e), ), with a 5-year OS AUC of 0.811 (Figure 4(f)). The distribution of the risk scores, OS statuses, and OS times of the 457 total included patients are shown in Figures 5(a) and 5(b). The heatmap in Figure 5(c) shows the expression distributions of the 13 ARDElncRNAs in the low-risk and high-risk groups, with the color changing from green to red, indicating rising trends from low expression to high expression levels.

3.6. Risk Factors Predictive of Overall Mortality in LUAD

Univariate Cox regression analysis was conducted to investigate the influences of clinicopathological factors, including the age at diagnosis, sex, T status, N status, M status, stage, surgical margin resection status, and risk score of the ARDElncRNA-based signature. Age at diagnosis and sex were not significant factors associated with OS; however, T2 (vs T1, HR 1.535, 95% CI 1.063-2.216, ), T3/4 (vs T1, HR 3.208, 95% CI 2.001-5.142, ), positive N status (vs negative N status, HR 2.693, 95% CI 1.985-3.653, ), positive M status (vs negative M status, HR 2.236, 95% CI 1.304-3.835, ), stage II (vs stage I, HR 2.765, 95% CI 1.903-4.017, ), stage III/IV (vs stage I, HR 3.629, 95% CI 2.524-5.216, ), positive surgical margin (vs negative surgical margin, HR 4-027, 95% CI 2.249-7.212, ), and high risk score (vs low risk score, HR 2.036, 95% CI 1.499-2.767, ) were significant factors associated with OS (Table 1 and Figure 6(a)). Multivariate Cox regression analysis indicated that only positive surgical margin (vs negative surgical margin, HR 3.428, 95% CI 1.808-6.498, ) and high risk score (vs low risk score, HR 1.823, 95% CI 1.315-2.526, ) were independent risk factors predictive of overall mortality in LUAD (Table 2 and Figure 6(b)). The flowchart of this study is shown in Figure 7.

4. Discussion

Autophagy plays a dual role in suppressing and promoting initiation or progression in different phases of cancer. Autophagy was found to mediate the secretion of immune-modulating factors that promote cellular proliferation and leads to an invasive cancer phenotype [23]. The upregulation of autophagy facilitates cancer survival under stress circumstances and increases cancer growth and aggressiveness; therefore, efforts to inhibit autophagy to improve cancer therapeutic effects have attracted great interest [24].

Growing evidence indicates a close correlation between autophagy and lung cancer. PAQR3 was demonstrated to suppress the tumor progression of NSCLC cells by modulating EGFR-regulated autophagy [25]. The downregulation of autophagy facilitated the anti-LUAD efficacy of Shh pathway suppression, thus highlighting a potential approach for LUAD therapy [26]. miR-150-mediated autophagy dysfunctions were found to induce ER stress and the DNA damage response and contribute to NSCLC development [27]. With regard to cisplatin resistance in NSCLC, the autophagy inhibition of cancer stem cells identified by CD133 expression could promote the efficacy of cisplatin against NSCLC [28].

This study initially identified 30 DEAGs between 497 LUAD tissues and 54 normal tissues. Gene functional enrichment analysis revealed that the main mechanisms were involved in carcinogenesis, including bladder cancer, EGFR tyrosine kinase inhibitor resistance, and platinum drug resistance. Moreover, the immunological pathways were also suggested to correlate with the DEAGs. Targeting autophagy was found to be an alternative and novel strategy in cancer immunology. Rocaglamide, a natural product, could enhance the natural killer cell-mediated lysis of NSCLC cells by targeting ULK1, which is required for autophagy initiation and autophagy inhibition [29]. SIRPαD1-Fc, as a CD47-targeting fusion protein, promoted macrophage-mediated phagocytosis and cytotoxicity by inhibiting autophagy, which highlighted a potential approach for NSCLC treatment involving simultaneously targeting CD47 and autophagy [30].

Considering the increasingly significant role of lncRNAs in cancer, this study screened 2997 DElncRNAs between cancer and normal samples, and 1183 ARDElncRNAs associated with 30 DEAGs were identified. The correlations between lncRNAs and autophagy in cancer biology have been widely investigated. The upregulation of the lncRNA GAS5 was found to enhance cisplatin sensitivity in NSCLC by inhibiting autophagy [31]. The overexpression of the lncRNA NBAT1 also inhibited autophagy by interacting with PSMD10 and suppressing ATG7 transcription in NSCLC cells, which led to reduced cell viability, clonogenicity, and chemoresistance [32]. The lncRNA MSTO2P promoted lung cancer cell proliferation and autophagy by upregulating EZH2 [33]. The lncRNA BLACAT1 promoted ATG7 expression through miR-17, facilitated autophagy, and promoted the chemoresistance of NSCLC cells through the miR-17/ATG7 signaling pathway [34]. In this study, 78 OS-associated ARDElncRNAs and 16 corresponding associated DEAGs were identified and provided more molecular targets to investigate the underlying mechanism on carcinogenesis and progression of LUAD.

Taking clinical and prognostic information into consideration, 78 of 1183 ARDElncRNAs were associated with OS in LUAD. A signature consisting of 13 ARDElncRNAs was developed from the training cohort and further validated in the validation group and total group. Univariate and multivariate Cox regression analyses were conducted to investigate the influences of clinicopathological factors and the risk score of the ARDElncRNA-based signature, and only positive surgical margin (vs negative surgical margin) and high score (vs low score) based on the ARDElncRNA-based signature were independent risk factors predictive of overall mortality in LUAD.

However, there were some limitations in our study. First, the ARDElncRNAs were identified based on the expression correlations between DElncRNAs and DEAGs. The underlying mechanisms and molecular correlations between ARDElncRNAs and autophagy need to be investigated. Moreover, the ARDElncRNA signature was developed and validated using a retrospective cohort from TCGA, and the predictive efficacy needs to be further proven in other prospective cohorts.

5. Conclusions

This study first investigated the correlation between autophagy-associated lncRNAs and prognosis in LUAD and built a robust predictive signature of 13 ARDElncRNAs to predict OS.

Abbreviations

TCGA:The Cancer Genome Atlas
DEAGs:Differentially expressed autophagy genes
DElncRNAs:Differentially expressed lncRNAs
OS:Overall survival
AUC:Areas under curve
ROC:Receiver operatoring characteristic
SCLC:Small cell lung cancer
NSCLC:Non-small-cell lung cancer
LUAD:Lung adenocarcinoma
KEGG:Kyoto Encyclopedia of Genes and Genomes
ARDElncRNAs:Autophagy-related DElncRNAs
FC:Fold change.

Data Availability

The authors declare that the data supporting the findings of this study is available within the article.

Conflicts of Interest

The authors declare that they have no competing interests.

Authors’ Contributions

XZ designed this study; MZ, WHS, and HYD collected the data; MZ and XZ performed the statistical analysis; MZ and XZ wrote the paper. All authors read and approved the final version of the paper.

Acknowledgments

This study was supported by the Chongqing Science and Technology Commission (No. cstc2019jcyj-msxmX0841).