Identification a Novel Autophagy Related Long Non-Coding RNA Prognostic Signature for Pancreatic Adenocarcinoma
- 1. Department of General Surgery, Heping Hospital, Changzhi Medical College, China
Abstract
Pancreatic adenocarcinoma (PAAD) is one of the most common malignant tumors in the digestive system. At present, the prognosis and 5-year survival time are still unsatisfied. Autophagic genes have been demonstrated as a crucial factor in pancreatic cancer progression, we studied the potential prognostic value of autophagy associated long non-coding RNA (lncRNA) in patients with pancreatic cancer. In our research, we summarized five autophagy-related lncRNAs based on The Cancer Genome Atlas (TCGA) pancreatic cancer patient’s data. According to the prognostic lncRNAs, we established a prognostic risk signature and then divided all the patients into low-risk or high-risk groups based on their risk scores. The overall survival (OS) time in the high-risk group is shorter than low risk group (HR=3.75, 95%CI: 2.45-5.73, p<0.001). The autophagy-related lncRNA signature was an independent prognostic predictor with an AUC value of 0.694 (1 year) and 0.703 (5 year). Nomogram was constructed to predict the patients’ survival probabilities based on the risk scores. Gene set enrichment analysis was performed to detect the signaling pathway involved in the different groups, which revealed the related genes were markedly enriched in multiple signaling pathways in high or low- risk group. Moreover, we examined these lncRNAs expression in HPNE cells and three pancreatic cancer cell lines including Mia- PaCa-2, CFPAC-1 and Panc-1. In addition, the biological function between the high and low risk groups was significantly different. We also analyzed the relationship between the autophagy-related lncRNAs signature and pancreatic cancer infiltration lymphocytes via CIBERSORT method in this study. To summarize, the 5-autophagy related lncRNAs we screened in this study has prognostic capability for PAAD and may play a crucial role in pancreatic cancer biology progress.
Keywords
• LncRNA; Autophagy; Prognostic signature; PAAD; TCGA
CITATION
He Y (2023) Identification a Novel Autophagy-Related Long Non-Coding RNA Prognostic Signature for Pancreatic Adenocarcinoma. JSM Biomar 6(1): 1015.
INTRODUCTION
Pancreatic Adenocarcinoma is an extremely aggressive digestive system tumor with insidious clinical manifestations and high degree of malignancy [1]. Due to early metastasis and local progression along with the lack of effective methods for early diagnosis, PAAD has been the 3rd leading cause of cancer related death in the United States [2]. Surgery is currently the only cure for PAAD, and even then, only 37 percent of PAAD patients survive more than five years [3]. Consequently, it is essential to find biomarkers for early diagnosis and accurate prediction of the risk degree for improving PAAD prognosis. Autophagy is an evolutionarily highly conserved intracellular degradation system intended to maintain cell homeostasis in response to different cellular stresses. Autophagy levels are usually at a low level under physiological conditions, whereas could be activated under oxidative stress, nutritional starvation, or multiple disease states [4,5]. Dysregulation of autophagy has been reported in malignant tumor, degenerative diseases of the nervous system, cardiovascular diseases, diabetes, and inflammatory disorders [6]. Autophagy can play a bidirectional regulatory role in tumors, which can either inhibit or promote the tumor progression based on the stage of tumor development. In pancreatic cancer, autophagy is involved in the growth and metabolism of PAAD. High levels of autophagy can both remove damaged cell components and provide metabolites for biosynthesis and energy production for tumor cells [7]. Recent study also showed that autophagy mediated immune escape in pancreatic cancer, could lead to immunotherapy failure [8]. Hence, it is vital to determine autophagy related biomarkers that could serve as effective the early diagnostic and prognostic biomarkers for PAAD patients.
The Human Genome Project revealed that there are 3 billion base pairs in the human genome, of which 1.5 percent encode proteins and 98.5 percent non- protein-coding genes, which were once considered junk genes. However, subsequent ENCODE projects have shown that about 75% of the human genome can be transcribed into RNAs, of which 74% are non-protein coding RNAs (ncRNAs). The long non-coding RNAs are a type of ncRNAs having more than 200 nucleotides with or without protein-coding capacity [9]. LncRNAs regulate important biological functions in cell growth via the form of RNA, including epigenetic regulation, transcriptional regulation, and post- transcriptional regulation [10,11].
Furthermore, some studies have shown that lncRNAs regulate autophagic functions in several cancers. For example, Wang et al., demonstrated that lncRNA-ATB promotes proliferation of hepatocellular carcinoma by activating autophagy [12]. Another study reveals that lncRNA HOTAIR regulates sunitinib-resistance of renal cancer by altering autophagy [13]. With new advances of microarray and gene detection technology, genome sequence had played a progressively important role in the exploration of biomarkers related to tumor diagnosis, treatment and prognosis [14]. Several autophagy related genes have been reported as biomarkers for cancers [15-17].
Therefore, we hypothesized that autophagy related lncRNAs may have the potential to be prognostic biomarkers for PAAD patients. We thoroughly summarized the relationship between autophagy related lncRNAs expression and clinicopathological features in 178 PAAD patients from TCGA database in this research and further constructed a 5-autophagy related lncRNAs prognostic risk model to estimate the PAAD patients’ prognosis individually and accurately.
In the present study, expression profiles and clinical data of 178 PAAD patients from TCGA were involved. The prognostic role of the 5-autophagy related lncRNAs signature was identified by multifaceted analysis. The relationships between the signature and immune cell type fractions, immune checkpoint regulators, mutation profile and functional analyses were further evaluated to explore underlying value of the signature.
MATERIALS AND METHODS
Collection of Data
The PAAD patient’s gene transcriptome raw data together with the matching clinical materials were downloaded from TCGA database (https://portal.gdc.cancer.gov/repository). A total of 182 tissues with gene expression profiles were collected in this study, including 178 PAAD tissues and 4 normal pancreas tissues. The 178 PAAD patients had complete follow-up time and clinical characteristics. In addition, 165 healthy controls from GTEx database were included for this study. According to the Genome Reference Consortium Human Build 38 (GRCh38) information in GENCODE website, we annotated all the lncRNAs and mRNAs that gathered from TCGA dataset (https://www.gencodegenes. org/human/). Finally, a total of 14,142 lncRNAs were recognized for all the patients’ transcriptome data sites.
Identification of the autophagy related lncRNAs
We downloaded the autophagy related genes from the Human Autophagy Database, which offering details about the mechanisms and the regulation of autophagy (HADb: http:// www.autophagy.lu/). Transcriptome data matrix and clinical information of autophagy related genes of 178 patients with PAAD were acquired from TCGA database. Pearson correlation was performed to generate the relationship between all the lncRNAs and autophagy-related genes. LncRNAs with correlation coefficient |R2| > 0.6 and P < 0.001 was identified as autophagy related lncRNAs.
Construction of the Autophagy Related lncRNAs Prognostic Risk Model for PAAD
Initially, the correlation between autophagy related lncRNA and the prognosis of PAAD patients was evaluated by univariate and multivariate cox regression analysis. The prognostic lncRNAs (P-value < 0.001) in univariate analysis were selected for the additional multivariate regression cox analysis to generate the prognostic risk model. The risk score of each patient was calculated depending on the following formula: Risk score = coef (lncRNAgene) × expr (lncRNAgene), where coef (lncRNAgene) and expr (lncRNAgene) separately represented the autophagy related lncRNAs survival correlation coefficient and expression level. Cox analysis was performed to establish a prognostic risk model for predicting PAAD patient’s survival. The median risk scores were used to divide the patients into high-risk and low-risk groups. The Kaplan–Meier survival curve was used to analyze the OS of patients in the two groups. Principal component analysis (PCA) and three-dimensional PCA analysis were conducted to reduce data dimension and separate the patient’s distribution based on the autophagy related genes and autophagy related lncRNAs expression profiles.
Furthermore, univariate and multivariate cox regression analyses were performed to detect whether the risk score was independent prognostic factor for PAAD. The receiver operating characteristic (ROC) curve was performed to compare the prognostic value between the lncRNAs prognostic signature and other clinicopathological values.
Set up a Nomogram
We used nomogram to predict prognosis, including factors such as grade, stage, and risk score. Subsequently, the nomogram was calibrated. These work was done through the “rms”, “foreign”, and “survival” package in R.
Cell culture
HPNE cell lines (Normal human pancreatic duct epithelial cell line) was obtained from Nanjing Medical University. Human pancreatic cancer cell lines (Mia-PaCa-2, Panc-1 and CFPAC-1) were acquired from ATCC: Global Bioresource Center. The pancreatic cancer cell lines and HPNE cells were cultured in the incubator with 37°C and 5% CO2 concentration. Cell culture medium are consisted of RPMI 1640 medium (Gibico, United States) or DMEM (Gibco, United States) with 10% fetal bovine serum (10%, FBS), 100 U/ml penicillin and 100 μg/ml streptomycin antibodies.
Total RNA isolation and quantitative real-time PCR (qRT-PCR)
To detect the expression level of autophagy related lncRNA, we used RNA Trizol reagent (Invitrogen, Carlsbad, CA, United States) to extract cell total RNA. Reverse transcription kit was purchased from TakaRa. We reverse transcribed the total RNA into cDNA after verifying the RNA quality. Primers used for qRT-PCR were synthesized from GenScript (Nanjing, China). Real-time fluorescent quantitative PCR was performed by SYBR Prime Script RT-PCR kit (Roche, Germany). The lncRNA expression level was analyzed using cycle threshold (CT) in the ΔΔCT, and the housekeeping gene GAPDH was selected as the internal parameters to standardize qRT-PCR data. All sequence of the primers used in this research are listed.
Gene Set Enrichment Analysis (GSEA)
We performed Gene set enrichment analysis in the PAAD patients via Hallmarks data sites to enrich the biological signaling pathways in high- or low- risk groups. In our research, we considered the gene sets with a false discovery rate (FDR) < 0.05 and a normalized P value < 0.05 were significantly different in the two subgroups.
Verification of Tumor Infiltration Immune Cells
We used CIBERSORT method to obtain the fraction of immune cell types, and Spearman was used to analyze the correlation between the autophagy related lncRNAs and these immune cells.
Predict chemotherapy responses
To evaluate the response to chemotherapy drugs, we applied public pharmacogenomics database Genomics of Drug Sensitivity in Cancer (GDSC) to predict the chemotherapy response. R package “pRRophetic” was used to calculate the half-maximal inhibitory concentration (IC50).
Statistical analysis
Autophagy related lncRNAs with expression levels P < 0.05 were regarded as statistically significant, which were further used to establish an autophagy related lncRNA-mRNA co-expression relationship via Cytoscape software (version 3.5.1; Cytoscape Consortium, USA). PCA analysis was used to detect whether the autophagy related genes and the 5-autophagy related lncRNAs expression profiles dimensionality have been reduced effectively and visually. GSEA was used to analyze the functional biological states. The Kaplan–Meier method was performed to compare OS time of PAAD patients in each group. Moreover, univariate and multivariate cox regression analyses were applied to identify important prognostic factors. ROC curve was expended to measure the predictive efficiency between the prognostic risk scores and other clinical parameters. The qRT-PCR experiments were analyzed by PRISM 7. Statistical analysis was performed using R software (version 4.0.2).
RESULTS
The expression and CNV status of autophagy related genes in PC
The overall idea of this study was shown in Figure 1. Firstly, we downloaded PAAD patients and healthy controls data from TCGA and GTEx. Figure 2A showed the patients characteristics from TCGA database. We recognized a total of 14142 lncRNAs and 19658 mRNA, which was obtained from the TCGA pancreatic cancer database). A list of 232 autophagy related genes were downloaded from the Human Autophagy Database (HADb: http://www.autophagy.lu/,). Waterfall diagram showed the autophagy related gene mutations in PAAD patients remarkably. TP53 (54%) and CDKN2A (17%) were the two genes with the highest mutation frequency in the 158 PC samples. Subsequently, we analyzed the expression differences of autophagy related genes and lncRNAs in PAAD patients and healthy controls [Figure 2B]. Heatmap showed differentially expressed autophagy related genes. The expression of PTK6, NRG3, TP63, IFNG, IL24, BIRC5, CXCR4, APOL1, CDKN2A, and ATG9B in tumor tissues was significantly higher than that in normal tissues [Figure 2C]. In addition, Figure 2D-E displayed the CNV (Copy Number Variations information) and location of autophagy related genes on the chromosome in PC, which showed the CNV of autophagy related genes were more in loss status.
Identified autophagy related lncRNAs and established an autophagy related lncRNAs signature for PAAD
We further performed a pearson correlation analysis between the lncRNAs and the autophagy related genes using |R| > 0.6 and P < 0.001 as the filter criterion to distinguish autophagy related lncRNAs. Ultimately, 492 lncRNAs were recognized and their expression profiles was listed in supplementary table 6. Depend on the autophagy related lncRNAs data, we used Univariate cox regression analysis and Kaplan-Meier (KM) method to screen prognostic related lncRNAs in 492 autophagy associated lncRNAs. We ranked the prognostic autophagy related lncRNAs in ascending order by their KM and univariate cox regression analysis P values (all less than 0.001). The results showed a total of 20 lncRNAs have prognostic value for PAAD patients [Table 1]. Univariate cox regression analysis also revealed that 19 of the prognostic lncRNAs were belong to protective factors, only one lncRNA (AC245041.2) was risk factor [Table 2]. Additionally, we analyzed co-expression relationship between all the prognostic lncRNAs and autophagy related genes based on the 177 PAAD patients’ data from TCGA. Circos plot displayed a strong positive correlation between the lncRNAs, suggesting a co-activation relationship or a role in similar biological processes [Figure 3A]. Consequently, Multivariate cox regression analysis results indicated five autophagy related lncRNAs were suitable candidates for constructing the prognostic risk model based on the lowest Akaike information criterion (AIC=770.28) [Table 3]. Among the screened autophagy related lncRNAs that were included in the prognostic signature, AC064836.3, AL022328.4, FLVCR1-DT and AC005332.6 were considered as protective factors (HR values 1), whereas AC245041.2 was considered as risk factors (HR values > 1). Then, overall survival analysis was performed depended on the expression of the selected 5 autophagy related lncRNAs. The outcomes showed these autophagy related lncRNAs were meaningfully associated with the overall survival of PAAD patients (P < 0.01; Figure 3 B-F).
Evaluation prognostic signature containing 5 autophagy-related lncRNAs for PAAD
The risk score for each PAAD patient in the TCGA dataset was calculated using the following formula for the autophagy-related lncRNA signature: risk score = 0.137527218542142 × expression level of AC005332.6)
(-0.248704379564644 × expression level of AC064836.3) + (-
0.624755831968904 × expression level of AL022328.4) + (-
0.741903940236085 × expression level of FLVCR1-DT) + (-
+ (0.217095145019759 × expression level of AC245041.2).
Then, PAAD patients were apportioned into high-risk (n = 88) and low-risk (n = 89) groups using the median risk score as the cut-off point. Moreover, Kaplan-Meier survival curve analysis showed that the overall survival of PAAD patients with high-risk scores was significantly shorter than those with low risk scores (HR: 3.75, 95% CI: 2.45-5.73, p< 0.001, Figure 4A). A principal components analysis (PCA) and three-dimensional PCA analysis based on the five-autophagy related lncRNAs showed two significantly different distribution patterns between high risk and low-risk groups [Figure 4 B-C]. Additionally, the 3-year survival rates were approximately 11.9% (95% CI: 0.0528-0.268) and 56.8% (95% CI: 0.446-0.725) for the high-risk and low-risk patients, respectively.Then rank all patients according to the risk scores calculated by autophagy- related lncRNAs prognostic risk model. The scatter dot plot demonstrated that the overall survival of the PAAD patients correlated with the risk scores, and patients in higher risk scores field revealed lower survival time [Figure 4 D-E].
Furthermore, heatmap exposed distinct differences in the levels of the 5 prognostic related lncRNAs in the high- and low risk PAAD patients. High-risk patients expressed higher levels of risk factors (AC245041.2), while higher levels of protective factors (AC064836.3, AL022328.4, FLVCR1-DT and AC005332.6) were found in low-risk patients [Figure 4 F].
The autophagy-related lncRNAs signature is an independent prognostic factor
Univariate and multivariate cox regression analyses were performed to verify whether autophagy-associated lncRNA prognostic risk score was an independent prognostic factor for PAAD patients. Univariate analysis results demonstrated that autophagy related lncRNAs prognostic risk score (HR: 1.510, 95%CI: 1.158-1.969, P: 0.002) were significantly related with OS. However, age, gender, AJCC stage and TNM stage have no obvious association with OS in these TCGA data [Figure 5A]. Multivariate analyses also indicated autophagy related lncRNAs prognostic risk score (HR: 1.499, 95%CI: 1.122-2.003, P: 0.006) were significantly associated with OS and could be an independent prognostic factor [Figure 5B]. All these data demonstrated that the autophagy- related lncRNA prognostic signature is an independent prognostic factor for PAAD patients. Additionally, the one-year ROC curve analysis demonstrated that the AUC value for the autophagy related lncRNAs prognostic signature was 0.694, which was higher than the AUC values for age (AUC=0.534), gender (AUC=0.597), grade (AUC=0.607), AJCC stage (AUC=0.450), T stage (AUC=0.504), N stage (AUC= 0.518) and M stage (AUC=0.467) [Figure 5C]. Furthermore, the five year ROC curve analysis showed the same results, risk score AUC value was also the highest during the other factors (AUC=0.703) [Figure 5D].
Establish the prediction nomogram and verify the prognostic lncRNAs expression in vitro
Nomograms has been reported as an effective clinical tools to accurately predict survival time for a patient by calculating the nomogram score based on the points assigned for each prognostic factor included in the nomogram [18].
Thus, nomogram was formed to precisely estimate the 1-, 3-, and 5-year survival probabilities by using risk score calculated from the autophagy-related lncRNA prognostic model and other clinicopathological factors, including age, gender, grade, T stage, M stage and N stage [Figure 6A]. The calibration curve analysis showed that the actual and the predicted 1- and 5-year survival times were in accordance with the reference line [Figure 6B]. These results demonstrated that the nomogram using the autophagy-related lncRNAs prognostic signature risk scores was reliable. Due to the lack of T1 and T4 stage patients in the TCGA database for pancreatic cancer, more data still be needed in the future to generate more accurate model. To verify the expression level of autophagy- related lncRNAs in PAAD cells, we used qRT-PCR analysis to detect normal pancreatic duct cells (HPNE) and pancreatic cancer cells (Mia- PaCa-2, Panc-1 and CFPAC-1). As showed in [Figure 6 C-G], the results revealed that AC245041.2 was obviously overexpressed in PAAD cells, the high expression was associated with poor survival, HR > 1, which indicated AC245041.2 may play a role as an oncogene in PAAD. Additionally, AC064836.3 was downregulated in PAAD cell lines compared with HPNE cells, which is consistent with the better survival based on TCGA data base. AC005332.6, FLVCR1-DT and AL022328.4 were overexpressed in Mia-PaCa-2, Panc-1 and CFPAC-1 cell lines, however, their high expression was related with better survival, HR < 1, so the internal mechanism still needed to be further studied.
Establish coexpression network and functional enrichment analysis
We further investigated the underlying roles of the 5 lncRNAs in PAAD via generating an autophagy related lncRNA and mRNA co-expression network by Cytoscape software. Moreover, the Graphical method and Sankey diagram unveiled the relationship between the mRNAs and the 5 screened risk or protective lncRNAs [Figure 7 A-B]. Kyoto Encyclopedia of Genes and Genomes (KEGG), pathway analysis confirmed that neuroactive ligand-receptor interaction, cAMP signaling pathway and insulin secretion were the top three most possible pathways for survival differences between high and low risk groups [Figure 7C]. Additionally, the Gene Ontology (GO) results showed that signal release was the top enrichment biological process, while transport vesicle and voltage- gated ion channel activity were the top enrichment cellular component and molecular function [Figure 7D]. Circos plot showed the association between different biological process [Figure 7 E].
Gene set enrichment analysis
GSEA was further performed to detect the signaling pathways enriched between the high- and low-risk groups. The results showed the altered genes in the high- risk PAAD patients belong to pathways related to cytokines and cancer. For regulating cytokines, IFN- α response, TGF-β signaling pathways, androgen response, cholesterol homeostasis and protein secretion were enriched. Additionally, MYC and Notch signaling pathways have been reported in various malignant tumors [Figure 8A]. This suggested that activation of pathways regulation cytokines and tumor growth function in the high-risk group may contribute to negative prognosis or worse survival outcomes. Pancreas beta cells and spermatogenesis were involved in the low-risk group [Figure 8 B]. Several spermatogenesis related genes were proved to play an essential role in cancers including breast cancer and testicular cancer, however, have never been reported in pancreatic cancer [19,20]. The results proposed that high prognostic signature risk score was associated with cytokines and cancer related signaling pathways, while low prognostic signature risk score was correlated with pancreas beta cells and spermatogenesis function [Figure 8C]. All these outcomes provided important clues for us to future examine the potential personalized therapies for PAAD patients with different risk scores.
Association of clinicopathological variables and comparison of the immune status
When we detected the relationships between multiple clinical factors and the autophagy related lncRNAs risk scores, the resulted showed the pancreatic cancer immune scores had significant difference with the risk scores [Figure 9A]. We further analyzed the differences of tumor microenvironment in the high- or low-risk groups. A total of 22 tumor infiltrating immune cells were screened via CIBERSORT method. The results demonstrated that naïve B cells, M0 macrophages, resting dendritic cells and plasma cells exhibited a higher expression in low risk group (P<0.05).Whereas the M2 phenotype macrophages had a higher expression in high risk group [Figure 9B]. Moreover, we also research the association between the autophagy related lncRNAs risk score and tumor-infiltrating immune cells. Spearman’s correlation analysis revealed that the risk score was positively related with 2 tumor infiltrating immune cells (M0 and M2 phenotype macrophages), however, was negatively correlated with naïve B cells, regulatory T cells, CD8+T cells, and plasma cells [Figure 9C]. M2 phenotype macrophages have been proved as a tumor-promoting factor in pancreatic cancer [21]. Hence, higher autophagy related lncRNAs risk scores may promote M2 macrophages infiltration in PAAD. In addition, we also analyzed the correlation between each lncRNA in the signature and tumor infiltrating immune cells [Figure 9D]. In brief, our findings indicated that the autophagy related risk scores could discriminate different characteristics of tumor immune cells in PAAD.
Immune checkpoint modulators and tumor mutation burden
It is well known that immune checkpoint modulators and tumor mutation burden play an important role in tumor progression. Therefore, we also evaluated the relationship between the signature and these two factors. The results showed that there were significant differences of the expression of immune checkpoint modulators in the high and low risk groups. CD44, CD276 and TNFSF9 were higher in the high-risk group than in the low-risk group, while others were higher in the low risk group [Figure 10A]. In addition, we found that the tumor mutation burden in the high-risk group was significantly higher than that in the low-risk group [Figure 10B]. Figure 10 C-D show the difference mutation frequency in the two groups. The order of mutation frequency in high-risk group was as follows: KRAS > TP53 > CDKN2A > SMAD4 > TTN > MUC16 > DAMTS12. And the order in low-risk group was as follows: KRAS > TP53 > SMAD4 > TTN > RNF43 > MUC16 > RYR1. Moreover, cancer stem cell like properties analysis showed that risk score was positively correlated with RNAss [Figure 10E].
Evaluation of immunotherapy response and chemo therapy response
We further assessed potential response to immunotherapy of each patient [Figure 11 A-C]. And the high-risk group had a lower potential for immune dysfunction and immune escape, which may indicate that the high-risk group responded better to immunotherapy. Our prediction of chemotherapy response found that seventeen drugs were more sensitive in the low-risk group, meaning that these drugs might be more suitable for patients in the low-risk group [Figure 11 D].
DISCUSSION
Pancreatic cancer is a malignant digestive system tumor with extremely high mortality rate [22]. Aggressive multimodal therapy was fully utilized in many malignancies, nevertheless, overall survival has not improved for pancreatic cancer, and treatment outcomes remain unsatisfactory [23, 24]. With the further study on tumor clinical management, the prognostic factors such as tumor size, tumor grade and CA19-9 level have been gradually clarified. Next generation sequencing biotechnology is being widely used to predict cancer recurrence and metastasis by detecting transcriptome expression levels. The role of autophagy in tumorigenesis has been reported for several cancers, including PAAD [25-29].
Autophagy related genes are involved in the regulation of autophagy level in vivo, and the abnormal expression of autophagy related genes can lead to a variety of diseases, including cancer. Previous studies have focused on the role of specific autophagy related genes in PAAD progression. In pancreatic cancer, autophagy process is a metabolic requirement, which could be used to inhibit the expression of MHC-I on the PAAD cells surface, thus blocking antigen presentation and achieving immune escape [8,30].
Recently, several lncRNAs have been proved as regulators in multiple cancers via directly or indirectly targeting autophagy related genes. For example, Fan et al. reported that silencing of lncRNA PRRT3-AS1 could suppress prostate cancer proliferation via promoting autophagy progress [31]. Hence, lncRNAs with co expression relationship to autophagy related genes may become potential diagnostic biomarkers and therapeutic targets for PAAD patients. Whereas there has been no autophagy related molecular analysis to recognize lncRNAs prognostic risk model for PAAD patients. Consequently, it is essential to set up a lncRNA signature for predicting the prognosis of PAAD patients. In our research, we used bioinformatics and statistical tools to systematically analyze the prognostic accuracy of autophagy related lncRNAs in PAAD. Firstly, we analyzed the differential expression of autophagy related genes and lncRNAs between the PAAD patients and healthy controls. Here, compared to previous similar studies, a special feature of our study is the inclusion of 165 additional healthy controls from the GTEx database, which makes our data more convincing.
According to the expression level of PAAD patient’s lncRNAs in TCGA database, 20 autophagy-associated lncRNAs that noticeably related with OS were screened by univariate cox regression analysis. Further, 5 autophagy-related lncRNAs, AC064836.3, AL022328.4, FLVCR1-DT, AC005332.6 and AC245041.2 were chosen to structure a prognostic risk model based on their =effects in the multivariate cox regression analysis. All the PAAD patients was divided into high- or low-risk group based on the expression level of the 5- autophgy related lncRNAs. In our research, PAAD patients in low-risk group had longer survival time than those in high-risk group. Furthermore, PCA and three-dimensional PCA analysis based on the five-autophagy related lncRNAs, clearly showed two distribution patterns between the high- and low-risk groups. ROC curve analysis validated the prognostic accuracy of the autophagy-related lncRNA prognostic signature in the PAAD patients. The risk score based on the autophagy related lncRNA prognostic risk model was an independent prognostic factor based on multi-variate cox regression analysis. In our study, the prognostic risk model of autophagy associated lncRNA was superior to other conventional clinical parameters in predicting prognosis. Nomogram has proven to be a resultful and dependable clinical tool for predicting survival in cancer patients. Thus, we generated a sturdy nomogram including the prognostic risk scores determined by the autophagy-related lncRNAs risk model to improve prognostic prediction of PAAD patients.
Moreover, calibration plots proved that the actual and predicted 1- and 5-year survival rates based on the nomogram were consistent. In general, the autophagy related lncRNA prognostic risk model precisely predicts survival outcomes of PAAD patients in our study.
For verifying the expression level of the selected autophagy related lncRNAs in PAAD cells, qRT-PCR analysis was performed in normal pancreatic duct cells (HPNE) and pancreatic cancer cells lines (Mia-PaCa-2, Panc-1 and CFPAC-1). The outcomes indicated that AC245041.2 and AC064836.3 were consistent with the results based on TCGA data base. We also evaluated the 5 lncRNAs, which were associated with autophagy related genes expression in PAAD patients and constructed the lncRNA-mRNA co-expression network. Additionally, GO and KEGG functional enrichment analyses showed that multiple biological processes were enriched. GSEA results also revealed distinct differences in the several signaling pathways between the high- and low risk groups. Cytokines- and cancer-related pathways were enriched in the high-risk group, whereas pancreas beta cells and spermatogenesis pathways were involved in the low-risk group. These data give us some clues for further research in PAAD treatment.
Moreover, we found differences in immune scores between the high and low risk groups. We then further evaluated the tumor-infiltrating immune cells and found our autophagy related risk scores could distinguish different characteristics of tumor infiltrating immune cells in PAAD. We found that M2 phenotype macrophages were more expressed in the high-risk group and were proportional to the risk score. As is known to all, tumor associated macrophages (TAMs) are important innate immune cells in the tumor microenvironment and play an important role in the occurrence and development of tumors [32,33], including M1 macrophages, which can cause inflammatory responses and improve anti-tumor immunity, and M2 macrophages, which have anti-inflammatory effects and play a role in tumor development and metastasis [34,35]. The M2 phenotype of macrophages has been shown to be a tumor-promoting factor in pancreatic cancer [21]. This result intrigued us that a higher risk score might promote M2 phenotype macrophage infiltration in PAAD.
In addition, analyses of immune checkpoint modulators and TMB found obvious differences between high and low risk groups. There were also significant differences between the two groups in response to immunotherapy and chemotherapy. This may be beneficial to the clinical selection of appropriate treatment for different patients. There are still some limits in our research. Initially, our findings need further validation in other independent research cohorts to verify the strength of the prognostic characteristics of autophagy associated lncRNAs. Next, this research was based on a single cohort of 178 patients from the public TCGA database and the data used for analysis was relatively insufficient. Moreover, we only detected the selected autophagy related lncRNAs expression in HPNE and PAAD cell lines. For getting more reliable results, the autophagy related lncRNA risk model should be validated in tumor tissue and in vivo experiments future. Finally, we found these five lncRNAs with prognostic value were only the first step. More importantly, we need to further study the internal mechanism of these lncRNAs affecting the PAAD patients’ prognosis.
In conclusion, our results revealed an autophagy related lncRNA prognostic risk model, which could accurately predict the survival outcomes of PAAD patients. Depend on this risk model we could recognized the PAAD patients into high or low-risk groups. Moreover, we also established and validated a prognostic nomogram via combining the autophagy related lncRNA prognostic risk model scores and other clinical characteristics. Our study also demonstrated that autophagy related risk scores could regulate the distribution of tumor immune cells in PAAD. The autophagy related lncRNAs risk model may provide us the clue for further research the mechanism of the autophagy related genes in regulating tumor growth and have potency to be prognostic and diagnostic biomarkers for PAAD therapy.
REFERENCES
1. SIEGEL RL, MILLER KD, FUCHS HE, Ahmedin Jemal. Cancer Statistics, CA: A Cancer J Clinicians. 2021; 71: 7-33.
2. CIVES M, STROSBERG JR. Gastroenteropancreatic Neuroendocrine Tumors. CA Cancer J Clin, 2018; 68: 471-487.
3. MIZRAHI JD, SURANA R, VALLE JW, Rachna T Shroff, Pancreatic cancer. Lancet. 2020; 395: 2008-2020.
4. LEVINE B, KROEMER G, Autophagy in the pathogenesis of disease. Cel. 2008; 132: 27-42.
5. KOCAK M, EZAZI ERDI S, JORBA G, I Maestro, J Farrés, V Kirkin, A Martinez, O Pless et al, Targeting autophagy in disease: established and new strategies. Autophagy. 2022; 18: 473-495.
6. LEVINE B, KROEMER G. Biological Functions of Autophagy Genes: A Disease Perspective. Cell. 2019; 176: 11-42.
7. YANG S, WANG X, CONTINO G, Marc Liesa, Ergun Sahin, at Haoqiang Ying, et al. Pancreatic cancers require autophagy for tumor growth. Genes Develop. 2011; 25: 717-729.
8. YAMAMOTO K, VENIDA A, YANO J, Douglas E Biancur , Miwako Kakiuchi , Suprit Gupta , et al. Autophagy promotes immune evasion of pancreatic cancer by degrading MHC-I. Nature. 2020; 581: 100 105.
9. SUMAN P, CHHICHHOLIYA Y, KAUR P, S Ghosh , A Munshi. Long non coding RNAs involved in different steps of cancer metastasis. Clin Transl Oncol. 2022; 24: 997-1013.
10. MISHRA S, VERMA SS, RAI V, Nikee Awasthee, Srinivas Chava, at Kam Man Hui, et al. Long non-coding RNAs are emerging targets of phytochemicals for cancer and other chronic diseases. Cell Mol Life Sci. 2019; 76: 1947-1966.
11. DAHARIYA S, PADDIBHATLA I, KUMAR S, Sanjeev Raghuwanshi, Adithya Pallepati, Ravi Kumar Gutti, et al. Long non-coding RNA: Classification, biogenesis and functions in blood cells. Mol Immunol. 2019; 112: 82-92.
12. WANG CZ, YAN GX, DONG DS, He Xin , Zhao-Yu Liu . LncRNA-ATB promotes autophagy by activating Yes-associated protein and inducing autophagy-related protein 5 expression in hepatocellular carcinoma. World J Gastroenterol. 2019; 25: 5310-5322.
13. LI D, LI C, CHEN Y, Lichen Teng , Yan Cao ,at Wentao Wang et al. LncRNA HOTAIR induces sunitinib resistance in renal cancer by acting as a competing endogenous RNA to regulate autophagy of renal cells. Cancer cell Int. 2020; 20: 338.
14. BOKERS S, URBAT A, DANIEL C, Kerstin Amann, Kenneth GC Smith, Marion Espéli, Lars Nitschke, et al. Siglec-G deficiency leads to more severe collagen-induced arthritis and earlier onset of lupus-like symptoms in MRL/lpr mice. J Immunol, 2014; 192: 2994-3002.
15. AN CH, KIM MS, YOO NJ, Sang Wook Park, Sug Hyung Lee . Mutational and expressional analyses of ATG5, an autophagy-related gene, in gastrointestinal cancers. Pathology, research and practice. 2011; 207: 433-7.
16. LEE JW, JEONG EG, LEE S H, Nam Jin Yoo, Sug Hyung Lee. Somatic mutations of BECN1, an autophagy-related gene, in human cancers. APMIS. 2007; 115: 750-756.
17. FU LL, YANG Y, XU HL, Cheng, X. Wen L. Ouyang, et al. Identification of novel caspase/autophagy-related gene switch to cell fate decisions in breast cancers. Cell proliferation. 2013; 46: 67-75.
18. BALACHANDRAN V P, GONEN M, SMITH J J, Ronald P DeMatteo. Nomograms in oncology: more than meets the eye. The Lancet Oncol. 2015; 16: e173-180.
19. OMISANJO OA, BIERMANN K, HARTMANN S, Lukas C Heukamp, Violetta Sonnack, Anne Hild, et al. DNMT1 and HDAC1 gene expression in impaired spermatogenesis and testicular cancer. Histochemistry Cell Biol. 2007; 127: 175-181.
20. SARVESTANI FM, SAFAEI A, TALEI A, Mehdi Dianatpour. The Evaluation of Cancer Testis Gene PIWIL2 Expression Levels as a New Prognostic Biomarker for Breast Cancer. Clinical laboratory. 2016; 62: 1469-1475.
21. YANG Y, GUO Z, CHEN W, Xiaofeng Wang , Meng Cao , Xuan Han , et al. M2 Macrophage-Derived Exosomes Promote Angiogenesis and Growth of Pancreatic Ductal Adenocarcinoma by Targeting E2F2. Mol Ther. 2021; 29: 1226-1238.
22. GUGENHEIM J, CROVETTO A, PETRUCCIANI N. Neoadjuvant therapy for pancreatic cancer. Updates Surg. 2022; 74: 35-42.
23. ANSARI D, DEL PINO BELLIDO C, BAUDEN M, Roland Andersson. Centrosomal Abnormalities in Pancreatic Cancer: Molecular Mechanisms and Clinical Implications. Anticancer Res. 2018; 38: 1241-1245.
24. KRATSCHMER C, LEVY M. Targeted Delivery of Auristatin-Modified Toxins to Pancreatic Cancer Using Aptamers. Mol Ther Nucleic acids. 2018; 10: 227-236.
25. MUKUBOU H, TSUJIMURA T, SASAKI R, Yonson Ku. The role of autophagy in the treatment of pancreatic cancer with gemcitabine and ionizing radiation. Int J Oncol. 2010; 37: 821-828.
26. DONADELLI M, DANDO I, ZANIBONI T . Gemcitabine/cannabinoid combination triggers autophagy in pancreatic cancer cells through a ROS-mediated mechanism. Cell Death Dis. 2011; 2: e152.
27. WANG P, ZHANG J, ZHANG L, Zhengfei Zhu, Jie Fan, Lianyu Chen, et al. MicroRNA 23b regulates autophagy associated with radioresistance of pancreatic cancer cells. Gastroenterol. 2013; 145: 1133-1143.
28. KINSEY C, BALAKRISHNAN V, O’DELL MR, Jing Li Huang , Laurel Newman , Christa L Whitney-Miller, et al. Plac8 links oncogenic mutations to regulation of autophagy and is critical to pancreatic cancer progression. Cell Rep. 2014; 7: 1143-1155.
29. PERERA RM, STOYKOVA S, NICOLAY BN, Kenneth N Ross, Julien Fitamant, Myriam Boukhali, et al. Transcriptional control of autophagy- lysosome function drives pancreatic cancer metabolism. Nature. 2015; 524: 361-365.
30. REYES-CASTELLANOS G, ABDEL HADI N, CARRIER A. Autophagy Contributes to Metabolic Reprogramming and Therapeutic Resistance in Pancreatic Tumors. Cells. 2022; 11: 426.
31. FAN L, LI H, WANG W. Long non-coding RNA PRRT3-AS1 silencing inhibits prostate cancer cell proliferation and promotes apoptosis and autophagy. Experimental Physiol. 2020; 105: 793-808.
32. BINNEMARS-POSTMA K, STORM G, PRAKASH J. Nanomedicine Strategies to Target Tumor-Associated Macrophages. Int J Mol Sci. 2017; 18: 979.
33. BRASSART-PASCO S, BREZILLON S, BRASSART B, Laurent Ramont, Jean-Baptiste Oudart, Jean Claude Monboisse, et al. Tumor Microenvironment: Extracellular Matrix Alterations Influence Tumor Progression. Front Oncol. 2020; 10: 397.
34. JEANNIN P, PAOLINI L, ADAM C, Yves Delneste. The roles of CSFs on the functional polarization of tumor-associated macrophages. FEBS J. 2018; 285: 680-699.
35. KAWAMURA K, KOMOHARA Y, TAKAISHI K, Hidetaka Katabuchi, Motohiro Takeya. Detection of M2 macrophages and colony stimulating factor 1 expression in serous and mucinous ovarian epithelial tumors. Pathol Int. 2009; 59: 300-305.