Neoplastic and healthy cell clusters are identified in glioma
To gain insights into the subpopulations and cellular diversity of gliomas, we first analyzed the scRNA-seq data of 6341 IDH-A glioma cells from 10 patients [9] by performing t-distributed stochastic neighbor embedding (t-SNE) based on the top 5% highly variable genes in expression. As shown in Fig. 1a, clusters 1 to 10 were separately composed of the cells from a single one of 10 patients, whereas clusters of 11–13 were mixed with the cells originating from distinct patients (Fig. 1b). Since tumor cells are more heterogeneous than normal cells, the cells in clusters 1–10 could be neoplastic, while the cells in clusters 11–13 might be healthy cells. We further validated the neoplastic and healthy cell clusters based on the expression profile of the known markers of mature human neurons astrocytes, oligodendrocytes, microglia/macrophages, and endothelial cells [25]. Previous study showed that the expression of neoplastic cell markers like EGFR can be used to identify neoplastic cells with high sensitivity and specificity [7]. We also found that neoplastic cell clusters of gliomas express significantly higher level of EGFR than that of healthy cell clusters (Fig. 1c). However, those neoplastic cell clusters barely expressed the myeloid-specific gene PTPRC (CD45), but this gene was broadly expressed in healthy cells.
To further explore the TME of non-neoplastic cells in glioma, we used the reported marker genes for different residential cell types to define the identity of each cell cluster [25]. As expected, we observed large expression differences among several cell types including neurons, oligodendrocytes, microglia/macrophages, and endothelial cells, in those three healthy cell clusters (Fig. 1d). We then annotated the cells with known marker genes and found that microglia/macrophages occupy the largest portion of healthy cells (Fig. 1d).
Macrophage and microglia cells are distinguished from non-neoplastic clusters
Previous study suggested that 237 human-mouse homologs are enough to discriminate discrete human glioma TAMs from macrophage/microglia cluster through the principal component analysis (PCA) based on scRNA-seq data [26]. We then explored those 237 genes in the macrophage/microglia cell cluster of glioma using PCA and found that macrophage and microglia cells are stratified into two subpopulations, respectively (Fig. 2a). Moreover, Muller et al. previously identified 66 potential maker genes that were strongly differently expressed between macrophage and microglia cells [26]. We also observed that those 66 marker genes of macrophage/microglia are respectively enriched in the two subpopulations identified by us, which further confirmed that one subpopulation is macrophages while the other is microglia (Fig. 2b).
To investigate the polarization heterogeneity of macrophages, we compared the expression level of specific markers for M1 (e.g. CD64, SOCS1 IDO, CD86, CD80, IL1R1, TLR2 and TLR4) and M2 (such as MRC1, TGM2, CCL22, CD163, TLR1, TLR8 and CLL22) macrophages [27]. Intriguingly, we found that the markers of M1/M2 co-exist in macrophages and no significant correlation is observed between the marker expression of M1 and M2 (Fig. 2c). Accordingly, our results further confirm previous finding that macrophages may have mixed M1/M2 phenotypes.
Expression correlation analysis reveals significant autocrine ligand-receptor pairs in glioma
Since the expression level of ligands and their receptors could reflect the level of cell-to-cell communication [28], we tried to identify the autocrine network among glioma neoplastic cells. We first collected 696 ligands and 653 cognate receptors (2419 pairs of chemokine/cytokine) from a public interaction dataset [20], and then compared their expression between neoplastic and non-neoplastic cells. By screening the expression of both receptors and ligands increased in glioma neoplastic cells with MAST [22], 186 ligand-receptor pairs were detected. Considering that co-expression of a ligand and its corresponding receptors is necessary for cell-cell communication through secreting signals, thus we subsequently calculated the Spearman’s correlation coefficients for the 186 ligand-receptor pairs using the TCGA Low-Grade Glioma (LGG) dataset [29]. A total of 16 pairs with significant Spearman’s correlation coefficients higher than 0.4 (P-value < 0.05) were identified (Fig. 3a, Additional file 1: Table S1). Furthermore, down-regulation of NOTCH1 and its ligand DLL1 has been reported to inhibit tumor proliferation in glioma cell lines, which may play the function of promoting tumor proliferation [30]. Our finding showed that NOTCH1 and DLL1 are both highly expressed in tumor cells (Fig. 3b). Additionally, previous research suggested that the secretion of NLGN3 could promote glioma growth [31], and we found that its binding partners are also highly expressed in tumor cells (Fig. 3c).
Specifically, GDNF (Glial Cell Derived Neurotrophic Factor) has been hypothesized to promote tumor growth and invasion in prostate cancer and its effect of resisting treatment could be related to the expression of its receptor GFRA1 [32]. We found that GDNF together with its two receptors in autocrine ligand-receptor pairs are highly expressed in glioma neoplastic cells (Additional file 1: Figure S1a). We also conducted gene functional enrichment analysis and found that those 16 ligand-receptor pairs are mainly involved in the pathways of tumor invasion, cell adhesion and cytoskeleton, cell growth and proliferation (Fig. 3d). Therefore, those 16 ligand-receptor pairs identified by us have important roles and may be tightly associated with the development of glioma.
The crosstalk between CSCs and macrophages is associated with the survival of glioma patients
To discover the crosstalk between CSCs and macrophages based on ligands and receptors, we first selected the stem-like cells with high expression of the 90 genes associated with stemness in glioma (Additional file 1: Table S2) [9]. Then we compared the expression of those ligands and receptors between stem-like and other differentiated tumor cells. The ligands and receptors that were with higher expression in macrophages than that of microglia were further identified. Next, we explored the cell-to-cell communications between stem-like cells and macrophages based on those highly expressed ligand-receptor pairs.
Interestingly, the crosstalk between stem-like cells and macrophages that showed ligands highly expressed in stem-like cells and receptors highly expressed in macrophages was mainly related to invasion and angiogenesis (Fig. 4a, Additional file 1: Figure S2). We found that the CSCs highly express some ligands encoding collagens I and III, such as COL1A1, COL1A2, COL3A1, COL6A1 and COL6A2. A previous study suggested that cancer cells and leukocytes could migrate instantly along collagen fibers [33]. Our results also implied similar mechanism that those ligand-receptor pairs may help glioma cells to migrate in glioma. Moreover, the receptors that highly expressed in macrophages, like CD44, ITGB1, ITGB3 and ITGA5, contained the integrin members. Integrins are known to directly bind the components of extracellular matrix (ECM) and could provide the traction necessary for cell motility and invasion [34, 35]. Intriguingly, we observed that the patients with high expression of COL1A1, ITGB1 and ITGB3 are prone to with significantly poor prognostic outcomes in TCGA LGG datasets (Fig. 4b c, f and Additional file 1: Figure S3; P-value < 0.05). Furthermore, our result shows that BGN (Biglycan) is highly expressed in CSCs, and BGN has been reported to promote tumor invasiveness via inducing integrin-β1 expression [36]. Additionally, laminins (e.g. LAMA4, LAMA5 and LAMC1) may take part in the interactions between CSCs and macrophages as well. We found that the glioma patients with higher expression of LAMA4 is significantly associated with shorter survival days (Fig. 4e, P-value <0.05).
Besides, we also identified a ligand-receptor pair that related to angiogenesis. Previous research suggested that PLGF (Placenta growth factor) could recruit macrophage lineage cells in inflammatory sites [37] [38]. We observed that CSCs remarkably express PLGF and its only receptor VEGFR1 is highly expressed in macrophages which may help macrophages enrich in tumor site. Thus our result implies that CSCs could recruit macrophages via PLGF-VEGFR1 system.
The ligands in macrophages and the receptors in CSCs are also correlated with patient survival
To further explore the interactions between CSCs and macrophages, we analyzed how macrophages communicated with CSCs via ligand-receptor pairs. For the ligands and receptors separately highly expressed in macrophages and CSCs, we identified 24 ligand-receptor pairs that showed crosstalk from macrophages to stem-like cells (Fig. 5a). We found that these 24 ligand-receptor pairs are mainly involved in the pathways of angiogenesis, integrin and apoptosis (Additional file 1: Figure S4).
It has been suggested that TNF expressed by macrophages can promote angiogenesis through IL-8 and VEGF to help tumor cells escape from immune system [39] [40]. Moreover, IL-8 has been shown to be strongly correlated with brain cancer grades [41, 42]. However, we found that IL-8 is mainly expressed by macrophages rather than tumor cells in glioma. Furthermore, we observed that the patients with lower expression of SDC1 is significantly associated with longer survival days (Fig. 5b, P-value <0.05). The high expression of SDC1 is crucial for promoting stem-like cell invasion and metastasis in cell-extracellular matrix adhesion [43]. C3 plays a key role in complement system and could support tumor growth by immunosuppression as well [44]. We observed that C3 is highly expressed in macrophages, and the glioma patients with higher expression of C3 tend to have significantly decreased survival days (Fig. 5c, P-value <0.05). Other ligand-receptor pairs many also influence tumor progress in different ways. For example, Kast et al. suggested that FN1 may play a key role in extracellular matrix components to increase tumor invasion and IL-18 could promote tumor growth [45]. Our result indicates that both FN1 and IL-18 are significantly associated with poor prognosis of glioma patients (Fig. 5d and e, P-value <0.05). In total, 34 ligand and receptor genes in the crosstalk between CSCs and macrophages are significantly associated with survival in the LGG dataset (P-value <0.05).
Prognostic model based on ligand-receptor interactions achieves high accuracy
To predict the prognosis of glioma patients, we built a model based on the ligand-receptor interactions between macrophages and CSCs using advanced machine learning algorithm XGBoost [46]. Firstly, we performed principle component analysis (PCA) based on the genes of ligand-receptor pairs in IDH-mutated and IDH-wildtype gliomas, and found that IDH-mutated and IDH-wildtype gliomas mixed together in PCA (Additional file 1: Figure S5), suggesting that these two types of gliomas have little difference in the expression of the ligand-receptor pairs. Then, we randomly chose 75% of the gene expression dataset of all the TCGA glioma samples as training set and the rest of 25% as test set. Then the training set was used to train a XGBoost model to predict the prognosis of glioma patients. We observed that this XGBoost model achieve a precision of 0.80 and a recall of 0.78 in the test set, which is much better than models built on randomly selected genes (Additional file 1: Figure S6).
To further validate the result of our XGBoost model, we tested the RNA-seq data of 325 glioma samples from CGGA (Chinese Glioma Genome Atlas) [19]. Strikingly, our model can significantly divide the patients into high-risk and low-risk groups as well, and the survival time of these two patient groups are remarkably different (P-value = 0), which further confirmed the performance of our robust model (Fig. 6a). The expression level of ITGB5 is the most important feature in the model followed by SEMA4F, IL18 and HMGB1 (Fig. 6b). Among the top 10 important genes, ITGB5, IL18 and CXCL1 are highly expressed in macrophages, while the receptors of SCD1, TNFRSF6B, ACVR1B are highly expressed in CSCs and the ligands of PODXL2, HMGB1, SEMA4F and APP are mainly expressed in CSCs (Fig. 6c). In addition, previous studies suggested that ITGB5 and IL18 are related to tumor invasiveness [40, 43], whereas CXCL1 and SEMA4F are associated with cell proliferation [47, 48]. Moreover, HMGB1 and SCD1 have been demonstrated to be closely correlated with drug resistance in glioma [49, 50]. However, little is known about the role of ACVR1B and PODXL2 on gliomas.
