H3.3 G34W Promotes Growth and Impedes Differentiation of Osteoblast-Like Mesenchymal Progenitors in Giant Cell Tumor of Bone

Giant cell tumor of bone (GCT) is a locally aggressive yet rarely metastasizing neoplasm that accounts for 20% of benign bone tumors (1). Tumors are composed of three major compartments: multinucleated osteoclast-like giant cells, macrophage-like monocytes, and stromal cells. Although the primary cause of morbidity in GCT is bone resorption by monocyte-derived giant osteoclasts, the actual neoplastic component is formed by cells from the poorly defined mesenchymal stromal compartment, which bears heterozygous somatic mutations in H3F3A, one of two genes encoding histone 3 variant H3.3 (2). These mutations lead to substitution of glycine 34 to tryptophan or, more rarely, to leucine (G34W/L, >92% of cases; ref. 2). Unlike H3.3 glycine 34-to-arginine or valine substitutions (H3.3 G34R/V) that arise in brain tumors and co-occur with partner TP53 and ATRX mutations (3, 4), H3.3 G34W/L mutations are the only recurrent molecular alteration in GCT, potentially driving this destructive bone disease (2, 5).

Although the almost universal occurrence of H3.3 G34W/L mutations in GCT (92%) strongly suggests a causal role, how these mutations drive tumorigenesis in the bone remains unknown. Moreover, the cellular processes that result in multinucleated giant cell recruitment have yet to be elucidated. Other H3 mutations, highly recurrent in a wide range of cancers, affect residues that are direct targets of post-translational modifications (PTM). These lysine to methionine substitutions on K27 (K27M) and K36 (K36M) of the H3 tail impair enzymatic activities of the corresponding methyltransferases (6). This results in genome-wide changes in the deposition and distribution of major antagonistic repressive (H3K27me3) or activating (H3K36me3, H3K36me2, H3K27ac) chromatin marks, and redistribution of the repressive polycomb receptor complex (PRC) 1 and 2 to promote aberrant transcription and tumorigenesis (7–13). In contrast, H3.3G34 mutations are less well understood. They exclusively affect noncanonical H3.3 histones, suggesting an underinvestigated H3.3-dependent role in their tumorigenesis. Indeed, deposition of H3.3 is replication-independent and follows a specific temporospatial pattern that may be relevant to GCT formation. H3.3 marks active euchromatin at enhancers and gene bodies, but also represents the main H3 variant in silent pericentromeric and telomeric heterochromatin (14). The presence of mutant H3.3G34 could therefore affect chromatin locally in any of these key genomic regions, contributing to tumorigenesis. Although H3.3G34 mutations occur on a residue that does not undergo PTMs, they affect the adjacent K36 residue, leading to decreased H3.3K36me3 (6) and potential reciprocal increase of H3.3K27me3 on the mutant H3.3 tail (15). Structural studies suggest that H3.3G34 mutants cause steric hindrance of the catalytic groove of SETD2, an H3K36 trimethyltransferase (16). This hypothesis is further supported by the finding that loss-of-function SETD2 mutations are mutually exclusive with H3.3G34R/V mutations in gliomas, suggesting a convergence of effects (17). H3.3G34R/V have also been suggested to impede H3.3K36 access by KDM2A, an H3K36me2 demethylase (18), or inhibit members of the KDM4 H3K9/K36 demethylase family, resulting in increased H3K36me3 and H3K9me3 at select loci (19).

Altogether, there is no consensus on the mechanism of action of H3.3G34R/V/W/L mutations. The cell type–specific mutation patterns indicate a strong context dependence that remains poorly understood, and that has so far hampered accurate tumor modeling. GCT models have been mainly restricted to overexpression systems that increase mutant histone expression beyond that observed in tumors (15, 20). However, accounting for dosage is key, because oncohistones contribute to less than 20% of the total H3 pool (6). To address these challenges, we generated H3.3 G34W (referred to hereafter as G34W) isogenic models using CRISPR–Cas9 in patient-derived GCT stromal cells as well as patient-derived orthotopic xenografts (PDOX). We profiled the epigenome, transcriptome, and secreted Golgi proteome of cellular models, and generated transcriptomic data of GCT primary tumors and PDOX at single-cell resolution. Our combined modeling and profiling efforts shed light on the mechanism by which G34W transforms cells to promote tumorigenicity and enables interactions with the tumor microenvironment (TME) of the bone to foster disease.

To determine the role of the H3.3 G34W oncohistone on GCT tumorigenesis in pure neoplastic stromal cells, we used CRISPR–Cas9 to target the H3F3A G34W mutation (referred to hereafter as G34W) in three immortalized (Im-GCT) stromal cell lines established from primary tumors (Supplementary Table S1; Supplementary Fig. S1A). We generated clones carrying loss-of-function insertions or deletions on the mutant allele (G34W-KO) or clones where the mutation was corrected to wild-type (G34-WT), and confirmed absence of oncohistone expression using a previously validated G34W-specific antibody (ref. 21; Fig. 1A; Supplementary Fig. S1B). G34W-KO and G34-WT clones showed significantly decreased proliferation (Fig. 1B; Supplementary Fig. S1C) and colony-forming ability in vitro (Fig. 1C; Supplementary Fig. S1D) in the isogenic cell lines. Because G34W-KO and G34-WT clones behaved similarly in in vitro assays, they are hereafter referred to as edited clones. Orthotopic tibial and subcutaneous tissue implantation of luciferase-tagged G34W and edited Im-GCT-4072 cells showed a weekly increase of bioluminescence and tumor formation only in G34W-mutant cells (Fig. 1D and E; Supplementary Fig. S1E–S1H). Notably, no tumors formed when two distinct edited clones were implanted in tibia and subcutaneously, with more than 1 year follow-up post-implantation (Fig. 1E; Supplementary Fig. S1E).

To evaluate whether reintroduction of H3.3 G34W in edited clones could rescue tumor formation, we overexpressed (O/E) H3.3 G34W and H3.3WT constructs in two Im-GCT-4072 G34W-edited clones (Supplementary Fig. S1I). Edited lines O/E G34W formed tumors, albeit with reduced penetrance and increased latency relative to parental G34W lines (Fig. 1E; Supplementary Fig. S1J). Histology was similar to parental and unedited lines (Fig. 1F and G; Supplementary Fig. S1K), reinforcing the oncogenic role of the G34W mutation in GCT.

Our xenograft model was aggressive, prone to disseminate locally and distally (Supplementary Fig. S1L and S1M). Orthotopic tibial tumors exhibited osteolytic properties seen in GCTs (Fig. 1F and G). Notably, murine osteoblast and osteoclast cells were recruited to G34W stromal cells adjacent to bone tissue, characteristic of a reactive bone-remodeling process (Fig. 1G). Eponymous multinucleated osteoclasts were evident within the neoplastic stromal cell mass extruding outside bone and were most abundant in regions of xenograft tumors containing differentiated stromal cells (Fig. 1F). Taken together, these results show that G34W is required for tumor formation and promotes recruitment of multinucleated osteoclasts in orthotopic xenograft tumors.

