Skip to main content

Immune cell gene expression signatures in diffuse glioma are associated with IDH mutation status, patient outcome and malignant cell state, and highlight the importance of specific cell subsets in glioma biology

Abstract

The tumor micro-environment (TME) plays an important role in various cancers, including gliomas. We estimated immune cell type-specific gene expression profiles in 3 large clinically annotated glioma datasets using CIBERSORTx and LM22/LM10 blood-based immune signatures and found that the proportions and estimated gene expression patterns of specific immune cells significantly varied according to IDH mutation status. When IDH-WT and IDH-MUT tumors were considered separately, cluster-of-cluster analyses of immune cell gene expression identified groups with distinct survival outcomes. We confirmed and extended these findings by applying a signature matrix derived from single-cell RNA-sequencing data derived from 19 glioma tumor samples to the bulk profiling data, validating findings from the LM22/LM10 results. To link immune cell signatures with outcomes in checkpoint therapy, we then showed a significant association of monocytic lineage cell gene expression clusters with patient survival and with mesenchymal gene expression scores. Integrating immune cell-based gene expression with previously described malignant cell states in glioma demonstrated that macrophage M0 abundance significantly correlated with mesenchymal state in IDH-WT gliomas, with evidence of a previously implicated role of the Oncostatin-M receptor and macrophages in the mesenchymal state. Among IDH-WT tumors that were enriched for the mesenchymal cell state, the estimated M0 macrophage expression signature coordinately also trended to a mesenchymal signature. We also examined IDH-MUT tumors stratified by 1p/19q status, showing that a mesenchymal gene expression signature the M0 macrophage fraction was enriched in IDH-MUT, non-codeleted tumors. Overall, these results highlight the biological and clinical significance of the immune cell environment related to IDH mutation status, patient prognosis and the mesenchymal state in diffuse gliomas.

Introduction

Gliomas are most aggressive and most common primary brain tumors that account for about 80% of all brain malignancies and include lower-grade gliomas (LGG) and glioblastomas (GBM) [46]. LGG typically range from grade 2 to 3 and have better survival, while GBM are categorized as high-grade tumor (grade 4) with poor survival [46]. Previous studies have highlighted the tumor heterogeneity within both LGG and GBM, where LGG were classified into three subtypes: IDH gene mutation with the chromosome 1p/19q codeletion, IDH mutation without 1p/19q codeletion, and IDH wildtype [7], while GBM were categorized into four molecular subtypes: proneural, neural, classical, and mesenchymal [64]. Recently, in the new classification by the World Health Organization (WHO), diffuse gliomas in adults can broadly be divided into three molecular subtypes based on their IDH mutation and 1p/19q codeletion status: (i) IDH wildtype glioma (GBM), (ii) IDH-mutant glioma without 1p/19q codeletion (IDH-mutant astrocytoma), and (iii) IDH mutant glioma with 1p/19q codeletion (IDH-mutant and 1p/19q co-deleted oligodendroglioma) [39].

Current treatment mainly involves neurosurgical resection followed by radiotherapy and temozolomide (TMZ) chemotherapy. Despite such therapies, GBM remains an incurable tumor with a median survival of only 15 months [68]. Recently, different immunotherapeutic strategies have been explored [14], such as immune checkpoint inhibitors and chimeric antigen receptor T cells used in several extracranial cancers [5, 18, 19]. Prior literature suggests that histopathologically similar tumors can respond differently to specific treatments [1, 11, 42]. To account for this heterogeneity, the tumor microenvironment (TME) can be crucial to select the immunotherapy responses [60]. Thus, the complex interplay between tumor cells and their tumor immune microenvironment can influence the outcome of immunotherapy and many other anti-cancer therapies to stratify patients [12, 19, 33]. Recent studies have noticed subtype-specific enrichment of immune cells [61] and their prognostic association in multiple tumors [4, 27, 28, 48, 74]. For example, higher levels of immune cell infiltration are associated with HER-2 positive and triple-negative breast cancers [18, 30]. The immunological heterogeneity of gliomas has also been investigated, revealing the predominant antitumor immune response in the mesenchymal subtype of GBM [13] and identification of immune-specific subtypes in diffuse LGG [72]. Recently, a comprehensive study has analyzed the tumor microenvironment of the brain and demonstrated significant enrichment of tumor-associated macrophages (TAMs) between glioma subtypes and brain metastasis [34]. Additional work has highlighted the role of macrophages in inducing mesenchymal-like state in glioblastoma, showing that macrophage-derived oncostatin M interacts with its receptor on glioma cells which in turn promotes mesenchymal state via STAT3 signaling [24]. Such complex interplay between tumor cells and the tumor immune microenvironment in glioma is therefore worthy of further exploration where deconvolution of bulk expression data in clinically annotated datasets could yield additional insights into this biology as well as its clinical relevance.

Deconvolution of gene expression data is an increasingly used tool to estimate immune cell type proportions, as well as immune cell specific gene expression profiles. CIBERSORTx is an analytical tool that has been proposed to impute gene expression profile from bulk data, providing an estimation of the abundances of member cell types in a mixed cell population [45, 57]. The LM22 and LM10 signature matrices is commonly used to estimate cell type abundance and gene expression of 22 and 10 immune cell types, respectively. Using this approach, we hypothesized that the glioma subtypes classified based on their IDH mutation and 1p/19q codeletion status may harbor distinct tumor immune microenvironments and that specific microenvironmental signatures would correlate to patient outcome within IDH-WT and IDH-MUT gliomas. We also generate a single cell RNAseq dataset of 19 glioma samples to generate a signature matrix highly relevant for deconvolution of bulk glioma expression data to extend key findings from the LM10 and LM22 signatures. To fully explore and evaluate these hypotheses, we use 3 independent glioma datasets (The Cancer Genome Atlas (TCGA) and two Chinese Glioma Genome Atlas (CGGA) datasets) that offered the opportunity for reproducibility and validation of findings. In addition to demonstrating the subtype specific prognostic impact of infiltrating immune cells and their association with the tumor component, our study uncovered a comprehensive framework for future studies designed to characterize the interplay between tumor cells and their surrounding immune cells.

Methods

Glioma datasets and immune cell identification

