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

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.

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/ kimsc 77/ 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 CIB-ERSORTx with its default parameters and B-mode batchcorrection [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-memrest, 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, AClike 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].

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 Coxregression 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.

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.

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).
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 immunebased 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.
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 DCactive 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)  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, immunebased 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 CIBER-SORTx 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 epithelialmesenchymal 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 MESlike cells. 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.
Unsupervised clustering of the checkpoint inhibitorsamples 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 epithelialmesenchymal 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 microenvironmenttargeted 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 cellstates 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, MESlike 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   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. 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.