High expression of stromal signatures correlated with macrophage infiltration, angiogenesis and poor prognosis in glioma microenvironment

Introduction

Glioma is the most deadly central nervous system (CNS) tumor. Patients with gliomas always show poor outcome, and their median overall survival remains from 14.6 to 16.8 months. The prognosis of gliomas patients is correlated with tumor subtype, age and sex (Mortazavi, Mortazavi & Paknahad, 2018). Based on tumor morphology and molecular alterations, the WHO updated a new classification system in 2016 (Molinaro et al., 2019). Gliomas are classified based on the absence or presence of mutations in IDH and the absence or presence of 1p/19q chromosomal co-deletion (Lim et al., 2018):IDH-mutant, 1p/19q-co-deleted present a benign prognosis; IDH-mutant, non-1p/19q-codeleted have intermediate survival outcomes; and IDH-wild-type have an unfavorable prognosis. The poor prognosis of gliomas has stimulated the discovery of new treatment strategies such as immunotherapy (Gan et al., 2017; Broekman et al., 2018). However, gliomas present significant resistance to immunotherapy. Resistance of innate immunity prevents the tumor from an immune response, resistance of adaptive immunity deactivates tumor-infiltrating immune cells (Jackson, Choi & Lim, 2019). Thus, in order to carry out further study of immunotherapy, it is vital to figure out the tumor microenvironment in gliomas. Previous study integrated several tumor microenvironment related genes based on The Cancer Genome Atlas (TCGA) database applied by ESTIMATE algorithmbased immune scores (Chen et al., 2018). Even so, the correlation between tumor-infiltrating immune cells and stromal signatures in gliomas are still unknown.

The Chinese Glioma Genome Atlas (CGGA) is a database which contains 2,000 glioma samples from Chinese cohorts. This database was finished a large scale update in July 2019, and 693 mRNA sequencings and matched clinical data were provided, which includes: histology, WHO grades, age, gender, chemotherapy, radiotherapy, overall survival. CIBERSORT is a computational method for quantifying 22 immune cell subtypes fractions from gene expression profile such as microarray or RNA-seq (Chen et al., 2018). Flow cytometry verified this computational method and several studies showed the landscape of tumor-infiltrating immune cells can be accurately determined by this method in various malignant tumors such as colon cancer, breast cancer, lung cancer and prostate cancer (Zhang et al., 2019; Ali et al., 2016; Liu et al., 2017).

In this study, we extracted 693 RNA-seq and matched clinical data from the CGGA database. Then, we analyzed the differences of 22 immune cell subtypes between normal brain tissue and gliomas, and revealed their prognosis value respectively. We also compared immune cell subtypes among different phenotypes such as: histology, WHO grades, age, gender, chemotherapy, radiotherapy, IDH and 1p/19q chromosomal co-deletion. Next, we identified five stromal signatures with high macrophages infiltration and angiogenesis which significantly contribute to prognosis by functional enrichment analysis and protein–protein interaction networks. And the stromal signatures were adopted to construct a classifier. The classifier showed powerful predictive ability (AUC = 0.748). Patients with a high risk score (RS) showed poor survival. Finally, we validated this classifier in the TCGA database and the result was consistent with the CGGA database. The results provide a comprehensive view of the tumor microenvironment and may stimulate new strategies of immunotherapy in Chinese glioma patients.

Methods

Result

The landscape of tumor infiltrating immune cell subtypes in gliomas

Based on the 693 RNA-seq from CGGA database, the different infiltration of 22 immune cell subtypes between normal brain tissue and gliomas were analyzed by t-test. The fraction of M0 macrophages, M2 macrophages, resting NK cells, plasma cells, CD8+ T cells was significantly higher in gliomas tissue than normal brain tissue (Fig. 1A). On the contrary, the fraction of naïve B cells, resting dendritic cells, Eosinophils, Neutrophils, M1 macrophages, activated mast cells, follicular helper T cells (Tfh), γδ T cells and resting CD4+ memory T cells was significantly lower in gliomas tissue than normal brain tissue (Fig. 1B).

