Genome-Wide Identification of Autophagy Prognostic Signature in Pancreatic Cancer

Abstract

Background:

Autophagy plays a vital role in cancer development. However, there is currently no comprehensive study regarding the effects of autophagy-related genes (ARGs) on pancreatic cancer prognosis. Thus, this study aimed to establish an autophagy-related signature for predicting the prognosis of patients with pancreatic cancer.

Methods:

We identified and validated differentially-expressed ARGs using data from The Cancer Genome Atlas (TCGA) database, Genotype-Tissue Expression project (GTEx) and Expression Omnibus (GEO) database. We performed Cox proportional hazards regression analysis on the differentially-expressed ARGs to develop an autophagy-related signature. We tested the expression of these genes through western blotting and verified their prognostic values through gene expression profiling and interactive analyses (GEPIA).

Results:

We identified a total of 21 differentially-expressed ARGs and screened 4 OS-related ARGs (TP63, RAB24, APOL1, and PTK6). Both the training and validation sets showed that the autophagy-related signature was more accurate than the Tumor Node Metastasis (TNM) staging system. Moreover, the western blotting result showed that the expression of TP63, APOL1, and PTK6 was high, whereas that of RAB24 was low in cancer tissues.

Conclusion:

This 4-ARG signature might potentially help in providing personalized therapy to patients with cancer.

Introduction

Pancreatic cancer is an aggressive malignancy with a 5-y survival rate of <8%. ¹ Although pancreatectomy is the preferred treatment approach, this is only possible for patients with small lesions and no metastasis. ² However, even with surgical resection, most patients eventually die due to recurrence. ³ The keys for the treatment of pancreatic cancer are early detection and treatment; however, most patients are treated at an advanced stage. Despite many recent advances, new therapies are still urgently required. In recent years, autophagy has received a lot of attention with regard to tumorigenesis and the development of cancers. ⁴

Autophagy is an essential physiological process that degrades cellular components through lysosomes and is critical for maintaining homeostasis of the internal environment. ⁵ At present, it is believed that autophagy plays a bidirectional role in the process of tumorigenesis. In the early stage of carcinogenesis, cells can maintain normal structural and metabolic stability by eliminating damaged organelles. ⁶ However, during tumor development, autophagy can promote tumor growth through the provision of nutrients. ⁷ The role of specific autophagy-related genes (ARGs) in the development and progression of cancer has been investigated in previous studies. For instance, the decrease in the levels of autophagy related 5 (ATG5) has been shown to be associated with tumor metastasis and shorter survival time. ⁸ Loss of ATG5 or ATG7 in the pancreas with activated mutations in KRAS can prevent the progression of precancerous lesions to invasive cancers. ⁹

However, to date, no large-scale ARG-expression profile screening for the identification of a prognostic signature for pancreatic cancer has been performed. Furthermore, few studies have established a hierarchical prognostic risk assessment system for cancer, based on an ARG signature. The purpose of this study was to establish an autophagy prognosis signature to effectively predict the individualized prognostic information of patients with cancer.

Materials and Methods

Data Collection