We downloaded publicly available transcriptomic datasets for gliomas from two independent resources: TCGA data portal and Chinese Glioma Genome Atlas (CGGA). From TCGA, we considered LGG and GBM while from CGGA, CGGA325 [3, 75], and CGGA693 [38, 67] were used. These encompass a total of 1721 high to low grade gliomas. Read count matrices from TCGA were normalized to TPM (Transcripts Per Kilobase Million) using a shiny R based COEX-seq package (https://github.com/kimsc77/COEX-seq) [32], while the CGGA datasets were already processed and ready to use for the downstream analysis. These bulk datasets were then subjected to deconvolution against a signation matrix by using CIBERSORTx with its default parameters and B-mode batch-correction [57]. For this study, we employed LM22 and LM10 signatures, where LM22 encompass signatures from 22 individual immune cell types while LM10 is a condensed matrix comprising of 10 broad immune cell types from LM22. Notably, In case of LM10, monocyte includes M0, M1 and M2 macrophages into a single unit while T-CD4 includes T-naïve, T-mem-activ, T-mem-rest, T-gam-del, T-help and T-reg cells. We also included an additional checkpoint inhibitor treated glioma dataset encompassing 29 tumors (GSE121810) and was deconvolved with LM10 signatures due to its small sample size.

Tumor cell state deconvolution

To deconvolve previously desctibed tumor cell states we used publicly available scRNA Seq datasets from GSE70630 and GSE89567 for IDH-MUT samples, and from GSM3828672 for IDH-WT samples. The cells were labeled by using a previously published method https://github.com/jlaffy/scalop with its scalop::sigScores and scalop::maxcol_strict functions to score cells using previously defined meta-modules. Cells from IDH-WT tumors were processed with MES-like, NPC-like, AC-like and OPC-like meta-modules [44], while those from IDH-MUT were processed with Oligo-like, Astro-like and Stem-like specific genes [63]. Cell having a minimum score difference of 0.3 between the two maximum scoring modules were selected to generate signature matrix by using CIBERSORTx [57]. These malignant cell state signatures from both IDH-WT and IDH-MUT tumors were then independently used to deconvolve the 3 bulk transcriptomic datasets with the default parameters and B-mode batch-correction.

Cell type-specific gene expression profiling

Unsupervised clustering of gene expression deconvolution data was performed by Seurat v3. For the dimensionality reduction, we used most variable 3000 genes with “Seurat::RunUMAP” function and identified cell type specific clusters [23]. We also explored subclusters of each cell types by using K-means clustering with nstart = 25, iter.max = 500, algorithm = "Hartigan-Wong" and silhouette method to compute the optimal number of clusters. These were visualized by using R based Rtsne package and ggplot function.

Identifying closely related clusters of cell type clusters

To identify similarly behaving clusters from all the cell types or cell states based on their sample memberships we employed a binary approach, where "1" denotes the presence of a sample in a cluster while "0" represents its absence. These clusters were then subjected to unsupervised clustering by using a R based ConsensusClusterPlus package with distance = pearson and maxK = 6, reps = 50 [70]. We used NbClust R package to determine the optimal number of clusters and the results were visualized by using a R based ComplexHeatmap package [20].

Functional annotation of the identified clusters

We used R based singscore package to score each sample against an immune signature, pathways, or biological processes [15]. The hallmark epithelial-mesenchymal transition (EMT) gene-set used in the analysis were downloaded from Molecular Signatures Database (MSigDB) [37]. We obtained Cytolytic (CYT) score by calculating the geometric mean of GZMA and PRF1 expression [51] in each dataset from each of the identified clusters. Similarly, T cell exhaustion (TCE) score was calculated by using PD-1, CTLA-4, TIM-3, LAG-3, CD160, 2B4, TIGIT, CD39, and BTLA genes [9]. We also attempted to estimate the STAT3-based gene scores by using STAT3, ELK3, DTX3L, AIDA, NEDD9, KLHDC8A, TWSG1, NAGA, MYO1C, SH3PXD2B, SLC35F5, HERC5, C5orf15, ZMYM6, TPM4, DAP, SNAP23, RHOJ, HMG20B, ZCCHC9, NAMPT, SLC43A3, BIRC2, BACE2, ITGB1, ITFG3, AGXT2L2, GNG12, PALLD, IGF2BP2, NUP37, CTNNA1, GMPPA, BRCA1, TMEM51, RPN2, FZD1, PTPN12, SHQ1, and NAA38 genes [58]. We used GENECODIS 4.0 beta for gene set enrichment analysis for the prognostic genes derived from M0 macrophages [8].

Survival analysis and statistical testing

Kaplan–Meier survival curves were computed by a R based “survival::survfit” function and visualized by “survminer::ggsurvplot” function. Multivariate Cox-regression was performed by “survival::coxph" function. To examine the statistical significance between the groups we either used chi-Square test for 2X2 matrix or Fisher’s exact test for other combinations. The evaluated significance levels were denoted by *** when p values < 0.0001; by ** when p values < 0.001; or by * when p values < 0.05. The overall workflow followed in this article is described in Additional file 1: Figure S1.

Single cell RNA-sequencing and data processing

Fresh tumor tissues were processed using Brain Tumor Dissociation Kit (P) and gentleMACS Dissociator (Miltenyi Biotec) following manufacturer’s protocol. Erythrocytes were removed using a density gradient separation medium (Lympholyte-H; Cedarlane Laboratories) and resulting cell suspensions were cleared from tissue debris using Debris Removal Solution (Miltenyi Biotec). Cells were counted and assessed with Luna-FL Automated Cell Counter (Logos Biosystems) and run on 10 × Chromium Next GEM Chip with target recovery of 6,000 cells per lane. Library generation, quality control, and sequencing were all performed following 10 × Genomics Single Cell 3’ Reagent Kit User Guide.

Read demultiplexing and alignment to the GRCh38 human reference genome was performed using the Cell Ranger Single Cell Software v2.0 (10× Genomics) with the cellranger’s mkfastq and count functions, respectively. Raw count matrices were filtered for the minimum number of genes detected per cell (greater than 200) and the percentage of mitochondrial unique molecular identifier (UMI) counts (less than 20%). Filtered barcodes were then merged and clustered using the Seurat v3 package in R as per the developers’ vignettes. Cell annotation and malignant cells identification were performed using scmap and CONICsmat. These annotated cells were used to generate signature matrix using CIBERSORTx with its default parameters and S-mode batch-correction.

Results

M2 macrophage abundance varies across glioma subtypes and correlates with patient survival

We first examined the composition of infiltrating immune cells in gliomas, we utilized CIBERSORTx [45] on 3 diffuse glioma datasets (TCGA, CGGA325 and CGGA693), using the validated LM22 signature matrix (Additional file 13: Table S1 and Additional file 2: Figure S2A). Similar to previous reports, we also found M2 macrophages were the most abundant cell type across all IDH mutant (IDH-MUT) and wildtype (IDH-WT) tumors [61, 65]. We also observed that specific immune cell distributions were associated with IDH status. For example, M0 and M1 macrophages, neutrophils and T-helper cells were more abundant in IDH-WT gliomas, while monocytes, resting dendritic cells, activated mast cells, and activated natural killer cells were significantly enriched in IDH-MUT tumors (Fig. 1A). However, within IDH-MUT tumors, monocytes were more represented in the non-codeleted tumors, while resting dendritic cells and T-helper cells were abundant in the1p/19q codeleted tumors (Fig. 1B). These findings are in line with previous studies [34, 43, 49, 69]. Next, we evaluated the association of tumor-infiltrating immune cell fractions with survival outcomes. Survival plots revealed that high proportion of M2 macrophages can predict poor prognosis across all datasets in both IDH wildtype (except TCGA) and mutant tumors (Fig. 1C). Previous reports have highlighted an association between immune infiltrate and DNA damage [61]. In line, we also found positive correlations for M2 macrophage abundance with various genomic alterations like number of segments, fraction of genome altered, and homologous recombination defects in TCGA datasets, where ‘‘fraction altered’’ represents the fraction of bases deviating from baseline ploidy (defined as above 0.1 or below -0.1 in log2 relative copy number (CN) space), while ‘‘number of segments’’ represents total number of segments in each sample’s copy number profile (Additional file 2: Fig. 2B). These results agree with prior studies and suggest altered tumor microenvironments in diffuse gliomas based on IDH-mutation status, as well as correlations of the microenvironment with genomic alterations in the neoplastic cells.

Fig. 1
figure 1

An overview of the infiltrated immune cells in glioma and its subtypes. A Differentially represented cell types between IDH-MUT (green) and IDH-WT (orange) tumors are depicted by using boxplots. Similarly, B boxplots representing relative abundance of the cell types between Chr1p/19q codeleted (cyan) and non-codeleted (dark green) samples in IDH-MUT tumors. C Forest plots displaying prognostic association of M2 macrophages in both IDH-WT and IDH-MUT tumors. The abbreviations used for each cell type is described in Additional file 13: Table_S1

Immune cell expression clusters are associated with IDH status, tumor grade, and patient survival

In addition to estimating the proportions of individual cell types, CIBERSORTx also estimates the gene expression profile of each cell type in each sample. Figure 2A shows the transcriptomic clustering of 22 cell types across all samples independently in each of the three cohorts, Interestingly, gene expression profiles 15 of the 22 immune cell types revealed a clear separation between IDH mutant and wildtype samples in all 3 datasets (Fig. 2A). We also analyzed expression profiles of the genes associated with macrophages and T cells from IDH-MUT and IDH-WT tumors. This set of analysis highlighted the expression-based cell-type specificity between the two subtypes. For example, in IDH-WT tumors, genes expressed in M0 macrophages were significantly enriched for mesenchymal pathways while in IDH mutants it was immunological pathways (Fig. 2B).

Fig. 2
figure 2

Gene expression clustering and the separation between IDH-MUT and IDH-WT tumors. A UMAPs denoting 15 independent gene expression clusters identified in all 3 datasets harboring a clear separation between IDH-MUT (green) and IDH-WT (orange) tumors. B Most significant pathways and related genes derived from selected cell types demonstrates comparative transcriptional profiles between IDH-WT and IDH-MUT tumors

Focusing further on M2 macrophages, we observed fine clusters that distinctly segregate by tumor IDH status across all three datasets (Additional file 2: Fig. S2C and D) and exhibit distinct survival characteristics (Additional file 2: Fig. S2E). Importantly Cox multivariate analyses, after adjusting for relevant covariates (IDH mutation status, tumor grade and patient age) indicated that these unsupervised cluster groups from deconvolved M2 macrophage gene expression data remained independent predictors of patient outcome (Additional file 2: Fig. S2F). Additional analysis using LM10 signatures, (where monocyte includes M0, M1 and M2 macrophages into a single unit while T-CD4 includes T-naïve, T-mem-activ, T-mem-rest, T-gam-del, T-help and T-reg cells) also showed distinct cluster representing each cell type. The majority of these cell types again showed a clear separation between IDH-MUT (green) and IDH-WT (orange) tumors (Additional file 3: Fig. S3A). For example, monocytes showed clean clusters that were distinctly segregated by their IDH status and distinct survival supports our results even with fewer cell types (Additional file 3: Fig. S3B, C and D).

For comparison, we next examined M2 macrophage expression-based clusters in non-glioma TCGA tumors but did not find such significant survival associations in most tumor types (Additional file 4: Figure S4), suggesting that the prognostic association of macrophage-based gene expression groups may be relatively specific to gliomas.

Tumor groups with distinct immune signature patterns reveal distinct clinical features

Based on the above analyses, multiple expression-based clusters were identified for each immune cell type. For instance, there were 3 clusters for M2 macrophages, 2 clusters for regulatory T cells, 4 clusters for monocytes, and so on. The unsupervised (k-means) clustering of each immune cell type-based clusters were next examined together to identify identical groups of immune-based clusters in glioma. Towards that end, we used cluster-of-cluster analysis (COCA), an unsupervised method to integrate multi-omics data into biologically relevant sub-classes [26]. By applying COCA to these immune based clusters for all 22 immune cell types, we found distinct clusters of glioma samples significantly grouped by their IDH status, as shown in Fig. 3A and B, supporting immune cell signatures underlie the major glioma clinical subtypes.

Fig. 3
figure 3

Cluster of cluster analysis of immune based clusters. A Heatmaps representing the hierarchical clustering of LM22 based clusters broadly differentiating IDH-WT tumors from IDH-MUT. B Histograms showing quantitative distribution IDH specific tumors across each cluster

We then separated the samples according to IDH mutation status (e.g., IDH-WT and IDH-MUT samples were analyzed separately) and performed COCA on these tumor groups in all 3 datasets. Among IDH-WT tumors (Additional file 5: Figure S5A), tumor clusters derived from this unsupervised analysis showed significant survival associations across all 3 datasets. Among IDH-MUT tumors (Additional file 5: Figure S5B), 2 of 3 datasets (TCGA and CGGA325) showed significant survival associations among these immune cell-based clusters, while the 3rd dataset (CGGA693) showed a statistical trend. These cluster memberships are listed in Additional file 14: Table S2.

The data represented in Additional file 5: Figure S5, included all 22 cell types and to further characterize important immune cell types, next we focused on those individual cell types whose gene-expression-based cluster groups were prognostic. Within IDH-WT tumors, we tested the gene expression cluster groups from each cell type to identify which ones showed a significant association with patient prognosis (Additional file 6: Figure S6). These included expression clusters derived from M0 macrophages, M2 macrophages, B-mem and DC-active cells were consistently prognostic across all 3 datasets in IDH-WT gliomas (Additional file 6: Figure S6A). Among IDH-MUT tumors, M0 macrophages, T-help and B-mem cell derived cluster groups were found to be consistently prognostic (Additional file 6: Figure S6B). We then applied these respective cluster groups to COCA analyses (Fig. 4) and found that samples defined by these COCA cluster groups remained prognostic in multivariate Cox analyses even after adjusting for clinical variables (patient age and tumor grade) across all 3 datasets for both IDH-WT (Figs. 4A and B) and IDH-MUT tumors (Figs. 4C and D). This results further point to the potential clinical relevance of immune cell biology in predicting patient outcome.

Fig. 4
figure 4

Tumor groups with distinct immune signatures based on selected cell types. A Heatmaps representing hierarchical clustering of the clusters identified from selected immune cell types from IDH-WT tumors. The above annotation bars representing the distributions of cluster assignments, tumor grade, MGMT promoter methylation, EGFR expression, TERT expression (surrogating the status of TERT promoter mutation) and Chr7 gain & Chr10 loss followed by their Kaplan–Meier curves denoting their survival differences and B forest plots displaying prognostic association of these immune based clusters for IDH-WT. Similarly, C heatmaps representing the hierarchical clustering of the clusters identified from selected immune cell types from IDH-MUT tumors. The above annotation bars representing the distributions of cluster assignments, tumor grade and Chr1p/19q loss with Kaplan–Meier curves below the heatmap denoting their survival differences and D forest plots displaying prognostic association of these immune based clusters for IDH-MUT tumors

We then tested for macrophage abundance among these immune based clusters. Based on previous reports, macrophages and microglia can be distinguished by ITGA4 expression, where elevated expression in ITGA4 represents macrophage marker [6, 34]. We found increased levels of ITGA4 in poor-surviving groups as compared to cluster groups with better survival across all 3 datasets for IDH-WT tumors and in 2/3 datasets for IDH-MUT tumors (Additional file 7: Figure S7A), suggesting a link between macrophage abundance, immune-based clusters, and the patient outcome in both IDH-WT and IDH mut gliomas.

We next explored cytolytic and T cell exhaustion scores and their survival associations with the identified immune cell-based clusters. With this set of analyses, we observed that higher levels of cytolytic and T cell exhaustion scores among the poor surviving groups (Additional file 7: Figure S7B and C) and are in line with the previous reports [22, 71], further highlighting the biological or clinical significance of tumor microenvironments and their prognosis in gliomas.

Malignant cell states interact with immune microenvironment signature

To further characterize the role of the TME and their potential interactions with the malignant cell states, we estimated the composition of tumor cell states in all the 3 glioma datasets. To better understand the composition of their malignant cell states we leveraged previously published single-cell resources along with their cell-annotations to derive the signature matrix. As previously established, the diverse malignant cells in IDH-WT tumors have been described as converging to a 4 cell-states: (AClike, MES-like, NPClike and OPClike) [44] while in case of IDH-MUT gliomas, 3 malignant cell-states have been proposed: Astro-like, Oligo-like and Stem-like [63]. Malignant cell-state proportions derived with these signatures were generated using CIBERSORTx and are depicted in Additional file 8: Figure S8A for IDH-WT and Additional file 8: Figure S8B for IDH-MUT tumors, respectively. With the malignant cell state deconvolution, we examined the interplay between these malignant cell-states and immune cell-based (COCA) survival groups. Towards this end, we tested for the distribution of each of the cell states among each of the immune derived survival clusters and found in the case of IDH-WT tumors, the MES-like component was higher in poor surviving groups (as defined in Fig. 4) for all three datasets (Fig. 5A). To further characterize this finding, we built on results from several recent reports which have underlined the role of macrophages to induce MES-like cell-state by implicating macrophage-derived oncostatin M (OSM) with its receptors and activating STAT3 signaling in glioma [24, 47, 66]. Along this line, we found a positive correlation for the proportion of MES-like state and estimated OSMR expression within M0 macrophages (Fig. 5B), but not within M1 macrophages, M2 macrophages or monocytes (Additional file 9: Figure S9), suggesting specificity within macrophage subtypes in a potential interaction with the MES-like malignant cell state. Moreover, the enrichment of an epithelial-mesenchymal transition signature as well as STAT3 signaling markers in poor-surviving groups further characterize these results (Fig. 5D and E), suggesting a role of M0 macrophage abundance in the MES-like cell state in glioma. We also included in-house scRNA RNA seq data from 19 gliomas encompassing 10 IDH-WT and 9 IDH-MUT tumors (Additional file 15: Table S3). With this set of analysis, we were able to identify 18 distinct clusters (Additional file 10: Fig. S10A and B), including the neoplastic and the immune cells (Additional file 10: Fig. S10C, D and Additional file 15: Table S3). These single cells derived cell types were then used to generate the signature matrix for deconvolution (Additional file 10: Fig. S10E and Additional file 15: Table S3). We also assessed the expression deconvolution data from macrophages and malignant cells and found a statistically significant correlation between the expressions of OSM and OSMR expressions from the macrophages and malignant cells respectively (Fig. 5C). This set of analysis support the connection between macrophages and MES-like cells.

Fig. 5
figure 5

Deconvolving malignant cell states and their interaction with immune based clusters. A Differential distributions of IDH-WT specific malignant cell-states between the two immune-based survival groups are depicted by boxplots. B Scatter plots representing a positive correlation between the proportions of M0 macrophages and MES-like component of IDH-WT tumors across all 3 datasets. C Scatter plots representing a significant positive correlation between the expressions of OSM from macrophages and OSMR from malignant cells in IDH-WT tumors identified by scRNA based deconvolution. D Boxplots depicting differential distribution of epithelial-mesenchymal transition markers between the two immune-based survival groups. Similarly, E boxplots depicting differential distribution of STAT3 signaling markers between the two immune-based survival groups

Additionally, we identified 368 MES-like corelating genes from M0 macrophages that were common in all 3 datasets that, in turn, showed enrichment for mesenchymal pathways, including extracellular matrix organization and angiogenesis along with the immunological pathways (Additional file 11: Figure S11A). In contrast, 95 genes from M2 macrophages which were common in all 3 datasets showed enrichment only for immunological pathways (Additional file 11: Figure S11B). We next analyzed IDH-MUT tumors and interestingly in 1p/19q non-codeleted tumors (compared to co-deleted tumors) we observed a similar enrichment for genes from mesenchymal pathways in M0 macrophages with overlapping correlated genes compared to the M0 gene expression signatures in the MES-like state of IDH-WT glioma (Additional file 12: Figure S12).

Extension of findings to checkpoint inhibitor-treated glioma patients

Next, we investigated infiltrating immune cell compositions in a predominantly IDH-WT 29-sample dataset collected from patients treated with pembrolizumab [10]. To explore the role of tumor-microenvironment in response to immunotherapy, we relied on the LM10 signatures comprising 10 broad immune cell types better suited for smaller datasets than the LM22. Notably, LM10 broadly covers all types of monocyte-derived cells (including M0, M1 and M2 macrophages) in a single unit. As expected from LM10-based deconvolution, this dataset demonstrated enrichment of the monocyte cell fraction across all the patients analyzed (Fig. 6A) and the predicted malignant cell state proportions are shown in Fig. 6B. Interestingly a higher fraction of monocytes was significantly associated with a shorter overall survival (HR = 5.0, P < 0.001) (Fig. 6C), suggesting that cells with monocytic lineage may play a role on checkpoint inhibitor-related patient outcome in glioma. We also repeated the LM10 based deconvolution for TCGA and CGGA datasets and observed higher fraction of monocytes was significantly associated with patient survival only in case of CGGA693 but not in CGGA325 and TCGA tumors (Fig. 6D). Such an inconsistency among the non-immunotherapy treated samples compared to the checkpoint inhibitor-treated tumors suggest the predictive value of monocytes.

Fig. 6
figure 6

Immunotherapy-treated glioma datasets. A Stacked bar plots representing relative proportion of 10 broader category of immune cell types from LM10 across samples treated with checkpoint inhibitor. B Stacked bar plots representing the relative proportions of 4 malignant cell-states across all immunotherapy treated samples. C Forest plot representing the association between their dichotomous immune cell fractions and their overall survival. D Monocyte-based gene expression clusters represented by tSNE followed by E Kaplan–Meier curves denoting distinct survival between the two clusters. F Boxplots depicting differential distribution of epithelial-mesenchymal transition markers between the two immune-based survival groups derived from monocytes from samples treated with checkpoint inhibitor. G Scatter plots representing significant positive correlations between the epithelial-mesenchymal transition and the T cell exhaustion scores both in patients treated with a checkpoint inhibitor, as well as larger cohorts. H Forest plot representing the insignificant or inconsistent association between their dichotomous monocytic fractions in non-immunotherapy or non-checkpoint inhibitor treated IDH-WT tumors (from TCGA, CGGA325 and CGGA693) and their overall survival

Unsupervised clustering of the checkpoint inhibitor-samples based on monocytic expression-based clusters (as above) revealed two clusters with distinct survival outcome (Fig. 6E and F). We also found epithelial-mesenchymal transition (EMT) markers to be significantly enriched in the poor surviving group (Fig. 6G). We also observed a positive correlation between the epithelial-mesenchymal transition and the T cell exhaustion scores both in patients treated with a checkpoint inhibitor, as well as larger cohorts (Fig. 6H). The finding of an association of epithelial-mesenchymal transition with T-cell exhaustion in glioma may be analogous to what has been reported in other cancer types [41]. These results together emphasize prognostic associations of monocyte/macrophage lineage cells in gliomas in the context of patients treated with a checkpoint inhibitor, and link mesenchymal gene expression in monocytic lineage cells with a T-cell exhaustion score in gliomas overall.

Discussion

Emerging evidence has highlighted the critical involvement of the immune microenvironment in tumor development and progression. These immune-related biomarkers may have a potential to predict patient outcome and therapy responsiveness [2, 31, 54, 73], and thus, have inspired numerous microenvironment-targeted therapies in multiple tumors [33, 34, 50]. In addition, various reports have discussed prognostic roles of immune related genes in gliomas, but the characterization of immune cell compartments has not been fully elucidated.

In this work, we explored the landscape of glioma infiltrating immune cells using CIBERSORTx deconvolution as summarized in Fig. 1. We found M2 macrophages constitute a large portion of glioma microenvironment. Macrophages are a large component of many tumor types and have been an attractive target for glioblastoma therapy [21, 36, 40, 55]. Deconvolution of the tumor bulk transcriptomic data enabled us to investigate immune cell distributions in glioma and their clinical outcome correlated at a much greater resolution. We found that specifically M0 macrophages were enriched in IDH-WT gliomas while monocytes were relatively more abundant in IDH-MUT gliomas. Additionally, higher levels of M2 macrophages were observed to be associated with poor prognosis. These differentially distributed monocyte/macrophage lineage cells and their association with patient outcome further highlight their prognostic roles in tumor progression. Targeting these specific cells to inhibit their tumor-promoting effect and reprogramming them into an anti-tumor phenotype could be a potential therapeutic approach for glioma.

Gene expression profiling of all the deconvolved immune cell types resulted into sample groups segregated by their IDH status, as discussed in Figs. 2 and 3. Further clustering of these sample groups resulted into distinct clusters which can predict patient outcome. These immune based clusters were broadly characterized by tumor grades, and in case of IDH-MUT tumors, also by their 1p/19q codeletion status. Based on WHO 2021 classification and other previous reports has already established that Chr7 gain & Chr10 loss, TERT mutation along with EGFR expression and MGMT promoter methylation are the strongest predictor of survival in glioblastomas [25, 35, 39]. Several reports suggest elevated levels of ITGA4 denoting macrophage abundance can also predict poor prognosis in gliomas [6, 34, 56]. We also observed similar associations between the Chr7 gain & Chr10 loss, expressions of EGFR, TERT (surrogating the status of TERT promoter mutation) and MGMT promoter methylation of the immune based clusters with the patient outcome (Fig. 4A and Additional file 5: Figure S5A). Similarly, the prognostic nature of these immune based clusters was also characterized by higher levels of ITGA4 expression as well as higher scores for cytolytic and T cell exhaustion, as reported previously [22, 51, 71] (Additional file 6: Figure S6), providing corroboration that CIBERSORTx-based deconvolution of the bulk datasets yields result that are consistent with prior investigations. These findings indicate that tumor microenvironment plays a major role in characterizing tumor subtypes [29, 53, 59]. These identified sample groups may exhibit a distinct immune context, prognosis, and immunotherapy benefit, which supports the idea that the immune environment is of vital importance in predicting patient prognosis and evaluating the response rate of checkpoint inhibitor immunotherapies. Later, clusters derived from a few selected cell types that involve M0 and M2 macrophages along with B-mem, DC-activ and T-help, were able to identify tumor subgroups that can also predict patient outcome. This set of analysis again emphasizes the prognostic role of M0 and M2 macrophages along with other immune cells. Further analyses may facilitate more opportunities to investigate other immune cells playing important roles in complex tissues.

Recent advances in tumor biology have revealed interactions between the tumor cells and their adjacent microenvironment, highlighting the importance of TME and their involvement in tumor growth and the development of metastasis [52, 62]. Towards this we extended our search by analyzing the proportion of malignant cell-states by using a previously published single cell data and their labels to derive a new signature matrix [44, 63]. To study the interplay of malignant cell-states and immune cell-based clusters, we checked for the distribution of malignant-states among the immune based clusters. With this set of analysis, we found in IDH-WT, MES-like component was found to be higher in poor surviving groups while in IDH-MUT, Oligo-like tumor cell-state was found to be higher in better surviving group (Additional file 7: Figure S7C). Based on recent literature, proportion of TCGA-MES subtype can be correlated with macrophage abundance [24, 66]. On the similar note, in IDH-WT we also found positive correlation between the M0 macrophages and MES-like proportions (r = 0.26, 0.45, and 0.38 in TCGA, CGGA325, and CGGA693 respectively in Fig. 5B) followed by enrichment for genes from mesenchymal pathways along with the immunological pathways in M0 macrophages while M2 macrophages demonstrated the enrichment only for immunological pathways. Together these results highlight the specificity within macrophage subtypes in context to the MES-like cell state in gliomas, as described in Additional file 9: Figure S9 and summarized in Fig. 7A. Additionally, IDH mut tumor also harbor a similar enrichment for genes from mesenchymal pathways in M0 macrophages over other cell types when comparing 1p/19q non-codeleted tumors with codeleted tumors (in Fig. 7B and Additional file 10: Figure S10). A deeper understanding of interactions like this will help not only to explore the cancer growth and metastasis but also to elucidate the mechanisms of action of classical drugs that have been discovered by empirical approaches. Such studies may highlight the importance of deconvolution approaches that can facilitate designing and development of novel cancer drugs.

Fig. 7
figure 7

Enrichment for genes from mesenchymal pathways in M0 macrophages while M2 macrophages demonstrated the enrichment for immunological pathways in A IDH-WT and B IDH-MUT tumors

We expanded our analysis our findings in an independent checkpoint inhibitor-treated glioma samples encompassing 29 high-grade glioma patients treated with pembrolizumab [10]. Here due to the small sample size we used LM10 signature matrix instead of LM22 for deconvolution, where monocytes are expected to encompass all types of macrophages along with monocytes itself. Results showed monocytes were the most predominant cell type and can predict poor prognosis in this glioma samples treated with checkpoint inhibitor (Fig. 6). Further clustering of the deconvolved gene expression from monocytes resulted into poor and better surviving groups where epithelial-mesenchymal transition markers were found to be enriched in the poor surviving group. As monocytes and macrophages share their lineages [17], these findings are consistent with those observed from TCGA and CGGA datasets. These results underline the prognostic potential of monocyte/macrophage lineage cells in gliomas, including patients treated with a checkpoint inhibitor.

Strengths of our study includes the use of clinically annotated genomically characterized glioma data, with the inclusion of 3 independent datasets, allowing the opportunity to examine reproducibility across cohorts. Our study has several limitations, which includes the use of the LM22 gene signature matrix, which was derived primary from cells in the blood and not from tumors. It is likely that variation in immune gene expression signatures exist within solid tumors. In addition, the distinction of M0, M1 and M2 macrophages may not best delineate the function of these cells in gliomas. That said, a recent report indicates that the M0 (unpolarized) state may best represent the majority of macrophages in gliomas [16]. Thus, while the M0/M1/M2 signatures in the LM22-derived matrix may not best represent glioma biology, our result suggests, at minimum that a subset of macrophages exist that may play specific roles in glioma biology, including relationships with the MES-like state in IDH-WT GBM, as well as associations with 1p/19q co-deletion status in IDH-MUT glioma. An explanatory model summarizes these main interactions (Fig. 7). Based on these findings, future analyses will examine signature matrices derived from glioma single cell RNA sequencing data to further characterize the interplay of macrophages and tumor cell characteristics. Even with these limitations, our findings suggest proof-of-principal to stimulate further and more detailed studies on the role of immune cell gene expression and biology in the clinical and biological aspects of both IDH-WT and IDH-MUT diffuse gliomas.

Conclusion

In conclusion, our analysis relies on a large sample set of gliomas to demonstrate immune cell gene expression signatures, derived by deconvolution with CIBERSORTx correlate with IDH mutation status in glioma and also with patient outcome within IDH-WT and IDH-MUT gliomas. The mesenchymal cell state in IDH-WT GBM showed association with M0 and M2 macrophages where it correlates with mesenchymal signatures in M0 macrophages, and immune signatures in M2 macrophages. Together these results highlighted the prominent association of monocytic lineage cells, specially the M0 macrophages, with MES-like state and the patient outcome thus provide insights for future investigation to better understand glioma biology and developing better immunotherapeutic approaches in gliomas.

References

  1. Asai M, Kato Y, Kawai S, Watanabe K, Yomota M, Okuma Y, Hosomi Y, Hishima T, Okamura T (2019) Management of cardiac tamponade during nivolumab of lung cancer with intrapericardial bleomycin: case report. Immunotherapy 11:467–472. https://doi.org/10.2217/imt-2019-0003

    CAS  Article  PubMed  Google Scholar 

  2. Ascierto ML, Kmieciak M, Idowu MO, Manjili R, Zhao Y, Grimes M, Dumur C, Wang E, Ramakrishnan V, Wang XY et al (2012) A signature of immune function genes associated with recurrence-free survival in breast cancer patients. Breast Cancer Res Treat 131:871–880. https://doi.org/10.1007/s10549-011-1470-x

    CAS  Article  PubMed  Google Scholar 

  3. Bao ZS, Chen HM, Yang MY, Zhang CB, Yu K, Ye WL, Hu BQ, Yan W, Zhang W, Akers J et al (2014) RNA-seq of 272 gliomas revealed a novel, recurrent PTPRZ1-MET fusion transcript in secondary glioblastomas. Genome Res 24:1765–1773. https://doi.org/10.1101/gr.165126.113

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  4. Batchu S (2021) Immunological landscape of Neuroblastoma and its clinical significance. Cancer Treat Res Commun 26:100274. https://doi.org/10.1016/j.ctarc.2020.100274

    Article  PubMed  Google Scholar 

  5. Bersanelli M, Buti S (2017) From targeting the tumor to targeting the immune system: transversal challenges in oncology with the inhibition of the PD-1/PD-L1 axis. World J Clin Oncol 8:37–53. https://doi.org/10.5306/wjco.v8.i1.37

    Article  PubMed  PubMed Central  Google Scholar 

  6. Bowman RL, Klemm F, Akkari L, Pyonteck SM, Sevenich L, Quail DF, Dhara S, Simpson K, Gardner EE, Iacobuzio-Donahue CA et al (2016) Macrophage ontogeny underlies differences in tumor-specific education in brain malignancies. Cell Rep 17:2445–2459. https://doi.org/10.1016/j.celrep.2016.10.052

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  7. Cancer Genome Atlas Research N, Brat DJ, Verhaak RG, Aldape KD, Yung WK, Salama SR, Cooper LA, Rheinbay E, Miller CR, Vitucci M et al (2015) Comprehensive, integrative genomic analysis of diffuse lower-grade gliomas. N Engl J Med 372:2481–2498. https://doi.org/10.1056/NEJMoa1402121

    CAS  Article  Google Scholar 

  8. Carmona-Saez P, Chagoyen M, Tirado F, Carazo JM, Pascual-Montano A (2007) GENECODIS: a web-based tool for finding significant concurrent annotations in gene lists. Genome Biol 8:R3. https://doi.org/10.1186/gb-2007-8-1-r3

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  9. Catakovic K, Klieser E, Neureiter D, Geisberger R (2017) T cell exhaustion: from pathophysiological basics to tumor immunotherapy. Cell Commun Signal 15:1. https://doi.org/10.1186/s12964-016-0160-z

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  10. Cloughesy TF, Mochizuki AY, Orpilla JR, Hugo W, Lee AH, Davidson TB, Wang AC, Ellingson BM, Rytlewski JA, Sanders CM et al (2019) Neoadjuvant anti-PD-1 immunotherapy promotes a survival benefit with intratumoral and systemic immune responses in recurrent glioblastoma. Nat Med 25:477–486. https://doi.org/10.1038/s41591-018-0337-7

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  11. Costa LB, Queiroz MA, Barbosa FG, Nunes RF, Marin JFG, Dzik C, Buchpiguel CA (2019) Pseudoprogression on PSMA PET imaging of a mCRPC patient under anti-PD1 treatment. Eur J Nucl Med Mol Imaging 46:1576–1577. https://doi.org/10.1007/s00259-019-04328-1

    Article  PubMed  Google Scholar 

  12. Devaud C, John LB, Westwood JA, Darcy PK, Kershaw MH (2013) Immune modulation of the tumor microenvironment for enhancing cancer immunotherapy. Oncoimmunology 2:e25961. https://doi.org/10.4161/onci.25961

    Article  PubMed  PubMed Central  Google Scholar 

  13. Doucette T, Rao G, Rao A, Shen L, Aldape K, Wei J, Dziurzynski K, Gilbert M, Heimberger AB (2013) Immune heterogeneity of glioblastoma subtypes: extrapolation from the cancer genome atlas. Cancer Immunol Res 1:112–122. https://doi.org/10.1158/2326-6066.CIR-13-0028

    CAS  Article  PubMed  Google Scholar 

  14. Fecci PE, Heimberger AB, Sampson JH (2014) Immunotherapy for primary brain tumors: no longer a matter of privilege. Clin Cancer Res 20:5620–5629. https://doi.org/10.1158/1078-0432.CCR-14-0832

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  15. Foroutan M, Bhuva DD, Lyu R, Horan K, Cursons J, Davis MJ (2018) Single sample scoring of molecular phenotypes. BMC Bioinformatics 19:404. https://doi.org/10.1186/s12859-018-2435-4

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  16. Gabrusiewicz K, Rodriguez B, Wei J, Hashimoto Y, Healy LM, Maiti SN, Thomas G, Zhou S, Wang Q, Elakkad A et al (2016) Glioblastoma-infiltrated innate immune cells resemble M0 macrophage phenotype. JCI Insight. https://doi.org/10.1172/jci.insight.85841

    Article  PubMed  PubMed Central  Google Scholar 

  17. Geissmann F, Manz MG, Jung S, Sieweke MH, Merad M, Ley K (2010) Development of monocytes, macrophages, and dendritic cells. Science 327:656–661. https://doi.org/10.1126/science.1178331

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  18. Glajcar A, Szpor J, Hodorowicz-Zaniewska D, Tyrak KE, Okon K (2019) The composition of T cell infiltrates varies in primary invasive breast cancer of different molecular subtypes as well as according to tumor size and nodal status. Virchows Arch 475:13–23. https://doi.org/10.1007/s00428-019-02568-y

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  19. Gnjatic S, Bronte V, Brunet LR, Butler MO, Disis ML, Galon J, Hakansson LG, Hanks BA, Karanikas V, Khleif SN et al (2017) Identifying baseline immune-related biomarkers to predict clinical outcome of immunotherapy. J Immunother Cancer 5:44. https://doi.org/10.1186/s40425-017-0243-4

    Article  PubMed  PubMed Central  Google Scholar 

  20. Gu Z, Eils R, Schlesner M (2016) Complex heatmaps reveal patterns and correlations in multidimensional genomic data. Bioinformatics 32:2847–2849. https://doi.org/10.1093/bioinformatics/btw313

    CAS  Article  PubMed  Google Scholar 

  21. Gutmann DH, McLellan MD, Hussain I, Wallis JW, Fulton LL, Fulton RS, Magrini V, Demeter R, Wylie T, Kandoth C et al (2013) Somatic neurofibromatosis type 1 (NF1) inactivation characterizes NF1-associated pilocytic astrocytoma. Genome Res 23:431–439. https://doi.org/10.1101/gr.142604.112

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  22. Haddad AF, Chen JS, Oh T, Pereira MP, Joshi RS, Aghi MK (2020) Higher cytolytic score correlates with an immunosuppressive tumor microenvironment and reduced survival in glioblastoma. Sci Rep 10:17580. https://doi.org/10.1038/s41598-020-73793-8

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  23. Hao Y, Hao S, Andersen-Nissen E, Mauck WM 3rd, Zheng S, Butler A, Lee MJ, Wilk AJ, Darby C, Zager M et al (2021) Integrated analysis of multimodal single-cell data. Cell 184(3573–3587):e3529. https://doi.org/10.1016/j.cell.2021.04.048

    CAS  Article  Google Scholar 

  24. Hara T, Chanoch-Myers R, Mathewson ND, Myskiw C, Atta L, Bussema L, Eichhorn SW, Greenwald AC, Kinker GS, Rodman C et al (2021) Interactions between cancer cells and immune cells drive transitions to mesenchymal-like states in glioblastoma. Cancer Cell 39(779–792):e711. https://doi.org/10.1016/j.ccell.2021.05.002

    CAS  Article  Google Scholar 

  25. Hegi ME, Diserens AC, Gorlia T, Hamou MF, de Tribolet N, Weller M, Kros JM, Hainfellner JA, Mason W, Mariani L et al (2005) MGMT gene silencing and benefit from temozolomide in glioblastoma. N Engl J Med 352:997–1003. https://doi.org/10.1056/NEJMoa043331

    CAS  Article  PubMed  Google Scholar 

  26. Hoadley KA, Yau C, Wolf DM, Cherniack AD, Tamborero D, Ng S, Leiserson MDM, Niu B, McLellan MD, Uzunangelov V et al (2014) Multiplatform analysis of 12 cancer types reveals molecular classification within and across tissues of origin. Cell 158:929–944. https://doi.org/10.1016/j.cell.2014.06.049

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  27. Huang L, Chen H, Xu Y, Chen J, Liu Z, Xu Q (2020) Correlation of tumor-infiltrating immune cells of melanoma with overall survival by immunogenomic analysis. Cancer Med 9:8444–8456. https://doi.org/10.1002/cam4.3466

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  28. Jiang J, Pan W, Xu Y, Ni C, Xue D, Chen Z, Chen W, Huang J (2020) Tumour-infiltrating immune cell-based subtyping and signature gene analysis in breast cancer based on gene expression profiles. J Cancer 11:1568–1583. https://doi.org/10.7150/jca.37637

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  29. Job S, Rapoud D, Dos Santos A, Gonzalez P, Desterke C, Pascal G, Elarouci N, Ayadi M, Adam R, Azoulay D et al (2020) Identification of four immune subtypes characterized by distinct composition and functions of tumor microenvironment in intrahepatic cholangiocarcinoma. Hepatology 72:965–981. https://doi.org/10.1002/hep.31092

    CAS  Article  PubMed  Google Scholar 

  30. Karn T, Jiang T, Hatzis C, Sanger N, El-Balat A, Rody A, Holtrich U, Becker S, Bianchini G, Pusztai L (2017) Association between genomic metrics and immune infiltration in triple-negative breast cancer. JAMA Oncol 3:1707–1711. https://doi.org/10.1001/jamaoncol.2017.2140

    Article  PubMed  PubMed Central  Google Scholar 

  31. Kim K, Jeon S, Kim TM, Jung CK (2018) Immune gene signature delineates a subclass of papillary thyroid cancer with unfavorable clinical outcomes. Cancers (Basel). https://doi.org/10.3390/cancers10120494

    Article  PubMed Central  Google Scholar 

  32. Kim SC, Yu D, Cho SB (2018) COEX-Seq: convert a variety of measurements of gene expression in RNA-Seq. Genomics Inform 16:e36. https://doi.org/10.5808/GI.2018.16.4.e36

    Article  PubMed  PubMed Central  Google Scholar 

  33. Klemm F, Joyce JA (2015) Microenvironmental regulation of therapeutic response in cancer. Trends Cell Biol 25:198–213. https://doi.org/10.1016/j.tcb.2014.11.006

    Article  PubMed  Google Scholar 

  34. Klemm F, Maas RR, Bowman RL, Kornete M, Soukup K, Nassiri S, Brouland JP, Iacobuzio-Donahue CA, Brennan C, Tabar V et al (2020) Interrogation of the microenvironmental landscape in brain tumors reveals disease-specific alterations of immune cells. Cell 181(1643–1660):e1617. https://doi.org/10.1016/j.cell.2020.05.007

    CAS  Article  Google Scholar 

  35. Li J, Liang R, Song C, Xiang Y, Liu Y (2018) Prognostic significance of epidermal growth factor receptor expression in glioma patients. Onco Targets Ther 11:731–742. https://doi.org/10.2147/OTT.S155160

    Article  PubMed  PubMed Central  Google Scholar 

  36. Liang P, Chai Y, Zhao H, Wang G (2020) Predictive analyses of prognostic-related immune genes and immune infiltrates for glioblastoma. Diagnostics (Basel). https://doi.org/10.3390/diagnostics10030177

    Article  PubMed Central  Google Scholar 

  37. Liberzon A, Birger C, Thorvaldsdottir H, Ghandi M, Mesirov JP, Tamayo P (2015) The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst 1:417–425. https://doi.org/10.1016/j.cels.2015.12.004

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  38. Liu X, Li Y, Qian Z, Sun Z, Xu K, Wang K, Liu S, Fan X, Li S, Zhang Z et al (2018) A radiomic signature as a non-invasive predictor of progression-free survival in patients with lower-grade gliomas. Neuroimage Clin 20:1070–1077. https://doi.org/10.1016/j.nicl.2018.10.014

    Article  PubMed  PubMed Central  Google Scholar 

  39. Louis DN, Perry A, Wesseling P, Brat DJ, Cree IA, Figarella-Branger D, Hawkins C, Ng HK, Pfister SM, Reifenberger G et al (2021) The 2021 WHO classification of tumors of the central nervous system: a summary. Neuro Oncol 23:1231–1251. https://doi.org/10.1093/neuonc/noab106

    CAS  Article  PubMed  Google Scholar 

  40. Lu J, Li H, Chen Z, Fan L, Feng S, Cai X, Wang H (2019) Identification of 3 subpopulations of tumor-infiltrating immune cells for malignant transformation of low-grade glioma. Cancer Cell Int 19:265. https://doi.org/10.1186/s12935-019-0972-1

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  41. Mak MP, Tong P, Diao L, Cardnell RJ, Gibbons DL, William WN, Skoulidis F, Parra ER, Rodriguez-Canales J, Wistuba II et al (2016) A patient-derived, pan-cancer EMT signature identifies global molecular alterations and immune target enrichment following epithelial-to-mesenchymal transition. Clin Cancer Res 22:609–620. https://doi.org/10.1158/1078-0432.CCR-15-0876

    CAS  Article  PubMed  Google Scholar 

  42. Masuhiro K, Shiroyama T, Nagatomo I, Kumanogoh A (2019) Unique case of pseudoprogression manifesting as lung cavitation after pembrolizumab treatment. J Thorac Oncol 14:e108–e109. https://doi.org/10.1016/j.jtho.2018.12.024

    Article  PubMed  Google Scholar 

  43. Najafi M, Hashemi Goradel N, Farhood B, Salehi E, Nashtaei MS, Khanlarkhani N, Khezri Z, Majidpoor J, Abouzaripour M, Habibi M et al (2019) Macrophage polarity in cancer: a review. J Cell Biochem 120:2756–2765. https://doi.org/10.1002/jcb.27646

    CAS  Article  PubMed  Google Scholar 

  44. Neftel C, Laffy J, Filbin MG, Hara T, Shore ME, Rahme GJ, Richman AR, Silverbush D, Shaw ML, Hebert CM et al (2019) An integrative model of cellular states, plasticity, and genetics for glioblastoma. Cell 178:835-849 e821. https://doi.org/10.1016/j.cell.2019.06.024

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  45. Newman AM, Steen CB, Liu CL, Gentles AJ, Chaudhuri AA, Scherer F, Khodadoust MS, Esfahani MS, Luca BA, Steiner D et al (2019) Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol 37:773–782. https://doi.org/10.1038/s41587-019-0114-2

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  46. Ostrom QT, Gittleman H, Truitt G, Boscia A, Kruchko C, Barnholtz-Sloan JS (2018) CBTRUS statistical report: primary brain and other central nervous system tumors diagnosed in the United States in 2011–2015. Neuro Oncol 20:iv1–iv86. https://doi.org/10.1093/neuonc/noy131

    Article  PubMed  PubMed Central  Google Scholar 

  47. Ozawa T, Riester M, Cheng YK, Huse JT, Squatrito M, Helmy K, Charles N, Michor F, Holland EC (2014) Most human non-GCIMP glioblastoma subtypes evolve from a common proneural-like precursor glioma. Cancer Cell 26:288–300. https://doi.org/10.1016/j.ccr.2014.06.005

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  48. Pages F, Galon J, Dieu-Nosjean MC, Tartour E, Sautes-Fridman C, Fridman WH (2010) Immune infiltration in human tumors: a prognostic factor that should not be ignored. Oncogene 29:1093–1102. https://doi.org/10.1038/onc.2009.416

    CAS  Article  PubMed  Google Scholar 

  49. Poon CC, Gordon PMK, Liu K, Yang R, Sarkar S, Mirzaei R, Ahmad ST, Hughes ML, Yong VW, Kelly JJP (2019) Differential microglia and macrophage profiles in human IDH-mutant and -wild type glioblastoma. Oncotarget 10:3129–3143. https://doi.org/10.18632/oncotarget.26863

    Article  PubMed  PubMed Central  Google Scholar 

  50. Quail DF, Joyce JA (2017) The microenvironmental landscape of brain tumors. Cancer Cell 31:326–341. https://doi.org/10.1016/j.ccell.2017.02.009

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  51. Rooney MS, Shukla SA, Wu CJ, Getz G, Hacohen N (2015) Molecular and genetic properties of tumors associated with local immune cytolytic activity. Cell 160:48–61. https://doi.org/10.1016/j.cell.2014.12.033

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  52. Seager RJ, Hajal C, Spill F, Kamm RD, Zaman MH (2017) Dynamic interplay between tumour, stroma and immune system can drive or prevent tumour progression. Converg Sci Phys Oncol. https://doi.org/10.1088/2057-1739/aa7e86

    Article  PubMed  PubMed Central  Google Scholar 

  53. Shen R, Li P, Li B, Zhang B, Feng L, Cheng S (2019) Identification of distinct immune subtypes in colorectal cancer based on the stromal compartment. Front Oncol 9:1497. https://doi.org/10.3389/fonc.2019.01497

    Article  PubMed  Google Scholar 

  54. Shen S, Wang G, Zhang R, Zhao Y, Yu H, Wei Y, Chen F (2019) Development and validation of an immune gene-set based Prognostic signature in ovarian cancer. EBioMedicine 40:318–326. https://doi.org/10.1016/j.ebiom.2018.12.054

    Article  PubMed  Google Scholar 

  55. Simmons GW, Pong WW, Emnett RJ, White CR, Gianino SM, Rodriguez FJ, Gutmann DH (2011) Neurofibromatosis-1 heterozygosity increases microglia in a spatially and temporally restricted pattern relevant to mouse optic glioma formation and growth. J Neuropathol Exp Neurol 70:51–62. https://doi.org/10.1097/NEN.0b013e3182032d37

    CAS  Article  PubMed  Google Scholar 

  56. Sorensen MD, Dahlrot RH, Boldt HB, Hansen S, Kristensen BW (2018) Tumour-associated microglia/macrophages predict poor prognosis in high-grade gliomas and correlate with an aggressive tumour subtype. Neuropathol Appl Neurobiol 44:185–206. https://doi.org/10.1111/nan.12428

    CAS  Article  PubMed  Google Scholar 

  57. Steen CB, Liu CL, Alizadeh AA, Newman AM (2020) Profiling cell type abundance and expression in bulk tissues with CIBERSORTx. Methods Mol Biol 2117:135–157. https://doi.org/10.1007/978-1-0716-0301-7_7

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  58. Tan MSY, Sandanaraj E, Chong YK, Lim SW, Koh LWH, Ng WH, Tan NS, Tan P, Ang BT, Tang C (2019) A STAT3-based gene signature stratifies glioma patients for targeted therapy. Nat Commun 10:3601. https://doi.org/10.1038/s41467-019-11614-x

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  59. Tang C, Ma J, Liu X, Liu Z (2020) Identification of four immune subtypes in bladder cancer based on immune gene sets. Front Oncol 10:544610. https://doi.org/10.3389/fonc.2020.544610

    Article  PubMed  PubMed Central  Google Scholar 

  60. Tang H, Qiao J, Fu YX (2016) Immunotherapy and tumor microenvironment. Cancer Lett 370:85–90. https://doi.org/10.1016/j.canlet.2015.10.009

    CAS  Article  PubMed  Google Scholar 

  61. Thorsson V, Gibbs DL, Brown SD, Wolf D, Bortone DS, Ou Yang TH, Porta-Pardo E, Gao GF, Plaisier CL, Eddy JA et al (2018) The immune landscape of cancer. Immunity 48:812-830 e814. https://doi.org/10.1016/j.immuni.2018.03.023

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  62. Ungefroren H, Sebens S, Seidl D, Lehnert H, Hass R (2011) Interaction of tumor cells with the microenvironment. Cell Commun Signal 9:18. https://doi.org/10.1186/1478-811X-9-18

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  63. Venteicher AS, Tirosh I, Hebert C, Yizhak K, Neftel C, Filbin MG, Hovestadt V, Escalante LE, Shaw ML, Rodman C et al (2017) Decoupling genetics, lineages, and microenvironment in IDH-mutant gliomas by single-cell RNA-seq. Science. https://doi.org/10.1126/science.aai8478

    Article  PubMed  PubMed Central  Google Scholar 

  64. Verhaak RG, Hoadley KA, Purdom E, Wang V, Qi Y, Wilkerson MD, Miller CR, Ding L, Golub T, Mesirov JP et al (2010) Integrated genomic analysis identifies clinically relevant subtypes of glioblastoma characterized by abnormalities in PDGFRA, IDH1, EGFR, and NF1. Cancer Cell 17:98–110. https://doi.org/10.1016/j.ccr.2009.12.020

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  65. Vidyarthi A, Agnihotri T, Khan N, Singh S, Tewari MK, Radotra BD, Chatterjee D, Agrewala JN (2019) Predominance of M2 macrophages in gliomas leads to the suppression of local and systemic immunity. Cancer Immunol Immunother 68:1995–2004. https://doi.org/10.1007/s00262-019-02423-8

    CAS  Article  PubMed  Google Scholar 

  66. Wang Q, Hu B, Hu X, Kim H, Squatrito M, Scarpace L, deCarvalho AC, Lyu S, Li P, Li Y et al (2017) Tumor evolution of glioma-intrinsic gene expression subtypes associates with immunological changes in the microenvironment. Cancer Cell 32:42-56 e46. https://doi.org/10.1016/j.ccell.2017.06.003

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  67. Wang Y, Qian T, You G, Peng X, Chen C, You Y, Yao K, Wu C, Ma J, Sha Z et al (2015) Localizing seizure-susceptible brain regions associated with low-grade gliomas using voxel-based lesion-symptom mapping. Neuro Oncol 17:282–288. https://doi.org/10.1093/neuonc/nou130

    Article  PubMed  PubMed Central  Google Scholar 

  68. Weller M, van den Bent M, Tonn JC, Stupp R, Preusser M, Cohen-Jonathan-Moyal E, Henriksson R, Le Rhun E, Balana C, Chinot O et al (2017) European Association for Neuro-Oncology (EANO) guideline on the diagnosis and treatment of adult astrocytic and oligodendroglial gliomas. Lancet Oncol 18:e315–e329. https://doi.org/10.1016/S1470-2045(17)30194-8

    Article  Google Scholar 

  69. White K, Connor K, Clerkin J, Murphy BM, Salvucci M, O’Farrell AC, Rehm M, O’Brien D, Prehn JHM, Niclou SP et al (2020) New hints towards a precision medicine strategy for IDH wild-type glioblastoma. Ann Oncol 31:1679–1692. https://doi.org/10.1016/j.annonc.2020.08.2336

    CAS  Article  PubMed  Google Scholar 

  70. Wilkerson MD, Hayes DN (2010) ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics 26:1572–1573. https://doi.org/10.1093/bioinformatics/btq170

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  71. Woroniecka KI, Rhodin KE, Chongsathidkiet P, Keith KA, Fecci PE (2018) T-cell dysfunction in glioblastoma: applying a new framework. Clin Cancer Res 24:3792–3802. https://doi.org/10.1158/1078-0432.CCR-18-0047

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  72. Wu F, Wang ZL, Wang KY, Li GZ, Chai RC, Liu YQ, Jiang HY, Zhai Y, Feng YM, Zhao Z et al (2020) Classification of diffuse lower-grade glioma based on immunological profiling. Mol Oncol 14:2081–2095. https://doi.org/10.1002/1878-0261.12707

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  73. Yang W, Lai Z, Li Y, Mu J, Yang M, Xie J, Xu J (2019) Immune signature profiling identified prognostic factors for gastric cancer. Chin J Cancer Res 31:463–470. https://doi.org/10.21147/j.issn.1000-9604.2019.03.08

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  74. Yao Y, Yan Z, Lian S, Wei L, Zhou C, Feng D, Zhang Y, Yang J, Li M, Chen Y (2020) Prognostic value of novel immune-related genomic biomarkers identified in head and neck squamous cell carcinoma. J Immunother Cancer. https://doi.org/10.1136/jitc-2019-000444

    Article  PubMed  PubMed Central  Google Scholar 

  75. Zhao Z, Meng F, Wang W, Wang Z, Zhang C, Jiang T (2017) Comprehensive RNA-seq transcriptomic profiling in the malignant progression of gliomas. Sci Data 4:170024. https://doi.org/10.1038/sdata.2017.24

    CAS  Article  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

This research is supported in part by the Intramural Research Program of the National Institutes of Health (NIH), National Cancer Institute (NCI), Center for Cancer Research (CCR). This work utilized the computational resources of the NIH HPC Biowulf cluster. The results shown here are based upon data generated by the TCGA Research Network: https://www.cancer.gov/tcga and Chinese Glioma Genome Atlas: http://www.cgga.org.cn/. We thank the NCI-Center for Cancer Research Single Cell Analysis Facility led by Dr. Michael Kelly for their help generating scRNA data and Dr. Eytan Ruppin for the helpful scientific discussions. Data generated for this study will be shared by the lead contact upon reasonable request.

Funding

Open Access funding provided by the National Institutes of Health (NIH).

Author information

Affiliations

Authors

Contributions

Bharati Mehani and Saleembhasha Asanigari have analyzed the data and organized all Figures. Bharati Mehani, Saleembhasha Asanigari, Hye-Jung Chung, Karen Dazelle, Arashdeep Singh, Sridhar Hannenhalli, and Kenneth Aldape have conceived the idea and contributed to data generation and interpretation. Bharati Mehani and Kenneth Aldape have written the first draft. Kenneth Aldape has designed the study and arranged for the required funds. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Kenneth Aldape.

Ethics declarations

Competing interests

The authors declare no competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Figure S1

. Flowchart to represent the workflow followed in this study.

Additional file 2: Figure S2

. A) Stacked bar plots representing the relative proportion of 22 immune cell types across all IDH-WT and IDH-MUT samples, where each color indicates each cell-type. B) Scatter plots representing a significant positive correlation between the proportion of macrophage M2 with number of segments or with fraction altered or with homologous recombination defects, where ‘‘fraction altered’’ represents the fraction of bases deviating from baseline ploidy (defined as above 0.1 or below -0.1 in log2 relative copy number (CN) space), while ‘‘number of segments’’ represents total number of segments in each sample’s copy number profile. C) M2 macrophage-based gene expression clusters represented by tSNE, followed by D) bar-plots showing IDH specific enrichment in each cluster, E) Kaplan–Meier curves denoting distinct survival between the clusters, and F) forest plots representing the survival differences corrected by IDH status, grade, and age.