A total of nine clinical parameters were analyzed, which includes: age, gender, histology, WHO grades, primary or recurrence, chemotherapy, radiotherapy, mutations in IDH and 1p/19q co-deletion. As a result, gender, primary or recurrence, chemotherapy showed no significant difference of immune cells. Monocytes were decreased in elderly patients (Fig. 2A). M0 macrophages, Tregs and activated dendritic cells were increased in high grade gliomas, whereas monocytes were decreased (Fig. 2B). The fraction of M0 macrophages, Tregs and γδ T cells was higher in glioblastoma (GBM) than astrocytoma (AOA), whereas monocytes and activated mast cells was lower (Fig. 2C). The fraction of activated mast cells, monocytes and resting CD4+ memory T cells was higher in IDH mutant than wildtype, while M0 macrophages, Tregs, γδ T cells and Tfh was lower (Fig. 2D). Tregs were increased in 1p/19q co-deletion, whereas naïve CD4+ T cells were decreased (Fig. 2E).

Prognostic value of tumor-infiltrating immune cells in gliomas

In this study, we divided all 22 immune cell subtypes into high/low infiltration group based on median cutoff. And we analyzed the correlation between 22 immune cell subtypes and patients prognosis by generating Kaplan–Meier survival curve. Result showed that high infiltration of M0 macrophages (P = 0.0104) and resting CD4+ memory T cells (P = 0.027) had an unfavorable prognosis. (Figs. 3A and 3B) Based on the previous studies, M2 macrophages showed pro-tumor function as well as M1 macrophages showed anti-tumor function in gliomas. We constructed a macrophage score: (M0 + M2)/M1. And results showed gliomas patients with high macrophage score correlated with poor outcome. (Fig. 3C).

Functional enrichment analysis of gene signature of macrophages

As the infiltration of macrophages play a significant role in glioma patients’ outcomes, we further investigated the potential mechanism and gene signature. First, we divided gliomas patients into high/low macrophage-infiltrating groups. Functional enrichment analysis was performed by Gene Set Enrichment Analysis Gene Set Enrichment Analysis (GSEA). GO categories include BP, MF, CC and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway were analyzed. As a result, a total of 5 GO terms of biological process, 2 GO terms of cellular component, 8 GO terms of molecular function and 4 KEGG pathway terms were identified to be significant (FDR < 0.1 and P < 0.001) (Fig. 4). Interestingly, the 3 enrichment GO terms with high macrophage infiltration were correlated with angiogenesis, so we obtained the 102 core gene signatures for next investigation.

Investigation of stromal signatures correlated with macrophage infiltration and angiogenesis

We built PPI networks to further investigate the 102 core gene signatures which correlated with high macrophage infiltration and angiogenesis by STRING tool. Cytoscape and Cytohubba were used to analyze hub genes (Fig. S3). And the gene set which contain 149 stromal signatures was provided by ESTIMATE and MSigDB. Then, we compared 102 hub genes with 149 stromal signatures. Eventually we generated five stromal signatures including HSPG2, FOXF1, KDR, COL3A1, SRPX2 (Fig. 5A). To investigate the prognostic value of five stromal signatures, we analyzed overall survival by Kaplan–Meier survival curves in the CGGA database. The result show that five stromal signatures were significantly correlated with poor outcome (P < 0.05) (Figs. 5B–5F).

Construction and verification of prognostic classifier based on stromal signatures

According to univariate COX regressions, all of five stromal signatures were significantly correlated with poor outcome in gliomas patients (HR > 1) (Fig. 6A). We calculated the RS of each sample based on the multivariate COX coefficient, and the low/high risk groups were defined according to the median cutoff RS. Then, K–M plot showed patients with high RS had lower survival than patient with low RS (P < 0.001) (Figs. 6B and 6E). Next, we generated univariate and multivariate analyses which enrolled clinical features and stromal signatures in the overall set. The result showed that stromal signature classifier was an independent factor for gliomas patients (Figs. 6C and 6D). Furthermore, we found all of five stromal signatures were correlated with WHO grades (P < 0.05) (Fig. 6F). As for histology, the expression of stromal signatures in GBM was significantly higher than AOA and oilgodendroglioma. The expression of COL3A1, HSPG2 and SRPX2 in patients with chemotherapy was higher than patients without chemotherapy. And the expression of COL3A1, KDR and SRPX2 was higher than patients with radiotherapy than patient without radiotherapy (Fig. S1). Principal component analysis showed significant different between groups of low and high RS (Fig. 6G). A receiver operating characteristic (ROC) curve was plotted to assess the effect of the classifier. The AUC of the classifier was 0.748. Although stromal classifier was less efficient than the WHO grade, it was much superior to histological type (Fig. 6H). To verify the classifier, we generated an independent cohort based on TCGA dataset. And the classifier also showed strong predictive ability in TCGA cohort (AUC = 0.726) (Fig. 6I).

