Discovering a specialized programmed-cell death patterns for prognostic model of pancreatic ductal carcinoma via machine learning

doi:10.21203/rs.3.rs-4670808/v1

Download PDF

Research Article

Discovering a specialized programmed-cell death patterns for prognostic model of pancreatic ductal carcinoma via machine learning

https://doi.org/10.21203/rs.3.rs-4670808/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

Substantial evidence implicates programmed cell death (PCD) in pancreatic ductal adenocarcinoma (PDAC) pathophysiology. Through advanced machine learning paradigms, our study identified 103 PCD-relevant hub genes. Employing a comprehensive panel of 167 algorithmic configurations, spanning 15 unique machine learning approaches, we analyzed the prognostic relevance of these PCD-linked features across diverse cohorts. Our systematic analysis yielded a groundbreaking prognostic indicator, the Cell Death Index (CDI), poised to markedly improve PDAC outcome predictions. Demonstrating notable accuracy in both prognosis and immunotherapy response forecasting, the CDI facilitated the development of an enhanced nomogram. Additionally, we pinpointed targeted therapeutic agents for PDAC patients classified according to specific CDI profiles, advancing personalized medicine strategies. MYOF, identified as a central hub gene, exhibited markedly heightened expression in PDAC tissues versus adjacent non-malignant tissues, as evidenced by quantitative PCR. Further probing revealed MYOF's critical role in mediating proliferation, viability, invasion, and migration in PDAC cells, underscoring its potential significance as a therapeutic target warranting further investigation.

programmed cell-death

machine learning

multi-omics analysis

deep learning

immunotherapy

Pancreatic ductal adenocarcinoma (PDAC) is characterized by its aggressive nature and poor prognosis, presenting substantial challenges in diagnosis, treatment, and overall prognosis. Over the past two decades, the global prevalence of this disease has notably risen (1). The lack of early diagnostic methods and the occurrence of atypical symptoms contribute to the rapid progression and inoperability of the majority of PDAC in clinical settings (2). Despite advancements in research and treatment, PDAC exhibits significantly lower survival rates compared to other malignancies (3, 4).

Programmed cell death (PCD), also known as regulated cell death, plays a crucial role in both pathological and physiological processes, including maintaining cell homeostasis, eliminating damaged or senescent cells, and tumorigenesis. These processes are biologically controlled by specific signaling pathways and molecular effectors (5–7). PCD is divided into two subgroups: apoptotic and non-apoptotic. Non-apoptotic PCD is further categorized into necroptosis, ferroptosis, pyroptosis, and entosis (8). Dysregulated PCD contributes to tumorigenesis, and uncontrolled PCD is a key feature of malignancy (9). The overexpression of the EGFR family in PDAC often leads to resistance to apoptosis (10). Additionally, Ras, a downstream target of EGFR, is identified as the most commonly mutated gene in PDAC. Mutated Ras also decreases apoptosis (11). In addition, PDAC develops uncontrolled resistance to immunotherapy and chemotherapy due to mutations in the PCD pathway (12, 13). Thus, it may be possible to override drug-resistant by controlling PCD.

In this study, we discovered PCD-related hub genes in the context of PDAC to understand its contribution to the disease. Then, we developed a predictive PCD-related signature using machine learning algorithms to facilitate personalized treatment options.

2.1 Datasets

A total of 178 pancreatic ductal adenocarcinoma (PDAC) samples were obtained from the Cancer Genome Atlas (TCGA) database (14), 455 PDAC samples were sourced from the International Cancer Genome Consortium (ICGC) database (15), and 91 PDAC samples were downloaded from the Gene Expression Omnibus (GEO) database (16) with the FPKM format in the bulk transcriptome level. Differentially expressed genes (DEGs) in PDAC were extracted from the Gene Expression Profiling Interactive Analysis (GEPIA) database (17), with the significance threshold at p-value < 0.05 and fold change (log2FC) > 1. A total of 788 genes that regulate the 5 kinds of PCD mutually were extracted, including 200 genes of apoptosis, 159 of pyroptosis, 110 of ferroptosis, 184 of entotic cell death, and 135 of necroptosis from the MSigDB database.

2.2 Identification of PCD-related hub genes in PDAC

The weighted gene co-expression network analysis (WGCNA) was utilized to investigate the hub genes in several PCD pathways in PDAC, which is a machine learning algorithm to construct gene modules based on similarly expressed genes (18). Initially, Pearson correlation among DEGs was calculated and then transformed into an adjacency matrix by indexing to generate a scale-free network. In such a network, a subset of hub genes exhibits high connectivity with other genes, whereas the majority of genes show limited connectivity, reflecting more closely the biological conditions. We then aimed to identify the appropriate exponent (soft-threshold power) by performing a linear regression analysis between the frequency of adjacency matrix values and the corresponding adjacency matrix values. This approach is believed to maximize overall connectivity and enhance the strength of the regression analysis. Additionally, the adjacency matrix was then transformed into topological overlap matrix (TOM) to minimize bias from other genes (19). Eventually, gene modules were constructed based on topological overlap matrix.

