Autophagy-related lncRNA and its related mechanism in colon adenocarcinoma

Feifei Tan,Zhongyin Zhou (✉)

Department of Gastroenterology,Renmin Hospital of Wuhan University,Hubei Key Laboratory of Digestive System,Wuhan 430060,China

Abstract Objective Colon cancer is a type of cancer with high morbidity and mortality,of which adenocarcinoma is the most common type.Numerous studies have found that long noncoding RNAs (lncRNAs) are related to the occurrence and development of colon cancer.Autophagy is a key metabolic process in the human body and has a role in affecting cancer growth.In this study,our aim was to explore the correlation between lncRNAs and colon adenocarcinoma (COAD) from the perspective of autophagy.Methods A series of bioinformatics methods were used to explore the correlation between lncRNA and COAD from the perspective of autophagy.Results Four autophagy-related lncRNAs related to the prognosis of COAD were identified:EB1-AS1,LINC02381,AC011462.4,and AC016876.1.These four lncRNAs may act as oncogenes involved in the occurrence and development of COAD.The prognostic model was established,and the accuracy of the model was verified by the receiver operating characteristic curve.The risk score of the model could independently predict the prognosis of patients and was preferable to other clinical indicators,with higher values indicating a worse prognosis of the patients.Gene Set Enrichment Analysis was performed for these four lncRNAs,which showed that the high expression group of these were enriched in the basal cell carcinoma pathway.To make it more convenient for clinicians to use,we constructed a nomogram based on age and risk score,which can be used to evaluate the one-,three-,and five-year survival rates of patients.Conclusion These results can help us understand the mechanism of action of lncRNA on COAD from the perspective of autophagy and may provide new directions for the diagnosis and treatment of COAD.The EB1-AS1 gene in this study is a potential candidate biological target for COAD treatment in the future.

Key words:colon adenocarcinoma (COAD);prognostic model;long noncoding RNA (lncRNA);EB1-AS1

According to Chenet al’s statistical analysis of cancer data from the National Central Cancer Registry of China from 2009 to 2011,colon cancer was ranked in the top five cancers in China regarding new incidence and mortality[1].Surgery is still the main treatment for colon cancer,though combinations of radiotherapy,chemotherapy,and neoadjuvant chemotherapy may be used in cases where surgery alone is unable to treat the cancer[2].With developments in biomedical research in recent years,molecularly targeted drugs have become an option for the non-surgical treatment of patients[3].However,the current median overall survival rate from colon cancer is only approximately two years[4].At present,the incidence and mortality of colon cancer are increasing in China;particularly concerning is the increased incidence among young people[5].Therefore,there is an urgent need to find new treatment targets for colon cancer to improve the prognosis of patients.

Unlike messenger RNA (mRNA),long noncoding RNA (lncRNA) does not participate in gene expression as a template for protein translation but can affect various biological activities in the human body by regulating protein synthesis and was found to be strongly correlated with cancer[6].There are almost 8000 cancerrelated lncRNAs in The Cancer Genome Atlas (TCGA) that can affect different stages of cancer,including cell proliferation,apoptosis,and metastasis[7].Autophagy is an important metabolic pathway for maintaining homeostasis in the human body and is responsible for the degradation of macromolecular substances,such as damaged organelles and long-lived macromolecular proteins[8].Studies have reported that autophagy can promote the occurrence,development,and metastasis of cancer and is related to the invasiveness of cancer cells[9].In addition to being regulated by the expression of autophagy-related genes,lncRNAs are also involved in autophagy regulation[10].The most common type of colon cancer is adenocarcinoma (colon adenocarcinoma,hereafter referred to a COAD) originating from the epithelium of the colonic mucosa.Therefore,understanding autophagy-related lncRNAs and their molecular mechanisms in COAD is beneficial for the treatment of colon cancer and may bring improvements to colon cancer therapy.To achieve this,we analyzed the expression level of autophagy-related genes in COAD tissues that were extracted from the TCGA and used a coexpression (Cox) analysis to obtain the related lncRNAs.The independent prognostic genes that were strongly correlated with COAD were then used to establish a clinically relevant prognosis model.To understand the usefulness of this prognosis model to clinicians,it was compared to other clinically relevant indicators.We also investigated the role of these lncRNAs in COAD using Gene Set Enrichment Analysis (GSEA).Finally,we used the risk score and patient age to prepare a nomogram to predict the survival of patients,whose accuracy we then verified using relevant biological methods.