We next characterized gene-expression changes associated with G34W using RNA sequencing (RNA-seq; Supplementary Table S1) and H3K27ac chromatin immunoprecipitation sequencing (ChIP-seq), a mark associated with transcription and active regulatory elements (22). Removal of G34W led to a global effect on the transcriptome, as shown by the segregation of samples by G34 status in all three isogenic lines using principal component analysis (PCA; Fig. 2A; Supplementary Fig. S2A). Despite their different genetic backgrounds, with a large expected individual variability, we identified common dysregulated pathways across the three GCT cell lines. Pathway analysis of upregulated genes in G34W revealed extracellular matrix (ECM) organization genes (e.g., COL6A1/3, EMILIN2, SOX9), a function associated with cartilage and bone-forming chondrocytes and osteoblasts, enriched in two isogenic lines (Fig. 2B and C; Supplementary Fig. S2B and S2C; Supplementary Tables S2 and S3). This transcriptional upregulation of ECM organization genes was associated with H3K27ac deposition at the promoters of these genes, indicative of epigenetic activation (Supplementary Fig. S2D and S2E). Notably, H3K27ac profiles displayed distinct genome-wide patterns between G34W and edited lines (Fig. 2D; Supplementary Fig. S2F). On the other hand, pathways related to actin filament-based processes and muscle contraction were strongly depleted in two of three G34W cell lines (e.g., ACTA2, CNN1, LMOD1, TAGLN; Fig. 2B and C; Supplementary Fig. S2B and S2C; Supplementary Tables S2 and S3). Intersection of commonly dysregulated genes across the three isogenic cell lines revealed 27 genes, all of which were downregulated in G34W-mutant cells. These included genes associated with muscle contraction (ADORA1, MYL1, PKP2, TNNT2; Fig. 2E and F; Supplementary Fig. S2G and S2H), suggesting that G34W in neoplastic stromal cells is potentially affecting muscle differentiation, a known lineage for mesenchymal progenitor cells of bone marrow origin (23). To further validate these findings, we obtained expression data from a previously reported dataset of primary G34W cell lines (n = 11 G34W vs. n = 6 WT; ref. 20), which confirmed downregulation in G34W of pathways related to muscle contraction (Supplementary Table S3). Altogether, these results indicate that G34W downregulates genes involved in muscle functions, possibly impairing further differentiation in a mesenchymal progenitor committed to myogenesis.

The G34W phenotypes we observed in driving tumor formation and altering the transcriptome are likely to originate at the epigenome. Thus, we first characterized the in cis effects of G34W on histone marks using histone mass spectrometry (hMS) on four primary G34W GCT cell lines (Supplementary Table S1), comparing patterns between mutant and WT peptides within each line. Mutant G34W contributed to only 2.6% of the total H3 pool despite accounting for 25.7% of total H3.3 (Fig. 3A). As expected (6, 15), we observed significant H3.3K36me3 decrease and concomitant H3.3K27me3 gain on the mutant peptide (Fig. 3B; Supplementary Fig. S3A and S3B), reflecting the in cis effects of G34W on the minority of H3.3-mutant histones within the overall H3 pool. Next, we compared isogenic G34W and edited lines and observed a modest but significant global decrease in H3K27me3 levels on both H3.3 and canonical H3.1/2 in G34W lines (Supplementary Fig. S3C), whereas no change was observed in overall H3K36me3 levels (Supplementary Fig. S3D). This suggests that, even though H3.3 G34W accounts for less than 3% of the total H3 pool, this mutant histone exerts an effect in trans, beyond mutant H3.3, that can be detected for the repressive H3K27me3 mark.

To further assess the effects of G34W genome-wide, we profiled H3.3, G34W, H3K36me3, and H3K27me3 deposition by ChIP-seq in isogenic G34W and edited lines from Im-GCT-4072 (Supplementary Table S1). G34W histones tracked with H3.3 in the parental line at H3.3 peaks and genome-wide (Fig. 3C; Supplementary Fig. S3E). Consistent with the expected distribution of H3.3 in euchromatin (14), G34W histones tracked with highly expressed genes with the largest H3.3 abundance at their transcriptional start and end sites (Supplementary Fig. S3F and S3G). Next, to identify G34W-enriched loci displaying in cis H3.3K36me3 loss, we compared total H3K36me3 deposited at H3.3 peaks between isogenic lines, as there are no available antibodies specific for H3.3K36me3. On the basis of the known effects of this mutation on H3.3K36me3, our expectation was to mainly observe areas with decreased H3K36me3 deposition in G34W compared with edited cells. Interestingly, we observed both gains and losses of H3K36me3 at H3.3 peaks when comparing the isogenic pairs (Fig. 3D–F). In G34W-mutant cells, areas which gained H3.3 deposition had concurrent H3K36me3 gain (Fig. 3E and G). Notably, areas that gained H3K36me3 also had H3K27me3 loss (Fig. 3E and G). Conversely, genomic regions that lost H3.3 deposition had concurrent H3K36me3 loss and H3K27me3 gain (Fig. 3F and G). Together, these results suggest a genome-wide H3.3 redistribution upon G34W removal (Fig. 3H; Supplementary Fig. S3H), implicating an H3.3 G34W–dependent chromatin remodeling process in GCT beyond initial sites of H3.3 deposition. This H3.3 redistribution may result from the spread of the H3K27me3 mark upon G34W-induced loss of the antagonistic H3.3K36me3 mark. Spread of H3K27me3 may, in turn, induce gene silencing and secondary eviction of H3.3 to other genomic areas.

Our investigation of H3K27me3 deposition in the isogenic pairs further identified global changes associated to G34W (Supplementary Fig. S4A–S4D). This mark was redistributed genome-wide: in G34W-mutant cells, H3K27me3 was lost predominantly in intergenic regions (P < 0.001; χ² test), whereas it was equally gained across genomic compartments, including promoter and intragenic regions (Fig. 4A; Supplementary Fig. S4E). In addition, G34W reexpression in edited lines also resulted in a modest shift of H3K27me3 deposition from intergenic to promoter and genic regions (Supplementary Fig. S4F). Redistribution of H3K27me3 was further confirmed by SUZ12 distribution, a core component of the PRC2 complex that catalyzes this repressive mark, as it mirrored H3K27me3 deposition in G34W-mutant cells (Fig. 4B; Supplementary Fig. S4G). Notably, genes with a chromatin state consistent with an in cis G34W effect (H3.3K36me3 loss/H3K27me3 gain) were enriched in pathways related to actin–myosin contractile functions (Supplementary Fig. S4H and S4I; Supplementary Tables S2 and S3), whereas genes associated with actin filament–based processes showed similar chromatin patterns (Fig. 4C and D; Supplementary Fig. S4J). These results are in keeping with our transcriptomic data, where downregulated genes between isogenic pairs are mainly involved in muscle function (Fig. 2B and E; Supplementary Fig. S2B and S2G; Supplementary Table S3).

