Similar pathology and gene expression in medulloblastomas from progenitors or stem cells
To initiate an oncogenic stimulus in CGNPs, we bred SmoM2 mice with Math1-Cre mice to generate Math1-Cre/SmoM2 (M-Smo) pups. To initiate an oncogenic stimulus earlier in brain development in pluripotent stem cells that give rise to CGNPs, we bred SmoM2 mice to Gfap- Cre mice to generate Gfap-Cre/SmoM2 (G-Smo) pups. Both M-Smo and G-Smo genotypes developed cerebellar tumors with 100% frequency and all tumors showed the small, round blue cell morphology typical of medulloblastoma and tendency to invade adjacent brain (Fig. 1a).
Microarray comparison of gene expression in samples from M-Smo and G-Smo tumors showed similar transcriptomic profiles with 64/41349 probes sets detecting statistically significant signals, representing 33 identified transcripts (Supplementary Table 1). One of these differentially expressed transcripts, Protamine 1 (Prm1) was included within the Gfap-Cre transgene, and thus expected to be differentially detected. Finding differential Prm1 provided an internal validation of the assay, while finding only 32 other differential transcripts demonstrated the high similarity between the tumors. Similarly, the same work flow applied to the previously published microarray datasets from G-Smo and M-Smo tumors (18) identified 66 out of 45105 probe sets, representing 54 differentially detected transcripts (Supplementary Table 2). Only 1 gene was differentially expressed in both studies. Both microarray comparisons demonstrated highly similar average gene expression profiles of G-Smo and M-Smo tumors, with minimal consistently observed differences.
Different, clinically-relevant behaviors of medulloblastomas from progenitors or stem cells
In contrast to the similarity in pathology and average gene expression of G-Smo and M- Smo tumors, we noted marked differences in the clinically relevant parameters of survival time and response to therapy. As in the prior study (18), the survival times of untreated G-Smo mice were typically shorter than those of M-Smo mice (Fig 1b). Moreover, while M-Smo mice showed significantly improved survival after radiation therapy (xRT), consistent with our prior studies (20), xRT did not extend the survival of G-Smo mice (Fig. 1c,d). The prognosis of G-Smo and M-Smo mice was therefore markedly different, despite their common oncogenic driver, tumor pathology and similarity in bulk transcriptomic studies.
scRNA-seq identifies difference between medulloblastomas from progenitors or stem cells
To examine the differences between G-Smo and M-Smo tumors with cellular resolution, we analyzed both tumor types using scRNA-seq. Transcriptomic analysis at single cell resolution allowed us to examine tumor sub-populations that might be overlooked in bulk transcriptomic studies. We harvested medulloblastomas from 5 M-Smo and 6 G-Smo at P15, dissociated the tumors and subjected the cells to bead-based Drop-seq analysis, as previously described (10). Putative cells identified by bead-specific barcodes were subjected to exclusion criteria described in Methods, to address the common problems of gene drop out, unintentional cell-cell multiplexing and premature cell lysis (21, 22). 5930 out of 11984 putative M-Smo cells, and 8699 out of 16489 G-Smo cells met criteria and were included in the analysis. To compare the two genotypes at similar sequencing depths, we randomly downsampled the G-Smo transcript counts to 60% of the original depth, as recommended in prior studies (23).
We subjected the scRNA-seq data from M-Smo and G-Smo tumors to a single principle component analysis (PCA) followed by Louvain clustering, as in our prior studies comparing M- Smo tumors with and without treatment (10). We identified 17 principle components that described >78% of the variance and used UMAP to place cells in a 2-dimensional space according to their distances in PCA space, with Louvain clusters color-coded (Fig. 2a). We noted that cells in the same clusters localized close together in the UMAP, supporting the validity of the cluster assignments. To determine the biological relevance of the clusters, we generated cluster-specific differential gene expression profiles (Supplementary Table 3) by comparing for each gene the expression by cells within each cluster to the expression by all cells outside the cluster using Wilcoxon rank-sum test. We then used cluster-specific gene expression patterns to infer the type of cells represented by each cluster.
Using these methods, we identified 8 clusters as different types of stromal cells typical of brain tissue, including astrocytes, oligodendrocytes, macrophages/microglia, vascular cells, fibroblasts and ciliated cells resembling ependymal or choroid plexus cells (Table 1; Fig. 2b). These 8 clusters localized as discrete single-cluster units on the UMAP projection. The other 14 clusters localized to a multi-cluster complex in which each cluster shared a border with other clusters. We identified these 14 clusters as tumor cells in a range of states that paralleled CGNP development, from proliferative cells at different phases of the cell cycle to non-proliferative cells at different stages of neural differentiation (Table 1). We identified proliferative clusters by expression of proliferation marker Mki67, Cyclin expression and SHH transcription factor Gli1, and further characterized proliferative cells as quiescent, cycling or M-phase enriched based on cluster-specific gene expression (Table 1). The non-proliferative clusters showed successive expression of early to later differentiation markers Barhl1, Cntn2, Rbfox3 and Grin2b (Table 1; Fig 2c), as in our prior study of M-Smo tumors (10). We included CGNs as the most differentiated cell type within this group.
Cluster
|
Cell type designation
|
distinctive markers
|
0
|
early differentiating CGNP-like tumor cells
|
Pde1c, Nhlh1/2
|
1
|
proliferative, quiescent tumor cells
|
Hes1, Ccnd1
|
2
|
proliferative, quiescent tumor cells
|
Srebf1, Ccnd2
|
3
|
proliferative, quiescent tumor cells
|
Srebf1, Ccnd1
|
4
|
proliferative, cycling tumor cells
|
Top2a, Lig1, Esco2
|
5
|
differentiating CGNP-like tumor cells
|
Mtss1, Cntn2
|
6
|
proliferative, cycling tumor cells
|
Hells, rrm2
|
7
|
proliferative, cycling tumor cells, M phase enriched
|
Cdc20, Cenp1
|
8
|
proliferative, cycling tumor cells, M phase enriched
|
Aspm, Ccnb1
|
9
|
late differentiating CGN-like
|
Pcp4, Car10
|
10
|
differentiating CGN-like tumors
|
Cntn2, Nhlh1
|
11
|
proliferative, cycling tumor cells
|
Hells, Lig1, Pclaf
|
12
|
proliferative, cycling tumor cells
|
Hells, Lig1, Gli1
|
13
|
CGNs and CGN-differentiated tumor cells
|
Gabra6, Vsnl1
|
14
|
oligodendrocytes
|
Ptpz1, Sox10, Fabp7, Olig1/2
|
15
|
M2 microglia/macrophages
|
Mrc1, C1qb, Aif
|
16
|
astrocytes
|
Aldoc, Aqp4, Fabp7
|
17
|
fibroblasts
|
Lum, Vtn
|
18
|
endothelial cells
|
Cldn5, Flt1
|
19
|
M1 microglia/macrophages
|
Selplg, Siglech, C1qa, Aif
|
20
|
ependymal/choroid plexus cells
|
Rsph1, Folr1
|
21
|
myelinating oligodendrocytes
|
Opalin, Plp
|
Table 1. Identification of clusters as specific types of tumor and stromal cells
To compare the populations within M-Smo and G-Smo tumors, we deconvoluted the UMAP by genotype (Fig 2d). For quantitative comparison, we determined the number of cells from each replicate animal in each cluster, normalized to the total number of cells from that animal, and then compared the cluster populations from replicate M-Smo and G-Smo mice (Fig. 2d,e). We found that most tumor cell clusters were similarly populated in M-Smo and G-Smo tumors. However, Clusters 1, 2 and 7, within the proliferative region, were significantly enriched in G-Smo tumors (p = 0.008 for each by Wilcoxon rank-sum test), while cluster 13, comprising CGNs at the differentiated pole, was significantly enriched in the M-Smo tumors (p = 0.023).
Statistically significant enrichment of smaller magnitude was also seen in fibroblast (p = 0.023) and differentiated oligodendrocyte clusters (p = 0.008) in M-Smo tumors.
The similarity in the populations of most clusters in G-Smo and M-Smo tumors was consistent with the similarity of these tumors in bulk transcriptomic studies. The differential representation of specific types of tumor and stromal cells in G-Smo and M-Smo tumors, however demonstrated differences in tumor subpopulations that could not be detected using bulk transcriptomic analysis. G-Smo tumors showed increased proliferation, while M-Smo tumors showed increased differentiation, and were specifically depleted in the cell types represented by clusters 1,2 and 7. The expression patterns of all genes detected in our studies can be plotted and compared in G-Smo and M-Smo UMAPs through our web-based application: https://gsmovmsmoviewer.shinyapps.io/GvMviewer/.
G-Smo tumors show larger populations of cells expressing stem cell markers
To characterize further clusters 1, 2 and 7 that were composed predominantly of G-Smo cells, we defined the set of genes differentially up-regulated in these clusters compared to all tumor cells in M-Smo mice, excluding stromal cell types (Supplementary Table 4). We excluded stromal cells in order to prevent stromal gene expression patterns from obscuring differences in tumor cell gene expression. We noted that Clusters 1,2 and 7 showed increased expression of genes associated with stem cell phenotype, including Nes, Vim, Olig1 and Olig2. Feature plots of each of these genes confirmed increased expression in G-Smo tumors, particularly in the region of Clusters 1,2 and 7 (Fig. 2f). We selected Olig2 for further study, because recent functional genetic studies have shown that Olig2+ tumor cells in medulloblastoma are cancer stem cells that play an important role in medulloblastoma initiation and recurrence (24).
Different temporal patterns of stem cell behavior in M-Smo and G-Smo tumors
To confirm differential expression of Olig2 at the protein level and to compare the temporal course of Olig2 expression efficiently between genotypes, we labeled tumor sections using immunohistochemistry. Our scRNA-seq data showed that both oligodendrocytes and tumor stem cells express Olig2, and that oligodendrocytes can be distinguished from stem cells by the expression of Sox10 (Fig. 2b,f). Labeling tumor sections with OLIG2 and SOX10 antibodies demonstrated both OLIG2+/SOX10+ cells that we considered to be oligodendrocytes and OLIG2+/SOX10- cells that we considered to be OLIG2-expressing tumor stem cells, equivalent to the Olig2+ cells of clusters 1, 2 and 7. Prior studies found that OLIG2+ tumor stem cells were more numerous in the early stages of medulloblastoma tumorigenesis and decreased as tumors enlarged (18). To compare OLIG2 dynamics over time in G-Smo and M-Smo tumors, we analyzed sections of both genotypes at P5 and P15.
We found that at P5 both M-Smo and G-Smo comprised similarly large fractions of OLIG2+/SOX10- cells (Fig 3b,c). The OLIG2+/SOX10- fraction decreased over time in both genotypes, but the decrease was more marked in M-Smo tumors (Fig. 3d,e). Prior studies showed that sub-curative cytotoxic treatment of medulloblastoma with radiation or chemotherapy induces proliferation of stem cells in the perivascular niche, and that these cells express OLIG2 (24, 25). We compared the dynamics of the OLIG2+/SOX10- population in G-Smo and M-Smo tumors recurring after treatment. We subjected P11 G-Smo and M-Smo mice to cytotoxic treatment with etoposide at a dose calibrated to produce regression followed by recurrence, and then we quantified OLIG2+/SOX10- cells 4 days later, at P15. Compared to untreated tumors at P15, we noted increased expression of OLIG2 in SOX10- cells, particularly in perivascular regions (Fig. 3f,g). Quantification of OLIG2+/SOX10- cells in tumor sections showed that these populations were similar at P5, but significantly smaller in M-Smo tumors at P15, and that etoposide induced a significant increase in M-Smo tumors to levels similar to untreated P15 G- Smo tumors (Fig 3h). The OLIG2+/SOX10- populations in M-Smo tumors were therefore not fixed, but rather were dynamic and varied inversely with tumor size, declining over time as tumors grew, and increasing after tumor shrinkage imposed by etoposide. The OLIG2+/SOX10- fractions in G-Smo tumors, in contrast, showed a smaller decline over time and varied less across all conditions tested.
Similar range of cell fates in G-Smo and M-Smo tumors
We examined whether the differences in stem cell populations of G-Smo and M-Smo tumors were accompanied by an expansion of tumor cell fates. We previously found by lineage tracing using the 3’ Yfp sequence of SmoM2 that SmoM2-induced tumorigenesis in M-Smo mice expanded the typically neuronally-restricted Atoh1 lineage to include astrocytic and oligodendrocytic progeny (10). Comparison of Yfp expression in G-Smo and M-Smo tumors showed that tumor lineage in both genotypes included neural progenitor-like tumor cells, differentiated neurons, astrocytes and oligodendrocytes (Fig. 4a). Cells within the fibroblast and macrophage/microglia clusters did not express Yfp, and while individual Yfp+ cells were noted in the endothelial and ependymal clusters, these rare cells did not indicate a significant trend as they were not observed in more than one replicate of either genotype. Glial cells were therefore the only stromal cell types that derived from the tumor lineage in either stem cell-derived or progenitor-derived tumors.
We compared the expression of specific neural markers to assess differences in tumor cell fates. In G-Smo tumors, the Gfap-Cre transgene is expected to activate SmoM2 in a lineage that is broader than the Atoh1 lineage activated by Math1-Cre in M-Smo tumors. The Atoh1 lineage comprises rhombic-lip derived CGNP and unipolar brush cell populations, marked by Barhl1 and Eomes (26) respectively. The set of neuronal cell types with SmoM2 activation in G- Smo tumors includes these rhombic-lip derived populations, NESTIN+/ATOH1- EGL cells that also generate CGNs, and ventricular-zone derived GABAergic interneurons and progenitors, marked by expression of Ascl1, Pax3 and Pax2 (10).
To determine if G-Smo tumors contained more cells resembling ventricular zone-derived progenitors, we compared expression of Barhl1, Eomes, Ascl1, Pax3 and Pax2 in G-Smo and M-Smo tumors. We found that both genotypes consisted predominantly of Barhl1+ cells. We detected a significant genotype-specific difference within the Atoh1 lineage, with a higher proportion of Barhl1+ cells in G-Smo tumors (p < 1.0 x 10-15, two-proportions z-test) and a higher proportion Eomes+ in M-Smo tumors (p < 1.0 x 10-15) (Fig. 4 b-c). However, we did not detect statistically significant differences in the populations of Ascl1, Pax3 and Pax2 expressing cells that define the lineage of GABA-ergic interneurons derive from the ventricular zone (Fig. 4d). Although Gfap-Cre activated SmoM2 in both the rhombic lip and ventricular zone lineages, cells showing ventricular zone lineage were not expanded in G-Smo tumors. Based on the similarity of cell types showing Yfp in G-Smo and M-Smo tumors, and the similar, rare expression of ventricular zone markers, we conclude that differences between G-Smo and M-Smo tumors derive from effects of the timing of oncogenic event on cells that progress through the CGN developmental trajectory, rather than from the recruitment of interneuron lineage cells for tumor growth.
Different stromal populations in M-Smo and G-Smo tumors
To determine if M-Smo and G-Smo tumors may interact differently with stromal cells in their microenvironments, we compared gene expression in endothelial cell, myeloid cells and fibroblasts. We selected these cell types because they were the most numerous cell types outside the tumor lineage. To identify tumor-specific changes in these populations, and to distinguish tumor-specific effects common to both tumor genotypes from effects specific to individual tumor genotypes, we combined the scRNA-seq data from G-Smo and M-Smo tumors with previously obtained data from WT cerebella at P7 (10). We obtained an initial grouping of cells from tumors and WT mice by cell type using the Harmony algorithm, which co-clustered tumor cells and their most similar normal progenitors (27). Using Harmony, we generated a UMAP combining G-Smo, M-Smo and WT cells, color-coded the clusters and analyzed cluster- specific gene expression profiles; of the proliferative cell types, CGNPs and medulloblastoma cells, grouped together in a set of Barhl1+ clusters, while interneuron progenitors formed a separate group distinguished by Pax3 and Pax2 (Supplementary Fig. 1, Supplementary Table 5). We identified each stromal cell type based on gene expression (Fig. 5a, Supplementary Table 5), and isolated the endothelial, macrophage/microglial and fibroblastic populations. We then subjected the cells of each isolated cell type to a new PCA to sub-cluster each cell type.
Endothelial cells show cancer-specific changes without reflecting developmental differences between tumors
Endothelial cells showed significant differences between WT and tumor, but did not show statistically significant differences between G-Smo and M-Smo tumors. Unsupervised analysis defined 2 clusters (Fig. 5b). Cluster E0 was populated mostly by cells from WT cerebella with both tumor genotypes contributing similar fractions, while cluster E1 was populated predominantly by cells from the tumors, with no statistically significant difference in G-Smo and M-Smo contributions (Fig. 5c,d). Both clusters showed widespread expression of the endothelial markers Pecam1 and Cldn5, confirming their endothelial cell identity (Fig. 5e). Each cluster showed cluster-specific gene expression (Fig. 5f,g, Supplementary Table 6). Cluster E1, which was predominantly populated by tumor-derived endothelial cells, showed increased expression of genes likely to contribute to malignancy, including the VEGF receptor Flt1, the p-Glycoprotein Abcb1a (aka Mdr1), and the CXCR4 ligand Cxcl12 (aka Sdf1), which has been shown to promote medulloblastoma growth and glioblastoma-endothelial interactions (28-31). scRNA-seq thus detected biologically relevant differences between endothelial populations in WT cerebella and tumors, with tumor-specific effects that were shared between tumor genotypes, consistent with their overall similarity.
Developmental differences alter subsequent tumor-associated myeloid populations
The cells with myeloid characteristics, in contrast to endothelial cells, differed significantly between M-Smo and G-Smo tumors. Unsupervised analysis grouped the myeloid-like cells into 5 clusters, M0-M4 (Fig. 6a). Projection of C1qb expression confirmed that clusters M0-M3 were populated by myeloid cells (Fig 6b). In contrast, cluster M4, which was the least populated, was C1qb- and expressed Cnn3 and Meis1 (Fig 6b); this marker pattern identified cluster M4 as choroid plexus epithelial cells that have been noted to cluster with myeloid cells in other scRNA- seq analyses (32).
To identify the types myeloid cells, we defined the sets of genes up-regulated by cells within each cluster, compared to cells in the other 4 clusters (Supplementary Table 7), and then projected these genes on a UMAP along with other known markers of phenotype (Fig 6c-f). Clusters M0 and M1, which localized to clusters on one side of the UMAP, expressed Cx3cr1 which distinguished them as microglia, while clusters M2 and M3, opposite in the UMAP, showed minimal Cx3cr1, indicating that they were macrophages (Fig. 6c). Cells of cluster M0 showed highly specific expression of Sparc (Fig. 6c), identifying these cells as mature, ramified microglia (33). Cluster M1 microglia specifically expressed Mrc1, Igf1 and Wfdc17 (Fig. 6d) which have all been linked to an M2-like, anti-inflammatory phenotype (34-36). Cluster M2 macrophages specifically expressed MHCII components, including Histocompatibility 2, Class II Antigen E alpha (H2-Ea) and the Invariant Polypeptide of Major Histocompatibility Complex, Class II Antigen-associated (Cd74) as well as Il1b and Ccr2 (Fig. 6e), consistent with a pro-inflammatory M1-like phenotype (37-39). Importantly Ccr2+ macrophages have previously been shown to exert an anti-tumor effect in a medulloblastoma (39). Cluster M3 macrophages specifically expressed Cd163 and Mrc1 (Fig. 6f), also consistent with an anti-inflammatory M2 phenotype (40, 41). Myeloid cells thus resolved into microglial and macrophage populations, each with M1- like and M2-like subsets (Table 2).
Cluster
|
Cell type designation
|
distinctive markers
|
M0
|
Mature microglia
|
Cx3cr1, Sparc
|
M1
|
M2 microglia
|
Igf1, Wfdc17, Mrc1
|
M2
|
M1 macrophages
|
Ccr2, Cd74, H2-Ea, Il1b
|
M3
|
M2 macrophages
|
Cd163, Mrc1
|
M4
|
choroid plexus epithelial cells
|
Cnn3, Meis1
|
Table 2. Identification of clusters with myeloid characteristics (M0-M4)
Each cluster was differently populated by cells from G-Smo, M-Smo and WT cerebella (Fig 6g,h). The M1-like microglial clusters M0 and the non-myeloid cluster M4 comprised cells from all 3 genotypes. The M2-like microglial cluster M1 was tumor-specific and comprised cells from both G-Smo and M-Smo tumors. The M1-like macrophage cluster M2 predominantly comprised cells from M-Smo tumors. The M2-like macrophage cluster M3 comprised cells from M-Smo tumors and WT cerebella. These analyses showed that M-Smo tumors included M1-like and M2-like microglia (clusters M0 and M1), and M1-like and M2-like macrophages (cluster M2 and M3). In contrast the myeloid populations of G-Smo tumors were more limited, and the macrophage populations of cluster M2 and M3 were significantly smaller.
To confirm differential marker expression at the protein level, we visualized the protein expression of the MHCII glycoprotein coded by H2-Ea (H2-EA) and the pan macrophage/microglial marker IBA1 using immunohistochemistry (Fig. 6i). Consistent with the scRNA-seq data, comparison of H2-EA and IBA1 immunostaining in sections of each genotype showed that H2-EA+ macrophages were more prevalent in M-Smo tumors (Fig. 6j). These results confirm that protein expression patterns match transcript data from scRNA-seq and demonstrate that immunohistochemical staining for MHCII components was sufficient to distinguish between G-Smo and M-Smo tumors.
Differential cytokine expression in G-Smo and M-Smo tumors.
To consider potential mechanisms for the differences in macrophage populations in the two tumor types, we compared cytokine expression. To analyze a complete list of known cytokines and chemokines, we used the set of 232 genes tagged with the Gene Ontology term “Cytokine Activity”. For each gene, we conducted differential expression testing between G-Smo and M-Smo tumors. Macrophage Migration Inhibitory Factor (Mif) was the only differentially expressed cytokine, and was higher in G-Smo tumors compared to M-Smo tumors (Fig. 6k). Intriguingly, tablis a ligand for CD74 (42), and the Cd74+ population was markedly lower in G- Smo tumors. Prior studies in melanoma have shown that intercellular communication through the MIF-CD74 interaction is immunosuppressive, and that blocking MIF-CD74 interaction increases tumor-associated M1 macrophages (43). Based on the inverse correlation of MIF and CD74 in our tumor model and on the prior studies of MIF-CD74 interaction in melanoma, we propose that MIF functions similarly in medulloblastoma to reduce the infiltration of CD74- expressing Cluster M2 macrophages, and acts more effectively in G-Smo tumors, which have higher Mif expression.
Developmental differences influence tumor fibroblast populations
The fibroblast populations in G-Smo, M-Smo and WT, similarly formed 3 clusters in unsupervised analysis, with significant differences between G-Smo and M-Smo tumors (Fig. 7a,b). Cluster F0 was WT-specific, Cluster F1 was WT and M-Smo-specific, and Cluster F2 was tumor-specific and populated with cells from both G-Smo and M-Smo tumors (Fig. 7b,c). Each cluster demonstrated different sets of cluster-specific gene expression patterns (Fig. 7d-f, Supplementary Table 8). We noted genotype specific effects on genes related to retinoid signaling, with Fabp5 expressed predominantly in WT fibroblasts and Rbp4 and Crabp2 expressed by fibroblasts from both WT cerebella and M-Smo tumors, but not by fibroblasts from G-Smo tumors (Fig. 7g). These differential patterns show that non-cell autonomous effects of the timing of oncogenesis occur in fibroblastic stroma as well as in tumor-associated myeloid cells.