Enrichment scores of PCD pathways including apoptosis, ferroptosis, necroptosis, entotic cell death, and pyroptosis were calculated by single-sample gene set enrichment analysis (ssGSEA) (20). Pearson correlation analyses were then performed between gene modules and these enrichment scores. The PCD-related gene module was selected for further investigation if it had a p-value less than 0.05 and a relevance score greater than 0.5. The genes with geneTraitSignificance greater than 0.4 and geneModuleMembership greater than 0.7 were identified. Ultimately, PCD-related hub genes were selected based on their interactions across apoptosis, ferroptosis, necroptosis, entotic cell death, and pyroptosis pathways.

2.3 Functional and pathway enrichment analysis

The biological functions and signaling pathways of the PCD-related hub genes were examined using Gene Ontology (GO), the Kyoto Encyclopedia of Genes and Genomes (KEGG), and Gene Set Enrichment Analysis (GSEA). These analyses were conducted with the R packages "clusterProfiler" and "org.Hs.eg.db".

2.4 Establishment of the PCD-related prognostic model in PDAC

The study identified PCD-related hub genes with prognostic relevance using univariate-cox regression analysis via the “survival” R package with the criteria of p < 0.05. Based on prognostic PCD-related hub genes, a total of 167 algorithmic combinations consisted of 15 machine learning methods were utilized within the cohorts in our research, including “Gradient Boosting Machine (GBM)”, “Random Survival Forest (RSF)”, “Least Absolute Shrinkage and Selection Operator (LASSO)”, “Ridge”, “Elastic network”, “Step-Cox”, “Support Vector Machine (SVM)”, “Support Vector Machine Recursive Feature Elimination (SVM-RFE)”, “Coxboost”, “Principle component analysis”, “partial least squares regression for COX (PLSR)”, XGBoost, CatBoost, AdaBoost and Deep learning. The average C-index of each algorithmic combination across the entire training, test, and validation sets was used as a criterion to assess the superiority of the combinations. The combination with the highest average C-index was chosen for subsequent model construction. Consequently, the Deep Learning-Random Survival Forest (RSF) combination was selected. Using the features identified via Deep Learning, we applied RSF to establish our PCD-related prognostic model by calculating the Cell Death Index (CDI). Patients were classified as high-risk or low-risk based on the median CDI cutoff. The predictive ability of the model was assessed using the "survival" and "timeROC" R packages, generating Kaplan-Meier survival analyses and time-dependent receiver operating characteristic (ROC) area under the curve (AUC) values over a period of 1 to 5 years in TCGA-PDAC, ICGC-PACA-AU, ICGC-PACA-CA, GSE28735 and GSE78229.

2.5 Construction of the Nomogram in PDAC

2.6 Assessment of gene mutation frequency in PDAC

For this study, single nucleotide variants (SNVs) in PDAC were downloaded from TCGA database and cleaned by R package "maftools". The variant classification, variant type, and gene mutation frequency in PDAC were analyzed, and differences across groups were calculated using the chi-square test. Survival differences between KRAS, TP53, and CDKN2A mutant and wild-type PDAC were identified by K-M analysis.

2.7 Assessment of tumor microenvironment and immune infiltration analysis

Several machine learning algorithms of immune infiltration, including "EPIC", "ESTIMATE", "MCPcounter", and "Xcell" were utilized to identify immune microenvironment via the corresponding methods in the R package “IOBR”, and differences across groups were calculated using the Wilcoxon test. The significance threshold was set at p-value < 0.05. To predict immunotherapy response, the Tumor Immune Dysfunction and Exclusion (TIDE) algorithm (an online analysis tool) was employed. The expression differences of immune checkpoint genes between high-risk and low-risk groups, including PDCD1 (PD-1), PD-L1 (CD274), CTLA4, CD47, BTLA, TIGIT, TNFRSF4, TNFRSF9, and VTCN1, were investigated using the Wilcoxon test.

2.8 Exploration of the drug sensitivity in PDAC patients

Gemcitabine, docetaxel, paclitaxel, and oxaliplatin are frequently utilized in the treatment of PDAC. This study aims to estimate the half-maximal inhibitory concentration (IC50) of these drugs for PDAC using data from the CellMiner database and the expression levels of PCD-related hub genes.

2.9 Validation of the significant prognostic PCD-related gene (MYOF)

Multi-variate Cox regression analysis was performed to identify independent risk factor, including UNX1, NTAN1, FNDC3B, VCAN, MYOF, ANO6, MXRA5, SRPX2, CORO1C, and LIMS, and MYOF was selected for further investigation. The expression of MYOF between PDAC and normal tissue in GSE28735, GSE62454, and ICGC-PACA-CA cohorts were conducted. Then, survival differences between PDAC with high MYOF expression and those with low MYOF expression were calculated using Kaplan-Meier analysis.

