# An Autophagy-Related Long Noncoding RNA Signature Contributes to Poor Prognosis in Colorectal Cancer

Authors: Jingsun Wei; Xiaoxu Ge; Yang Tang; Yucheng Qian; Wei Lu; Kai Jiang; Yimin Fang; Maxwell Hwang; Dongliang Fu; Qian Xiao; Kefeng Ding
Journal: Journal of Oncology
 (2020)

---

## Abstract

Purpose. Colorectal cancer is one of the most common malignant primary tumors, prone to metastasis, and associated with a poor prognosis. As autophagy is closely related to the development and treatment of colorectal cancer, we investigated the potential prognostic value of long noncoding RNA (lncRNA) associated with autophagy in colorectal cancer. Methods. In this study, we acquired information on the expression of lncRNAs in colorectal cancer from the Cancer Genome Atlas (TCGA) database and found that 860 lncRNAs were associated with autophagy-related genes. Subsequently, univariate Cox regression analysis was used to investigate 32 autophagy-related lncRNAs linked to colon cancer prognosis. Subsequently, eight of the 32 autophagy-related lncRNAs (i.e., long intergenic nonprotein coding RNA 1503 [LINC01503], ZEB1 antisense RNA 1 [ZEB1-AS1], AC087481.3, AC008760.1, AC073896.3, AL138756.1, AL022323.1, and TNFRSF10A-AS1) were selected through multivariate Cox regression analysis. Based on these autophagy-related lncRNAs, a risk signature was constructed, and the patients were divided into high- and low-risk groups. Results. The high-risk group’s overall survival time was significantly shorter than that of the low-risk group p<0.0001. Receiver operating characteristic curve analysis was performed to further confirm the validity of the model (area under the curve: 0.689). Moreover, multivariate regression suggested that the risk score was a significant prognostic risk factor in colorectal cancer. Gene set enrichment analysis showed that these gene sets are significantly enriched in cancer-related pathways, such as Kirsten rat sarcoma viral oncogene homolog (KRAS) signaling. Conclusion. The risk signature of eight autophagy-related lncRNAs has prognostic potential for colorectal cancer. These autophagy-related lncRNAs may play a vital role in the biology of colorectal cancer.


---

## Body

## 1. Introduction

Colorectal cancer, a disease caused by the interaction of genetic and environmental factors, is one of the most common malignancies worldwide [1]. Statistics show that, in 2017, nearly 140,000 individuals in the United States of America were diagnosed with colon cancer, resulting in 50,000 deaths [2]. In China, colorectal cancer is one of the most common tumors, and its incidence is gradually increasing [3]. Surgery combined with postoperative chemotherapy is the current main treatment for colorectal cancer [4]. The lack of early definitive diagnosis and metastasis is the major cause of death due to colon cancer. In 2014, the 5-year survival rate of patients with local disease was 90.3%, while that of patients with local and distant metastases was 70.4 and 12.5%, respectively [5]. Therefore, screening of the colorectal cancer-related risk signature to evaluate the prognosis is of great significance.Autophagy is a “self-eating” phenomenon in cells, a physiological process in which membranes enclose organelles and proteins in cells and degrade them [6]. The effect of autophagy on cells is dual; it can promote and inhibit the formation and development of tumors and play varied roles in different tumors. In the early stage of tumor development, autophagy plays a role in inhibiting tumor growth. In the later stages, autophagy helps tumor cells grow in a low vascularization environment [7, 8]. In recent years, many studies found targets in autophagy-related pathways to treat cancer [9]. Besides, studies have confirmed that autophagy-related protein Beclin1 is upregulated in colorectal cancer [10]. This finding indicated that autophagy plays a vital role in this disease.Long noncoding RNA (LncRNA) has a variety of regulatory modes. It mainly interacts with proteins, DNA, or RNA to regulate gene expression at the epigenetic, transcriptional, and posttranscriptional levels, and participates in many biological activities, such as transcriptional activation, transcriptional interference, and nuclear transport [11]. LncRNA can regulate the function and activity of autophagy-related DNA, RNA, or protein or affect autophagy-related stress factors and energy receptors, thereby participating in the regulation of cell autophagy [12, 13]. For example, HAGLROS is a 699 bp lncRNA that can inhibit cell autophagy through the mammalian target of the rapamycin (mTOR) signaling pathway, thereby promoting tumorigenesis and progression [14, 15]. Previous studies have shown that lncRNA can be used as a potential biomarker for COAD prognosis [16]. However, the autophagy-related lncRNAs are rarely studied, the precise mechanisms of autophagy-related lncRNA in colorectal cancer remain unknown; hence, further investigations are warranted.At present, research on the relationship between autophagy-related lncRNA and tumors has received extensive attention. In this study, autophagy-related lncRNAs may have latent value in the prognosis of patients with colorectal cancer, and they can be potential therapeutic targets. The present study’s objective was to establish autophagy-related lncRNA markers and provide clues for understanding their role in colorectal cancer.


## 2. Methods

### 2.1. Acquiring Information of Patients with Colorectal Cancer

