Revealing a cancer-associated fibroblast-based risk signature for pancreatic adenocarcinoma through single-cell and bulk RNA-seq analysis

Abstract

Purpose: Proliferation of stromal connective tissue is a hallmark of pancreatic adenocarcinoma (PAAD). The engagement of activated cancer-associated fibroblasts (CAFs) contributes to the progression of PAAD through their involvement in tumor fibrogenesis. However, the prognostic significance of CAF-based risk signature in PAAD has not been explored.

Methods: The single-cell RNA sequencing (scRNA-seq) data sourced from GSE155698 within the Gene Expression Omnibus (GEO) database was supplemented by bulk RNA sequencing data from The Cancer Genome Atlas (TCGA) and microarray data retrieved from the GEO database. The scRNA-seq data underwent processing via the Seurat package to identify distinct CAF clusters utilizing specific CAF markers. Differential gene expression analysis between normal and tumor samples was conducted within the TCGA-PAAD cohort. Univariate Cox regression analysis pinpointed genes associated with CAF clusters, identifying prognostic CAF-related genes. These genes were utilized in LASSO regression to craft a predictive risk signature. Subsequently, integrating clinicopathological traits and the risk signature, a nomogram model was constructed.

Results: Our scRNA-seq analysis unveiled four distinct CAF clusters in PAAD, with two linked to PAAD prognosis. Among 207 identified DEGs, 148 exhibited significant correlation with these CAF clusters, forming the basis of a seven-gene risk signature. This signature emerged as an independent predictor in multivariate analysis for PAAD and demonstrated predictive efficacy in immunotherapeutic outcomes. Additionally, a novel nomogram, integrating age and the CAF-based risk signature, exhibited robust predictability and reliability in prognosticating PAAD. Moreover, the risk signature displayed substantial correlations with stromal and immune scores, as well as specific immune cell types.

Conclusions: The prognosis of PAAD can be accurately predicted using the CAF-based risk signature, and a thorough analysis of the PAAD CAF signature may aid in deciphering the patient's immunotherapy response and presenting fresh cancer treatment options.

Keywords: cancer-associated fibroblasts; immunotherapy; nomogram; pancreatic adenocarcinoma; risk signature.

Figure 1

Profiling of CAF subpopulations. (A) UMAP plot displaying the distribution of four distinct CAF subpopulations post-clustering. (B) UMAP plot illustrating the expression of CAF marker genes (ACTA2, FAP, PDGFRB, NOTCH3, DCN, and COL1A2). (C) Dot plot showcasing the top 5 marker gene expressions across the four CAF clusters. (D) Relative proportions and cell numbers within each sample for the four CAF clusters. (E) UMAP plot delineating the distribution between malignant and non-malignant cells. (F) KEGG enrichment analysis of DEGs observed across the four CAF clusters.

Figure 2

Tumor-related pathway characteristics in CAF clusters. (A) Heatmap illustrating GSVA scores for ten tumor-related pathways enriched in both malignant and non-malignant CAF cells. (B) Comparison of malignant and non-malignant cell proportions across different CAF clusters. Comparative analysis of GSVA scores for ten tumor-related pathways between malignant and non-malignant cells within (C) CAF_1, (D) CAF_2, (E) CAF_3, and (F) CAF_4 clusters. *P < 0.05; **P < 0.01; ***P < 0.001; ns, not statistically significant.

Figure 3

Identification of CAF-associated hub genes with prognostic significance. (A) Volcano plot illustrating DEGs between tumor and normal tissues in the TCGA-PAAD cohort. (B) Volcano plot showcasing prognosis-related genes identified through univariate Cox regression analysis. (C) Functional enrichment analyses encompassing GO (BP, CC, and MF) and KEGG analyses of CAF-related DEGs. (D) Trajectory plot depicting each independent variable with lambda in the LASSO model for PAAD. (E) LASSO coefficient profiles highlighting the seven genes in PAAD. The plot shows coefficient profiles against the log (lambda) sequence. (F) Multivariate Cox coefficients for each gene in the risk signature. (G) Kaplan-Meier curves illustrating the risk model constructed using the seven genes in the TCGA-PAAD cohort. (H) ROC curves displaying the risk model constructed with the seven genes in the TCGA-PAAD cohort. Kaplan-Meier curves of the risk model constructed with the seven genes in the validation datasets (I) GSE78229 and (J) GSE85916.

Figure 4

Creation of a nomogram using CAF-related gene signature for PAAD prognostication. (A) Univariate and (B) multivariate Cox regression analyses involving risk score and clinicopathological characteristics. (C) Development of a nomogram model amalgamating age and risk score. (D) Calibration plots showcasing the prediction accuracy for 1-, 2-, and 3-year survival probabilities. (E) Decision curve analysis illustrating the nomogram’s utility. (F) Time-ROC curve analysis comparing the predictive performance of the nomogram against other factors.

Figure 5

Genetic profile of the seven genes in the risk signature. (A) Waterfall diagram illustrating SNV mutations of the seven pivotal genes. (B) Enrichment heatmap displaying key pathways associated with SNV data in PAAD. (C) CNV alterations in the seven crucial genes, showcasing instances of gain, loss, and absence of alterations. (D) Heatmap visualizing correlations between the seven pivotal genes and Aneuploidy Score, Homologous Recombination Defects, Fraction Altered, Number of Segments, and Nonsilent Mutation Rate. (E) Heatmap revealing gene-pathway correlations. (F) Heatmap illustrating enrichment scores for pathways. *P < 0.05, **P < 0.01, ***P < 0.001.

Figure 6

Risk signature response to immunotherapy in IMvigor210 and GSE78220 cohorts. (A) Prognostic differences among IMvigor210 cohort subgroups based on the risk score. (B) Variations in risk scores within IMvigor210 cohort responses to immunotherapy. (C) Distribution of immunotherapy responses among risk score groups in the IMvigor210 cohort. (D) Prognostic differences among subgroups of early-stage patients in the IMvigor210 cohort based on the risk score. (E) Prognostic differences among subgroups of advanced-stage patients in the IMvigor210 cohort based on the risk score. (F) Prognostic differences among subgroups of the GSE78220 cohort based on the risk score. (G) Variations in risk scores among GSE78220 cohort responses to immunotherapy. (H) Distribution of immunotherapy responses among risk score groups in the GSE78220 cohort.