Additional file 3: Figure S3

. Gene expression clustering and the separation between IDH-MUT and IDH-WT tumors. A) UMAPs denoting 10 independent gene expression clusters identified in all 3 datasets harboring a clear separation between IDH-MUT (green) and IDH-WT (orange) tumors. B) Bar-plots showing IDH specific enrichment in each cluster, C) Kaplan–Meier curves denoting distinct survival between the clusters, and D) forest plots representing the survival differences corrected by IDH status, grade, and age.

Additional file 4: Figure S4

. Deconvolved M2 macrophages gene expression profiles in non-glioma tumors. Unsupervised clustering of M2 macrophages gene expression profiles. A) tSNE plot representing the samples colored by their cluster groups. B) Kaplan–Meier curves estimating survival probability for each unsupervised cluster.

Additional file 5: Figure S5

. Tumor groups with distinct immune signatures based on LM22 clusters. Heatmaps representing the hierarchical clustering of LM22 clusters. The above annotation bars representing the distributions of cluster assignments, tumor grade, MGMT promoter methylation, EGFR expression, TERT expression (surrogating the status of TERT promoter mutation) and Chr7 gain & Chr10 loss or Chr1p/19q loss with Kaplan–Meier curves below the heatmap denoting their survival differences between these immune-based clusters for A) IDH-WT and B) IDH-MUT tumors.