2.10 Samples and real time quantitative PCR

Tumor tissue and paired adjacent tissues from patients diagnosed with pancreatic adenocarcinoma after surgery were procured from the First Affiliated Hospital of Chongqing Medical University and stored at -80°C. Reverse transcription-PCR (qRT-PCR) was used for quantitative analysis of gene expression in both pancreatic cancer and adjacent tissues after total RNA extraction. Trizol reagent was utilized to extract RNA from the tissues, following the manufacturer's instructions. The extracted RNA was then reverse transcribed into complementary DNA (cDNA), using RT primers and a reverse transcription reaction mix. The reaction mixture comprised of RNAse inhibitors (Sangon Biotech, Shanghai, China), MMLV RT enzyme (P7040L, Enzymatics, USA), buffer solution (B7040L, Enzymatics, USA), and dNTPs (7DN1, HyTest Ltd, Finland). The cDNA samples underwent qRT-PCR analysis utilizing qPCR Master Mix on a Gene Amp PCR System 9700. Gene-specific primers were designed for both the target gene MYOF and the reference gene GAPDH (Supplementary Table 1). The sequences for the forward and reverse primers were provided. The expression of the targeted genes was relatively quantified using the 2-ΔΔCt method.

2.11 Cell culture and regents

The PANC-1 cell line was obtained from the Cell Bank of the Type Culture Collection in Shanghai. Standard protocols were used for cell culture, with DMEM from Gibco supplemented with 10% FBS and 1% penicillin/streptomycin. The cells were maintained at 37°C in a humidified incubator with a 5% CO2 concentration. To achieve lentiviral overexpression, we followed the manufacturer's guidelines for transducing MYOF-overexpressing lentiviruses from GeneCopoeia in Guangzhou, China. For knockdown, we employed MYOF-specific short hairpin RNAs (siRNAs) obtained from Shanghai Genechem Co.

2.12 CCK-8, Wound Healing Assay and Transwell Invasion Assay

We examined the proliferative effects of PANC-1 cells with and without MYOF downregulation using the Cell Counting Kit-8 (CCK-8) assay. PANC-1 cells, including both control and si-MYOF cells, were placed individually into 96-well plates, each group seeded at a density of 8 × 10^3 cells per well. The cells were then incubated until adherence was achieved, followed by the addition of 10 µL of CCK-8 reagent to each well at 0, 24, 48, and 72-hour intervals. The CCK-8 reagent was applied to all wells in both the control and si-MYOF groups. Subsequently, the cells were incubated for two hours after the introduction of the CCK-8 reagent. The spectrophotometer then assessed the absorbance of each well at a wavelength of 450 nm following a two-hour incubation period.

The Transwell assay is used to assess the invasion ability of PANC-1 cells with diverse levels of MYOF expression. Matrigel-coated Transwell chambers from BD Sciences in Sparks, MD, USA are employed, and each of the two groups in the upper chamber is seeded with 5 × 10^4 cells. The groups consist of PANC-1 cells exhibiting downregulation of MYOF expression, PANC-1 cells that lack MYOF downregulation, PANC-1 cells overexpressing MYOF, and PANC-1 cells without MYOF overexpression. After a 24-hour incubation period, cells that penetrate the Matrigel and descend to the lower part of the chamber are fixed with a 4% formaldehyde solution. Following fixation, the fixed cells are treated with a 0.2% crystal violet solution for approximately 20 minutes. To remove excess staining solution and non-invading cells from the upper portion of the Transwell chamber, carefully wipe the inner surface with a cotton swab. Count the number of invading cells that are stained on the underside of the Transwell chamber by using an inverted microscope.

The quantity of invading cells serves as a gauge for the invasion capacity of PANC-1 cells under various MYOF expression conditions. PANC-1 cells with downregulated or unaltered MYOF expression are cultured in separate wells of a 6-well plate, achieving complete confluency to create a nearly full monolayer before evaluation. Once the cells reach the desired level of confluence, a wound is created by manually scraping the cell monolayer using a sterile pipette tip that is 200 µL in size. This procedure produces a gap, or "wound," in the cell layer that is imaged and monitored at two different time points: 0 and 24 hours after the wound creation. Images are taken using an inverted microscope. The wound area at 0 and 24 hours is then analyzed using ImageJ to determine the extent of wound closure. The decrease in wound area indicates that cells have migrated into the gap and have the capability to effectively promote wound healing.

2.13 Statistical analysis

