An Immunity-Related Gene Model Predicts Prognosis in Cholangiocarcinoma

Abstract

The prognosis of patients with cholangiocarcinoma (CCA) is closely related to both immune cell infiltration and mRNA expression. Therefore, we aimed at conducting multi-immune-related gene analyses to improve the prediction of CCA recurrence. Immune-related genes were selected from the Gene Expression Omnibus (GEO), The Cancer Genome Atlas (TCGA), and the Immunology Database and Analysis Portal (ImmPort). The least absolute shrinkage and selection operator (LASSO) regression model was used to establish the multi-gene model that was significantly correlated with the recurrence-free survival (RFS) in two test series. Furthermore, compared with single genes, clinical characteristics, tumor immune dysfunction and exclusion (TIDE), and tumor inflammation signature (TIS), the 8-immune-related differentially expressed genes (8-IRDEGs) signature had a better prediction value. Moreover, the high-risk subgroup had a lower density of B-cell, plasma, B-cell naïve, CD8+ T-cell, CD8+ T-cell naïve, and CD8+ T-cell memory infiltration, as well as more severe immunosuppression and higher mutation counts. In conclusion, the 8-IRDEGs signature was a promising biomarker for distinguishing the prognosis and the molecular and immune features of CCA, and could be beneficial to the individualized immunotherapy for CCA patients.

Keywords: LASSO; TCI; cholangiocarcinoma; immunity; prognosis.

Figure 1

Identification of immune-related differentially expressed genes in CCA from the dataset. (A) Volcano plots of DEGs in the GSE76297 dataset. (B) Volcano plots of DEGs in the GSE26566 dataset. (C) Volcano plots of DEGs in the GSE119336 dataset. (D) Volcano plots of DEGs in the GSE89749 dataset. (E) Volcano plots of DEGs in the TCGA dataset [x-axis: log2(FC); y-axis: −log10(FDR) for each gene. Genes with FDR < 0.01 and FC >1.5 or <−1.5 were considered as DEGs in TCGA. Blue: downregulated genes; Gray: non-differential genes; Red: upregulated genes]. (F) Overlapping analyses of DEGs in GSE76297, GSE26566, GSE119336, and GSE89749 groups; DEGs shared within 2 datasets or more were regarded as credible DEGs in each Venn diagram. (G) Overlapping analysis of GEO, TCGA, and ImmPort datasets.

Figure 2

Construction of an 8-IRDEGs signature from the TCGA cohort. (A) Tenfold cross-validation for tuning parameter selection in the LASSO model. The dotted vertical lines are drawn at the optimal values by minimum criteria (lambda. min, left vertical dotted line) and 1-SE criteria (lambda.1se, right vertical dotted line). (B) LASSO model at optimal lambda value; 8 mRNAs with non-zero coefficients were selected.

Figure 3

Correlations between the prognostic signature-derived risk score and infiltration abundances of multiple immune cells. (A) B cells, (B) B-cell plasma, (C) B-cell naïve, (D) CD4+ T-cell memory, (E) CD4+ T cell (Th1), (F) CD4+ T cell (Th2), (G) CD8+ T cell, (H) CD8+ T-cell naïve, and (I) CD8+ T-cell memory (Pearson correlation analysis).

Figure 4

Evaluation of the 8-IRDEGs signature for relapse in the TCGA cohort. (A) Distribution of the risk score derived from the signature. Patients are ranked according to the corresponding risk score. (B) Survival status of CCA patients in different risk subgroups. (C) The Kaplan–Meier survival curve of recurrence-free for patients between two different groups. (D) Time-dependent ROC curve at 1, 3, and 5 years. (E) Comparison of prognostic accuracy between the signature and single mRNAs. (F) Comparison of prognostic accuracy between the signature and clinical characteristics. p-values were calculated using the log-rank test. HR, hazard ratio; AUC, area under the ROC curve; RFS, recurrence-free survival.

Figure 5

Mutation analysis and GSEA. (A) Significantly mutated genes in the mutated CCA samples of the low-risk subgroup. (B) Significantly mutated genes in the mutated CCA samples of the high-risk subgroup. [Samples (columns) are arranged to emphasize mutual exclusivity among mutations. The right panel shows the mutation percentage, and the top panel shows the overall number of mutations. The color coding indicates the mutation type.] (C) Gene sets enriched in the low-risk subgroup (p < 0.05). (D) Gene sets enriched in the high-risk subgroup (p < 0.05). (E) KEGG pathway in the low-risk subgroup (p < 0.05). (F) KEGG pathway in the high-risk subgroup (p < 0.1).

Figure 6

The Kaplan–Meier survival analysis for immune functions. (A) The difference of immune functions between the high-risk and low-risk subgroups. (B) The KM curve of B cells. (C) The KM curve of Mast cells. (D) The KM curve of T helper cells. (E) The KM curve of CD8+T cells. (F) The KM curve of the checkpoint. (G) The KM curve of Macrophages. (H) The KM curve of Neutrophils. (I) The KM curve of NK cells. (J) The KM curve of Tfh cells. (K) The KM curve of Th2 cells. (L) The KM curve of TIL. (M) The KM curve of Treg (*p < 0.05, ** p < 0.01).

Figure 7

Evaluation of the risk score formula for relapse prediction in the TCGA cohort. (A) Scatter plot for the distribution of risk score and relapse status of individual patients. (B) Survival status of CCA patients in the two 8-IRDEGs signature subgroups. (C) The Kaplan–Meier survival curve of recurrence-free for patients between two different groups. (D) Time-dependent ROC curve at 1, 3, and 5 years. (E) Comparison of prognostic accuracy between the signature and single mRNAs. (F) Comparison of prognostic accuracy between the signature and clinical characteristics. p-values were calculated using the log-rank test. HR, hazard ratio; AUC, area under the ROC curve; RFS, recurrence-free survival.