Human TME components are depleted in GBM PDOXs upon in vivo passaging

To characterize the recapitulation of TME in preclinical models in vivo we studied PDOXs derived from aggressive high-grade gliomas. Up to date, we have derived a cohort of 46 PDOX models by intracranial implantation of patient-derived organoids (Fig. 1A, Additional file 1: Table S1). Our protocol is based on short-term cultures of patient tumor tissue fragments that form 3D organoids ex vivo [16]. Organoids are further implanted intracranially to form orthotopic xenografts in the brain. While use of NOD/SCID and NSG mice allows for higher engraftment rate, well-established models can be recapitulated also in nude mice with less immunodeficient background. We have previously shown that this protocol allows for an efficient recapitulation of histopathological and molecular features of patient tumors [18, 19, 49] (Additional file 2: Fig S1A-B), without a loss of cancer stem-like properties [23]. Our current cohort comprises 42 IDH wild-type GBMs and 4 IDH mutant high-grade astrocytomas. In-depth assessment of genetic and DNA methylation features in 39 models revealed diverse GBM profiles regularly observed in patients (Fig. 1B, Additional file 1: Table S1). While human tumor cells can be expanded in vivo by serial passaging, we hypothesized that human components of the TME are depleted upon orthotopic xenografting and are replaced by the equivalent mouse counterparts. As expected, flow cytometry confirmed depletion of these TME components already in the first PDOX generation (Fig. 1C). These results were consistent with our prior profiling of human and mouse cells in PDOXs [18, 19, 22]. They also underscore the challenges associated with preserving long-term human TME components within preclinical models.

Fig. 1

Composition of the mouse-derived TME in GBM PDOXs. A Schematic of the preclinical modeling of GBM tumors in PDOXs. Created with Biorender.com. See PDOX characteristics in Additional file 2: Fig S1 and Additional file 1: Table S1. B Oncoplot of glioma-specific somatic mutations, gene amplifications, and deep deletions in the PDOX cohort. Longitudinal PDOXs are highlighted with color. MGMT promoter methylation status in PDOX models in depicted. C Flow cytometric analysis showing depletion of human CD31+ endothelial cells and CD45+ immune cells upon xenografting. Examples are shown for 3 GBM patient tumors (single viable cells) and respective PDOXs models at the first passage (single viable human cells, characterized as GFPneg population in GFP+ NOD/SCID mice). D Top: UMAP projection of scRNA-seq data showing the overall gene expression profile of TME cell types. scRNA-seq data combined the biological groups: nude mouse normal brain (Nu-NB), PDOXs (9 models), C57BL6/N mouse normal brain (BL6-NB), GL261 tumor (3 collection time points: early, middle, late). Bottom: Proportions of TME cell types across different tumors and normal brains. Statistical difference between PDOXs (n = 9) and GL261 (n = 3) was evaluated with two-tailed Student’s t test with Bonferroni correction (***p < 0.001, *p < 0.1); Cell types are color-coded; OPCs: oligodendrocyte progenitor cells

scRNA-seq analyses identify major TME components in GBM PDOXs