R version 4.3.1 software was utilized for all static analysis, with statistical significance defined as a P value of less than 0.05. The following algorithms were utilized: GBM), RSF, LASSO, Ridge, Elastic network, Step-Cox, Support Vector Machine (SVM), SVM-RFE, Coxboost, PCA, PLSR, XGBoost, CatBoost, AdaBoost and Deep learning. The algorithms were executed using various R packages, including "randomForestSRC", "e1071", "gbm", "glmnet", "xgboost", "catboost", "adabag", and "h2o". The implemented machine learning algorithms in R were superpc, survivalsvm, gbm, PLSR, COXBoost, Ridge, Lasso, ElasticNet, stepwise Cox, and RSF. The packages "superpc", "survivalsvm", "gbm", "plsRcox", "CoxBoost", "glmnet", and "randomForestSRC" were utilized as they provide the necessary functions and implementations for the aforementioned algorithms.

3.1 Identification of PCD-related hub genes in PDAC

A total of 9,222 DEGs were identified between PDAC and corresponding normal tissue from the GTEx database. The DEGs were discovered based on a fold change greater than 1 (|log2FC| > 1) and a statistically significant p-value less than 0.05 (P < 0.05). For further analysis, 178 cases of PDAC were selected, and Pearson correlations among the DEGs were calculated. We then picked a soft-threshold power of β = 8 to create an adjacency matrix and a topological overlap matrix (TOM) (Fig. 1a). Using this TOM, DEGs were divided into 19 gene modules: black (514 genes), blue (793 genes), brown (735 genes), cyan (197 genes), green (624 genes), green-yellow (281 genes), gray (427 genes), gray-60 (83 genes), light cyan (164 genes), light green (59 genes), magenta (392 genes), midnight blue (177 genes), pink (449 genes), purple (327 genes), red (549 genes), salmon (278 genes), tan (281 genes), turquoise (1,195 genes), and yellow (650 genes) (Supplementary Table 2), as presented in Fig. 1(b-c). To assess the relevance between the above gene modules and enrichment scores of PCD, we applied Pearson's analysis (Fig. 1d). Thus, the brown gene module (735 genes) significantly correlated with PCD pathways, including apoptosis (cor = 0.75, p < 2.2*10^-16), entotic cell-death (cor = 0.64, p < 2.2*10^-16), ferroptosis (cor = 0.61, p < 2.2*10^-16), necroptosis (cor = 0.60, p < 2.2*10^-16), and pyroptosis (cor = 0.65, p < 2.2*10^-16), as shown in Fig. 1e. And their hub genes were extracted based on geneTraitSignificance greater than 0.4 and geneModuleMembership greater than 0.7. The 103 PCD-related hub genes were selected by interaction among apoptosis, entotic cell death, ferroptosis, necroptosis, and pyroptosis (Fig. 1f and Supplementary Table 3).

3.2 Functional and pathway enrichment analysis

KEGG analysis revealed significant enrichment of PCD-related hub genes in focal adhesion, the PI3K-Akt pathway, the Hippo pathway, and the AGE-RAGE pathway (Fig. 2a). GO analysis indicated significant enrichment in extracellular matrix organization, extracellular structure organization, and external encapsulating structure organization (Fig. 2b). Additionally, GSEA analysis demonstrated that PCD-related hub genes were significantly enriched in collagen degradation, collagen formation, and extracellular matrix organization (Fig. 2c). These findings suggest that PCD-related hub genes regulate cell death through interactions with extracellular matrix components.

3.3 Integrated development of a prognostic model in PDAC

The study identified 51 PCD-related hub genes with prognostic relevance using univariate-cox regression analysis via the “survival” R package with the criteria of p < 0.05 (Supplementary Table 4). These genes were then used in our integration program (ten-fold cross-validation) to create the prognostic models within the TCGA-PAAD cohort as training set, with two datasets (PACA-CA and GSE28735) as validation sets, which involved in 15 different machine learning algorithms (167 combinations). The average C-index of each algorithmic combination across the entire training, test, and validation sets was used as a criterion to assess the superiority of the combinations (Fig. 3a). The combination with the highest average C-index was chosen for subsequent model construction. According to the rank of average C-index, the Deep Learning-Random Survival Forest (RSF) combination was selected. Using the features identified via Deep Learning (RUNX1, NTAN1, FNDC3B, VCAN, MYOF, ANO6, MXRA5, SRPX2, CORO1C, and LIMS1), we applied RSF to establish our PCD-related prognostic model by calculating the Cell Death Index (CDI) (Fig. 3b). Patients were classified as high-risk or low-risk based on the median CDI cutoff (Fig. 3e). The predictive ability of the model was assessed using the "survival" and "timeROC" R packages, generating Kaplan-Meier survival analyses and time-dependent receiver operating characteristic (ROC) area under the curve (AUC) values over a period of 1 to 5 years. The Kaplan-Meier survival analysis demonstrated that patients with low CDI had significantly better overall survival (OS) than those with high CDI (P < 0.05) (Fig. 3c and Fig. 4a, c, e, g). The time-related ROC curves were exhibited with the AUCs of 1-, 3-, and 5-year in each cohort, respectively (Fig. 3d and Fig. 4b, d, f, h). Functional and pathway enrichment analysis of DEGs between high-CDI and low-CDI group by GSEA had indicated that these genes enriched in developmental biology, immune system, and infectious disease etc.(Fig. 3g).