To assess whether PRC2 redistribution (gain and loss of H3K27me3) in G34W lines is accompanied by secondary chromatin remodeling, we profiled broad intergenic marks with known mutually exclusive deposition patterns relative to H3K27me3, namely H3K36me2 (9) and H3K9me3 (24). An interesting interplay of these marks in intergenic regions was apparent: loci that gained H3K27me3 in G34W lines consistently exhibited loss of H3K36me2 (Fig. 4D and E), whereas those that lost H3K27me3 gained either H3K9me3 or H3K36me2 (Fig. 4E and F). The TNNI1/TNNT2 locus exemplifies H3K27me3 gain at the expense of H3K36me2 in G34W lines, specifically repressing the TNNI1 and TNNT2 genes encoding subunits of the troponin complex which regulates myofibril contractility (Fig. 4D). In contrast, the locus comprising the negative regulator of muscle differentiation, BMP2, displays H3K36me2 gain at the expense of H3K27me3 (Fig. 4F). The LARGE1 locus, which encodes for a glycosyltransferase important for muscle function (25), shows spread of a heterochromatin H3K9me3 domain into territories occupied by H3K27me3 and H3K36me2 in edited lines, concurrent with repression in G34W (Fig. 4F). Similarly, the spread of H3K9me3 into the MYH cluster likely results in the repression of sarcomeric myosin proteins encoded in the locus (Supplementary Fig. S4K).

Collectively, these results indicate that deposition of G34W parallels that of wild-type H3.3, and that G34W in cis effects may promote PRC2/H3K27me3 silencing of H3K36me3-depleted nucleosomal substrates, leading to subsequent H3.3 eviction, redistribution, and gene silencing. Recruitment of PRC2 to genic regions may occur at the expense of weaker intergenic substrates where loss of H3K27me3 is replaced by the spread of neighboring antagonistic marks, such as H3K36me2 and H3K9me3. These epigenetic changes are coupled to transcriptional changes affecting mesenchymal cell state, with muscle contractile genes likely constituting early targets of G34W in mesenchymal progenitors at a differentiation stage yet to be determined, permissive to the generation of neoplastic stromal cells in GCT.