Due to the depletion of human TME upon xenografting, we hypothesized that growth of human GBM cells in mouse brains is supported by diverse cell types of mouse origin forming the TME in GBM PDOX tumors. To assess TME composition in an unbiased manner, we performed scRNA-seq on mouse-derived cells. PDOXs were derived by intracranial implantation of GBM organoids into nude mice, which have the least immunocompromised background compared to NOD/SCID and NSG strains. We selected nine genetically and phenotypically diverse models, which were derived from treatment-naïve and recurrent IDH wild-type GBMs (Fig. 1B, Additional file 1: Table S1, Additional file 2: Fig S1A) [18]. Four models represented longitudinal tumors derived from GBM patients prior and after standard-of-care treatment, which included radiotherapy and TMZ (LIH0347: T347/T470, LIH0192: T192/T233) [18]. All tumors were collected at endpoint, when tumors were fully developed. Tumor tissues were microdissected, following the previously evaluated MRI and histopathological features of each model, to ensure minimal contamination of healthy mouse brain. Mouse-derived cells of the TME were purified and processed by Drop-seq. In total we obtained 15,366 cells from nine PDOXs. The data were combined with the Pires-Afonso et al. dataset [24] of the TME of the GL261 syngeneic orthotopic GBM mouse model derived in C57BL6/N (BL6/N) wild-type mice (3 time points, 2492 cells in total). Normal brain controls were included for both mouse strains (1692 cells/nude brain, 1972 cells/BL6/N brain). Unsupervised clustering and uniform manifold approximation and projection (UMAP) analysis based on 21,522 cells and 24,067 genes in total, revealed nine major cellular clusters present in all samples (Fig. 1D, Additional file 1: Table S3). Cell clusters were identified based on the expression of cell type-specific markers (Additional file 2: Fig S2A) and included well-known components of the normal brain and GBM TME, such as astrocytes, endothelial and ependymal cells, pericytes, oligodendrocytes, and oligodendrocyte progenitor cells (OPCs), as well as immune cells (Fig. 1D). All major cellular subpopulations were present in PDOXs, GL261, and normal brain controls. Similar to patients, myeloid cells constituted the major immune component in the PDOX and GL261 models. As expected, T lymphocytes were largely depleted and few functional B and NK cells were detected in PDOXs (Additional file 2: Fig S2B), while the majority of infiltrated lymphocytes in the GL261 TME were T and NK cells, with increased proportions upon tumor development (3.3–13.5%, Fig. 1D). GL261 also displayed higher proportions of oligodendrocytes (15–23%) than PDOXs (0.15–3.3%). In accordance with our previous report [19], PDOXs with stronger angiogenic features (P13, T16, P3) had higher proportions of endothelial cells than more invasive PDOXs (Fig. 1D). No correlation between histopathological features and the abundance of myeloid cells was observed (Fig. 1D, Additional file 2: Fig S1B).

The TME composition exhibited a patient-specific trend, e.g., longitudinal models (LIH0347: T347/T470, LIH0192: T192/T233) showed similar cellular proportions, where PDOXs derived from the LIH0192 patient showed a high percentage of myeloid cells, whereas PDOXs of the LIH0347 patient were particularly abundant in astrocytes. This suggests a potential influence of the genetic background of tumor cells on the TME composition, as has been suggested in human GBMs [13]. Since mesenchymal GBMs were described to contain the highest proportions of myeloid cells, we examined transcriptomic heterogeneity of human GBM tumor cells in PDOXs. PDOXs with high myeloid content did not show an increased abundance of mesenchymal-like GBM tumor cells (Additional file 2: Fig S2C). We also did not observe major differences in the TME composition between PDOXs derived from treatment-naïve and recurrent GBMs. This can be explained by the fact that these models show similar genetic profiles at recurrence (Additional file 1: Table S1) [18] and do not display transcriptomic evolution toward the mesenchymal subtype (Additional file 2: Fig S2C).

We further assessed the ontogeny of cells with active cell cycle programs. Assessment of the expression of key cell cycle marker genes (e.g., Top2a, Mki67, Additional file 2: Fig S2D) revealed that cells with an active cell cycle gene expression program were in general clustered together in the “cycling cells” cluster, with an exception of cycling myeloid cells that were also identified in the original “myeloid cells” cluster. The “cycling cells” cluster was identified in all conditions ranging from 0.5 to 1.6% for normal brains and 1.7 to 7.5% for PDOXs and GL261 tumors (Fig. 1D) and was composed of different cellular entities (Additional file 2: Fig S2E). Importantly, cells exhibiting active cell cycle gene expression profiles constituted a minor proportion within the cells in each cell type (Additional file 2: Fig S2F). The highest proportion of cycling cells was observed for pericytes (28%), ependymal cells (19.5%), and oligodendrocytes (14.5%). Cycling myeloid cells were investigated in the follow-up analyses.

In summary, the mouse-derived TME in PDOXs is composed of cellular types relevant to human GBM [4].

TME subpopulations in PDOXs show transcriptional adaptation toward GBM-specific phenotypic states