3.5 Construction of Nomogram for predicting survival in PDAC

In this study, integrating various clinical information (sex, age, maximum tumor dimension, residual tumor, TNM stage, tobacco, alcohol, diabetes, and chronic pancreatitis), both univariate and multivariate Cox regression analysis suggested CDI could serve as an independent prognostic factor in PDAC patients (Figure S1a-b). Subsequently, a clinical nomogram incorporating CDI and the aforementioned clinical features (TNM stages and residual tumor) was constructed for predicting overall patient survival (OS) at 1, 2, and 3 years using the "rms" R package in TCGA-PDAC cohort (Fig. 5a). Our nomogram exhibited a AUC value for 3-year (0.944) (Fig. 5b) and high C-index (0.92) (Fig. 5c), indicating its superior performance. And the calibration curve described the predictive reliability of the nomogram model at 1-, 2-, and 3-year intervals (Fig. 5d). This indicates the high reliability and precision of the nomogram model.

3.6 Assessment of gene mutation frequency in PDAC

The research involved 84 patients in the low-risk group and 82 patients in the high-risk group. Missense mutations were the most commonly variant classification, and SNPs were the most frequent variant type in PDAC (Figure S2a-b). TP53, KRAS, and CDKN2A were frequently mutated across the groups (Fig. 6a-b). In the high-risk group, 61 patients had KRAS mutations, 57 had TP53 mutations, and 20 had CDKN2A mutations. In contrast, in the low-risk group, 39 patients had KRAS and TP53 mutations, and 8 had CDKN2A mutations. We found that the low-risk group had a significantly higher mutation frequency than the high-risk group (P < 0.05) (Table S2). Additionally, patients with wild-type genes had significantly better overall survival (OS) compared to those with mutations (P < 0.05), except for CDKN2A (Fig. 6c-e).

3.7 Immune infiltration and immunotherapy response

We applied several machine learning algorithms, including CIBERSORT, Quantiseq, and Xcell, to predict immune infiltration between the high-CDI and low-CDI groups. The CIBERSORT algorithm indicated that the low-CDI group had a higher percentage of naive B cells and CD8 + T cells (Fig. 7a). Additionally, Quantiseq analysis revealed a higher percentage of dendritic cells in the low-CDI group (Fig. 7b). Furthermore, Xcell analysis showed that the low-CDI group exhibited a higher percentage of naive CD4 + T cells, M2 macrophages, memory B cells, and NKT cells (Fig. 7c). The TIDE algorithm was used to evaluate the effectiveness of immunotherapy in PDAC, which showed that the low-CDI group had lower exclusion scores and the high-CDI group had low dysfunction scores. This suggests that the low-CDI group displays a relatively lower abundance of suppressive immune cells but a higher level of T cell dysfunction (Fig. 7d). The analysis of expression differences in immune checkpoint genes between the groups showed no significant differences in PDCD1, CD274, and CTLA4 between the two groups. However, the high-CDI group exhibited higher CD47 expression levels, while the low-CDI group had higher expression levels of BTLA and TNFRSF4 (Fig. 7e). Furthermore, low CDI had significantly better overall survival (OS) than those with high CDI (P < 0.05) in IMvigor210CoreBiologies cohort (Fig. 7f). Besides, CR/PR group had lower CDI than SD/PD group in in IMvigor210CoreBiologies (Fig. 7g).

3.8 Investigation of the potential therapeutic drugs for PDAC patients

To investigate the clinical applicability of PCD-related hub genes, we assessed the IC50 values of three common drugs for PDAC: docetaxel, paclitaxel, and oxaliplatin. The IC50 of docetaxel was significantly negatively associated with ANO6 and CORO1C. Additionally, the IC50 of paclitaxel was significantly negatively correlated with NTAN1 and CORO1C. Furthermore, the IC50 of oxaliplatin was significantly negatively associated with FNDC3B, MYOF, CORO1C, and SRPX2 (Fig. 8).

3.9 Validation of the significant prognostic PCD-related gene (MYOF) in Vitro

MYOF expression was significantly increased in PDAC compared to corresponding normal tissues in the GSE28735, GSE62452, and PACA-CA datasets (Fig. 9a). Multivariate Cox regression analysis identified MYOF as an independent risk factor (Fig. 9b). PDAC patients with high levels of MYOF expression had a worse prognosis compared to those with lower levels (Fig. 9c). Consequently, further investigation focused on MYOF. The mRNA expression of MYOF was notably greater in PDAC cells as opposed to adjacent tissues (P < 0.05) (Fig. 9d). A CCK-8 assay was conducted to evaluate cell viability, revealing that interfering with MYOF expression significantly reduced the viability of the PANC-1 cell line (Fig. 9g). Additionally, a wound healing assay performed on the PANC-1 cell line demonstrated that overexpression of MYOF enhanced cell migration, whereas the si-MYOF group exhibited an inhibitory effect (Fig. 9h). Furthermore, a transwell migration assay demonstrated that overexpression of MYOF boosted the migration of the PANC-1 cell line, whereas suppression of MYOF resulted in a significant reduction in cell migration (Fig. 9i).

