Computational immune infiltration analysis of pediatric high-grade gliomas (pHGGs) reveals differences in immunosuppression and prognosis by tumor location

Immunotherapy for cancer has moved from pre-clinical hypothesis to successful clinical application in the past 15 years. However, not all cancers have shown response rates in clinical trials for these new agents. igh-grade gliomas, in particular, have proved exceedingly refractory to immunotherapy. In adult patients, there has been much investigation into these failures


INTRODUCTION
Pediatric high-grade gliomas (pHGGs) are aggressive brain tumors in children with poor prognoses and limited therapeutic options. A frequent mutation in pHGG subtypes is amino acid substitutions in histone tails, specifically on histone H3.1 and H3.3. Lysine-to-methionine (H3.1/3.3-K27M) mutations occur in brainstem and midline tumors almost exclusively, and indicate the worst prognosis among pHGGs. Hemispheric tumors arise in the cerebral cortex and are often H3-WT but sometimes feature H3.3-G34R/V mutations, which have a worse prognosis than H3-WT but significantly better than H3.1/3.3-K27M [1]. Radiotherapy is the standard of care for brainstem tumors, while hemispheric tumors may add chemotherapy or targeted therapy in combination with radiotherapy depending on detected mutations [2]. There have been no significant advances in pHGG therapy and these cancers are in desperate need of inventive and efficacious modalities. Clinical trials have recently begun investigating immuno-modulating therapies for pHGG, including vaccines (NCT01130077, NCT03334305, NCT03615404), immune checkpoint blockade (NCT03690869), and cytokine therapy (NCT03330197). For these interventions to work properly, there must be cytotoxic immune cells present either in the tumor, or in the peripheral blood that can traffic to the tumor site, to be stimulated and become more active. These trials are not designed to account for mutational and anatomical differences among pHGG patients, which may play a role in the efficacy of immunotherapies if immune infiltration differs by these factors. At present, there have been limited investigations on the immune status of pHGGs that include hemispheric and brainstem tumors, and how different immune cell subtypes may contribute to patient prognosis.
Within this report, we use the computational methods CIBERSORT, MCP Counter, and CIBERSORTx to investigate a pHGG patient dataset, which includes both H3-WT hemispheric and H3-K27M brainstem gliomas. We find that distributions of immune cells differ between these tumor locations and that positive patient prognoses can be predicted by immune cell types. The presence of regulatory T-cells, memory B-cells, eosinophils, and dendritic cells are positively prognostic for hemispheric tumors, but not brainstem tumors. We further investigate immunosuppressive factors in the RNA-Seq data and find brainstem tumors possess a gene signature consistent with a more immunosuppressive microenvironment. We also determine genetic correlates of patient survival that suggest cytokines and growth factors are influential in the progression of pHGGs and may represent targets for developmental therapeutics.

RNA-Seq dataset and CIBERSORT/CIBERSORTx/MCP counter analysis
Raw RNA-Seq from dataset published by Mackay et al [1] was analyzed by the CIBERSORT algorithm with the standard LM22 matrix [3,4]. Patient demographic characteristics of the population analyzed are summarized in Table 1. The dataset was used for MCP Counter analysis [5] and CIBERSORTx analysis. Individual patients were segmented into H3-WT hemispheric tumors and H3-K27M brainstem tumors for further analysis and matched to survival data. Midline and G34R/V hemispheric tumors were not used for survival analyses due to a lack of statistical power. CIBERSORT values per patient and immune cell type were classified as significant (P < .05), nonsignificant (P > .05), and undetectable (P-value could not be computed). The p-value from CIBERSORT output is a p-value for the global deconvolution of each sample. MCP Counter scores are presented with a cut-off value of 10, excluding < 2% of all data points. CIBERSORTx in B-mode (batch correction) was applied to the normalized bulk RNA-seq gene expression dataset with the LM22 gene signature to estimate the relative fractions of 22 immune cell types.