As all relevant cell types were detected in the TME, we further tested to what extent murine TME cells in PDOXs are instructed by human GBM cells and acquire GBM-specific molecular profiles. To investigate the transcriptomic changes of each cell type, we compared identified cell population in GBM tumors with the corresponding cells in the naïve nude mouse brain. We detected pronounced transcriptomic differences across all the populations of the TME that were generally stronger than the changes observed in GL261 tumors when compared to normal BL6/N brain (Fig. 2A, Additional file 2: Fig S3A), thus indicating effective crosstalk between human GBM tumor cells and mouse-derived TME. The more prominent changes may result from the different tumor tissue structures between PDOXs and GL261. While PDOXs well recapitulate the cellular tumor niche with large areas of tumor-TME crosstalk, GL261 creates circumscribed tumors mostly representing the angiogenic niche with limited infiltration of brain-derived components of the TME, which remain high in the surrounding normal brain structures (Additional file 2: Fig S1A). Transcriptomic adaptation was most pronounced in myeloid cells, endothelial cells, astrocytes, and OPCs (Additional file 1: Table S4). Myeloid cells in PDOXs displayed transcriptomic programs linked to cell migration, inflammation, and cytokine production (Fig. 2B). Furthermore, key “homeostatic” Mg genes including P2ry12, Tmem119, and Gpr34 [50, 51] were downregulated (Fig. 2B). In parallel, myeloid cells overexpressed GBM-specific TAM markers such as Spp1 (Osteopontin), Fn1, Cst7, and Ch25h pointing toward reciprocal crosstalk with GBM cells and transition to TAMs. We confirmed myeloid cell adaptation across the PDOX and GL261 models by qPCR of FACS-sorted CD11b+ cells (Additional file 2: Fig S3B).

Fig. 2

Transcriptomic adaptation of GBM-educated TME subpopulations in PDOXs. A Differentially expressed genes (DEGs) between the TME of PDOX and GL261 versus corresponding normal brains in identified cell types (FDR ≤ 0.01, |log2FC|≥ 1, Wilcoxon rank sum test with Benjamini–Hochberg correction). B Top four gene ontology terms characterizing DEGs in PDOX versus Nu-NB. C Gene expression levels of exemplary DEGs for distinct cell types in four biological groups: nude mouse brain (Nu-NB), PDOXs (9 models combined), C57BL6/N mouse normal brain (BL6-NB), GL261 tumor (3 time points combined). D Expression levels of exemplary markers in distinct cell types detected and annotated in human GBM tumors by Darmanis et al. [35]

Interestingly, other cell types within the PDOX TME also activated biological processes linked to phenotypic states in GBM. For example, OPCs activated programs of tissue inflammation and regeneration (e.g., Pdgfra, Cspg4, Cspg5, Cacng4, Fig. 2B,C, Additional file 2: Fig S3C), while astrocytes expressed genes linked to metabolic processes and cellular shape, suggesting ongoing reactive gliosis (e.g., upregulated Gfap and Vim, downregulated Slc1a2 and Slc1a3, Fig. 2B,C, Additional file 2: Fig S3D). In agreement with our previous study [19], endothelial cells displayed an activated and proliferative phenotype associated with angiogenesis (Fig. 2B,C, Additional file 2: Fig S1B). We detected similar profiles in corresponding TME subpopulations of human GBM tumors (n = 3589 cells from 4 IDHwt GBMs [35], Fig. 2D). Altogether, these results point toward GBM-specific transcriptomic adaptations of myeloid cells and the main TME components in PDOXs.

GBM-educated myeloid cells in PDOXs are largely of microglial origin