PDAC is a deadly cancer with a poor prognosis, characterized by aggressive behavior and a high propensity for metastasis (2). However, PDAC treatment yielded limited success (21). Emerging evidence indicates that PCD significantly impacts PDAC. The propensity of PDAC to evade cell death is a plausible explanation for its unfavorable prognosis (22).

In previous investigations, in order to predict over survival of PDAC, several PCD-related predictive signature have been reported. By combining single-cell sequencing and transcriptome analysis, a necroptosis-related signature was constructed comprising POLR3GL, COL17A1, DDIT4, PDE4C, CLDN1, HMGA2, CENPF, and EPS8 to predict PDAC prognosis (23). In addition, 14 ferroptosis-related genes were selected to establish a prognostic model in PDAC (24). Interestingly, when considering apoptosis, pyroptosis, and necroptosis (PAN) together, it was found that a PAN-related prognostic signature effectively served as an indicator for predictive risk classification in gastric cancer (25). The multifaceted nature of PCD within tumor microenvironments defies simplistic prognostic models centered around a singular pathway. Given this complexity, there arises an imperative need for comprehensive genomic profiling to pinpoint key regulatory genes orchestrating various forms of programmed cell demise. This approach promises to enhance our predictive capabilities regarding tumor behavior and treatment response, thereby facilitating personalized oncological care.

To our knowledge, this study is the first to identify PCD-related hub genes across five PCD patterns—apoptosis, entotic cell death, ferroptosis, necroptosis, and pyroptosis—in PDAC. We discovered 103 genes as PCD-related hub genes in both internal and external cohorts. Our research has shown a significant correlation between PCD-related hub genes and overall survival in PDAC patients. These genes were incorporated into our integration program to create a prognostic signature using 15 distinct machine learning algorithms (167 combinations) in the TGCA-PAAD, GSE28735, and PACA-CA datasets. Among the 167 algorithm combinations, the combination of deep learning and RSF algorithms, which achieved the highest average C-index, was selected for further investigation. The model demonstrated outstanding predictive performance in both Kaplan-Meier survival analysis and time-dependent ROC analysis. Subsequent research revealed that the predictive model served as an independent risk factor in PDAC. Moreover, mutations in TP53, KRAS, and CDKN2A were more prevalent in the high-CDI group, indicating a poor prognosis for these patients. In particular, the low-CDI group was more sensitive to immunotherapy due to relatively lower levels of suppressive immune cells and higher levels of effector immune cells. Interestingly, the study suggests that patients exhibiting increased expression levels of ANO6, NTAN1, FNDC3B, MYOF, CORO1C, and SRPX2 genes respond better to chemotherapeutic agents such as docetaxel, paclitaxel, and oxaliplatin, which are widely used in the treatment of PDAC.

MYOF expression levels exhibited a significant increase in PDAC tissues compared to corresponding normal tissues in the GSE28735, GSE62452, and PACA-CA datasets. Moreover, patients with significantly elevated MYOF expression levels had a worse prognosis than those with lower levels. To increase the reliability of the predictive markers related to hub genes, we conducted expression and in vitro functional studies on the MYOF gene. Our study found a significant increase in MYOF expression in PDAC compared to normal tissue. MYOF was associated with key processes in PDAC cells, including proliferation, viability, invasion, and migration. Myoferlin (MYOF), a key factor in cell membrane repair, is overexpressed in PDAC and correlates with poor outcomes. Notably, MYOF localization in lysosomes is essential for lysosomal maintenance and provides an initial defense against membrane damage (26). Inhibiting MYOF can effectively restrain the proliferation and migration of PDAC cells (27). Additionally, depletion of MYOF induces mitochondrial fission, leading to decreased cell proliferation, reduced ATP production, and increased autophagy (28).

However, our study has several limitations. Further analysis is necessary to elucidate the underlying mechanisms of PCD-related hub genes. Additionally, the predictive signature's accuracy should be validated in a large cohort of PDAC patients. Moreover, further research is needed to explore the association between myoferlin (MYOF), regulated cell death, and the tumor microenvironment.

The aim of our study was to discover PCD-related hub genes in PDAC using bioinformatics and machine learning algorithms. Additionally, we developed a predictive algorithm to forecast patient survival and guide specific management strategies using 15 machine learning algorithms (167 cross-validations). Experimental findings underscore the involvement of MYOF in pivotal processes such as proliferation, viability, invasion, and migration in PDAC cells, suggesting its potential as a promising target for further investigation.

Consent for publication