Additional file 6: Figure S6

. Forest plots displaying prognostic association of the clusters from selected cell types which were consistently significant in A) IDH-WT tumors that involve i) M0 Macrophages, ii) M2 macrophages, iii) Dendritic activated cells and iv) B memory cells. Similarly, forest plots representing the prognostic clusters from B) IDH-MUT tumors which were significantly consistent across all datasets involving i) M0 Macrophages, ii) T helper and iii) B memory cells.

Additional file 7: Figure S7

. Box plots demonstrating distribution of A) ITGA4 expression, B) cytolytic scores, and C) T cell exhaustion scores between the two immune-based survival groups.

Additional file 8: Figure S8

. A) Stacked bar plots representing the relative proportion of 4 malignant cell states across all IDH-WT samples. B) Stacked bar plots representing the relative proportion of 3 malignant cell states in IDH-MUT tumors. C) Dodged boxplots depicting differential representation of 3 IDH-MUT specific malignant cell states between the two immune-based survival groups derived from 3 datasets.

Additional file 9: Figure S9

. Scatter plots representing the correlation between proportions of A) M1 macrophages, B) M2 macrophages, and C) monocytes with MES-like component of all IDH-WT tumors from each dataset.

Additional file 10: Figure S10

. A) UMAP representing unsupervised clustering of 89,926 cells from 19 tumors revealed 18 distinct clusters. B) Color-coded UMAP representing sample-wise composition of each cluster. C) Copy number-based identification of malignant cells and are marked in red and green. D) Color-coded UMAP representing distinct cell types identified from. E) IDH-WT tumors specific single cell derived signature matrix.