The abundance of myeloid cells was further confirmed by Iba1+ staining in PDOX tumors at the endpoint of tumor growth (Fig. 3A, Additional file 2: Fig S1B). While normal brains of nude mice showed Iba1+ myeloid cells with morphology corresponding to “surveilling” ramified Mg, GBM tumors in PDOXs displayed TAMs with different morphologies. The cellular tumor was in general occupied by Iba1+ TAMs showing amoeboid or hyper-ramified morphology. Tumors displaying less-invasive growth showed a gradient of TAM phenotypes, from ramified and hyper-ramified phenotypes at the invasive front toward amoeboid TAMs in the tumor center. Myeloid cells with macrophagocytic morphology were especially present in areas of pseudopalisading necrosis (PDOXs P13, T16). In these models, we observed a notable accumulation of TAMs at the tumor border. A pronounced accumulation of TAMs was also detectable at a much sharper delineated GL261 tumor border (Fig. 3A). Tumor showing invasive growth without a distinct tumor border showed more uniform, diffuse infiltration and activation of myeloid cells toward amoeboid states. Notably, certain resting and hyper-ramified morphologies remained detectable, particularly at the invasive front, in areas marked by lower tumor cell density (e.g., PDOXs P8 Fig. 3A; PDOXs T101, T192 Additional file 2: Fig S1B).This was in agreement with increased presence of amoeboid Iba1+ cells concomitant with increased tumor cell density upon tumor growth over time in mice (Fig. 3B).

Fig. 3

