Prognostic value and immune relevancy of a combined autophagy-, apoptosis- and necrosis-related gene signature in glioblastoma

Abstract

Background: Glioblastoma (GBM) is considered the most malignant and devastating intracranial tumor without effective treatment. Autophagy, apoptosis, and necrosis, three classically known cell death pathways, can provide novel clinical and immunological insights, which may assist in designing personalized therapeutics. In this study, we developed and validated an effective signature based on autophagy-, apoptosis- and necrosis-related genes for prognostic implications in GBM patients.

Methods: Variations in the expression of genes involved in autophagy, apoptosis and necrosis were explored in 518 GBM patients from The Cancer Genome Atlas (TCGA) database. Univariate Cox analysis, least absolute shrinkage and selection operator (LASSO) analysis, and multivariate Cox analysis were performed to construct a combined prognostic signature. Kaplan-Meier survival, receiver-operating characteristic (ROC) curves and Cox regression analyses based on overall survival (OS) and progression-free survival (PFS) were conducted to estimate the independent prognostic performance of the gene signature. The Chinese Glioma Genome Atlas (CGGA) dataset was used for external validation. Finally, we investigated the differences in the immune microenvironment between different prognostic groups and predicted potential compounds targeting each group.

Results: A 16-gene cell death index (CDI) was established. Patients were clustered into either the high risk or the low risk groups according to the CDI score, and those in the low risk group presented significantly longer OS and PFS than the high CDI group. ROC curves demonstrated outstanding performance of the gene signature in both the training and validation groups. Furthermore, immune cell analysis identified higher infiltration of neutrophils, macrophages, Treg, T helper cells, and aDCs, and lower infiltration of B cells in the high CDI group. Interestingly, this group also showed lower expression levels of immune checkpoint molecules PDCD1 and CD200, and higher expression levels of PDCD1LG2, CD86, CD48 and IDO1.

Conclusion: Our study proposes that the CDI signature can be utilized as a prognostic predictor and may guide patients' selection for preferential use of immunotherapy in GBM.

Keywords: CGGA; Cell death index; Glioblastoma (GBM); Immune infiltration; Prognostic; TCGA.

Fig. 1

A) Combined cell death index (CDI) was generated, which included the highest expression of genes involved in autophagy, apoptosis, and necrosis. B) The interaction between CDI signature genes in GBM. The circle size represented the effect of each signature gene on the prognosis, and the range of values calculated by Log-rank test was p < 0.001, p < 0.01, p < 0.05 and P < 0.1, respectively. The Autophagy, apoptosis and necrosis signature gene was marked with blue, yellow and red respectively. Green dots in the circle represent protective factors of prognosis and black dots in the circle represent risk factors of prognosis. The lines linking signature genes showed their interactions, and thickness showed the correlation strength between genes. Negative correlation was marked with blue and positive correlation with red

Fig. 2

A) Volcano plot showing the differential expression of genes between high CDI and low CDI patients (p < 0.05, |log2 fold-change|> 1). B) The differential analysis of copy number variations between two groups was visualized by Circos plot, which revealed that compared with the low CDI group, 225 (35.3%) genes were significantly amplified, and 309 (48.4%) were significantly deleted in the high CDI group. Red dots represented amplifications and blue dots represented deletions. C) Left panel: Circos plots of each risk group revealing the amplifications and deletions of chromosomes. Right panel: Boxplots inhibited more burdens of copy number amplifications and deletions in high CDI group. D) Waterfall plots of 15 most frequently mutated DEGs which were altered in 126 GBM samples. E) Waterfall plots showed the top 10 mutated in high risk and low risk group. F) The proportion of mutation status of PTEN, ATRX, TP53 and EGFR in the two groups. G) Violin plots of CDI risk score in individual samples of GBM patients, stratified by IDH, MGMT promoter and TERT mutation status

Fig. 3

Distribution of the CDI in TCGA training cohort. A) Classification of patients into different risk groups based on CDI. B) Distribution of patients’ survival time and status. C) Kaplan–Meier survival analysis suggested that patients’ OS were significantly different between high CDI and low CDI group. D) The prognostic performance of CDI demonstrated by ROC curves for predicting the 1.5-, 3.0-, and 4.5- year OS rates. E) Kaplan–Meier survival analysis suggested that patients’ PFS were significantly different between high CDI and low CDI group. F) The prognostic performance of CDI demonstrated by ROC curves for predicting the 1.5-, 3.0-, and 4.5- year PFS rates

Fig. 4

Distribution of the CDI in CGGA validation cohort. A) Classification of patients into different risk groups based on the CDI. B) distribution of patients’ survival time and status. C) Kaplan–Meier survival analysis suggested that patients’ OS were significantly different between high CDI and low CDI group. D) The prognostic performance of CDI demonstrated by ROC curves for predicting the 1.5-, 3.0-, and 4.5- year OS rates

Fig. 5

Volcano plot showing the differential expression of cytokines (n = 265 genes) between patients in the high and low significant gene expression of individual cell death pathway (p < 0.05): A) apoptosis, B) autophagy, C) necrosis, D) CDI. E) Gene Ontology (GO) functional enrichment of genes with higher expression in high CDI group

Fig. 6

A) Heatmap for immune responses based on XCELL, TIMER algorithms, CIBERSORT, quanTIseq, MCPcounter and ssGSEA among high and low risk group. B) Correlations between CDI riskscore and the immune cells using Spearman analysis. Negative correlation was marked with blue and positive correlation with red. Immune cell distribution and immune function in patients with the high and low significant gene expression of individual cell death pathway: C) apoptosis, D) autophagy, E) necrosis, F) CDI

Fig. 7

Comparison of immune checkpoint blockade–related genes expression levels in patients with the high and low significant gene expression of individual cell death pathway: A) apoptosis; B) autophagy, C) necrosis, D) CDI. E) Association analyses between immune checkpoint inhibitors CD86, CD40, PDCD1, CD48, CD200, IDO1 and PDCD1LG2 and CDI. F) Correlation plot showing the association between CDI risk model and CD86, PDCD1LG2, CD48, CD40, IDO1, CD200 and PDCD1

Fig. 8

Functional enrichment analysis of highly expressed genes (log2 fold-change > 1) in the A) high risk group; and B) low risk group. Pathways enriched in high and low CDI groups: Molecular Function (GO: MF) and Biological Process (GO: BP) of the (C, D) high CDI group and (E, F) the low CDI group. G) Heatmap illustrated the enrichment scores of 20 differentially enriched molecular pathways evaluated by GSVA analysis between low CDI and high CDI patients. Red represented high enrichment scores, and blue represented low enrichment scores. H) PPI network of hubba DEGs obtained from the DEGs network