Materials and methods

Data download and preprocessing

Construction of the lncRNA and mRNA matrix:tissues from COAD patients were obtained from TCGA to construct the lncRNA and mRNA matrix.Since the data used in this study were obtained exclusively from the TCGA database and strictly followed the TCGA publication guidelines (http://cancergenome.nih.gov/abouttcga/policies/publicationguidelines),it was not necessary to obtain the approval of the ethics committee.Autophagy-related genes were obtained from the Human Autophagy Database (HADb),and the expression levels of these genes were extracted from the mRNA matrix.The lncRNA matrix and the expression level of autophagy genes were used for the joint analysis.Autophagy-related lncRNAs were found by setting CorFilter ＞ 0.4 andP＜ 0.001,and the expression level of lncRNAs in the matrix was extracted by Cox analysis.

Building a prognostic signature

The obtained autophagy-related lncRNAs were combined with the clinical data obtained from the TCGA database [including age (divided into ＞ and ＜ 65 years old groups),T (the extent of the primary tumor),M (whether there is distant metastasis),and N (the involvement of local lymph nodes) staging,grading,and other related data].Currently,clinical prognosis is mainly based on the stage and grade of tumor cells and the patient’s age.Univariate cox (unicox) analysis was performed to obtain autophagic lncRNAs related to the prognosis of COAD,and then,multivariate cox (multicox) analysis was performed with the obtained genes to obtain independent prognostic lncRNAs that were strongly related to COAD.These lncRNAs were found to be related to specific clinical indicators by the Clinical correlation analysis,and the specific formula was as follows:

where n represents the number of prognostic lncRNAs,the regression coefficient is β,and EXP is the expression value.

Evaluation of the prognostic signature

By calculating the risk score of all samples in this study,the median expression level was obtained.With the median as the boundary,samples with higher expression levels were defined as the high expression level group,and those with lower expression levels were the low expression level group.The receiver operating characteristic (ROC) curve was drawn to evaluate whether the signature was representative of the groups.

An independent prognostic analysis was used to assess whether the prognostic signature could be used as a prognostic factor independent of the other clinical indicators.In order to better understand how these lncRNAs affect autophagy,Cox analysis of these genes and the mRNA was performed.In addition,the“survival”package in R4.0.3 was used to draw a survival curve to evaluate the impact of a single prognostic lncRNA and risk score on patient survival.

Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis

In order to better understand the ways in which these lncRNAs acted,KEGG enrichment analysis was performed for both groups using GSEA 4.1.0 (set atP＜ 0.05,|ES| ＞ 0.6,gene size ≥ 100,and false discovery rate ＜ 25%) to identify the biological pathways that were enriched in the high and low groups.

Construction and evaluation of the nomogram

In order to better serve clinicians,we included the age and risk score to draw nomograms and evaluated them.First,the ROC curve was used to evaluate their representativeness,and then a c-index was used to evaluate their predictive ability.Finally,a calibration curve for one,three,and five-year survival rates was established.

Results

The entire research process is shown in Fig.1.

Fig.1 The entire research process

Screening of autophagy-related lncRNA

Four-hundred and thirty-seven tissues were obtained from the TCGA,of which 398 were COAD tissues and 39 were adjacent tissues.A total of 257 autophagy-related genes were obtained from the HADb,and their respective expression levels were obtained.A total of 922 autophagyrelated lncRNAs were yielded from the Cox analysis.

Cox regression analysis and clinical correlation analysis

The unicox analysis showed that 18 lncRNAs were correlated with the prognosis of COAD,and anything withP＜ 0.05 was considered meaningful.Finally,four lncRNAs were obtained by the multicox analysis,namely EB1-AS1,LINC02381,AC011462.4,and AC016876.1.The clinical correlation analysis showed that these four lncRNAs were related to the classification and staging of COAD (Fig.2),and the higher their expression,the higher the COAD grade and staging level,meaning that these four lncRNAs may act as oncogenes in the occurrence and development of COAD.

Fig.2 Clinical correlation analysis,integrating prognostic genes with tumor staging and grading (ns,none significance;＊,P ＜ 0.05;＊＊,P ＜ 0.01;＊＊＊,P ＜ 0.001)

Construction and evaluation of risk signatures

The prognostic model of autophagy-related lncRNAs is given in Fig.3 (with the left showing the Independent prognostic analysis and significance denoted byP＜ 0.05).Both the high and low expression groups hadP-values ＜ than 0.05,indicating that the score could be used to independently judge the prognosis of patients.Age was also an independent prognostic risk factor (P＜ 0.05 in both expression groups).The ROC curve on the right showed that the area under the curve (AUC) of the risk score was 0.678,higher than the tumor cell grade and T/N/M stage,which indicated that the accuracy of the risk score in predicting survival was higher than other clinical features,thus verifying the accuracy of the model.

Fig.3 Receiver operating characteristic (ROC) curve evaluation prognostic model.Independent prognostic analysis of high and low expression groups,which includes age,gender,risk score,tumor stage,and tumor grade

Construction of the core lncRNA co-expression network

Fig.4 shows the co-expression network diagram of autophagy-related lncRNAs (left panel),where the diamond block was the prognostic gene and the purple ellipse was the autophagy mRNA.It was evident that the four genes interacted with the autophagy mRNA.The Sankey diagram in the right panel shows the interaction between them more visually.In addition,it showed that these four lncRNAs belonged to the risk group,which was consistent with the data presented in Fig.2.

Fig.5 further validates the results displayed in Fig.4b which includes the survival analysis,risk curve,and heat map,indicating that these genes play a negative role in cancer prognosis.

Fig.4 Co-expression network diagram and Sankey diagram of autophagy-related long noncoding RNAs (lncRNAs).The diamond blocks are prognostic genes,and the purple ellipses are autophagy-related messenger RNAs (mRNAs) showing the autophagy mRNAs with which these four prognostic genes interact.The Sankey diagram (right panel) shows the interaction between these more intuitively.All four genes are in the high-risk group

Fig.5 Correlation between these four prognosis genes and the prognosis of colon adenocarcinoma (COAD) based on survival analysis.The risk curve and heat map further validate the results in Fig.4b,indicating that these genes play a negative role in cancer prognosis

KEGG enrichment analysis

According to the KEGG analysis of COAD patients (Fig.6),the high expression group of these lncRNAs was mainly enriched in the basal cell carcinoma pathway,indicating that the basal cell carcinoma pathway may be the carcinogenic target of these lncRNAs.This provides new insights for molecular research of COAD.The pathways for the enrichment of the low expression group were abundant.During the transcription process,it mainly affected DNA replication,pyrimidine metabolism,DNA mismatch repair,nucleotide excision repair,pentose phosphate pathway,RNA metabolism,purine metabolism,and spliceosomes.In the translation process,it mainly affected the cell cycle,proteasomes,protein export and transportation,ribosomes,and AA-tRNA.

Fig.6 Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis of the four long noncoding RNAs (lncRNAs).High gene expression was related to the basal cell carcinoma pathway,and low expression was enriched in a variety of pathways,including gene transcription and translation

Construction and evaluation of a nomogram

To provide a quantitative method for predicting the probability of survival time,we used information from all samples to construct a nomogram that integrated the risk score and patient age (Fig.7).In the previous steps,we observed that age was related to prognosis,which Aquinaet al[11]also found in COAD patients.The older the age,the worse the prognosis.Therefore,in the nomogram,we integrated the age and risk score to provide patients with a comprehensive score to better predict the one-year,three-year,and five-year survival rates.The nomogram showed that the risk score made the biggest contribution to the nomogram.The accuracy of the model was verified by the ROC curve (AUC=0.708) and C-index=0.691.

Fig.7 A nomogram drawn from the risk score and age factors.The patient’s value was marked on each axis,and the above variables which include risk score and age were added.The results show that the risk score contributed the most to the nomogram.The receiver operating characteristic (ROC) curve (area under the curve;AUC=0.708) and C-index (=0.691) verified the accuracy of the model

Discussion

Globally,COAD has high morbidity and mortality rates[12].Current screening methods for colon cancer have shortcomings[13].The gold standard is colonoscopy,but bowel preparation and contraindications make it unacceptable for many patients.The commonly used intestinal tumor biomarker,carcinoembryonic antigen,lacks sensitivity and specificity.The fecal immunochemical test has a high false positive rate[13].Therefore,there is an urgent need to develop screening tests that are relevant for prognosis and suitable for all patients to improve the detection rate,prognosis,and five-year survival rate of patients with colon cancer.Thus,lncRNAs have attracted attention.Firstly,lncRNA detection is performed through taking plasma samples,which are easy to obtain and inexpensive[14].Secondly,numerous studies have shown that lncRNAs are closely related to the occurrence,development,and metastasis of COAD[15].The lncRNA activated by TGF-β can promote the epithelial-mesenchymal transition process by inhibiting the expression of E-cad,thereby affecting the occurrence of colon tumors.The lncRNA CASC1 can regulate the miR-4310/LGR5/Wnt/β-catenin signal transduction pathway to promote the proliferation and metastasis of COAD[16].HOTAIR may be related to distant metastasis and short survival rates from colon cancer by inhibiting the transcription of the tumor suppressor gene miR-34a[17].Autophagy is a major metabolic pathway in the human body,and many studies have found that it plays an important role in the occurrence and development of colon cancer,although the precise mechanism is not yet clear[18].

Therefore,this study was conducted to explore the relationship between lncRNAs and COAD from the perspective of autophagy.First,relevant genetic information and clinical information in TCGA and HADb were integrated,and four independent prognostic lncRNAs that strongly correlated with COAD were obtained by Cox analysis.These lncRNAs (EB1-AS1,LINC02381,AC011462.4,and AC016876.1) may act like oncogenes.These lncRNAs were used to construct a prognostic signature,and a series of biological processes was used to verify the signatures.Finally,a nomogram was made using risk scores and age,which was converted into specific numbers to predict the one-,three-,and five-year survival of patients.

In this study,EB1-AS1 was the most central among the four lncRNAs,that is,it had a greater prognostic correlation than the other lncRNAs.There have been many studies in recent years on the relationship between EB1-AS1 and COAD,but these have identified differences in the specific ways of action.For example,ZEB1-AS1 can promote cancer by binding to Mir-181A-5p and inhibiting the microRNA (miRNA)-induced β-catenin inhibitory pathway[19];ZEB1-AS1 can also inhibit miR-101 to promote the proliferation and metastasis of cancer cells.The expression level was positively correlated with the histological grade and T stage of the cancer,that is,the higher the expression level,the worse the prognosis of the patients[20].Current studies mainly focus on the interaction between ZEB1-AS1 and miRNA.MicroRNA is a kind of non-coding RNA that generally acts as a tumor suppressor gene,probably mainly by inhibiting transcription or mediating degradation[21].In this study,it was found that the autophagy effect of ZEB1-AS1,namely the interaction between lncRNA and mRNA,may play a very important role in the occurrence and development of COAD,which is a novel result from this study.

Studies have found that LINC02381 may play an inhibitory role in COAD by regulating the PI3K-Akt signaling pathway[22].In this study,the expression level of this lncRNA was inversely related to survival,indicating that it may act as an oncogene.In another study of autophagy-related lncRNAs,this gene was also considered to act as an oncogene[23].This indicates that lncRNA may affect the growth of COAD cells through a variety of ways.Whether it is inhibited or enhanced in cancer cells may be tissue-specific.As far as we are aware,no published studies currently exist for the other two lncRNAs (AC011462.4 and AC016876.1).

Based on the GSEA,we know that high expression of these four genes may play a role in promoting tumor cell recurrence and metastasis through the basal cell carcinoma pathway.A bioinformatics analysis previously showed that the basal cell carcinoma pathway may be the oncogenic target of lncRNAs[24].However,there is a lack of relevant experiments to confirm this,so this is a direction we can consider in the future.The lncRNAs from our study differ from those in previous studies[23,25-26],but we believe our findings may be more credible.The reasons for this are:(1) our study only identified four lncRNAs,which was less than in other studies,indicating that our study found more core lncRNAs;(2) our study identified one core lncRNA,the involvement of which in COAD has been confirmed in many experiments and is likely to be a future biological target of COAD treatment;(3) half of the genes we identified (EB1-AS1 and LINC02381) have been experimentally proven to be related to COAD.

Conclusion

In summary,four independent prognostic lncRNAs related to COAD were found in this study,and among these,EB1-AS1 is very likely to be a new biological target for COAD treatment.Moreover,these four lncRNAs were used to construct a prognostic signature that was superior to the prognostic indicators currently used in clinical practice.Finally,the possible carcinogenic pathways of these four lncRNAs were determined through the enrichment analysis.The EB1-AS1 gene and basal cell carcinoma pathway were specifically daintified and will be the focus of the future research direction of our team,and relevant experiments will be carried out to verify their roles in COAD.

Conflicts of interest

The authors indicated no potential conflicts of interest.