Ontogeny of GBM-educated myeloid cells. A Representative Iba1 staining in PDOXs depicting myeloid cells in invasive (PDOX P8), intermediate (PDOX P3), and angiogenic (PDOX P13) tumor growth, normal nude brain, and GL261 tumor. Tumor core and tumor border zones are highlighted. Arrows indicate examples of ramified (green), hyper-ramified (red), and amoeboid (black) myeloid cells. Scale bar: 50 µm. Sections were co-stained with hematoxylin. See more examples in Additional file 2: Fig S1B. B Representative Iba1 staining representing myeloid cells in PDOX P3 at different stages of tumor growth. Inserts represent sections of the entire mouse brains. Sections were co-stained with hematoxylin to visualize tumor cell density. C UMAP projection of reference-based mapping of myeloid cells from TME of GBM PDOXs and GL261 tumors and respective normal brains. Three myeloid cell entities were identified: microglia (Mg), peripheral monocyte-derived cells (Mo), and border-associated macrophages (BAMs). Inserts show expression of marker genes: pan-myeloid: Itgam (CD11b), Mg: P2ry12, Mo: Ly6c2, BAMs: Mrc1 (CD206). The color gradient represents expression levels. D Proportions of myeloid cell subpopulations in nude mouse normal brain (Nu-NB), PDOXs (9 models), C57BL6/N mouse normal brain (BL6-NB), GL261 tumors (3 collection time points: early, middle, late). Statistical difference between PDOXs (n = 9) and GL261 (n = 3) was evaluated with two-tailed Student’s t test with Bonferroni correction (***p < 0.001). E Flow cytometry-based quantification of CD45+CD11b+Ly6G−Ly6C−CD206− Mg, CD45+CD11b+Ly6G−Ly6C+CD206− Mo and CD45+CD11b+Ly6G−Ly6C−CD206+ BAMs in the Nu-NB and PDOX P3, P8, and P13 TME. For PDOXs P3 and P13 invasive zone and distant normal brain areas were also collected (n ≥ 3 mice/condition, mean ± SEM, one-way ANOVA with Tukey's HSD correction, ****p < 0.0001, *p < 0.05). See gating strategy in Additional file 2: Fig S4D

Due to pronounced heterogeneity within the TME, we hypothesized that myeloid cells in PDOXs are composed of cells of different ontogeny and phenotypic states. To examine the ontogeny of myeloid cells, we combined our dataset with previously published scRNA-seq data of myeloid cells in GL261 tumors, which assigned the ontogeny of TAMs to Mg, Mo, and BAMs based on transcriptomic profiles [7, 8, 24] (Fig. 3C,D, Additional file 2: Fig S4A-C). Referencing PDOX myeloid cells to Ochocka et al. [7] dataset confirmed a high abundance of Mg (81–100%) and a low proportion of Mo (0–19%) and BAMs (0–8%). P13 PDOX with pronounced angiogenic features showed a higher proportion of Mo and BAMs (7 and 5% respectively), although rare Mo were also present in invasive T101 PDOX. In contrast, GL261 tumors contained notably more Mo (Fig. 3D). Flow cytometry confirmed high proportions of Mg in the PDOX TME compared to Mo and BAMs (Fig. 3E, Additional file 2: Fig S4D-E). While BAM proportions showed a trend toward decreased proportions in the tumor cores, we observed higher levels, although not significant, of Mo in the tumor cores compared to other brain areas in PDOXs and normal brain. It was accompanied by increased levels of neutrophils and lymphocytes in tumor cores of PDOX P13 (Additional file 2: Fig S4F) suggesting that peripheral immune cell infiltration in PDOXs is very low and is limited to the tumor regions with highly disrupted blood–brain barrier.

Of note, we revealed major differences across GL261 datasets, with Ochocka et al. [7] and Pires-Afonso et al. [24] datasets containing 26–32% Mo, whereas Pombo-Antunes et al. [8] dataset carrying more than 78% Mo, according to the different adopted cell isolation strategies (Additional file 2: Fig S4C). In fact, within the first two studies, myeloid cells were isolated at also early stages of tumor development from a larger part of the tumor-containing hemisphere, whereas in Pombo-Antunes et al., cells were extracted at the endpoint and specifically from the tumor center. This highlights that Mg/Mo proportions depend on the sampling approach, suggesting diverse spatial locations in the tumor and adjacent brain regions.

Mg-derived TAMs display heterogeneous transcriptional programs

We further hypothesized that cells of Mg, Mo, and BAMs origin adapt their transcriptome toward GBM-specific states in PDOXs. To interrogate the phenotypic heterogeneity of myeloid cells in normal brain and tumors of different histopathological features, we next took advantage of our unique dataset containing normal brain and TME from PDOXs and GL261. To avoid batch effects arising from different scRNA-seq technologies, we performed reference-free analysis of our in-house Drop-seq dataset on the myeloid compartment, containing Mg, Mo, and BAMs. The analysis stratified myeloid cells into nine phenotypic states: Mg formed seven phenotypic clusters (CL0-6), whereas CL7 and CL8 displayed transcriptional profiles of Mo and BAMs, respectively (Fig. 4A–C, Additional file 1: Table S5). PDOXs showed pronounced transitions toward heterogeneous Mg states, although with variable proportions (Fig. 4B). Homeostatic Mg (Ho-Mg), highly enriched in the normal brain, grouped into two clusters (CL0-1), with CL1 showing lower expression levels of homeostatic genes (e.g., P2ry12, Tmem119, Gpr34, Fig. 4C). Importantly, this was not the result of Mg activation via enzymatic digestion, since the markers of enzymatically activated Mg (e.g., Erg1, Fos [52]) were expressed by a subset of Ho-Mg in CL0 (Additional file 2: Fig S5A). Five phenotypic states were observed to be enriched in Mg-derived TAMs (Mg-TAMs, CL2-6) (Fig. 4A–C, Additional file 1: Table S5). These included classical pro-tumorigenic Mg-TAMs, which were present at the highest levels in CL2 and CL3, and were high for, e.g., Spp1, Cst7, Cxcl13, and Apoe (Fig. 4C). Among these two groups, CL3 presented higher cytokine expression levels (Ccl3, Ccl4), suggesting stronger secretory properties and education by GBM when compared with CL2. Subset of CL3 cells showed also high transcriptional activity (e.g., Rpl17, Rps12, Additional file 2: Fig S5A). As expected, CL3 Mg-TAMs showed higher Ptprc (CD45) expression and lower levels of homeostatic Mg genes (Additional file 2: Fig S5A). This is reminiscent of the decrease of homeostatic genes in reactive Mg, known to occur in the GBM TME, but also under inflammatory and neurodegenerative conditions [53,54,55,56]. Mg transition toward CD45high and CCR2+ TAM states in the tumor core was further detected by flow cytometry (Fig. 4D). While Mg in distant brain areas (PDOXs P3, P13) resembled normal brain characteristics, the invasive niche (PDOX P3) showed partial activation of Mg toward Mg-TAMs.

Fig. 4

GBM-driven activation of Ho-Mg toward heterogeneous Mg-TAMs. A UMAP plot showing clusters of myeloid cells in PDOXs, GL261, and normal brain controls revealing nine distinct clusters (CL). CL0-6 represent microglia (Mg), CL7 peripheral monocyte-derived cells (Mo), CL8 border-associated macrophages (BAMs). B Proportions of cells assigned to nine clusters of myeloid cells in nude mouse normal brain (Nu-NB), PDOXs (9 models), C57BL6/N mouse normal brain (BL6-NB), GL261 tumor (3 collection time points: early, middle late). CL0-1: homeostatic Mg (Ho-Mg), CL2-6: Mg-derived tumor-associated macrophages (Mg-TAMs), CL7: Mo, CL8: BAMs. C Discriminative marker genes for each myeloid state (row z-scores of the expression levels). D Representative flow cytometry graphs and quantification of CD45high and CCR2+ cells in CD45+CD11b+Ly6G−Ly6C−CD206− Mg in PDOX P3, P8 and P13. For PDOXs P3 and P13 invasive zone and distant normal brain areas were also collected (n ≥ 3 mice/condition, mean ± SEM, one-way ANOVA with Tukey’s HSD correction, ****p < 0.0001, ***p < 0.001). E Relative transcription factor (TF) activity of regulons identified by SCENIC in myeloid clusters. Regulons with RSS < 0.05 are shown

Additional subpopulations included Mg-TAMs displaying astrocytic features (CL4, e.g., Sparcl1, Gfap), and expression of endothelial cell markers (CL5, e.g., Pecam1, Cldn5) at similar levels than the original astrocyte and endothelial cell clusters (Additional file 2: Fig S5A). We excluded contamination by other TME subpopulations within these clusters as these cells expressed myeloid cell markers, including Hexb, Csf1r, Itgam (CD11b), and Ptprc (Additional file 2: Fig S5A) and showed low scores for potential doublets (Additional file 2: Fig S5B). CL4 and CL5 cells also displayed expression of pro-tumorigenic markers, but at more variable levels than CL3 Mg-TAMs. However, we cannot exclude the possibility of detection of the mRNA from astrocytic and endothelial cells due to potential phagocytosis of these cells by Mg-TAMs. We also identified a separate Mg-TAM cluster corresponding to cells with activated cell cycle programs (CL6). While GFAP+ astrocytic-like, CD31+ endothelial-like and Ki67+ cycling Iba1+ myeloid cells were detected by immunohistochemistry and flow cytometry (Additional file 2: Fig S5C-G), further investigation is required to validate the presence of specific markers at the protein level in the spatial context.

To understand the interdependence between identified Mg states, we further performed the trajectory analysis of Mg cells. TSCAN-based trajectory analysis revealed a general transition from Ho-Mg to Mg-TAMs via CL2, with more profound differences found for CL5 endothelial-like and CL6 cycling Mg-TAMs (Additional file 2: Fig S5H-I). To further examine the transcriptomic differences between identified myeloid states, we next sought to reveal transcriptomic regulators by conducting SCENIC analyses [30] (Fig. 4E). We identified a high number of regulons for Ho-Mg states including Maf, Nr3c1, and Sox4. While CL2 Mg-TAMs showed again transitory features between Ho-TAMs and CL3, CL3 Mg-TAMs appeared regulated by transcription factors such as Hif1a, Stat1, Mafb, and Irf8 [57] suggesting a role of hypoxia in Mg state transitions toward pro-tumorigenic states. Astrocytic-like Mg-TAMs (CL4) showed high activity for Thra, Sox9, and Sox2, which are known to regulate astrocytic states, while cycling Mg-TAMs (CL6) regulons were enriched in classical cell cycle regulators (e.g., Tfdp1, Atf1). Importantly, although certain regulatory networks were shared between Mg-TAM states, Mo and BAMs, we also detected specific regulons unique for Mo (e.g., Fosl2, Irf1, Cebpb) or shared between Mo and BAMs (e.g., Bach1, Prdm1 and Klf4). These data further highlight the factual differences between transcriptomic states of Mg and show the impact of TME niches in shaping Mg heterogeneity in GBM.

Rare Mo and BAMs undergo phenotypic adaptation toward TAM features in GBM TME in PDOXs

We further aimed to investigate whether Mo and BAMs undergo phenotypic adaptation toward TAMs in GBM tumors developed in PDOXs. Since both cell types were rare in PDOXs, we reverted to the reference-based analysis, which allowed us to discriminate between three Mo states described initially by Ochocka et al. [7] (Additional file