Discussion

First of all, based on the RNA-seq data from CGGA and the CIBERSORT, we calculated 22 subtypes of tumor-infiltrating immune cells in gliomas. We comprehensively analyzed the tumor-infiltrating immune cells present in gliomas and revealed the prognosis value.

The landscape of tumor-infiltrating immune cells in gliomas as follows:

1. B cell lines: Naïve B cells can differentiate into antibody-secreting plasma cells after being stimulated by antigen (Liao et al., 2016). Our study found that naïve B cell was a significant decrease in gliomas, meanwhile, the plasma cells were significant increased. It was illustrated that naïve B cell may be differentiated into plasma cells by gliomas antigens.

2. T cell lines: Tregs are correlated with unfavorable prognosis in the several kinds of tumor microenvironment (e.g., ovarian cancer, breast cancer, kidney cancer and pancreatic cancer) (Mirzaei, Sarkar & Yong, 2017). Our results showed that Tregs were significantly increased in GBM and high WHO grade, which demonstrated Tregs may promote the development of gliomas (Banissi et al., 2009; Vandenberk & Van Gool, 2012). Previous study showed that tumor-infiltrating Tfh cells have protective roles by suppressing lymphoid tumor-promoting effects in breast cancer and colorectal cancer (Woroniecka et al., 2018). Our results also shown that Tfh cells was significantly decrease in gliomas. Previous study showed neutrophils and eosinophils are implicated in almost every stage of oncogenesis and display both anti- and pro-tumor properties (Curran, Evans & Bertics, 2011; Schneider, Kwan & Boockvar, 2018).

3. Myeloid cell lines: Our results showed that neutrophils and eosinophils were significantly decreased in gliomas, which indicated neutrophils and eosinophils may have negative correlation with gliomas. Mast cells play an important role in the growth of tumors (Põlajeva et al., 2011). However, the contribution of mast cells in the microenvironment of solid malignancies remains controversial. Previous studies have illustrated that gliomas are associated with a profound accumulation of mast cells, and STAT5 play an important role on the recruitment of mast cells to gliomas (Põlajeva et al., 2014). Our results showed that resting mast cells were increased in gliomas, whereas activated mast cells were significantly decreased in gliomas. In addition, the infiltration of activated mast cells of GBM was prominently lower than AOA. Consistent with previous study, gliomas promoted mast cells recruitment, but gliomas may also block resting mast cell activation. And glioblstoma may show stronger ability of preventing mast cell activation than AOA. Dendritic cell vaccinations have emerged as newly strategies in the treatment of gliomas (Garg et al., 2016; Wen et al., 2019). Adjuvant dendritic cell-based immunotherapy for gliomas patients can induce long-term survival (Li et al., 2018). We found resting dendritic cells were significantly decreased in glioma, which was consistent with previous studies that dendritic cells recruitment correlated with favorable outcome. Previous study showed monocytes played an important role in cancer development and progression (Põlajeva et al., 2011). And different monocytes displayed both anti- and pro-tumor function (Prosniak et al., 2013). Our results showed monocytes were significant decreased in high grade glioma and GBM, which indicated stimulating monocytes recruitment may prevent progression of high grade gliomas. What’s more, our results also showed monocytes were significantly lower in patients with radiotherapy than patients without radiotherapy, which may contribute to the resistance of gliomas radiotherapy.