The autophagy gene list was downloaded from the Human Autophagy Database (HADB, https://autophagy.lu/clustering/index.html). RNA-seq expression profiles were obtained from The Cancer Genome Atlas (TCGA) database and Genotype-Tissue Expression project (GTEx). GSE78229 was download from the Gene Expression Omnibus (GEO) database and was used as the validation set. All datasets are open to the public.

Screening of Differentially-Expressed ARGs

To further explore the differences in ARGs between cancer and normal tissues, we used 178 cancer tissues from TCGA and 165 normal pancreatic tissues from the GTEx website for differential analysis. The expression profiles were quantified as raw read counts. We used the “sva” package in R (version 3.6) to remove batch effects and DESeq2 for normalization and differential analysis. Identification of differentially-expressed ARGs was based on the cut-off criteria of adjusted P < 0.05 and |log₂FC| >1.5.

Functional Enrichment Analysis

Differentially-expressed ARGs were analyzed by Gene Ontology (GO) and the Kyoto Encyclopaedia of Genes and Genomes (KEGG) using the “clusterProfiler” package in R. Adjusted P < 0.05 was used as the cut-off criterion.

Development and Validation of the Autophagy-Related Signature

To discover ARGs associated with the overall survival (OS) of patients with cancer, univariate and multivariate Cox proportional hazards regression analyses were performed to construct the autophagy-related signature.

The autophagy-related signature formula was as follows:

$\text{Risk Score = }\sum \text{β n ×}Exp\text{n,}$

β n: coefficient of genes by multivariate Cox regression

Expn: expression level of genes by multivariate Cox regression

Patients were divided into 2 groups according to the median risk score. Survival analysis was performed using the Kaplan-Meier test (log-rank test) for the comparison of high-and low-risk groups. The area under the receiver operator characteristics (ROC) curve was used to evaluate the accuracy of the autophagy-related signature.

Western Blotting

The tissues were lysed by chilled radio-immunoprecipitation assay (RIPA) and the bicinchoninic acid (BCA) method was used to determine the protein concentration of each lysate. Proteins (80 µg) were separated using 10% sodium dodecyl sulfate-polyacrylamide gel electrophoresis. The procedures were performed according to the manufacturer’s protocol and references. ¹⁰ The primary antibodies used were as follows: TP63 (rabbit polyclonal antibody, 1:1000, Cat#12143-1-AP, Proteintech), APOL1 (rabbit polyclonal antibody, 1:1000, Cat#11486-2-AP, Proteintech), RAB24 (rabbit polyclonal antibody, 1:500, Cat#11445-1-AP, Proteintech), PTK6 (rabbit polyclonal antibody, 1:1000, Cat#18697-1-AP, Proteintech,), and GAPDH (murine monoclonal antibody, 1:2000, Cat#TA-08, ZSGB-BIO). The following secondary antibodies were used: HRP-labeled goat anti-rabbit antibody or murine IgG (1:2000, ZSGB-BIO, China).

Statistical Analysis

Statistical analyses and simulations were performed using SPSS 23.0 (IBM, Chicago, IL, USA) along with R software 3.6 and GraphPad Prism 8.0. Survival curves were estimated using the Kaplan-Meier test (log-rank test). The Cox proportional hazards model was applied for multivariate survival analysis. All statistical tests were 2-tailed, and P-value < 0.05 was considered statistically significant.

Results

Identification of Differentially-Expressed ARGs

We downloaded the mRNA expression profiles of 178 cancer and 165 normal pancreatic tissues from TGCA and GTEx. We also obtained a total of 231 ARGs from the HADb. Using adjusted P < 0.05 and | log₂FC |> 1.5 as the cut-off, we identified 13 upregulated and 8 downregulated genes. The box-plot shows the expression pattern of 21 differentially-expressed ARGs between non-tumor and cancer tissues (Figure 1). The full names and expression levels of these are shown (Table 1).

Figure 1.

Expression of 21 differentially-expressed autophagy-related genes between normal and cancer tissues.

Table 1.

The Expression Levels of Differentially-Expressed ARGs.

Gene	Full name	Log2 FC	P. adj
NRG3	Neuregulin 3	3.5111	2.17E-36
IFNG	Tnterferon gamma	2.6758	1.68E-26
TP63	Tumor protein p63	2.4991	1.00E-27
IL24	Interleukin 24	2.3703	6.70E-44
BIRC5	Baculoviral IAP repeat containing 5	2.3686	4.11E-50
PTK6	protein tyrosine kinase 6	1.7634	5.52E-44
APOL1	Apolipoprotein L1	1.7152	7.74E-56
CXCR4	C-X-C motif chemokine receptor 4	1.6670	4.46E-54
FAM215A	Family with sequence similarity 215 member A	1.6575	2.67E-15
ATG9B	Autophagy related 9B	1.5434	1.73E-22
CDKN2A	Cyclin dependent kinase inhibitor 2A	1.5373	4.96E-08
NLRC4	NLR family CARD domain containing 4	1.5151	5.53E-32
CCR2	C-C motif chemokine receptor 2	1.5026	6.90E-16
BNIP3	BCL2 interacting protein 3	−1.5364	2.26E-56
DAPK2	Death associated protein kinase 2	−1.5471	5.09E-47
PIK3C3	Phosphatidylinositol 3-kinase catalytic subunit type 3	−1.5668	2.26E-56
TM9SF1	Transmembrane 9 superfamily member 1	−1.5751	2.26E-56
RAB24	Member RAS oncogene family	−1.6029	8.70E-56
TMEM74	Transmembrane protein 74	−1.9033	6.08E-41
NRG2	Neuregulin 2	−2.0083	9.64E-50
SPNS1	Sphingolipid transporter 1 (putative)	−2.1907	2.26E-56

Functional Enrichment Analysis

GO analysis of 21 differentially-expressed ARGs revealed that these genes were involved in several essential biological processes (BPs), cellular components, and molecular functions (Figure 2A). We also found that these differentially-expressed genes were significantly enriched in autophagy and in processes utilizing autophagic mechanisms. While the autophagosome was enriched in cellular components (CC), the receptor-ligand and receptor regulator activities were significantly enriched in the molecular function (MF). Moreover, we found that the KEGG pathway was significantly associated with ARGs, including viral protein interactions with cytokines and cytokine receptors, cytokine-cytokine receptor interactions, platinum drug resistance, and the ErbB signaling pathway (Figure 2B).

Figure 2.

Functional enrichment of 21 differentially-expressed autophagy-related genes. (A) Gene Ontology (GO) and (B) Kyoto Encyclopaedia of Genes and Genomes (KEGG) analyses of differentially-expressed autophagy-related genes.

Development and Validation of the Autophagy-Related Signature

To further screen for ARGs significantly associated with OS in patients with cancer, we performed univariate and multivariate proportional hazards Cox regression analyses on differentially-expressed ARGs (Figure 3). Eventually, we identified a total of 4 independent prognostic ARGs (TP63, RAB24, APOL1, and PTK6) that were used to develop the autophagy-related signature. Applying the coefficient of multivariate Cox regression analysis as a weighting factor, the risk score was given by the following formula: 0.442 × TP63 expression value + 1.3721 × APOL1 expression value + 0.769 × PTK6 expression value − 0.863 × RAB24 expression value. Using a median risk score as the cut-off criterion, we calculated the risk score for each patient and plotted the distribution of the risk score, survival status, heatmaps, and OS time (Figure 4). Kaplan-Meier survival curves revealed that the OS of the high-risk group was much lower than that of the low-risk group (P < 0.001). We also found that the risk score, survival status, heatmap, and OS time in the validation set were consistent with those of the training set.

Figure 3.

Forest map of univariate (A) and multivariate (B) analyses.

Figure 4.

Distribution of autophagy-related signature risk score. From top to bottom are survival status, the heatmap, and Kaplan-Meier survival curve in the training (A) and validation (B) sets.

Furthermore, we used the X-title software (version 3.6.1) to choose the best quartile cut-off values of the signature risk score. We noticed that in the training group, the risk score values for the prediction of the 3-and 5-y OS were 0.749 and 0.807 (Figure 5A), while in the validation group, these values were 0.695 and 0.73, respectively (Figure 5B). The prediction ability of the autophagy-related signature was demonstrated to be more accurate than that of the TNM staging system in both the training and validation groups.

Figure 5.

Comparison of the areas under the curve (AUC) of the autophagy-related signature risk score and TNM stage system in the training (A) and validation (B) sets.

Prognostic Values of the 4 Dysregulated ARGs

When we examined the expression of TP63, RAB24, APOL1, and PTK6 in tissues using western blotting, we found that the expression of TP63, APOL1, and PTK6 was high, whereas that of RAB24 was low in cancer tissues (Figure 6A and B). We further explored the prognostic values of the 4 ARGs in cancer using GEPIA, which is an interactive tool based on a network database. ¹¹ Our results showed that the high expression of TP63, APOL1, and PTK6 was significantly correlated with shortened OS (P = 0.027, P = 0.017, and P = 9.4e-05, respectively). We additionally observed that the low expression of RAB24 in cancer tissues was significantly associated with poor prognosis (P = 0.012) (Figure 6C-F).The above results suggested that these 4 ARGs might play important roles in the occurrence, development, and prognosis of cancer.

Figure 6.

A, Expression of TP63, RAB24, APOL1, and PTK6 in tissues by western blotting (n = 4 for each group). B, Relative values of TP63, RAB24, APOL1, and PTK6 in the tissues mentioned in (A). C-F, Prognostic values of dysregulated TP63, RAB24, APOL1, and PTK6 in cancer tissues, as analyzed by GEPIA. Data are presented as mean ± SD. *, P < 0.05; **, P < 0.01; ***, P < 0.001.

Discussion

Pancreatic cancer is one of the most fatal malignancies. Due to ineffective early detection and treatment, patients are already in the middle and late stages when they are diagnosed, with a meagre 5-y survival rate. The TNM staging system is usually used for the classification and selection of treatment for patients with cancer. However, due to the heterogeneity of tumors, even at the same stage, the therapeutic effect might be different. In recent years, high-throughput sequencing has been widely used in the diagnosis and treatment of cancer. Moreover, there have been many studies on the mechanism of autophagy in cancer; however, these studies have only focused on a single autophagy gene.

In the present study, we conducted a comprehensive analysis of the ARGs in cancer and identified a total of 21 differentially-expressed genes between cancer and normal tissues. Considering that these differentially-expressed genes are involved in the process of cancer, we analyzed them using GO and KEGG pathway functional enrichment analyses, which demonstrated that they were accumulated in tumor-related signaling pathways such as platinum drug resistance, the ErbB signaling pathway, and the chemokine signaling pathway. Hence, these signaling pathways could explain the molecular mechanism of ARGs in cancer.

After univariate and multivariate analysis, we identified 4 ARGs (TP63, RAB24, APOL1, and PTK6) and developed an autophagy-related signature for patients with cancer. Using the median risk score signature as the cut-off value, the patients with cancer were divided into high-and low-risk groups. Survival curves indicated that the autophagy-related signature could adequately distinguish patients with cancer risk stratification. Compared with the TNM stage system, this autophagy-related signature was demonstrated to have more predictive ability in both the training and validation groups.

In the present study, we first established an autophagy-related signature including 4 ARGs (TP63, RAB24, APOL1, and PTK6), which were associated with the OS of patients with cancer. Among these, TP63 belongs to the TP53 family of tumor suppressor genes, ¹² and these genes induce cell cycle arrest and apoptosis. ^{13 -15} In particular, TP63, which encodes 2 subtypes of TAP63 and carcinogenic deltaNp63 subtypes, ¹⁴ has been involved in many cancers such as bladder, uterine, and breast cancer, as well as squamous cell carcinoma of the head and neck. ^{16 -19} Somerville and Xu ²⁰ found that the deltaNp63 subtype-driven enhancer reprogramming enhances tumor growth and the invasive ability of cells, promoting epithelial differentiation and metastasis of cancer cells. Moreover, they showed that the continuous expression of TP63 is important for the growth of cancer epithelial cells. RAB24 is an atypical member of the Rab GTPase family. ²¹ Amaya and Militello ²² found that RAB24 is necessary for normal cell division and might be involved in regulating chromosome segregation and cytokinesis. However, there have been few reports on the roles of RAB24 in cancer. APOL1, which encodes a secreted high-density lipoprotein, has been shown to be involved in the formation of most cholesterol lipids in the plasma, playing a vital role in host defense and the maintenance of intracellular homeostasis. ^{23 -25} The altered function of APOL1 has been linked to chronic kidney disease and cancer. ²⁴ Liu et al ²⁵ found that APOL1 could be used as a diagnostic marker for cancer; however, its mechanism in cancer remains unclear. PTK6 is a non-receptor intracellular tyrosine kinase and has been found to be over-expressed in several cancers such as lung, bladder, ovarian, cervical, gastric, head and neck cancers, as well as B-cell and T-cell lymphomas. ^26,27 More specifically, it has been reported to function in prolonging the S-phase and enhancing the active apoptosis induced by gemcitabine, ²⁸ suggesting that cancer subtypes with high expression of PTK6 are more sensitive to gemcitabine. Moreover, PTK6 has been shown to regulate the migration and invasion of cancer cells through ERK signals, and therefore, it might be a therapeutic target for cancer. ²⁹

However, there were some limitations in this study. First, the sample size included in the study was relatively small. Therefore, future studies will need larger sample sizes to verify the model. Second, we only used western blotting to detect the expression of these genes, but no functional experiments were performed to further reveal the potential mechanism of the autophagy-related signature. Therefore, more functional experimental studies are needed to evaluate the accuracy of the model in the future.

Conclusion

In summary, based on an in-depth analysis of ARGs, we identified 4 OS-related ARGs in patients with cancer, indicating that genes of the autophagy pathway might also serve as potential biomarkers and targets of therapeutic intervention in pancreatic cancer. We developed an autophagy-related signature based on differentially-expressed ARGs to predict the OS of patients with cancer and revealed the association of cancer tissues with specific signaling pathways. Meanwhile, both the training and validation sets revealed that the signature was superior to TNM staging in terms of predicting survival. The novel 4 autophagy-related signature might help the development of personalized therapy for patients with cancer.