The studies involving humans were approved by The First Affiliated Hospital of Chongqing Medical University Ethics committee. The all participants provided their written informed consent to participate in this study. All experiments were performed in accordance with relevant guidelines and regulations in the first affiliated hospital of Chongqing medical university. All methods were carried out in accordance with relevant guidelines and regulations. All experimental protocols were approved by The First Affiliated Hospital of Chongqing Medical University Ethics committee.

Availability of data and material

Publicly available datasets were analyzed in this study and all algorithms could be found in Supplementary Table 4. This data can be found here: TCGA: GDC Data Portal Homepage (cancer.gov), GEO: Home - GEO - NCBI (nih.gov), the ICGC: ICGC Data Portal | Retirement Notice, GTEx: https://xenabrowser.net/datapages/, Genecard database: https://www.genecard.org, MSigDB: GSEA | MSigDB (gsea-msigdb.org), GEPIA: GEPIA (Gene Expression Profiling Interactive Analysis) (cancer-pku.cn), TIDE: Tumor Immune Dysfunction and Exclusion (TIDE) (harvard.edu).

Competing interests

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Funding

This study was funded by the National Natural Science Foundation of China under grants No. 81871261.

Authors’ contributions

WZW and CY conceived the review, acquired data, and drafted the manuscript. HK and JSM undertook the initial research. WZW and HK was involved in writing and reviewing the manuscript. All authors read and approved the final manuscript.

Acknowledgements

We thank GEO, TCGA, ICGC, and GTEx databases for providing invaluable data for statistical analyses.

Statement

All the supplementary file would get published along with the main manuscript file. During the preparation of this work the author(s) used OpenAI in order to revise and refine the article. After using this tool, the author(s) reviewed and edited the content as needed and take(s) full responsibility for the content of the publication.

Siegel RL, Miller KD, Wagle NS, Jemal A. Cancer statistics, 2023. CA Cancer J Clin. 2023;73(1):17-48.
Zhou Z, Zhang Y, Li J, Weng S, Li J, Chen S, et al. Crosstalk between regulated cell death and immunity in redox dyshomeostasis for pancreatic cancer. Cell Signal. 2023;109:110774.
Kleeff J, Korc M, Apte M, La Vecchia C, Johnson CD, Biankin AV, et al. Pancreatic cancer. Nat Rev Dis Primers. 2016;2:16022.
Wood LD, Canto MI, Jaffee EM, Simeone DM. Pancreatic Cancer: Pathogenesis, Screening, Diagnosis, and Treatment. Gastroenterology. 2022;163(2):386-402.e1.
Peng F, Liao M, Qin R, Zhu S, Peng C, Fu L, et al. Regulated cell death (RCD) in cancer: key pathways and targeted therapies. Signal Transduct Target Ther. 2022;7(1):286.
Legrand AJ, Konstantinou M, Goode EF, Meier P. The Diversification of Cell Death and Immunity: Memento Mori. Mol Cell. 2019;76(2):232-42.
Tang D, Kang R, Berghe TV, Vandenabeele P, Kroemer G. The molecular machinery of regulated cell death. Cell Res. 2019;29(5):347-64.
Galluzzi L, Vitale I, Aaronson SA, Abrams JM, Adam D, Agostinis P, et al. Molecular mechanisms of cell death: recommendations of the Nomenclature Committee on Cell Death 2018. Cell Death Differ. 2018;25(3):486-541.
Hanahan D, Weinberg RA. Hallmarks of cancer: the next generation. Cell. 2011;144(5):646-74.
Yarden Y, Sliwkowski MX. Untangling the ErbB signalling network. Nat Rev Mol Cell Biol. 2001;2(2):127-37.
Deramaudt T, Rustgi AK. Mutant KRAS in the initiation of pancreatic cancer. Biochim Biophys Acta. 2005;1756(2):97-101.
Zhu H, Chen K, Chen Y, Liu J, Zhang X, Zhou Y, et al. RNA-binding protein ZCCHC4 promotes human cancer chemoresistance by disrupting DNA-damage-induced apoptosis. Signal Transduct Target Ther. 2022;7(1):240.
Mahadevan KK, LeBleu VS, Ramirez EV, Chen Y, Li B, Sockwell AM, et al. Elimination of oncogenic KRAS in genetic mouse models eradicates pancreatic cancer by inducing FAS-dependent apoptosis by CD8(+) T cells. Dev Cell. 2023;58(17):1562-77.e8.
Hu Z, Yuan J, Long M, Jiang J, Zhang Y, Zhang T, et al. The Cancer Surfaceome Atlas integrates genomic, functional and drug response data to identify actionable targets. Nat Cancer. 2021;2(12):1406-22.
Milella M, Luchini C, Lawlor RT, Johns AL, Casolino R, Yoshino T, et al. ICGC-ARGO precision medicine: familial matters in pancreatic cancer. Lancet Oncol. 2022;23(1):25-6.
Edgar R, Domrachev M, Lash AE. Gene Expression Omnibus: NCBI gene expression and hybridization array data repository. Nucleic Acids Res. 2002;30(1):207-10.
Li J, Miao B, Wang S, Dong W, Xu H, Si C, et al. Hiplot: a comprehensive and easy-to-use web service for boosting publication-ready biomedical data visualization. Brief Bioinform. 2022;23(4).
Liu C, Zhou Y, Zhou Y, Tang X, Tang L, Wang J. Identification of crucial genes for predicting the risk of atherosclerosis with system lupus erythematosus based on comprehensive bioinformatics analysis and machine learning. Comput Biol Med. 2023;152:106388.
Yip AM, Horvath S, editors. The Generalized Topological Overlap Matrix for Detecting Modules in Gene Networks. Biocomp; 2006: Citeseer.
Peng Z, Ye M, Ding H, Feng Z, Hu K. Spatial transcriptomics atlas reveals the crosstalk between cancer-associated fibroblasts and tumor microenvironment components in colorectal cancer. J Transl Med. 2022;20(1):302.
Sarantis P, Koustas E, Papadimitropoulou A, Papavassiliou AG, Karamouzis MV. Pancreatic ductal adenocarcinoma: Treatment hurdles, tumor microenvironment and immunotherapy. World J Gastrointest Oncol. 2020;12(2):173-81.
Chen X, Zeh HJ, Kang R, Kroemer G, Tang D. Cell death in pancreatic cancer: from pathogenesis to therapy. Nat Rev Gastroenterol Hepatol. 2021;18(11):804-23.
Chen L, Zhang X, Zhang Q, Zhang T, Xie J, Wei W, et al. A necroptosis related prognostic model of pancreatic cancer based on single cell sequencing analysis and transcriptome analysis. Front Immunol. 2022;13:1022420.
Jiang P, Yang F, Zou C, Bao T, Wu M, Yang D, et al. The construction and analysis of a ferroptosis-related gene prognostic signature for pancreatic cancer. Aging (Albany NY). 2021;13(7):10396-414.
Huo J, Xie W, Fan X, Sun P. Pyroptosis, apoptosis, and necroptosis molecular subtype derived prognostic signature universal applicable for gastric cancer-A large sample and multicenter retrospective analysis. Comput Biol Med. 2022;149:106037.
Gupta S, Yano J, Mercier V, Htwe HH, Shin HR, Rademaker G, et al. Lysosomal retargeting of Myoferlin mitigates membrane stress to enable pancreatic cancer growth. Nat Cell Biol. 2021;23(3):232-42.
Gu H, Zhang T, Li Y, He Y, Guan T, Kan W, et al. Discovery of 1,5-diaryl-1,2,4-triazole derivatives as myoferlin inhibitors and their antitumor effects in pancreatic cancer. Future Med Chem. 2022;14(20):1425-40.
Rademaker G, Hennequière V, Brohée L, Nokin MJ, Lovinfosse P, Durieux F, et al. Myoferlin controls mitochondrial structure and activity in pancreatic ductal adenocarcinoma, and affects tumor aggressiveness. Oncogene. 2018;37(32):4398-412.