Next, we further investigated tumor-infiltrating macrophages in gliomas. Based on their function, M0 macrophages can be differentiated into two categories: classically activated macrophages (M1 macrophages) and alternatively activated macrophages (M2 macrophages) by stimulation of various cytokines (Chen et al., 2019; Jiang et al., 2019; Lailler et al., 2019). Previous studies show gliomas contain a plenty of macrophages, and tumor-associated macrophages correlate with progression and angiogenesis of gliomas by several ways, including: microRNA-macrophage feedback loop (Bao & Li, 2019; Liu et al., 2019), extracellular lipid loading (Offer et al., 2019) and cytokines like OPN, IRGM, IL-6 (Hernandez-SanMiguel et al., 2019; Wang et al., 2018; Xu et al., 2019). In addition, previous study successfully reeducated the pro-tumor M2 macrophages toward anti-tumor M1 macrophages by a dual-targeting biomimetic treatment strategy which provides a good method for the pharmacotherapy of gliomas (Zhao et al., 2018). Our results showed that M0 and M2 were significantly increased in gliomas, while M1 were significant decreased. The results confirmed pro-tumor function of M0 and M2 as well as anti-tumor function of M1 in gliomas. And the infiltration of M0 macrophages showed correlated with decreased survival (P < 0.05). In addition, we found M0 macrophages were significantly lower in IDH mutant than IDH wildtype. This result may explain glioma patients with IDH-mutant always have favorable survival outcomes. After functional enrichment analysis, 3 gene sets includes 102 genes about angiogenesis were detected in high macrophage-infiltrating group, which was consistent with previous study.

Among 102 core genes, we found five stromal signatures correlated with high macrophages infiltration and angiogenesis indicated poor prognosis (P < 0.05). The five stromal signatures including HSPG2, FOXF1, KDR, COL3A1 and SRPX2. Sun et al. (2019) constructed the multicellular gene network between gliomas and macrophages, they found macrophage-related gene signature had good prognostic value for predicting resistance to targeted therapeutics and survival of glioma patients. We found the relationship between stromal signatures and macorphages signatures in gliomas. Based on precious studies, high expression of perlecan/HSPG2 in gliomas lead to tumor promotion through the transformation of brain extracellular matrix into tumor microenvironment (Hu et al., 2016). In addition, high expression of COL3A1 contributes to glioma cells proliferation and migration (Shin et al., 2017). And up-regulation of KDR promotes proangiogenic myeloid cells, result in low-grade to high-grade transition (Huang et al., 2017). Moreover, glioma with high tumor associated macrophages has increased expression of SRPX2 (Hung et al., 2016). We also found these five stromal signatures were correlated with clinical features like grades, histology, chemotherapy and radiotherapy. The expression of stromal signatures was significantly higher in patients with higher grade and worse pathological features like GBM. Then, we constructed a prognosis classifier based on HSPG2, FOXF1, KDR, COL3A1 and SRPX2. And the classifier showed strong predictive ability in both CGGA cohort and TCGA cohort. We confirmed the high infiltration of M2 macrophages in protein and histology level via CPTAC and the human protein atlas. Eventually, based on our investigation, up-regulation of five stromal signatures may lead to macrophage infiltration and angiogenesis in gliomas patients. And our classifier also provides strong predictive ability in prognosis. Our study may provide a novel immunotherapeutic strategy and a more clearly overview of tumor microenvironment in gliomas.

However, there still exited several limitations in our research. First, the calculation results from public database may show bias. Although we have verified the results in two independent cohorts, we should continue to deeper research. Second, the macrophage infiltration and angiogenesis function of five stromal signatures in gliomas need to further confirmation in vitro and in vivo.

Conclusions

Based on CGGA RNA-seq datasets, the landscape of 22 tumor-infiltrating immune cells in gliomas was analyzed by CIBERSORT. High infiltration of macrophages was correlated with poor outcome. After functional enrichment analysis, three gene sets includes 102 core genes about angiogenesis were detected in high macrophage-infiltrating group. Next, we found five stromal signatures indicated poor prognosis which including HSPG2, FOXF1, KDR, COL3A1, SRPX2. Five stromal signatures were adopted to construct a classifier. The classifier showed powerful predictive ability. Finally, we validated this classifier in TCGA and the result was consistent with CGGA. Our investigation of tumor microenvironment in gliomas may stimulate the new strategy in immunotherapy. The five stromal signatures correlated with poor prognosis also provide a strong predator of glioma patient outcomes.

Supplemental Information