Additional file 11: Figure S11

. Venn diagram demonstrating the number of MES-like correlating genes from A) M0 macrophages and B) M2 macrophages that were common in all 3 IDH-WT datasets followed by significantly enriched overlapping pathways.

Additional file 12: Figure S12

. Venn diagram highlighting the number of common genes upregulated in non-codeleted tumors compared to 1p/19q codeleted IDH-MUT tumors in all 3 datasets from A) M0 macrophages, B) M2 macrophages and C) Tcell CD8 cells, followed by significantly enriched overlapping pathways.

Additional file 13: Table S1

. Table listing abbreviations used for all the 22 cell types from LM22 and 10 cell types from LM10 utilized for deconvolution.

Additional file 14: Table S2

. Table listing cluster memberships for each sample derived from cluster of cluster analysis of immune based clusters.

Additional file 15: Table S3

. Table listing inhouse samples with scRNA Seq, with cell type annotation and scRNA derived signature matrix.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Mehani, B., Asanigari, S., Chung, HJ. et al. Immune cell gene expression signatures in diffuse glioma are associated with IDH mutation status, patient outcome and malignant cell state, and highlight the importance of specific cell subsets in glioma biology. acta neuropathol commun 10, 19 (2022). https://doi.org/10.1186/s40478-022-01323-w

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s40478-022-01323-w

Keywords

  • Glioma
  • Tumor microenvironment
  • Malignant cell-state
  • Deconvolution
  • CIBERSORTx
  • Prognosis