No competing interests reported.

Download PDF

Editorial decision: Revision requested
11 Jul, 2024
Editor assigned by journal
06 Jul, 2024
Submission checks completed at journal
06 Jul, 2024
First submitted to journal
01 Jul, 2024

You are reading this latest preprint version

Discovering a specialized programmed-cell death patterns for prognostic model of pancreatic ductal carcinoma via machine learning

Status:

Version 1

Abstract

Figures

1. Introduction

2. Methods

2.1 Datasets

2.2 Identification of PCD-related hub genes in PDAC

2.3 Functional and pathway enrichment analysis

2.4 Establishment of the PCD-related prognostic model in PDAC

2.5 Construction of the Nomogram in PDAC

2.6 Assessment of gene mutation frequency in PDAC

2.7 Assessment of tumor microenvironment and immune infiltration analysis

2.8 Exploration of the drug sensitivity in PDAC patients

2.9 Validation of the significant prognostic PCD-related gene (MYOF)

2.10 Samples and real time quantitative PCR

2.11 Cell culture and regents

2.12 CCK-8, Wound Healing Assay and Transwell Invasion Assay

2.13 Statistical analysis

3 Result

3.1 Identification of PCD-related hub genes in PDAC

3.2 Functional and pathway enrichment analysis

3.3 Integrated development of a prognostic model in PDAC

3.5 Construction of Nomogram for predicting survival in PDAC

3.6 Assessment of gene mutation frequency in PDAC

3.7 Immune infiltration and immunotherapy response

3.8 Investigation of the potential therapeutic drugs for PDAC patients

3.9 Validation of the significant prognostic PCD-related gene (MYOF) in Vitro

4 Discussion

5 Conclusion

Declarations

References

Additional Declarations

Supplementary Files

Status:

Version 1