Statistical analysis and graphing software
Statistical tests such as Log-Rank survival, Wilcoxon survival, Spearman correlation, ANOVA, and unpaired t-tests were performed with GraphPad Prism 8. Alpha was set at 0.05 and data shown in Figures 2A and 6A were corrected for multiple hypothesis testing whereas all other figures are depicting independent comparisons and do not require multiple hypothesis testing. All graphs were made in GraphPad Prism 8. For CIBERSORTx data, graphs were made in R (v.3.6.0) with R packages ggplot2 and ggpubr. Pvalues were calculated by a two-sided Wilcoxon rank-sum test.

RESULTS
We first examined distributions of immune cells that could be detected by the CIBERSORT platform and their differences between pHGG tumor location. Significant differences were found between WT hemispheric pHGGs and K27M brainstem pHGGs for regulatory T-cells, activated dendritic cells, and eosinophils ( Figure 1A, *P < .05). Addi-TA B L E 1 Description of patient population examined by CIBERSORT. Age, gender, and survival of patients studied are shown for the entire group and by histone mutation status and tumor location. These data are compiled from Mackay et al  tionally, G34R/V hemispheric pHGGs displayed fewer regulatory T-cells and NK cells than WT hemispheric pHGGs ( Figure 1A, **P < .01). As sample sizes were highest for WT hemispheric and K27M brainstem patients, we further analyzed cell type distribution by normalizing significance to total cell number. We found more detectable amounts of CD8, NK, M1 macrophages, and activated mast cells in hemispheric tumors, but more detectable amounts of activated dendritic cells (DCs) and neutrophils in brainstem tumors ( Figure 1B). Of the detectable total for each cell type, hemispheric tumors had more significant amounts of regulatory T-cells and eosinophils. We next compared if the significant immune infiltrate of each cell type held prognostic value for patient sur-vival outcomes. Our previously published work examined NK and CD8 + T cells [6], therefore, here we examined additional cell types. In the lymphoid compartment, we found that memory B-cells (Figure 2A), CD4 + regulatory T-cells ( Figure 2B, P = .01), and activated DCs ( Figure 2C) were positively prognostic in hemispheric pHGGs when patients had presence of each cell type. Notably, this did not hold true for brainstem DIPG patients, who showed no survival benefit for these cell types.
We next examined if cell types in the myeloid compartment were prognostic and varied by pHGG tumor location. We again found that brainstem tumors never benefit from immune infiltrate, but positive associations were found in hemispheric tumors for eosinophils ( Figure 3A,  Figure 3C). Interestingly, neutrophils were negatively prognostic for both hemispheric ( Figure 3D, P = .03) and brainstem ( Figure 3D, P = .01) locations. We compared hemispheric to brainstem tumors using the average patient survival for each type, aiming to profile which immune cell types could predict long-term survivors, or the "long-tail" seen in immunotherapy regimens [7]. Here we found significant differences by tumor location for NK cells, regulatory T-cells, dendritic cells, memory B-cells, eosinophils, monocytes, and M1 macrophages ( Figure 4A and B, *P < .05, ns = not significant). Cytotoxic CD8 Tcells and M2 macrophages could predict small numbers of long-term survivors, but the differences were nonsignificant, and detectable presence of these cells pushed sur-vival below the average for hemispheric pHGGs. However, it is important to note that macrophage and microglia gene expression profiles in brain as assessed by CIBERSORT may overlap and M1 values may be inclusive of other cell types present in the brain.
Given the complete lack of survival benefit in brainstem tumors across several immune cell types, we hypothesized the local tumor microenvironment may be immunosuppressive and lacking in inflammatory signals. We investigated the RNA-Seq data used for CIBERSORT and plotted immunosuppressive genes segregated by tumor location. Brainstem tumors uniformly expressed more immunosuppressive genes, with significant differences ( Figure 5A them to immune suppression in brain tumors [8,9,10,11,12]. We next sought to examine if we could detect secretory cytokines from immune cells within the bulk RNA-Seq data and if they differed by tumor location. Using patients with significant relative NK cell infiltrate as a model, we found hemispheric tumors expressed significantly more TGFβ1, but less IFNG and GZMB, than brainstem tumors ( Figure 5B, *P < .05). TGFβ family members are well known to suppress NK function [13], but NK cells are able to confer survival benefit in hemispheric tumors [14], suggesting NK activating signals are expressed highly enough in hemispheric pHGG to compensate. For tumors with NK infiltrate, the immunomodulatory genes GZMB and SLAMF6 correlate significantly with survival in hemispheric tumors ( Figure 5C, P < .05), and no immunosuppressive genes correlated in brainstem tumors. When examining all patients together, we found that IL10, FGL2, VEGFB, and VEGFC were significantly correlated with hemispheric pHGG survival ( Figure 5D, P < .05). In brainstem tumors, IL10 and IDO2 were significantly correlated ( Figure 5E, P < .05), suggesting that IL10 may be a common immunomodulator across pHGG subtypes.
We concluded our study by using alternate computational analyses of the same pHGG RNA-Seq data and performed MCP Counter analysis as well as CIBERSORTx. The MCP Counter method claims to offer superior intersample analytical capability versus CIBERSORT's intrasample leukocyte resolution [5]. We first plotted MCP Counter scores for pHGG tumor types, which directly correlate to proportion of the indicated cell type in the tumor sample. We found K27M brainstem tumors to have significantly more monocytes, while WT hemispheric tumors had significantly more fibroblasts and significantly less neutrophils compared to G34R/V hemispheric and WT midline tumors, respectively ( Figure 6A, *P < .05). We next probed survival data by taking the top 20% and bottom 20% of MCP Counter scores for each immune cell type and plotted WT hemispheric and K27M brainstem tumor survival. Like our CIBERSORT findings, we found NK cells and CD8 T-cells to benefit WT hemispheric survival but have no or detrimental effects for K27M brainstem survival ( Figure 6B Intriguingly, neutrophils were positively prognostic for WT hemispheric tumors but a negative factor for K27M brainstem tumors. CIBERSORTx is a relatively newly designed informatics platform and was built based on single-cell RNA-seq data and is expected to be more accurate for RNAseq data deconvolution relative to CIBERSORT which uses bulk expression array data.  NK cells, activated DCs, eosinophils, and neutrophils than WT tumors.

DISCUSSION
A handful of studies have investigated neuro-oncology immune infiltrates using a combination of computational and live tissue methods, as well as cohorts of adult glioblastoma, pediatric gliomas, or a mixture of the two. Tang et al recently published an immune risk score (IRS) based upon CIBERSORT data in cohorts of adult glioblastoma [15]. In their analysis, they found low numbers of activated NK cells correlated with poor patient prognosis, matching our observations in hemispheric pHGG [14]. However, they also found that significant infiltration of memory B-cells, activated dendritic cells, and M1 macrophages were negatively prognostic, the opposite of our observations. They also did not report on eosinophils, neutrophils, or regulatory T-cells, possibly because these datasets did not report significant infiltration of these cell types. For single-gene correlations, they found IDO and GZMB to be negatively prognostic with regards to their IRS score. This was again opposite of our observations, however, we correlated expression directly with patient survival, not risk score. Other reports have shown distinct phenotypes of immune cells present in adult gliomas compared to pediatric [16] that may explain these observational differences.
Bockmayr et al developed their own immune signature algorithm to analyze a dataset of over 1000 samples that included both adult glioma and pHGG [17]. Their analysis found that H3-WT gliomas had a significant enrichment in endothelial gene signatures compared to H3-mutated pHGGs, suggesting increased vascularization of H3-WT tumors. This hypothesis is being investigated using mouse models of DIPG and pHGG [18]. They also found that tumors rich in antigen-presenting cells (APCs), such as dendritic cells and helper T-cells, had a favorable prognosis if the tumor also contained cytolytic cells (CD8 T-cells and NK). H3-WT tumors in this cohort contain 6-7X more adult gliomas than pHGGs, and the authors did not separate these cases in their analyses. However, by examining H3-G34R/V pHGGs compared to H3-K27M pHGGs, we can make partial conclusions based upon tumor location in pediatric patients from this data. H3-G34R/V tumors presented much higher proinflammatory signaling compared to H3-K27M, and tumors with a DIPG diagnosis followed the same trend when compared to anaplastic astrocytoma and glioblastoma. However, we cannot rule out the contribution of the H3-G34R/V mutation in these observed phenotypes compared to H3-WT hemispheric pHGG, as we also observed significant differences in immune cell infiltration comparing G34R/V vs. WT hemispheric pHGGs. It is also important to note the differences in algorithms across platforms and the consequences on data generated. This was apparent in our current study since using CIBER-SORT, the presence of neutrophils (P < .05) related with poor survival, but in MCP Counter having high (top 20%) neutrophils correlated with improved survival.
Lieberman et al used both gene expression as well as IHC and functional assays to assess immune infiltrate in pediatric tumors exclusively, allowing direct comparisons of DIPG and hemispheric pHGG [19]. This type of hybrid approach of examining cell histological features as well as genomic data will be important in confirming our macrophage data from CIBERSORT analyses given overlaps in gene expression signatures for microglia and macrophage populations. Regardless, our data were concurrent with that of Lieberman et al in finding that DIPG express lower amounts of TBFβ1 but higher amounts of VEGFα compared to pHGG. Their overarching hypothesis was that lack of immunosurveillance in DIPG was responsible for the low number of immune infiltrates, particularly T-cells. However, we did not find that increased presence of dendritic cells correlated with improved survival in DIPG patients. Furthermore, it has been shown that DIPG patients produce tumor-specific T-cells for the K27M antigen [20], suggesting that tumor microenvironment and trafficking of cytotoxic leukocytes plays a larger role. Recent investigations suggest roles for resident macrophages [21] and chemoattractant gradients [22] in immune reactivity of gliomas, and further studies are needed to probe cytokine secretion capabilities of various pHGG tumor types and their associated resident macrophages.
Collectively, these data combined with our findings suggest brainstem pHGGs possess a harshly immunosuppressive microenvironment lacking in inflammatory signals [16,17,19], potentially explaining why immune infiltrate in these tumors is never positively prognostic compared to hemispheric pHGGs. It should be noted that local neuroinflammation caused by infused CAR-T therapy was shown to be fatal in mouse models [23], indicating that caution must be used when attempting to stimulate cytokines in brainstem pHGGs. Another cogent hypothesis is that vascular differences exist between these tumor locations [17,18,24], preventing the influx of immune cells to the tumor site. Our data showing mast cell and eosinophil associations with tumors are surprising. Little is known regarding the presence of these immune cells in healthy brain, however reports indicate they are mobilized to the brain during injury, and have been noted in adult glioblastoma patients [25,26]. Going forward, immunotherapeutic modalities for pHGG will need to consider tumor location when designing new interventions. We suggest that hemispheric pHGGs may respond well to vaccines, checkpoint blockade, and macrophage depletion, while brainstem DIPG will likely benefit more from carefully titrated adoptive cell therapies, epigenetic modulation, and new surgical delivery techniques.

C O N F L I C T O F I N T E R E S T
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

E T H I C A L S TAT E M E N T
This study did not does not contain any human subjects research or animal research performed by any of the authors. Data sets that were probed are publicly available.

A U T H O R C O N T R I B U T I O N S
CB analyzed and plotted CIBERSORT/MCP Counter output and RNA-Seq data, MF, RW, and SZ performed CIBER-SORT/CIBERSORTx/MCP Counter analysis, and RW and LW edited the manuscript. CB and JC conceived the hypotheses and wrote and edited the manuscript.

F U N D I N G I N F O R M AT I O N
This work was supported by the Schissler Foundation Fellowship and a University of Texas MD Anderson Center for Cancer Epigentics Scholar Award (to CB) and from the (i) National Institutes of Health: R21 NS093387 (to JC); R61/R33 NS111058 (to JC) and the (ii) Brain Tumor SPORE: P50 CA127001 (Developmental Research Project to JC).

A C K N O W L E D G M E N T S
The authors thank Dr. Michael Curran for suggestion of CIBERSORT/MCP Counter and Dr. Chris Jones for providing raw RNA-Seq data for analysis.