The RNA sequence transcriptome data and relevant clinical information of patients with colon adenocarcinoma (COAD) and rectum adenocarcinoma (READ) were acquired from the Cancer Genome Atlas (TCGA;https://cancergenome.nih.gov/) database.


### 2.2. Selection of lncRNA and Autophagy Genes

The TCGA database was used to obtain the profiles of lncRNAs in patients with colorectal cancer. The autophagy genes were obtained from the Human Autophagy Database (http://autophagy.lu/clustering/index.html). All expression data were normalized using the limma package (version 3.22.7). Pearson’s correlation test was performed to study the relationship between the autophagy-related genes (ATGs) and lncRNAs. A correlation coefficient  |R2| >0.4 and p<0.001 confirmed a correlation of the lncRNA with an ATG.


### 2.3. Establishment of the Risk Signature

We evaluated the prognostic value of autophagy-related lncRNAs based on the univariate and multivariate Cox regression analyses. We initially identified lncRNAs withp<0.05 through univariate analysis and subsequently incorporated those lncRNAs into a multivariate Cox regression analysis to establish risk scores. The following formula was used to calculate the risk score for each patient:(1)risk score=βgene1×exprgene1+βgene2×exprgene2+…+βgenen×exprgenen.The risk model was constructed to predict the survival of patients with colon cancer through Cox analysis. The risk scores were assigned by linear combination of the expression of the lncRNAs, which was weighted according to the regression coefficient (β). The Akaike information criterion was utilized to optimize the risk model. The subgroups (low- and high-risk groups) were established based on the median risk score. The expression level of lncRNA is defined as exp rgenen.


### 2.4. Gene Set Enrichment Analysis (GSEA)

The gene expression data were explored by GSEA. The GSEA software was downloaded to analyze the gene sets and deduce the function of the data. This approach can be utilized to study whether there was a statistically significant difference between the high- and low-risk groups. In the present study, we verified whether differentially expressed genes between the two subgroups were enriched during the autophagy process.


### 2.5. Cell Culture and Real-Time PCR

All the cells we used were cultured in RPMI‐1640 (Gibco, USA) with 10% fetal bovine serum (FBS, Life Technologies, Carlsbad, CA), 100°U/mL penicillin, and 100 mg/mL streptomycin. The cells were cultured at 37°C with 5% CO2. All cells were acquired from the American Type Culture Collection (ATCC, Rockville, MD).TRIzol Reagent (Invitrogen, USA) was performed to extract the RNA of the cell lines. Then, the cDNAs were acquired by the reverse transcription (PrimeScript™ RT Master Mix, Takara, Japan). Subsequently, qPCR was used using SYBR® Premix Ex Taq™ GC (Takara, Japan). And, the GAPDH was chosen as the internal control. The primers of the eight autophagy lncRNAs are listed in Table S1. The method of ΔCt was applied to analyze the data.


### 2.6. Statistical Analysis

The autophagy-lncRNA co-expression network was constructed using the CYTOSCAPE software (version 3.5.1;https://cytoscape.org/) [17]. Pearson’s correlation, Cox regression, and Kaplan-Meier curve analyses were performed using the R software (version 3.6; https://www.r-project.org/). The functional annotation indicated that the risk model was effe of the two subgroups was accomplished using the GSEA software (version 4.0.3; http://www.broadinstitute.org/gsea/index.jsp). p values < 0.05 denoted statistically significant differences.


## 2.1. Acquiring Information of Patients with Colorectal Cancer

The RNA sequence transcriptome data and relevant clinical information of patients with colon adenocarcinoma (COAD) and rectum adenocarcinoma (READ) were acquired from the Cancer Genome Atlas (TCGA;https://cancergenome.nih.gov/) database.


## 2.2. Selection of lncRNA and Autophagy Genes

The TCGA database was used to obtain the profiles of lncRNAs in patients with colorectal cancer. The autophagy genes were obtained from the Human Autophagy Database (http://autophagy.lu/clustering/index.html). All expression data were normalized using the limma package (version 3.22.7). Pearson’s correlation test was performed to study the relationship between the autophagy-related genes (ATGs) and lncRNAs. A correlation coefficient  |R2| >0.4 and p<0.001 confirmed a correlation of the lncRNA with an ATG.


## 2.3. Establishment of the Risk Signature

We evaluated the prognostic value of autophagy-related lncRNAs based on the univariate and multivariate Cox regression analyses. We initially identified lncRNAs withp<0.05 through univariate analysis and subsequently incorporated those lncRNAs into a multivariate Cox regression analysis to establish risk scores. The following formula was used to calculate the risk score for each patient:(1)risk score=βgene1×exprgene1+βgene2×exprgene2+…+βgenen×exprgenen.The risk model was constructed to predict the survival of patients with colon cancer through Cox analysis. The risk scores were assigned by linear combination of the expression of the lncRNAs, which was weighted according to the regression coefficient (β). The Akaike information criterion was utilized to optimize the risk model. The subgroups (low- and high-risk groups) were established based on the median risk score. The expression level of lncRNA is defined as exp rgenen.


## 2.4. Gene Set Enrichment Analysis (GSEA)

The gene expression data were explored by GSEA. The GSEA software was downloaded to analyze the gene sets and deduce the function of the data. This approach can be utilized to study whether there was a statistically significant difference between the high- and low-risk groups. In the present study, we verified whether differentially expressed genes between the two subgroups were enriched during the autophagy process.


## 2.5. Cell Culture and Real-Time PCR

All the cells we used were cultured in RPMI‐1640 (Gibco, USA) with 10% fetal bovine serum (FBS, Life Technologies, Carlsbad, CA), 100°U/mL penicillin, and 100 mg/mL streptomycin. The cells were cultured at 37°C with 5% CO2. All cells were acquired from the American Type Culture Collection (ATCC, Rockville, MD).TRIzol Reagent (Invitrogen, USA) was performed to extract the RNA of the cell lines. Then, the cDNAs were acquired by the reverse transcription (PrimeScript™ RT Master Mix, Takara, Japan). Subsequently, qPCR was used using SYBR® Premix Ex Taq™ GC (Takara, Japan). And, the GAPDH was chosen as the internal control. The primers of the eight autophagy lncRNAs are listed in Table S1. The method of ΔCt was applied to analyze the data.


## 2.6. Statistical Analysis

The autophagy-lncRNA co-expression network was constructed using the CYTOSCAPE software (version 3.5.1;https://cytoscape.org/) [17]. Pearson’s correlation, Cox regression, and Kaplan-Meier curve analyses were performed using the R software (version 3.6; https://www.r-project.org/). The functional annotation indicated that the risk model was effe of the two subgroups was accomplished using the GSEA software (version 4.0.3; http://www.broadinstitute.org/gsea/index.jsp). p values < 0.05 denoted statistically significant differences.


## 3. Results

### 3.1. Establishment of the Autophagy lncRNAs Co-Expression Network

A total of 257 ATGs were acquired from the Human Autophagy Database. Subsequently, two expression matrices of lncRNAs and ATGs in colorectal cancer were obtained through the TCGA database. We obtained co-expression results by performing a correlation analysis of lncRNAs and ATGs. The expression network of autophagy lncRNAs was established to discern the lncRNAs associated with autophagy. As a result, 860 autophagy-related lncRNAs were identified based on the filtering criteria of a correlation coefficient <0.4 andp<0.001. According to the value of the correlation coefficient, the top 10 lncRNAs are presented in Table 1.Table 1
Top 10 negative and positive co-expression of autophagy genes and related lncRNAs based on the value of the correlation coefficient.
ARGgenelncRNACor.p valueNegative-corMLST8EBLN3P−0.5452412192.66E − 45CAPN1PSMA3-AS1−0.5325150976.47E − 43CAPN1EBLN3P−0.5290227752.81E − 42MLST8OIP5-AS1−0.5281872363.98E − 42CAPN1AC125257.1−0.5175550713.08E − 40GABARAPAC020663.2−0.5149001248.92E − 40CAPN1OIP5-AS1−0.5141924011.18E − 39MLST8PSMA3-AS1−0.5083774631.17E − 38GAPDHNORAD−0.5071328741.89E − 38CAPNS1AL049840.3−0.5063564292.56E − 38Positive-corBIRC6Z68871.10.8573589512.19E − 165WDFY3AL157392.30.8452480673.62E − 156WDFY3AC097376.20.8448227767.37E − 156NAMPTHIF1A-AS20.8431250641.24E − 154WDFY3AC004492.10.8419500398.55E − 154WDFY3SCARNA90.8415461271.66E − 153WDFY3AC021078.10.8408786414.92E − 153WDFY3USP12-AS10.8408615695.05E − 153PTENAL163051.20.8378041566.92E − 151WDFY3SDCBP2-AS10.8342703571.80E − 148


### 3.2. Construction of a Risk Model Based on Autophagy-Related lncRNAs in Patients with Colorectal Cancer

After constructing the co-expression network, according to the 860 autophagy-related lncRNAs, we further assessed the prognostic value of lncRNAs through univariate Cox regression analysis. A total of 548 patients with colorectal cancer in the TCGA database were used for the analysis. The significance filtering criterion was set atp<0.05, and 32 autophagy-related lncRNAs were identified (Table 2).Table 2
The information of 32 selected lncRNAs significantly related to overall survival in colorectal cancer.
lncRNAKMBSEHRHR.95 LHR.95Hp valueLINC015030.0440.2330.0761.2621.0871.4662.25E − 03AC026471.40.0070.1300.0401.1391.0531.2321.15E − 03AC007128.10.0350.3160.1601.3721.0031.8764.79E − 02AC011462.40.0340.4180.1231.5191.1941.9326.71E − 04AL353194.10.0380.1480.0751.1591.0001.3435.00E − 02AC099850.30.024−0.0620.0250.9400.8950.9861.17E − 02AC004148.20.0040.1440.0731.1551.0021.3324.67E − 02AC068580.30.0020.3360.1251.4001.0961.7887.11E − 03ZEB1-AS10.0010.7660.1722.1501.5343.0148.74E − 06AC087481.30.0380.1770.0741.1941.0331.3811.67E − 02AL354836.10.0050.1140.0381.1211.0411.2062.39E − 03AC005261.30.0170.1760.0771.1921.0261.3862.20E − 02AC027796.40.0210.3810.0921.4641.2221.7533.42E − 05AC083843.20.0180.0690.0351.0721.0011.1474.53E − 02LINC023810.0440.1890.0791.2081.0351.4101.64E − 02AC068580.10.0330.2610.1041.2981.0581.5931.23E − 02LINC010110.0300.3500.1441.4201.0701.8831.51E − 02AC010973.20.0080.5950.1281.8131.4122.3303.19E − 06AC087741.10.0290.4210.1431.5241.1522.0153.13E − 03AC107375.10.0070.2350.1011.2651.0381.5421.98E − 02AC008760.10.0320.4500.1191.5681.2431.9801.52E − 04AC073896.30.024−0.5170.2210.5960.3870.9191.91E − 02SNHG70.0120.0350.0141.0361.0071.0661.38E − 02AL138756.10.0270.1980.0741.2191.0541.4107.54E − 03AL022323.10.004−0.3320.1690.7170.5150.9984.88E − 02AL512413.10.0100.3490.1281.4181.1021.8246.57E − 03AL161729.40.0490.3210.0971.3781.1411.6658.94E − 04AL162586.10.0000.4020.1071.4951.2111.8451.79E − 04TNFRSF10A-AS10.003−0.1380.0570.8710.7800.9731.43E − 02LINC001740.0130.2490.0951.2831.0651.5448.62E − 03AL022328.20.0470.2030.0971.2241.0131.4803.61E − 02AC020558.20.0130.3780.1591.4591.0681.9911.75E − 02Subsequently, multivariate Cox regression analysis was performed to ulteriorly identify the independent prognostic factors among the 32 lncRNAs. When the Akaike information criterion value reached the optimal value, eight of the 32 lncRNAs were identified as independent prognostic factors: long intergenic nonprotein coding RNA 1503 (LINC01503), ZEB1 antisense RNA 1 (ZEB1-AS1), AC087481.3, AC008760.1, AC073896.3, AL138756.1, AL022323.1, and TNFRSF10A-AS1 (Table3).Table 3
The coefficients and hazard ratio (HR) value of 8 autophagy-related lncRNAs estimated by multivariate Cox regression.
lncRNACoef.HRLINC015030.2219716511.248535983ZEB1-AS10.5218792311.68519154AC087481.30.2006657791.222216213AC008760.10.365005531.440521974AC073896.3−0.5894650540.554623899AL138756.10.1611738451.1748892AL022323.1−0.359889480.697753438TNFRSF10A-AS1−0.1419951460.867625468Furthermore, we constructed a network using these eight lncRNAs to exhibit the relationship between the prognostic lncRNAs and ATGs (Figure1). Next, we verified the relationship between these eight autophagy-related lncRNAs and patients’ prognosis with colorectal cancer through a survival curve (Figure 2). All eight lncRNAs were significantly associated with the prognosis of colorectal cancer p<0.05. According to the risk score acquired from the multivariate Cox regression analysis, we produced a Sankey diagram to exhibit the association among ATGs, autophagy-related lncRNAs, and risk types of lncRNAs (Figure 3). The results showed that three autophagy-related lncRNAs were protective factors (AC073896.3, AL022323.1, and TNFRSF10A-AS1), while the remaining five autophagy-related lncRNAs were risk factors (LINC01503, ZEB1-AS1, AC087481.3, AC008760.1, and AL138756.1).Figure 1
Co-expression network of autophagy genes and related prognostic lncRNAs. Red nodes represent autophagy-related lncRNAs. Blue nodes represent autophagy genes. LncRNA, long noncoding RNA.Figure 2
Kaplan–Meier survival curve of the eight selected autophagy-related lncRNAs in colorectal cancer in TCGA database. The Kaplan–Meier survival analysis indicated that patients with high expression of AC008760.1, AC087481.3, AL138756.1, LINC01503, and ZEB1-AS1 had worse prognosis. In contrast, those with high expression of AC073896.3, AL022323.1, and TNFRSF10A-AS1 had better prognosis. LncRNA, long noncoding RNA; TCGA, the Cancer Genome Atlas; LINC01503, long intergenic nonprotein coding RNA 1503; AS1, antisense RNA 1.
(a)(b)(c)(d)(e)(f)(g)(h)Figure 3
The relationships among autophagy-related genes, autophagy-related lncRNAs, and risk types in the Sankey diagram. lncRNA, long noncoding RNA.


### 3.3. Prognostic Effects of the Autophagy-Related lncRNA Risk Model on Colorectal Cancer

Subsequently, the autophagy-related lncRNA risk model was established based on the risk score. The patients with colorectal cancer were divided into two groups, namely, high- and low-risk groups (Figures4(a)–4(c)). Consequently, the risk score was a significant predictor of overall survival (OS) in patients with colorectal cancer, with the low-risk group exhibiting a longer OS than that observed in the high-risk group p<0.05. In addition, the log-rank test was used to obtain the Kaplan–Meier survival curve. The results also suggested that the high-risk group had a worse prognosis (Figure 5(a)). The receiver operating curve was used to measure the risk model’s effectiveness, and the area under the curve (0.689) indicated that the risk model was effective (Figure 5(b)). We also used the previously published lncRNAs signature [16] compared with our signature. The results showed that the previously published lncRNA signature was significantly related to OS (Figure S1A); however, the area under the curve was 0.436 (Figure S1B), indicating our signature’s superiority based on autophagy-related lncRNAs. In addition, univariate and multivariate Cox regression analyses were performed based on the risk score. Both regression analyses demonstrated that the risk score was a significant prognostic risk factor (univariate regression: p<0.05, hazard ratio: 1.384; multivariate regression: p<0.05, hazard ratio: 1.329) (Figures 6(a) and 6(b)). Collectively, the results suggested that the risk score obtained from our signature can be used as an independent prognostic factor to estimate the OS of patients with colorectal cancer.Figure 4
Risk analysis of autophagy-related lncRNAs in patients with colorectal cancer in TCGA database. (a) The risk score in the high- and low-risk groups for the risk signature in colorectal cancer. (b) Survival time of patients with colorectal cancer in high- and low-risk groups. (c) Expression of the eight selected autophagy-related lncRNAs in colorectal cancer. Red color represents increased expression. Blue color represents decreased expression. LncRNA, long noncoding RNA; TCGA, the Cancer Genome Atlas.
(a)(b)(c)Figure 5
Prognostic value of the risk signature based on coefficients and HR. (a) Survival analysis of the high- and low-risk groups according to the risk model in TCGA database. (b) The accuracy of the model was verified using multi-indicator ROC analysis. HR, hazard ratio; TCGA, the Cancer Genome Atlas; ROC, receiver operating characteristic.
(a)(b)Figure 6
Analysis of the risk signature as an independent prognostic indicator through (a) univariate and (b) multivariate cox regression analyses.
(a)(b)


### 3.4. The Risk Model Was Closely Related to Clinicopathological Features in Colorectal Cancer

After detecting the significant association between the risk score and the prognosis of colorectal cancer, we further investigated the relationship between the risk score and clinical traits in the clinical information TCGA database. The clinical characteristics included age, sex, stage status, and tumor-node-metastasis status. As shown in Table4, there was a strong relationship between the risk model of the eight lncRNAs and clinicopathological features, such as stage status and tumor-node-metastasis stage p<0.05. The risk score in the late stages (stage III-IV) tends to be higher compared with that obtained for the early stages (stage I–II), and it was the same for the T (T1–2, T3–4), N (N0, N1–2), and M (M0, M1) stages. These results suggested that the risk score may be closely related to the progression of colorectal cancer. Interestingly, we also found that males had a higher risk score compared with females; this observation is consistent with the epidemiology of colorectal cancer.Table 4
The association between the clinicopathological characteristics and risk score in colorectal cancer.
ClinicalGroupnMeanSDTpAge≤652021.2271.32−0.9460.345>652711.3361.119GenderFemale2221.1590.912−2.2700.024Male2511.4041.413StageStage I-II2741.0720.981−4.4340Stage III-IV1991.5881.415TT1-2941.0550.707−2.9610.003T3-43791.3471.298MM03731.1821.096−3.1580.002M11001.6891.501NN02831.0860.985−4.2410N1-21901.5921.432


### 3.5. GSEA

The GSEA was used to explore further functional annotation. The results showed that the differentially expressed genes were significantly enriched in some tumor-related pathways between the high- and low-risk groups. Twelve gene sets were upregulated in the high-risk group and two gene sets were significantly enriched at nominalp value <0.05 and false discovery rate <25%. Figure 7 shows that Kirsten rat sarcoma viral oncogene homolog (KRAS) signaling and myogenic gene sets were significantly enriched in the high-risk group (Figures 7(a) and 7(b)). The G2°M checkpoint and E2F target gene sets were significantly enriched in the low-risk group (Figures 7(c) and 7(d)). The results based on the autophagy-related lncRNAs showed that the KRAS signaling pathway might play an important role in autophagy in colorectal cancer.Figure 7
Gene set enrichment analysis showing the enriched hallmark pathways in the two groups based on TCGA database. (a), (b) KRAS and myogenic gene sets were significantly enriched in the high-risk group. (c), (d) G2°M checkpoint and E2F targets gene sets were significantly enriched in the low-risk group. TCGA, the Cancer Genome Atlas; KRAS, Kirsten rat sarcoma viral oncogene homolog.
(a)(b)(c)(d)


### 3.6. Expression of the lncRNAs in Cell Lines of Colorectal Cancer

The Q-PCR assay was established to study the eight autophagy-related lncRNAs in cell lines of colorectal cancer (CRC). Seven cell lines of CRC were used in this study: SW620, DLD1, HCT116, LOVO, SW480, HT29, and RKO. The value ofΔCt was used to confirm the expression of the eight lncRNAs in cell lines. A higher value of ΔCt means a lower level of lncRNAs expression. The results indicated that five risk factors (LINC01503, ZEB1-AS1, AC087481.3, AC008760.1, and AL138756.1) were highly expressed in colorectal cancer, and two protective factors (AC073896.3 and AL022323.1) were low (Figure 8(a)). Figure 8(b) shows the mean value of ΔCt of the each autophagy-related lncRNAs in all cell lines: AC073896.3 = 15.03, AL022323.1 = 18.8, TNFRSF10A-AS1 = 9.2, LINC01503 = 10.65, ZEB1-AS1 = 9.28, AC087481.3 = 10.31, AC008760.1 = 13.16, and AL138756.1 = 11.23. The results were generally consistent with our work.Figure 8
The expression levels of the eight autophagy-related lncRNAs in CRC cell lines. (a) The expression of the 8 lncRNAs in each cell lines. (b) The total expression levels of the eight autophagy-related lncRNAs in 7 cell lines. LncRNA, long noncoding RNA; CRC, colorectal cancer.
(a)(b)


## 3.1. Establishment of the Autophagy lncRNAs Co-Expression Network

A total of 257 ATGs were acquired from the Human Autophagy Database. Subsequently, two expression matrices of lncRNAs and ATGs in colorectal cancer were obtained through the TCGA database. We obtained co-expression results by performing a correlation analysis of lncRNAs and ATGs. The expression network of autophagy lncRNAs was established to discern the lncRNAs associated with autophagy. As a result, 860 autophagy-related lncRNAs were identified based on the filtering criteria of a correlation coefficient <0.4 andp<0.001. According to the value of the correlation coefficient, the top 10 lncRNAs are presented in Table 1.Table 1
Top 10 negative and positive co-expression of autophagy genes and related lncRNAs based on the value of the correlation coefficient.
ARGgenelncRNACor.p valueNegative-corMLST8EBLN3P−0.5452412192.66E − 45CAPN1PSMA3-AS1−0.5325150976.47E − 43CAPN1EBLN3P−0.5290227752.81E − 42MLST8OIP5-AS1−0.5281872363.98E − 42CAPN1AC125257.1−0.5175550713.08E − 40GABARAPAC020663.2−0.5149001248.92E − 40CAPN1OIP5-AS1−0.5141924011.18E − 39MLST8PSMA3-AS1−0.5083774631.17E − 38GAPDHNORAD−0.5071328741.89E − 38CAPNS1AL049840.3−0.5063564292.56E − 38Positive-corBIRC6Z68871.10.8573589512.19E − 165WDFY3AL157392.30.8452480673.62E − 156WDFY3AC097376.20.8448227767.37E − 156NAMPTHIF1A-AS20.8431250641.24E − 154WDFY3AC004492.10.8419500398.55E − 154WDFY3SCARNA90.8415461271.66E − 153WDFY3AC021078.10.8408786414.92E − 153WDFY3USP12-AS10.8408615695.05E − 153PTENAL163051.20.8378041566.92E − 151WDFY3SDCBP2-AS10.8342703571.80E − 148


## 3.2. Construction of a Risk Model Based on Autophagy-Related lncRNAs in Patients with Colorectal Cancer

After constructing the co-expression network, according to the 860 autophagy-related lncRNAs, we further assessed the prognostic value of lncRNAs through univariate Cox regression analysis. A total of 548 patients with colorectal cancer in the TCGA database were used for the analysis. The significance filtering criterion was set atp<0.05, and 32 autophagy-related lncRNAs were identified (Table 2).Table 2
The information of 32 selected lncRNAs significantly related to overall survival in colorectal cancer.
lncRNAKMBSEHRHR.95 LHR.95Hp valueLINC015030.0440.2330.0761.2621.0871.4662.25E − 03AC026471.40.0070.1300.0401.1391.0531.2321.15E − 03AC007128.10.0350.3160.1601.3721.0031.8764.79E − 02AC011462.40.0340.4180.1231.5191.1941.9326.71E − 04AL353194.10.0380.1480.0751.1591.0001.3435.00E − 02AC099850.30.024−0.0620.0250.9400.8950.9861.17E − 02AC004148.20.0040.1440.0731.1551.0021.3324.67E − 02AC068580.30.0020.3360.1251.4001.0961.7887.11E − 03ZEB1-AS10.0010.7660.1722.1501.5343.0148.74E − 06AC087481.30.0380.1770.0741.1941.0331.3811.67E − 02AL354836.10.0050.1140.0381.1211.0411.2062.39E − 03AC005261.30.0170.1760.0771.1921.0261.3862.20E − 02AC027796.40.0210.3810.0921.4641.2221.7533.42E − 05AC083843.20.0180.0690.0351.0721.0011.1474.53E − 02LINC023810.0440.1890.0791.2081.0351.4101.64E − 02AC068580.10.0330.2610.1041.2981.0581.5931.23E − 02LINC010110.0300.3500.1441.4201.0701.8831.51E − 02AC010973.20.0080.5950.1281.8131.4122.3303.19E − 06AC087741.10.0290.4210.1431.5241.1522.0153.13E − 03AC107375.10.0070.2350.1011.2651.0381.5421.98E − 02AC008760.10.0320.4500.1191.5681.2431.9801.52E − 04AC073896.30.024−0.5170.2210.5960.3870.9191.91E − 02SNHG70.0120.0350.0141.0361.0071.0661.38E − 02AL138756.10.0270.1980.0741.2191.0541.4107.54E − 03AL022323.10.004−0.3320.1690.7170.5150.9984.88E − 02AL512413.10.0100.3490.1281.4181.1021.8246.57E − 03AL161729.40.0490.3210.0971.3781.1411.6658.94E − 04AL162586.10.0000.4020.1071.4951.2111.8451.79E − 04TNFRSF10A-AS10.003−0.1380.0570.8710.7800.9731.43E − 02LINC001740.0130.2490.0951.2831.0651.5448.62E − 03AL022328.20.0470.2030.0971.2241.0131.4803.61E − 02AC020558.20.0130.3780.1591.4591.0681.9911.75E − 02Subsequently, multivariate Cox regression analysis was performed to ulteriorly identify the independent prognostic factors among the 32 lncRNAs. When the Akaike information criterion value reached the optimal value, eight of the 32 lncRNAs were identified as independent prognostic factors: long intergenic nonprotein coding RNA 1503 (LINC01503), ZEB1 antisense RNA 1 (ZEB1-AS1), AC087481.3, AC008760.1, AC073896.3, AL138756.1, AL022323.1, and TNFRSF10A-AS1 (Table3).Table 3
The coefficients and hazard ratio (HR) value of 8 autophagy-related lncRNAs estimated by multivariate Cox regression.
lncRNACoef.HRLINC015030.2219716511.248535983ZEB1-AS10.5218792311.68519154AC087481.30.2006657791.222216213AC008760.10.365005531.440521974AC073896.3−0.5894650540.554623899AL138756.10.1611738451.1748892AL022323.1−0.359889480.697753438TNFRSF10A-AS1−0.1419951460.867625468Furthermore, we constructed a network using these eight lncRNAs to exhibit the relationship between the prognostic lncRNAs and ATGs (Figure1). Next, we verified the relationship between these eight autophagy-related lncRNAs and patients’ prognosis with colorectal cancer through a survival curve (Figure 2). All eight lncRNAs were significantly associated with the prognosis of colorectal cancer p<0.05. According to the risk score acquired from the multivariate Cox regression analysis, we produced a Sankey diagram to exhibit the association among ATGs, autophagy-related lncRNAs, and risk types of lncRNAs (Figure 3). The results showed that three autophagy-related lncRNAs were protective factors (AC073896.3, AL022323.1, and TNFRSF10A-AS1), while the remaining five autophagy-related lncRNAs were risk factors (LINC01503, ZEB1-AS1, AC087481.3, AC008760.1, and AL138756.1).Figure 1
Co-expression network of autophagy genes and related prognostic lncRNAs. Red nodes represent autophagy-related lncRNAs. Blue nodes represent autophagy genes. LncRNA, long noncoding RNA.Figure 2
Kaplan–Meier survival curve of the eight selected autophagy-related lncRNAs in colorectal cancer in TCGA database. The Kaplan–Meier survival analysis indicated that patients with high expression of AC008760.1, AC087481.3, AL138756.1, LINC01503, and ZEB1-AS1 had worse prognosis. In contrast, those with high expression of AC073896.3, AL022323.1, and TNFRSF10A-AS1 had better prognosis. LncRNA, long noncoding RNA; TCGA, the Cancer Genome Atlas; LINC01503, long intergenic nonprotein coding RNA 1503; AS1, antisense RNA 1.
(a)(b)(c)(d)(e)(f)(g)(h)Figure 3
The relationships among autophagy-related genes, autophagy-related lncRNAs, and risk types in the Sankey diagram. lncRNA, long noncoding RNA.


## 3.3. Prognostic Effects of the Autophagy-Related lncRNA Risk Model on Colorectal Cancer

Subsequently, the autophagy-related lncRNA risk model was established based on the risk score. The patients with colorectal cancer were divided into two groups, namely, high- and low-risk groups (Figures4(a)–4(c)). Consequently, the risk score was a significant predictor of overall survival (OS) in patients with colorectal cancer, with the low-risk group exhibiting a longer OS than that observed in the high-risk group p<0.05. In addition, the log-rank test was used to obtain the Kaplan–Meier survival curve. The results also suggested that the high-risk group had a worse prognosis (Figure 5(a)). The receiver operating curve was used to measure the risk model’s effectiveness, and the area under the curve (0.689) indicated that the risk model was effective (Figure 5(b)). We also used the previously published lncRNAs signature [16] compared with our signature. The results showed that the previously published lncRNA signature was significantly related to OS (Figure S1A); however, the area under the curve was 0.436 (Figure S1B), indicating our signature’s superiority based on autophagy-related lncRNAs. In addition, univariate and multivariate Cox regression analyses were performed based on the risk score. Both regression analyses demonstrated that the risk score was a significant prognostic risk factor (univariate regression: p<0.05, hazard ratio: 1.384; multivariate regression: p<0.05, hazard ratio: 1.329) (Figures 6(a) and 6(b)). Collectively, the results suggested that the risk score obtained from our signature can be used as an independent prognostic factor to estimate the OS of patients with colorectal cancer.Figure 4
Risk analysis of autophagy-related lncRNAs in patients with colorectal cancer in TCGA database. (a) The risk score in the high- and low-risk groups for the risk signature in colorectal cancer. (b) Survival time of patients with colorectal cancer in high- and low-risk groups. (c) Expression of the eight selected autophagy-related lncRNAs in colorectal cancer. Red color represents increased expression. Blue color represents decreased expression. LncRNA, long noncoding RNA; TCGA, the Cancer Genome Atlas.
(a)(b)(c)Figure 5
Prognostic value of the risk signature based on coefficients and HR. (a) Survival analysis of the high- and low-risk groups according to the risk model in TCGA database. (b) The accuracy of the model was verified using multi-indicator ROC analysis. HR, hazard ratio; TCGA, the Cancer Genome Atlas; ROC, receiver operating characteristic.
(a)(b)Figure 6
Analysis of the risk signature as an independent prognostic indicator through (a) univariate and (b) multivariate cox regression analyses.
(a)(b)


## 3.4. The Risk Model Was Closely Related to Clinicopathological Features in Colorectal Cancer

After detecting the significant association between the risk score and the prognosis of colorectal cancer, we further investigated the relationship between the risk score and clinical traits in the clinical information TCGA database. The clinical characteristics included age, sex, stage status, and tumor-node-metastasis status. As shown in Table4, there was a strong relationship between the risk model of the eight lncRNAs and clinicopathological features, such as stage status and tumor-node-metastasis stage p<0.05. The risk score in the late stages (stage III-IV) tends to be higher compared with that obtained for the early stages (stage I–II), and it was the same for the T (T1–2, T3–4), N (N0, N1–2), and M (M0, M1) stages. These results suggested that the risk score may be closely related to the progression of colorectal cancer. Interestingly, we also found that males had a higher risk score compared with females; this observation is consistent with the epidemiology of colorectal cancer.Table 4
The association between the clinicopathological characteristics and risk score in colorectal cancer.
ClinicalGroupnMeanSDTpAge≤652021.2271.32−0.9460.345>652711.3361.119GenderFemale2221.1590.912−2.2700.024Male2511.4041.413StageStage I-II2741.0720.981−4.4340Stage III-IV1991.5881.415TT1-2941.0550.707−2.9610.003T3-43791.3471.298MM03731.1821.096−3.1580.002M11001.6891.501NN02831.0860.985−4.2410N1-21901.5921.432


## 3.5. GSEA

The GSEA was used to explore further functional annotation. The results showed that the differentially expressed genes were significantly enriched in some tumor-related pathways between the high- and low-risk groups. Twelve gene sets were upregulated in the high-risk group and two gene sets were significantly enriched at nominalp value <0.05 and false discovery rate <25%. Figure 7 shows that Kirsten rat sarcoma viral oncogene homolog (KRAS) signaling and myogenic gene sets were significantly enriched in the high-risk group (Figures 7(a) and 7(b)). The G2°M checkpoint and E2F target gene sets were significantly enriched in the low-risk group (Figures 7(c) and 7(d)). The results based on the autophagy-related lncRNAs showed that the KRAS signaling pathway might play an important role in autophagy in colorectal cancer.Figure 7
Gene set enrichment analysis showing the enriched hallmark pathways in the two groups based on TCGA database. (a), (b) KRAS and myogenic gene sets were significantly enriched in the high-risk group. (c), (d) G2°M checkpoint and E2F targets gene sets were significantly enriched in the low-risk group. TCGA, the Cancer Genome Atlas; KRAS, Kirsten rat sarcoma viral oncogene homolog.
(a)(b)(c)(d)


## 3.6. Expression of the lncRNAs in Cell Lines of Colorectal Cancer

The Q-PCR assay was established to study the eight autophagy-related lncRNAs in cell lines of colorectal cancer (CRC). Seven cell lines of CRC were used in this study: SW620, DLD1, HCT116, LOVO, SW480, HT29, and RKO. The value ofΔCt was used to confirm the expression of the eight lncRNAs in cell lines. A higher value of ΔCt means a lower level of lncRNAs expression. The results indicated that five risk factors (LINC01503, ZEB1-AS1, AC087481.3, AC008760.1, and AL138756.1) were highly expressed in colorectal cancer, and two protective factors (AC073896.3 and AL022323.1) were low (Figure 8(a)). Figure 8(b) shows the mean value of ΔCt of the each autophagy-related lncRNAs in all cell lines: AC073896.3 = 15.03, AL022323.1 = 18.8, TNFRSF10A-AS1 = 9.2, LINC01503 = 10.65, ZEB1-AS1 = 9.28, AC087481.3 = 10.31, AC008760.1 = 13.16, and AL138756.1 = 11.23. The results were generally consistent with our work.Figure 8
The expression levels of the eight autophagy-related lncRNAs in CRC cell lines. (a) The expression of the 8 lncRNAs in each cell lines. (b) The total expression levels of the eight autophagy-related lncRNAs in 7 cell lines. LncRNA, long noncoding RNA; CRC, colorectal cancer.
(a)(b)


## 4. Discussion

Colorectal cancer is a common type of malignant tumor. Owing to its rapid progression and easy metastasis, it has become the third leading cause of tumor-related mortality worldwide [18].At present, colorectal cancer treatment involves surgery, chemotherapy, radiotherapy, biological therapy, etc. However, the effectiveness and prognosis are unsatisfactory. The 5-year survival rate of patients with colorectal cancer is only 42.7% [19]. Based on the widely utilized high-throughput biological technology, genomic alterations are accurately detected to predict the risk of metastasis and prognosis in patients with cancer [20, 21]. With the development of sequencing technology, investigators found that lncRNA plays an important role in transcription, posttranscriptional regulation, chromatin modification, splicing regulation, genomic imprinting, cell cycle regulation, and epigenetic regulation [22, 23]. This indicates that lncRNAs may be useful in predicting the prognosis of colorectal cancer. In previous studies, Zhou et al. found an effective lncRNA signature for predicting the risk of recurrence in colon tumors [24]. In addition, a risk model composed of mRNA and lncRNA could also be applied to detect the early recurrence of colon cancer [25]. These studies further prove the important role of lncRNA in tumors. Autophagy has been a hot topic in tumor research. However, there are few studies on the autophagy-related lncRNAs for patients with colorectal cancer. Therefore, it is necessary to construct the risk signature based on autophagy-related lncRNAs.In the present study, the TCGA dataset was used to research the value of autophagy-related lncRNAs to predict prognosis of colorectal cancer. Firstly, we selected 860 autophagy-related lncRNAs through the co-expression of lncRNA and autophagy genes. Subsequently, eight autophagy-related lncRNAs were identified through multivariate Cox regression analysis. The risk signature was determined according to these autophagy-related lncRNAs, which can classify patients with colorectal cancer into high- and low-risk groups. Our results indicated that OS in the high-risk group was shorter than that noted in the low-risk group. Receiver operating curve analysis was used to validate the accuracy of the risk signature. The area under the curve of the risk score and stage was 0.689 and 0.695, respectively, indicating the risk signature’s accuracy. Moreover, we further validated the risk model via univariate and multivariate Cox regression analyses. Moreover, we also performed Q-PCR assay to analyze the expression levels of the eight autophagy-related lncRNAs. Finally, we concluded that the risk signature based on these eight lncRNAs was an independent factor for colorectal cancer and was significantly associated with OS.Autophagy is a highly conserved regulatory mechanism in eukaryotic cell evolution and an important degradation pathway for almost all cells from yeast to mammals [26]. More than 30 ATGs (e.g., ATG5, ATG12, and ATG16) are closely related to autophagy. The ATG12 binding process and LC3 modification process play important roles during autophagy. ATG12 is conjugated with ATG5 with the assistance of ATG7 to form the ATG5–ATG12 complex. LC3-I is finally transformed into LC3-II with the assistance of this complex, promoting autophagosome formation [27–29]. As mentioned above, autophagy plays different roles in tumors based on the type, stage, and genetic background of cancer [30]. It has been found that autophagy promotes the growth of cancer cells in a variety of tumor cells [31]. Under a stress state, the degree of autophagy increases to enhance cell adaptability, which is conducive to cell survival [32]. Yang et al. found that autophagy inhibition led to significant regression of tumors, which inversely proved that the degree of autophagy increased in pancreatic cancer cells [33]. In contrast to autophagy’s potential role in promoting tumor proliferation, autophagy exerts an inhibitory effect on tumor growth. At present, the mechanism through which autophagy is induced to inhibit tumor development is not fully understood. Thorburn et al. found that autophagy prevented early steps in developing epithelial neoplasms [34]. It was also reported that autophagy could prevent tumorigenesis by limiting necrosis and chronic inflammation [35]. Our results showed the ATG12 was co-expressed with ZEB1-AS1, AC073896.3, AL138756.1, TNFRSF10A-AS1, and AC008760.1. In our risk model, AC073896.3 and TNFRSF10A-AS1 were protective factors, whereas ZEB1-AS1, AL138756.1, and AC008760.1 were risk factors. The results of this study indicated that the identified ATGs might exert different effects on colorectal cancer through different lncRNAs. Furthermore, this also reflected the complex functions of autophagy in tumors.As described in Section3, among the eight autophagy-related lncRNAs, three lncRNAs were protective factors (AC073896.3, AL022323.1, and TNFRSF10A-AS1) and the remaining five lncRNAs were risk factors (LINC01503, ZEB1-AS1, AC087481.3, AC008760.1, and AL138756.1). We also performed Q-PCR assay to explore the expression of these lncRNAs in cell lines. The expression of 7 lncRNAs was consistent with the previous results: 2 protective factors were low expression and five risk factors were high in CRCs. However, the protective factor, TNFRSF10A-AS1, is expressed highly in CRCs. We thought it might be due to the feedback upregulation during the development of the tumor. At present, only LINC01503 and ZEB1-AS1 have been studied extensively. Xie et al. found that the expression of LINC01503 was increased in squamous cell carcinoma, and its overexpression can contribute to the malignant phenotypes of squamous cell carcinoma cells in terms of cell proliferation, migration, and invasion [36]. Moreover, LINC01503 can promote proliferation and metastasis in gastric tumor cells through the WNT signaling pathway [37]. A previous study indicated that LINC01503 promoted the proliferation and metastasis of colorectal tumor cells by regulating miR-4492/FOXK1 signaling [38]. ZEB1-AS1 is one of the most characteristic oncogene regulatory factors in tumor-related lncRNAs and is often overexpressed in a variety of cancers [39]. Li et al. suggested that ZEB1-AS1 was upregulated in hepatocellular carcinoma (HCC) and promoted metastasis; hence, it can be used to predict the prognosis of this disease [40]. In addition to HCC, ZEB1-AS1 can also promote metastasis in numerous cancer types, such as prostate cancer, glioma, bladder cancer, and colorectal cancer [41–44]. In colorectal cancer, the OS and recurrence-free survival rates of patients with increased ZEB1-AS1 expression were lower, indicating that ZEB1-AS1 is a promising biomarker in predicting the clinical outcome [44]. In this study, LINC01503 and ZEB1-AS1 were associated with autophagy and may play a role in colorectal cancer prognosis. Moreover, the other autophagy-related lncRNAs we found may also help predict the prognosis of colorectal cancer, although most of them have been poorly studied.Studies have shown that lncRNA is a novel regulator of autophagy. The lncRNA maternally expressed 3 (MEG3) is significantly downregulated in glioma tissues and cell lines, and its over expression can significantly inhibit cell proliferation and promote the apoptosis and autophagy of glioma cells [45]. The specific mechanism of lncRNA for the regulation of autophagy can be divided into three categories: (1) lncRNA acts as a competitive endogenous RNA-binding miRNA to regulate miRNA expression, thereby affecting the process of autophagy [46]; (2) lncRNA can also affect the expression of ATG genes and (3) lncRNA promotes tumor progression by regulating autophagy-mediated apoptosis through signaling pathways, such as the AKT/mTOR signaling pathway [47]. Recent studies revealed the relationship between signaling pathways and autophagy. The mTOR signaling pathway is crucial in regulating autophagy. When cellular nutrients and growth factors are abundant, mTOR1 is activated and phosphorylates key autophagy-related proteins, leading to the inhibition of autophagy; on the contrary, inhibition of mTOR1 leads to the induction of autophagy [48]. RAS signal transduction can also regulate autophagy through two main cellular pathways. Activation of RAS leads to increased stimulation of the PI3°K pathway and upregulation of mTOR1, resulting in inhibition of the autophagy system [49]. Furthermore, RAS activation can lead to increased autophagy by reducing the signal transduction of the mitogen-activated protein kinase/extracellular signal-regulated kinase (MAPK/ERK) pathway [49]. In colorectal cancer, studies found that autophagy was induced in HT-29 colon cancer cells following the stimulation of the MAPK/ERK signaling pathway [50]. Our results of the GSEA demonstrated that KRAS signaling gene sets were significantly enriched in the high-risk group, indicating that the selected eight autophagy-related lncRNAs were associated with KRAS signaling. This pathway may regulate autophagy, affecting the prognosis of patients with colorectal cancer.


## 5. Conclusion

We established a co-expression network of autophagy and related lncRNAs. Subsequently, the risk signature based on the eight autophagy-related lncRNAs was determined, and the risk model can predict the prognosis of patients with colorectal cancer. Owing to the various functions of lncRNAs, these autophagy-related lncRNAs may be promising targets for treating colorectal cancer. The results of this study provide new insight into exploring the prognosis of colorectal cancer.



---

*Source: 4728947-2020-10-23.xml*