GCT stromal cells bearing the G34W mutation are ill-defined and have been alternately described as possessing properties of mesenchymal stem cells, fibroblasts, or pre-osteoblasts based on in vitro differentiation assays and histologic studies (26). To better define all cellular components and stromal cell heterogeneity, we performed single-cell RNA-seq (scRNA-seq) on four G34W GCTs (Supplementary Table S1). Clustering analysis revealed 17 transcriptionally distinct clusters (Supplementary Fig. S5A), which we labeled as lymphoid, myeloid, endothelial, and putative stromal cells (Supplementary Fig. S5B) based on similarity to cell types within the Human Primary Cell Atlas (HPCA; Supplementary Fig. S5C) and known markers (Supplementary Fig. S5D; Supplementary Table S4). In contrast to other clusters which contained cells from all tumors, putative stromal cell clusters segregated by patient (Supplementary Fig. S5B), indicating increased variation within this compartment likely reflecting divergent clonal evolution. Putative stromal cell clusters were significantly enriched for the G34W mutation (Supplementary Fig. S5E) and showed similar pathway enrichment as isogenic G34W stromal cells (e.g., ECM organization; (Supplementary Fig. S5F; Supplementary Tables S3 and S4). Importantly, stromal clusters across patient tumors displayed the highest enrichment score for the transcriptomic G34W signature derived from our independent isogenic model (Fig. 5A). We therefore conclude that these clusters represent the neoplastic stromal compartment, and that the transcriptional changes induced by G34W in our cellular models are recapitulated in primary tumors.

As expected, the mononuclear stromal component was detected in all samples. This prompted us to use Harmony data integration to reduce the impact of technical and biological intertumoral stromal cell heterogeneity, and to query potential neoplastic stromal subpopulations shared by all samples. Clustering analysis identified four clusters within the stromal compartment (S1A, S1B, S2, and S3; Fig. 5B; Supplementary S5G–S5I), with osteopontin (SPP1) and smooth muscle actin (ACTA2) as the most discriminative markers between stromal clusters. On the basis of differential gene expression analysis and expression modules of genes correlated with SPP1 or ACTA2 (Fig. 5B; Supplementary Fig. S5H; Supplementary Table S4), clusters S1A and S1B likely represent a continuum of the same population, leading us to consider three populations: S1 (S1A, S1B), S2, and S3. S1 cells expressed genes associated with osteoblast and chondrocyte functions (e.g., SPP1, IBSP, MMP13), whereas S3 cells expressed markers of myofibroblasts (e.g., ACTA2, TAGLN, POSTN) and S2, with fewer differentially expressed genes, was considered an intermediate between both states (Fig. 5C; Supplementary Fig. S5I). We confirmed expression of both osteopontin (SPP1) and alpha smooth muscle actin (ACTA2/α-SMA) in patient tumors (Fig. 5D). In isogenic lines, G34W-mediated effects converge on epigenetic/transcriptional downregulation of muscle contractile genes and upregulation of ECM pathways. The S3 population is defined by similar pathways suggesting that cells within this subgroup are myofibroblast-like and committed to myogenesis, but are restricted by G34W from further differentiation. Indeed, lineage inference analysis identified two possible trajectories across the four populations: from S1A to S1B cells, or from S1A to S3 cells (Fig. 5E; Supplementary Fig. S5J). This suggests a continuum of progenitor stromal states, with G34W-mutant cells either transitioning between early S1A to S1B states or following a trajectory from an early S1A toward a possible myofibroblast-like progenitor S3 state.

To identify the cell of origin of the neoplastic stromal cells, we next mapped the transcriptional profiles of S1 to S3 populations to two recently published scRNA-seq murine bone marrow stroma datasets (refs. 27, 28; Fig. 5F and G). The control endothelial cell cluster mapped, as expected, to reference EC–arterial/arteriolar (27) and V1 arterial (28) cells. Strikingly, S1, S2, and S3 populations strongly mapped to the osteoprogenitor subset in the Baryawno dataset (27). In the Tikhonova dataset (28), they all mapped to the osteo-lineage O1 subset of COL2.3⁺ osteoblasts, which likely also represents an osteoprogenitor cell type because we show that O1 maps robustly to the osteoprogenitor subset and vice versa (Supplementary Fig. S5K). S1 and S2 populations mapped more strongly than S3 cells to the less-differentiated perivascular 4 (P4 Spp1^hi/Ibsp^hi) or “preosteoblasts” reference subsets (refs. 27, 28; Fig. 5E and F), consistent with their earlier position in the neoplastic trajectory. These findings indicate that GCT neoplastic stromal cells resemble osteoprogenitors or a closely related cell type, with some cells (S3) progressing to display contractile features (Supplementary Fig. S5L). We further verified that the three stromal subpopulations are also present in four PDOX tumors using single-nuclei RNA-seq (snRNA-seq; Supplementary Fig. S6A). Using the annotated clusters from patient tumors (Fig. 5B), we classified cells from PDOX tumors and found that they were formed predominantly by S3 cells (Supplementary Fig. S6B and S6C). Consistent with this result, tibial xenograft tumors displayed high expression of α-SMA/ACTA2, a marker of S3 cells, but low osteopontin expression (marker of S1 cells; Supplementary Fig. S6D). Enrichment of the S3 subtype in the xenograft models may explain why they exhibit a more malignant phenotype and pathology compared with relatively benign GCT.

To validate the myofibroblast differentiation potential of GCT stromal cells, we induced differentiation using TGFβ1, l-ascorbic acid, and platelet derived growth factor (PDGF-AB) in isogenic GCT lines and control human mesenchymal stromal cells (hMSC). Upon induction, GCT cells expressed high levels of myofibroblast markers including stress fiber-associated protein calponin 1 (CNN1; Fig. 5H and I) and α-SMA (Supplementary Fig. S6E). Notably, although both G34W and edited cells were able to differentiate into myofibroblasts, edited cells already exhibited higher expression of calponin 1 at baseline relative to G34W lines (Fig. 5H and I). This is consistent with transcriptomic data showing enrichment of actin-myosin contractile genes following G34W editing (Fig. 2B and C; Supplementary Fig. S2B and S2C). Together, these findings reinforce that neoplastic stromal cells transition from progenitor S1 to S3 states, where epigenetic reprogramming by G34W appears to impede differentiation into a terminal myofibroblast-like state.

In GCT, the formation of giant multinucleated histone WT osteoclasts from monocytic precursors (29) is an essential contributor to the pathogenesis of these tumors and is thought to be mediated by G34W-dependent secretion of factors such as RANKL (TNFSF11; ref. 30). In contrast to previous findings (20), neither RANKL nor its decoy receptor OPG (TNFRSF11B) were differentially expressed between G34W and edited clones (Fig. 6A; Supplementary Fig. S7A), suggesting that RANK/RANKL signaling, although required for osteoclastogenesis, is independent of G34W. We thus investigated G34W-mediated aberrant TME interactions that contribute to the giant cell phenotype. Using additional markers to build upon HPCA-based classification, we first elucidated cell types present in the GCT myeloid compartment (Fig. 6B; Supplementary Figs. S5C, S5D and S7B; Supplementary Table S4). We show that the putative osteoclast cluster is enriched for osteoclast differentiation and bone resorption pathways (Supplementary Fig. S7C). To investigate stromal cell ligands that promote osteoclastogenesis, we first determined the cell type giving rise to osteoclasts in GCT, using a combination of trajectory inference methods (Fig. 6C; Supplementary Fig. S7D–S7G). We identified a bifurcation event where monocytes can either differentiate into macrophages or give rise to osteoclasts through a highly cycling, preosteoclast intermediate. These findings confirm the proposed monocytic origin for osteoclasts in GCT (30) and describe a previously unappreciated fate for monocytes toward macrophage differentiation in GCT.

To identify ligands secreted by G34W stromal cells, we profiled the secreted proteome by isolating Golgi apparatus of isogenic GCT lines and performing mass spectrometry (MS; Fig. 6D; Supplementary Table S1), and identified significantly differentially secreted proteins between G34W and edited lines (Supplementary Fig. S7H; Supplementary Table S5). Ligand–receptor gene interactions were inferred between stromal cells and the myeloid cell compartment based on curated protein interaction databases (ref. 31; Fig. 6E; Supplementary Fig. S7I; Supplementary Table S6). Filtering these predicted ligand–receptor interactions for proteins differentially secreted by G34W lines, we observed stromal cell–specific expression of ECM ligands such as collagens (COL6A1/3, COL5A2) and proteoglycans (BGN) that are predicted to interact with specific integrin receptors on monocyte/macrophage (ITGB2) and osteoclast (ITGAV) cells (Fig. 6E; Supplementary Fig. S7I). Expression of COL6A1/3 and BGN was further confirmed in PDOX tumors by snRNA-seq (Supplementary Fig. S7J), and by IHC in primary GCT samples as well as PDOX tumors (Fig. 6F; Supplementary Fig. S7K). Interestingly, when intersecting the Golgi secretome and differentially expressed genes in S1 to S3 populations, the S3 cells were enriched for expression of G34W-specific secreted proteins (Fig. 6G; Supplementary Fig. S7L), potentially accounting for the overrepresentation of this subtype in our aggressive PDOX model. Indeed, 6 of the 24 genes specifically expressed by S3 cells are found in the ECM secretome, suggesting that S3 stromal cells are responsible for G34W-driven ECM remodeling, in line with our finding of increased osteoclasts in proximity to the differentiated stromal compartment of a xenograft tumor (Fig. 1F). Taken together, our results indicate that G34W stromal cells resemble osteoblast progenitors stalled in differentiation at a myofibroblast-like progenitor population (S3). These S3 cells secrete ECM remodeling proteins, promoting bone destruction by acting locally on myeloid cells in the TME. Therefore, G34W drives both components of GCT pathogenesis: it sustains the neoplastic transformation of the mononuclear stromal cells, as shown by our in vivo orthotopic model, and simultaneously enables the recruitment and formation of pathologic giant osteoclasts that largely contribute to the morbidity of this tumor.

We comprehensively characterized the epigenetic, transcriptomic (including at the single-cell level) and Golgi proteomics profiles of G34W GCT using primary tumors and isogenic models of tumor-derived cell lines. G34W has been associated with tumorigenesis in previous studies using xenotransplantation of the full tumor (32, 33) or transient siRNA knockdown of G34W in GCT-derived cell lines (34). We show that G34W removal is sufficient to prevent tumor formation in vivo, indicating that this mutation is necessary for GCT tumorigenesis. Although G34W occurs on a residue that does not undergo PTM, we show that tumorigenesis is likely the consequence of a global epigenetic remodeling process initiated by the in cis effects of this mutation in the stromal cell of origin (Fig. 7A). Indeed, we posit that the earliest event following H3.3K36me3 loss on G34W histones is the deposition of H3K27me3 in active genic regions normally enriched for H3.3 as well as G34W (representing 25% of total H3.3). H3K36me3 loss in genic regions may create improved substrates for the PRC2 complex and result in redistribution of H3K27me3 from lower affinity intergenic regions to promoter and genic regions, in turn promoting gene silencing and eviction of H3.3. Secondary and more widespread changes as a result of H3.3 redistribution may then induce transcriptional changes in the absence of direct epigenetic remodeling.

In intergenic regions, H3K27me3 loss is replaced by H3K36me2 or H3K9me3 in areas vacated by PRC2 (Fig. 7A). These compensatory mechanisms to H3K27me3 loss are likely to have distinct effects. While H3K36me2 enables repressive DNA methylation to be deposited through DNMT3A recruitment (35), H3K36me2-marked regions remain permissive to active regulatory states marked by H3K27ac (e.g., BMP2; ref. 36). In contrast, the constitutive heterochromatin H3K9me3 mark is strongly repressive (37), and its imposition into facultative heterochromatin at a few key loci (e.g., LARGE1, MYH) in G34W-mutant cells may result in undue silencing of genes associated with contractile functions.

The net effect of G34W-mediated epigenetic remodeling affects mesenchymal cell lineage commitment. Indeed, our transcriptomic and epigenetic data converge on a stalled osteoblast-like progenitor population, capable mainly of differentiation along a neoplastic trajectory toward a myofibroblastic cell state. Moreover, at single-cell resolution, G34W GCT stromal cells comprise distinct cell populations related by a neoplastic lineage trajectory with SPP1⁺ (S1) tumor cells at the origin leading to ACTA2⁺ myofibroblast progenitor (S3) tumor cells (Fig. 7B). A terminally differentiated myofibroblast cell state is blocked by G34W, as suggested by lower expression of myofibroblast markers such as calponin 1 in G34W-mutant cells. This block can be overcome in vitro using strong differentiation factors (e.g., TGFβ1), an observation which may have implications for GCT treatment, especially in recurring tumors or those where complete surgical resection is not feasible.

Notably, our data suggest that G34W induces the secretion of ECM-remodeling factors in the neoplastic S3 myofibroblast progenitor cells, including collagens and proteoglycans which interact directly with integrin receptors on myeloid cells (Fig. 7B). The presence of neoplastic myofibroblasts has been reported in primary and metastatic GCT using IHC and ultrastructural studies (38, 39). Moreover, scRNA-seq studies of several diseases where stroma plays a major role in pathogenesis have identified ACTA2+ stromal cells expressing ECM ligands (COL6A3, BGN) consistent with the S3 population (40–42). In GCT tumorigenesis, the initiating factor is the oncogenic G34W mutation which leads to a persistent activated progenitor state incapable of further transitioning toward a terminally differentiated myofibroblast state. The presence of myeloid cells in the bone microenvironment niche enables the recruitment and syncytia of osteoclasts through a G34W-mediated ECM-remodeling process enriched in S3 ACTA2⁺ cells. Notably, secreted collagen VI (COL6A1/3) maintains mechanical stiffness within the ECM, a function implicated in the activation of the mechanosensitive Ca²⁺ channel TRPV4, which is expressed on the plasma membrane of large osteoclasts and regulates terminal osteoclast differentiation (43, 44). We previously showed that H3 wild-type giant cell lesions of the jaw (GCLJ), a disease that is histologically and radiologically similar to GCT, carry TRPV4 gain-of-function mutations (45), suggesting possible convergence of effects between G34W-mutant GCT and TRPV4-mutant GCLJ on osteoclastogenesis and pathologic features.

In conclusion, we show that G34W is necessary to drive the two major pathologic features of GCT: the destruction of bone and the maintenance of proliferating osteoprogenitors. G34W in a neoplastic stromal ACTA2⁺ population mediates the secretion of factors that recruit osteoclasts within the bone TME, resulting in destruction of bone. Few transcriptional changes may be initially needed to maintain this aberrant progenitor state, and these are likely due to G34W's in cis effects, initiating the global epigenetic remodeling events we observed. This chromatin remodeling, whether a direct effect or indirectly mediated by the mutation, reflects a neoplastic maintenance of a progenitor state, which promotes GCT formation and helps maintain the neoplastic progenitor state. As these changes are potentially reversible, future therapies targeting the epigenome may be of benefit in GCTs.

Additional materials and methods can be found in the Supplementary Information.

Raw and processed sequencing data for bulk RNA-seq, scRNA-seq, snRNA-seq, and bulk ChIP-seq data were deposited into the Gene Expression Omnibus (GEO) under accession code GSE149211.

Further information and requests for resources/reagents should be directed to and will be fulfilled by the lead contact, Nada Jabado.

Protocols for this study involving collection of patient samples were approved by the Research Ethics and Review Board of McGill University (Monteral, Quebec, Canada) and affiliated Hospitals Research Institutes. Written informed consent was obtained from all research participants.

Fresh tumor specimens were collected in DMEM, washed twice in PBS, minced physically with surgical blades, and dissociated in Collagenase-Dispase (Sigma-Aldrich) with incubation at 37°C in 5% CO₂ and agitation for 30 minutes. The resulting cell pellet was treated with NH₄Cl to remove red blood cells, washed twice in DMEM (with 4.5 g/L glucose, l-glutamine, phenol red, by Wisent) containing 10% FBS and penicillin/streptomycin, plated in T-75 flasks in the same medium, and passaged at 75% confluency. DNA from various passages was assessed for H3F3A G34W mutation allele frequency by droplet digital PCR (details available on request). We thank Dr. David Allis for generously sharing the primary GCT-2611 cell line. Three lines were immortalized using hTERT retrovirus (see Supplementary Information for details). All lines were tested monthly for Mycoplasma contamination (MycoAlert Mycoplasma Detection kit by Lonza), and STR fingerprinting was regularly performed.

CRISPR–Cas9 editing (46) was performed using pSpCas9(BB)-2A-Puro (PX459 V2.0) plasmid, a gift from Feng Zhang (Addgene #48139, RRID: Addgene_62988). A single-guide RNA (sgRNA; IDT) targeting H3F3A G34W was cloned into the plasmid, and a single-stranded donor oligonucleotide (ssODN; IDT) template with H3F3A WT sequence (sequences in Supplementary Information) was transfected using Lipofectamine 3000 (Thermo Fisher Scientific). Transfected cells were selected using puromycin for 72 hours and then sorted into 96-well plates. Single-cell clones were expanded and screened for editing events (insertions, deletions, or repair to wild-type) by Sanger sequencing, and confirmed through targeted deep sequencing (Illumina MiSeq) and immunoblotting using a H3.3 G34W–specific antibody (1:500 dilution, RevMab Biosciences 31–1145–00, RRID:AB_2716434, see Supplementary Information). We refer to H3F3A^+/− and H3F3A^+/+ clones as “edited” and not “corrected” because H3F3A^+/− clones are hemizygous and the initial presence of G34W may have caused irreversible transcriptional/epigenetic changes.

Cells were seeded as 1,500 cells/well into 96-well plates (five technical replicates per line), and nuclei were stained with NucLight Rapid Red Cell Reagent (1:500 dilution, Essen Bioscience). Cells were imaged using IncuCyte ZOOM System real-time instrumentation (Essen Bioscience) every 2 hours (10× magnification), and images analyzed after 120 hours with IncuCyte ZOOMTM 2015A software.

Cells were seeded as 250 cells/well in 6-well plates. After 3 weeks, colonies were fixed with 100% methanol for 20 minutes, and stained with 0.5% crystal violet for 30 minutes. Colonies (containing >50 cells) were manually counted under an inverted microscope.

Cells were seeded as 6,000 cells/chamber of 8-well chamber slides in standard medium. At approximately 70% confluency, differentiation was induced with media containing 5 ng/mL TGFβ1 (R&D Systems 240-B-002/CF), 30 μmol/L l-Ascorbic acid 2-phosphate (Sigma A8960), and 5 ng/mL PDGF-AB (PeproTech 100–00AB; ref. 47). Fresh differentiation media were added twice weekly for 3 weeks. Immunofluorescence staining was performed using primary antibodies anti-α-smooth muscle actin (1:200, Abcam ab5694, RRID:AB_2223021) and anti-calponin 1 (1:200, Abcam ab46794, RRID:AB_2291941; see Supplementary Information for details).

The complete workflow for histone extraction, LC/MS, and data analysis was described in detail (refs. 48, 49; see Supplementary Information for brief description).

Cells were lysed using a 7 mL Dounce stainless tissue grinder, and a 2 mol/L sucrose solution was added to homogenates to obtain final 1.7 mol/L sucrose concentration. Golgi membranes were isolated by isopycnic centrifugation using discontinuous sucrose gradients (50), and processed for MS analyses as described previously (51). See Supplementary Information for full description of workflow.

All mice were housed, bred, and subjected to listed procedures according to the McGill University Health Center Animal Care Committee (Protocol #2010-5684) and in compliance with Canadian Council on Animal Care guidelines. Mice were monitored weekly and euthanized when xenograft tumor volume reached 2,000 mm³ or immediately at clinical endpoint when recommended by veterinary/biological services staff.

pSMAL-GFP-Luc lentiviral vector was a gift from Dr. Kolja Eppert (McGill University). Ninety-six hours after transduction with virus, GFP-positive cells were sorted using FACS.

1.2–1.5 × 10⁶ cells were prepared as a single-cell suspension in 200 μL PBS-50% Matrigel and implanted subcutaneously into the left flank of 8- to 12-week-old male immunodeficient NSG mice (NOD.Cg-Prkdc^scid Il2rg^tm1Wjl/SzJ, RRID:IMSR_JAX:005557). Tumor formation was monitored weekly using bioluminescence imaging and palpation starting 4 weeks postimplantation.

Seventy percent ethanol was used to clean the leg of 8- to 12-week-old male immunodeficient NRG mice (NOD.Cg-Rag1t^m1Mom Il2rg^tm1Wjl/SzJ, RRID:IMSR_JAX:007799). A predrilled hole was created in the tibial plateau using rotating movements with a 25G needle. Three microliters of PBS-20% Matrigel containing 1.8–2 × 10⁵ cells were implanted into the tibial diaphysis with a Hamilton syringe (27G needle). Tumor formation was monitored weekly using bioluminescence imaging and tumor volume measured with digital calipers.

Bioluminescence imaging was performed as described previously (52) using IVIS 100 system (Caliper). Anesthetized mice were placed in imager 7 minutes post–intraperitoneal injection of d-luciferin (150 mg/kg).

Xenograft tumors were fixed in 10% buffered formalin for 48 hours. Bone samples were decalcified using 10% EDTA for 3 to 4 weeks (solution changed weekly), followed by embedding in paraffin wax and sectioning at a thickness of 5 μm. Automated IHC was performed with Ventana Discovery Ultra. Slides were deparaffinized and rehydrated. Antigen retrieval was done using EDTA buffer, and slides were incubated with primary antibodies: H3.3 G34W (1:100, RevMAb Biosciences 31-1145-00, RRID:AB_2716434), Ki-67 (1:300, Abcam ab15580, RRID:AB_443209), alpha smooth muscle actin (1:200, Abcam ab5694, RRID:AB_2223021), osteopontin (1:500, Abcam ab8448, RRID:AB_306566), Col VI (1:50, Abcam ab6588, RRID:AB_305585), Biglycan (1:500, Abcam ab49701, RRID:AB_1523212). After washing, secondary antibody (anti-mouse or rabbit/mouse HRP) was added and DAB kit chromogen used to detect signal.

Slides were deparaffinized and rehydrated through graded ethanols to distilled water. Slides were placed in prewarmed TRAP Staining Solution Mix (4% pararosaniline, 4% sodium nitrite, 0.1 acetate buffer pH 11, naphthol AS-TR Phosphate, sodium tartare, pH 5.0), and incubated at 37°C until control was developed, then rinsed in distilled water. Slides were counterstained with 0.1% Fast Green for 1 minute, air-dried, and mounted. Osteoclasts stain red-violet, whereas background stains green using TRAP staining.

ChIP-seq was performed as described previously (ref. 8; see Supplementary Information for details). ChIP reactions for histone modifications were performed on Diagenode SX-8G IP-Star Compact by incubating 2 × 10⁶ cells of sonicated cell lysate with following antibodies: anti-H3K27me3 (1:40, Cell Signaling Technology 9733, RRID:AB_2616029), anti-H3K27ac (1:80, Diagenode C15410196, RRID:AB_2637079), anti-H3K36me3 (1:100, Active Motif 61021, RRID:AB_2614986), anti-H3.3 (1:66, Millipore 09–838, RRID:AB_10845793), anti–H3.3 G34W (1:66, RevMAb Biosciences 31-1145-00, RRID:AB_2716434), anti-H3K36me2 (1:50, Cell Signaling Technology 2901; anti-H3K9me3 (1:66, Abcam ab8898, RRID:AB_306848). SUZ12 ChIPs were performed manually using anti-SUZ12 (1:150, Cell Signaling Technology 3737, RRID:AB_2196850) antibody.

RNA-seq was performed as described previously (ref. 53; see Supplementary Information for details).

Description of statistical details for each experiment can be found in figure legends. Statistical significance was always adjusted for multiple testing and considered to be attained when P < 0.05.

The ChIP-seq pipeline from C3G's GenPipes toolset (RRID: SCR_016376; ref. 54; v3.1.0) was used. See the Supplementary Information for a complete treatment of sample selection, data processing, data visualization, H3.3 peak enrichment analysis, comparison of H3.3 and H3.3 G34W deposition patterns.

The abundance of each mark was quantified by counting the number of reads in genomic 10 kb bins and normalized to RPKM with the “multiBamSummary bins” functionality of deeptools (RRID:SCR_016366; ref. 55; v3.1.3). Differential enrichment of each mark at each bin was calculated as the log₂ fold-change (LFC) between average abundance in each condition. Thresholds of LFC ±1 (corresponding to a 2-fold-change) and ±0.58 (1.5-fold change) were used to determine large changes in H3K27me3 and SUZ12 in Im-GCT-4072 lines, respectively. For Im-GCT-3504 lines, we used a threshold of LFC ±0.58. For edited lines derived from Im-GCT-4072–overexpressing H3.3 G34W or H3.3WT, we used a threshold of LFC ±0.32 (1.25-fold change). To avoid overestimating LFCs, we ignored bins with low counts in both conditions (i.e., below-median average abundance in either condition).

Bins were classified as promoter, intragenic, or intergenic regions based on overlap with promoter (i.e., ±1.25 kb of a transcription start site), intergenic (i.e., covered in the reference), and intergenic regions (i.e., complement of promoter and intergenic regions) using the UCSC table browser's (RRID: SCR_00578) “whole” annotation for Ensembl's ensGene reference (GRCh37; n = 60,234 genes). The significance of differences in proportion of bins overlapping these regions that also gain or lose a chromatin mark was determined using a χ² test.

Enrichment of chromatin marks at promoter and gene bodies was obtained by counting reads over promoter and intragenic regions, respectively. Similar to differential gene expression, differential enrichment analysis was computed using DESeq2 (RRID: SCR_015687; ref. 56; v1.18.1), using unit-normalized library depth as size factors (i.e., total number of uniquely mapped reads). Genes with (i) sufficient coverage (baseMean > 25), (ii) LFCs exceeding a threshold of ±1 (2-fold) for H3K27me3 and H3K27ac and ±0.58 (1.5-fold) for H3K36me3, and (iii) P < 0.05, were considered as significantly differentially enriched. For Im-GCT-3504 lines, we used a threshold of ±0.58 (1.5-fold) for changes in H3K27me3. To ensure more accurate estimates of LFCs in the low-information setting (i.e., low counts and high dispersion), we enabled shrinkage of LFCs toward zero (56). P values were adjusted for multiple testing using Benjamini–Hochberg.

Bulk RNA-seq was analyzed as described previously (13). Please refer to Supplementary Information for complete information.

Please refer to Supplementary Information for complete information on single-cell RNA-seq data analysis.

Low-quality cells were excluded on the basis of outlier mitochondrial content (indicative of cellular stress or damage) and number of genes expressed using the R package Seurat (RRID: SCR_016341; refs. 57, 58; v3.0.0). The upper threshold on mitochondrial content varied from >15% to 20% of total genes detected. Because GCT samples contain multinucleated osteoclasts, where a high number of UMI counts is a biological feature of this cell type, we applied a more lenient upper threshold of 3 times the interquartile range. Filtered cells were then normalized together using SCTransform (v0.2.0; ref. 59; with parameter variable.features; n = 3,000), and mitochondrial content percent was regressed. To reduce the impact of ambient multispecies mRNAs in xenografts, we applied the SoupX (v1.4.5; ref. 60) algorithm.

To robustly identify cell types associated with each cell cluster, we first derived markers for each cluster through differential expression (DE) analysis using a Wilcoxon rank sum test configured to only test genes expressed in >25% of cells and with log-fold change ≥0.25. We then combined the three following approaches: (i) comparing cluster markers to known cell type markers from the literature (30, 61), (ii) classifying cell types from reference datasets based on correlation of gene expression profiles, and (iii) pathway enrichment analysis of cluster-specific markers. We used SingleR package (v1.0.5; ref. 62) with default parameters to assign cell types to each cluster based on their Spearman correlation to gene expression profiles of known cell types in the Human Primary Cell Atlas (HPCA) reference. Because the stromal and osteoclast cell types were not present in the reference, we performed pathway enrichment analysis on cluster-specific gene markers using G:Profiler (RRID: SCR_006809) with term size set to 1000 and deriving pathways from GO:BP only (63). See Supplementary Information for additional details on stromal and osteoclast cell identification.

To examine heterogeneity within stromal cell populations, stromal cells were reclustered (i.e., in the absence of nonstromal cells), and highly variable genes recalculated. PCA revealed that PC1 was defined by SPP1-correlated genes, and PC2 by ACTA2-correlated genes. To determine whether these genes, which represent major sources of variation within the data, defined specific subsets of cells, correlation analysis was performed using the correlateGenes function from R package Scran (RRID: SCR_016944; v1.12.1; ref. 64). This revealed distinct gene sets differentially enriched in stromal cells, notably “SPP1 module” and “ACTA2 module.”

Technical (i.e., sequencing batch) and biological (i.e., divergent tumor evolution and genetic backgrounds) variance can obscure similarities between neoplastic cells from different tumors. To validate that the gene modules identify distinct stromal cell subsets or states, we used the Harmony data integration method (v0.99.6; ref. 65; parameters max.iter.harmony = 15, theta = 2) to correct these sources of variation on the PCA embeddings. Although unsupervised clustering identified distinct S1A and S1B clusters, we collectively defined these stromal subtypes as S1 because they displayed similar enrichment for the “SPP1 module” and modest differential gene expression (Supplementary Fig. S5H). We confirmed that overcorrection did not occur, as the same nonstromal populations were recovered before and after Harmony integration. Both expression of the gene modules derived before correction and cluster-specific markers obtained by differential expression analysis after correction were in agreement.

To confirm the presence of the stromal cell subtypes S1 to S3 in xenograft tumors, we created a SingleR reference from Harmony-integrated patient GCT cells using the “aggregateReference” function (with power = 0). We then classified cell clusters from xenograft tumors using these references, as described earlier.

To determine cell trajectory within each lineage, we employed two complementary trajectory inference methods, Slingshot (RRID: SCR_017012; v1.3.1; ref. 66) and velocyto.R (RRID: SCR_018167; https://github.com/velocyto-team/velocyto.R, v0.6; ref. 67). To calibrate parameters for trajectory inference in stromal cells, we first analyzed a well-defined lineage, myeloid cells. To identify lineages, we applied unsupervised Slingshot analysis on a three-dimensional (3-D) diffusion map embedding. A 3-D diffusion map was constructed using the DiffusionMap function from the R package destiny (v2.14.0; ref. 68), with only high-quality cells included. High-quality cells were defined here as cells having < 5% mitochondrial gene content and genes with >10 counts in at least 100 cells. We then used Slingshot (with default parameters and without specifying a start cluster) to identify lineages in the diffusion map embedding (Supplementary Fig. S7E and S7F). Three-dimensional visualizations were created using R package plotly (v4.9.1; ref. 69). Because Slingshot cannot confidently establish the start of the lineage, we employed velocyto. Briefly, spliced and unspliced reads were counted using the Python implementation of velocyto (v0.17.13). RNA velocity of each gene passing the expression threshold (exonic and intronic reads with a minimum maximum-cluster average of 0.1 and 0.05, respectively, were retained) was estimated using the R package velocyto.R (v0.6) using a gene-relative model with parameters kCells = 20 and fit.quantile = 0.05 and projected onto UMAP embeddings. Consistent with Slingshot results and literature, velocyto predicted a bifurcated lineage with a mix of monocytes and preosteoclasts as a root. Slingshot and velocyto were then applied, using the same parameters, to reconstruct trajectories in stromal cells.

To better determine a putative cell-of-origin for stromal cells, we classified them using SingleR with custom references derived from two murine bone marrow scRNA-seq datasets from Baryawno and colleagues 2019 (27) and Tikhonova and colleagues 2019 (28), which better cover the diversity of cells found in in the bone milieu compared with HPCA. Mouse genes were humanized to their human orthologs using Ensembl Biomart (v2.42.0; refs. 70, 71). Genes without a human match and irrelevant genes (i.e., ribosomal and housekeeping-related) were removed to decrease the risk of spurious correlation. The humanized expression matrices were then converted into a SingleR reference using the function aggregateReference (with parameter power = 0). To ensure that humanization did not remove important genes, we confirmed that there were no errors when classifying clusters from the mouse datasets using their corresponding humanized references. To detect overlapping cell types across the two references, we classified cells from one reference using the other. The fact that neoplastic stromal cells match a cell type that was uniquely matched between the two references further supports the putative identity of the cell of origin (Supplementary Fig. S5K). Each reference was then used to classify the stromal cells from each of the stromal cell subtype populations.

Cell-to-cell interaction analysis in the GCT TME was performed as described previously (31) and described in the Supplementary Information.

Data processing, protein enrichment analysis, and identification of secreted factors by stromal cells are described in the Supplementary Information.

No disclosures were reported.

S. Khazaei: Conceptualization, data curation, formal analysis, writing-original draft. N. De Jay: Conceptualization, data curation, formal analysis, writing-original draft. S. Deshmukh: Conceptualization, data curation, formal analysis, writing-original draft. L.D. Hendrikse: Conceptualization, data curation, formal analysis, writing-original draft. W. Jawhar: Data curation, formal analysis. C.C.L. Chen: Conceptualization, data curation, formal analysis. L.G. Mikael: Resources, data curation, writing-review and editing. D. Faury: Methodology. D.M. Marchione: Data curation, methodology. J. Lanoix: Data curation, formal analysis, methodology. E. Bonneil: Data curation, formal analysis. T. Ishii: Data curation, methodology. S.U. Jain: Conceptualization, data curation, formal analysis. K. Rossokhata: Data curation, formal analysis. T.S. Sihota: Data curation, formal analysis. R. Eveleigh: Software, formal analysis. V. Lisi: Software, formal analysis. A.S. Harutyunyan: Data curation, formal analysis, methodology. S. Jung: Formal analysis. J. Karamchandani: Formal analysis. B.C. Dickson: Data curation, formal analysis. R. Turcotte: Investigation. J.S. Wunder: Data curation, formal analysis. P. Thibault: Data curation, formal analysis. P.W. Lewis: Conceptualization, formal analysis. B.A. Garcia: Data curation, formal analysis. S.C. Mack: Conceptualization, formal analysis. M.D. Taylor: Conceptualization, formal analysis. L. Garzia: Conceptualization, formal analysis. C.L. Kleinman: Conceptualization, formal analysis, supervision, funding acquisition, investigation, writing-original draft, project administration, writing-review and editing. N. Jabado: Conceptualization, resources, formal analysis, supervision, funding acquisition, visualization, writing-original draft, writing-review and editing.

This work was supported by funding from A Large-Scale Applied Research Project grant from Genome Quebec, Genome Canada, the Government of Canada, and the Ministère de l'Économie, de la Science et de l'Innovation du Québec, with the support of the Ontario Institute for Cancer Research through funding provided by the Government of Ontario to N. Jabado, M.D. Taylor, and C.L. Kleinman (13513); Fondation Charles Bruneau (to N. Jabado), NIH grant P01-CA196539 (to N. Jabado, P.W. Lewis, and B.A. Garcia; R01CA148699 and R01CA159859, to M.D. Taylor); the Canadian Institutes for Health Research (CIHR grant MOP-286756 and FDN-154307, to N. Jabado; and PJT-156086, to C.L. Kleinman); the Canadian Cancer Society (CCSRI grant 705182) and the Fonds de Recherche du Québec en Santé (FRQS) salary award (to C.L. Kleinman); NSERC (RGPIN-2016-04911) to C.L. Kleinman; CFI Leaders Opportunity Fund 33902 to C.L. Kleinman, Genome Canada Science Technology Innovation Centre, Compute Canada Resource Allocation Project (WST-164-AB); data analyses were enabled by compute and storage resources provided by Compute Canada and Calcul Québec. N. Jabado is a member of the Penny Cole Laboratory and the recipient of a Chercheur Boursier, Chaire de Recherche Award from the FRQS. This work was performed within the context of the International CHildhood Astrocytoma INtegrated Genomic and Epigenomic (ICHANGE) consortium, and the Stand Up To Cancer (SU2C) Canada Cancer Stem Cell Dream Team Research Funding (SU2C-AACR-DT-19-15, to M.D. Taylor, N. Jabado), provided by the Government of Canada through Genome Canada and the Canadian Institute of Health Research, with supplemental support from the Ontario Institute for Cancer Research, through funding provided by the Government of Ontario. Stand Up To Cancer Canada is a Canadian Registered Charity (Reg. # 80550 6730 RR0001). Research Funding is administered by the American Association for Cancer Research International - Canada, the Scientific Partner of SU2C Canada. M.D. Taylor is supported by The Paediatric Brain Tumour Foundation, The Terry Fox Research Institute, The Canadian Institutes of Health Research, The Cure Search Foundation, b.r.a.i.n.child, Meagan's Walk, Genome Canada, Genome BC, Genome Quebec, the Ontario Research Fund, Worldwide Cancer Research, V-Foundation for Cancer Research, Canadian Cancer Society Research Institute Impact grant and the Garron Family Chair in Childhood Cancer Research at the Hospital for Sick Children and the University of Toronto. N. De Jay is a recipient of a fellowship from FRQS and Réseau de Médecine Génétique Appliquée. We also acknowledge support from the We Love You Connie Foundation (to N. Jabado). S.C. Mack is supported by Cancer Prevention Research Institute of Texas (CPRIT) Scholar in Cancer Research award (RR170023), and Alex's Lemonade Stand Foundation (ALSF) A award. L. Garzia is supported by the Nicole et François Angers Sarcoma Research Chair from the MGH Foundation, C17 Ewing Sarcoma Grant, and the Terry Fox Research Institute Grant (number 1087).