Histone tail analysis reveals H3K36me2 and H4K16ac as epigenetic signatures of diffuse intrinsic pontine glioma

Diffuse intrinsic pontine glioma (DIPG) is an aggressive pediatric brainstem tumor. Most DIPGs harbor a histone H3 mutation, which alters histone post-translational modification (PTM) states and transcription. Here, we employed quantitative proteomic analysis to elucidate the impact of the H3.3K27M mutation, as well as radiation and bromodomain inhibition (BRDi) with JQ1, on DIPG PTM profiles. We performed targeted mass spectrometry on H3.3K27M mutant and wild-type tissues (n = 12) and cell lines (n = 7). We found 29.2 and 26.4% of total H3.3K27 peptides were H3.3K27M in mutant DIPG tumor cell lines and tissue specimens, respectively. Significant differences in modification states were observed in H3.3K27M specimens, including at H3K27, H3K36, and H4K16. In addition, H3.3K27me1 and H4K16ac were the most significantly distinct modifications in H3.3K27M mutant tumors, relative to wild-type. Further, H3.3K36me2 was the most abundant co-occurring modification on the H3.3K27M mutant peptide in DIPG tissue, while H4K16ac was the most acetylated residue. Radiation treatment caused changes in PTM abundance in vitro, including increased H3K9me3. JQ1 treatment resulted in increased mono- and di-methylation of H3.1K27, H3.3K27, H3.3K36 and H4K20 in vitro. Taken together, our findings provide insight into the effects of the H3K27M mutation on histone modification states and response to treatment, and suggest that H3K36me2 and H4K16ac may represent unique tumor epigenetic signatures for targeted DIPG therapy.

For example, histone H3K27 di-and tri-methylation (H3K27me2/3), transcriptionally silencing marks, are reduced in H3K27M DIPG, de-repressing transcription [13][14][15]. H3K27me2/3 loss is accompanied by increased H3K27ac and the formation of heterotypic H3K27M/ H3K27ac nucleosomes that colocalize with bromodomain proteins 2 and 4 (BRD2/4) and RNA polymerase II, suggesting a potential role of BRD proteins in DIPG pathogenesis [7,16]. We demonstrated pre-clinical efficacy of bromodomain inhibition (BRDi) with JQ1, a bromodomain and extra terminal domain (BET) inhibitor, resulting in reduction of H3K27ac and restoration of H3K27me3 [7]. These findings suggest that changes in epigenetic regulation of gene expression in DIPG may be reflected by histone codes, and detection of these signatures could shed light on mechanisms of tumorigenesis. Here, we employ targeted mass spectrometry via multiple reaction monitoring (MRM) to characterize 93 distinct PTM states in H3.3K27M mutant DIPG specimens, and investigate the effects of radiation and bromodomain inhibition on tumor epigenetic signatures.

Tissue specimens
Pediatric DIPG tumor tissue (n = 9, H3.3K27M mutant) was collected postmortem (IRB #Pro00001339). H3.3K27 wild-type (WT) adult glioma tumor tissue specimens (n = 2), and one non-tumor frontal lobe tissue specimen were collected during the course of treatment (IRB# STU00095863). All diagnoses were made radiographically by a neuroradiologist, then confirmed by neuropathological evaluation upon tissue collection. Patient identifiers were removed and a single sequential numerical identifier assigned to each specimen (Table 1). Tissues were stored at − 80°C prior to proteomic analysis (IRB #STU00202063).

Gene expression profiling
Total cell RNA was extracted for gene expression profiling based on manufacturer protocol (Qiagen). All samples were adjusted to a concentration of 20 ng/μl RNA prior to submission for gene expression analysis using the NanoString nCounter Plexset Analysis System DNA damage repair 96 gene codeset (NanoString Technologies, Seattle, WA). Gene expression values were determined and converted to fold change in expression level for each treatment condition compared to DMSO controls in order to determine statistically significant differences in gene expression levels between groups (NanoString). Functional pathways and gene ontology analysis was then performed on fold change expression values (Inginuity Pathways Analysis, Qiagen).

Histone preparation for targeted mass spectrometry
Nuclei were isolated in buffer (15 mM Tris-HCl pH 7.5, 60 mM KCl, 15 mM NaCl, 5 mM MgCl 2 , 1 mM CaCl 2 , 250 mM sucrose; 1 mM DTT, 1:100 HALT protease inhibitor (Thermo Scientific)), and 10 mM sodium butyrate added immediately prior to use containing 0.3% NP-40. Nuclei were centrifuged at 600 x g for five minutes at 4°C and washed twice with nuclear isolation buffer. Histones were extracted by via five volumes of H 2 SO 4 and incubated at room temperature one hour. Debris was removed by centrifugation at 4000 x g for five minutes. Trichloroacetic acid was added to the supernatant at a final concentration of 20% (v/v) and incubated overnight at 4°C. Precipitated histones were pelleted at 10, 000 x g for five minutes, washed once with ice-cold acetone with 0.1% HCl and twice with acetone at 15,000 x g for five minutes. Residual acetone was dried in a fume hood, and samples stored at − 20°C until further processing. Histones were derivatized with a single round of propionylation, digested with 1 μg trypsin (Promega, Madison, WI), and derivatized again with a single round of propionylation according to Garcia et al. [20].

Targeted LC-MS/MS of histone peptides
Histone peptides were resuspended in Loading Buffer (0.1% TFA in water) and analyzed on a Thermo TSQ Quantiva (Thermo Scientific) in-line with a Dionex nano-LC (Sunnyvale, CA). Peptides were loaded onto a trapping column (3 cm × 150 μm, packed with ProntoSIL C18-AQ, 3 μm, 200 Å resin (New Objective, Woburn, MA)) for 10 min at 2.5 μL/min with 100% Loading Buffer. Elution from the trapping column and separation on a PicoChip analytical capillary column (10 cm × 75 μm packed with ProntoSIL C18-AQ, 3 μm, 200 Å resin (New Objective)) was achieved by decreasing the percentage of Solvent A (0.1% formic acid in water) and increasing the percentage of Solvent B (0.1% formic acid in 95% acetonitrile) from 1 to 35% at a flow rate of 0.30 μL/min over 45 min. The peptides were introduced into the triple quadrupole mass spectrometer by electrospray from an emitter with a 10 μm tip. The instrument settings were as follows: collision gas pressure of 1.5 mTorr; Q1 peak width of 0.7 (FWHM); cycle time of 3 s; skimmer offset of 10 V; electrospray voltage of 2.5 kV. All injections were performed in technical triplicate. Targeted methods were performed as previously published [21]. A list of monitored peptide transitions is provided (Supplementary material 1).

Quantitation of histone modification states
Raw data were imported into Skyline with Savitzky-Golay smoothing for analysis. Peak area assignments were manually confirmed, and total peak areas used to determine the percent relative abundance of histone modifications. The percent relative peptide abundance was determined by dividing the area under the curve for a given peptide of interest by the sum of the areas under the curve for all peptides with the same sequence, then multiplying by 100. Heatmaps were generated with R using gplots package, and show modifications that were above the limit of detection in all samples.

Statistical analysis
Differences in peptide abundance were analyzed using independent-sample t test for between-group comparisons, or one-way ANOVA and Tukey's post hoc test for multi-group analyses using SPSS (IBM, Version 25), and expressed as mean ± standard deviation (SD). Pearson Correlation was performed to analyze the relationship between H3.3K27M mutation and detected histone codes. Multivariable analysis was performed using linear regression with Stepwise entry (p = 0.05) and removal (p = 0.10). All p values were two-sided, with p < 0.05 considered statistically significant.
We also quantified the relative abundance of H3.3K27 WT and H3.3K27M peptides in a given tissue sample. Similar to cell lines, 26.44% of H3.3K27 peptides in our mutant tumor specimens were H3.3K27M (SD = 10.7, range 5.2-44.5%, n = 9), with no mutant peptides detected in WT tissues. Lastly, we compared DIPG cell and tumor tissue epiproteomic profiles. Overall, modifications on the H3.1/ H3.3K27 N-terminal tail were similar in DIPG cell lines and tumor tissues. Specifically, there was no significant difference in abundance of the majority of H3.1/ H3.3K27 peptides quantified, with differential fold change (FC) in abundance > 2 or < 0.5 for only six peptide modification states between DIPG cells and tissues ( Table 2, Supplementary material 6). Significant differences in peptide modification states between groups were observed on H1.4K25, H3.3K36, H3K79, H3K9, and H4K20 peptides. Importantly, no difference in relative H3.3K27M peptide abundance was detected between DIPG tissue and cell lines (mean = 26.4 and 29.2, respectively, p = 0.686). Of note, unmodified peptides and modifications with < 5% relative abundance were excluded from this comparative analysis to ensure biological significance.

Radiation treatment alters DIPG histone codes in vitro
As radiation treatment (RT) is the standard therapy for DIPG, we examined the effects of RT on the glioma histone code in vitro. We treated H3.3K27M mutant DIPG cells (SF8628, DIPG007), H3.3K27 WT glioma cells (U87) and normal astrocytes (NHA) with 9 Gy radiation, and performed targeted histone mass spectrometry at 24 and 48 h post-treatment. Cell histone codes clustered first by origin, then by treatment condition (Fig. 5a). Comparative analysis revealed distinct differences in relative peptide abundance after RT, most notably in acetylated peptides (Supplementary material 8). In DIPG007 cells, the relative abundance of H2AK5ac decreased with RT at 48 h (Fig. 5b). In contrast, H2AK5ac abundance increased in NHA at 24 and 48 h, with no significant change in SF8628 or U87 (Fig. 5b). At 24 and 48 h post-RT, all lines except DIPG007 showed a significant increase in H4K16ac (Fig. 5c, Supplementary material 7). H3K9me3 also significantly increased in NHA, U87, and DIPG007, and trended towards a significant increase in SF8628, after RT (Fig. 5d, Supplementary  material 8).
To determine the potential biological impact of observed differences in histone modifications after RT relative to controls, we compared cell gene expression profiles in these same treatment groups. On functional pathways analysis of genes identified as statistically significantly differentially expressed by RT treated cells relative to controls (fold change in expression > 2 or < − 2, adjusted p-value < 0.05), we identified DNA double-strand break repair as a top canonical pathway activated in RT treated DIPG cell lines (p-value 1.10 × 10 − 21 , Supplementary material 9). The top implicated molecular and cellular function in RT treated cells, as well as the top network of molecular interaction, was DNA replication, recombination and repair (p-value 2.06 × 10 − 6 ),

Bromodomain inhibition alters DIPG histone codes in vitro
We previously reported pre-clinical efficacy of JQ1, a BET/Bromodomain inhibitor (BRDi), in DIPG [7]. We observed that BRD proteins are enriched at acetylated residues on histone H3 in DIPG cells and animal models, and that JQ1 treatment results in decreased global H3K27ac, with restoration of global H3K27 methylation [7]. To further elucidate the epigenetic effects of JQ1 in DIPG, we examined its effects on histone codes in SF8628, DIPG007, NHA, and U87 cells. As observed with RT, PTM profiles clustered by cell line, then treatment condition (Supplementary material 10). When analyzed individually, cell profiles cluster by treatment condition (Fig. 6a, Supplementary material 11). In DIPG007 cells, H3.3K27me1 decreased by 1.5% at 48 h after JQ1 treatment compared to vehicle control, representing the only significant change in modification state observed at that time point (Supplementary material 11). A small increase in H3.3K27M abundance was also detected in DIPG007 at 48 h (31.4% vs 26.9%, p = 0.067, Supplementary material 11). Surprisingly, in DIPG007, H3.3K27ac abundance was not significantly different from control at 48 h post-treatment, though a marginal decrease was observed at 24 h (1.65 and 1.95%, respectively, p = 0.065, Supplementary material 11). These results were validated in DIPG007 and SF8628 cells via western blot for H3K27M and H3K27ac (Supplementary material 12).
To determine the potential biological impact of observed differences in histone modifications after BRD inhibition relative to controls, we also compared cell gene expression profiles in these treatment groups. On functional pathways analysis of genes identified as statistically significantly differentially expressed by JQ1 treated cells relative to controls (fold change in expression > 2 or < − 2, adjusted p-value < 0.05), we identified nucleotide excision repair as the top canonical pathway activated in JQ1 treated DIPG cell lines (p-value 2.51 × 10 − 48 , Supplementary material 9). Other canonical pathways implicated by JQ1 treated cell gene expression profiles included DNA mismatch repair via nonhomologous end joining (NEMJ, p-value 2.53 × 10 − 32 ), while G2/M DNA damage check point regulation was the top toxicology list match (p-value 2.68 × 10 − 6 ). Top upstream regulator effect networks implicated ANXA7 activation (resulting in DNA binding and conformational modification, score 9.615), and TP73 inhibition (resulting in apoptosis, score − 15.254) in JQ1 treated DIPG cells relative to controls.

Discussion
DIPG is an incurable pediatric brain tumor, with a high rate of somatic missense mutations in genes encoding histone H3.1 and H3.3 proteins. Here, we present comprehensive characterization of histone N-terminal tail modification states in pediatric glioma specimens, including H3.3K27M mutant DIPG tissue, using novel, high-throughput proteomic analysis, revealing the effects of radiotherapy and BRDi on the tumor epigenetic landscape. The greatest differences in PTMs between H3K27M mutant and WT specimens were observed at H3K27, H3K36, and H4K16ac, while H3K36me2 was the highest co-occurring modification on the H3.3K27M peptide, and H4K16 was the most acetylated residue. These results suggest that H3K36me2 and H4K16ac may serve as epigenetic signatures in DIPG, and that targeting these modification states may represent a potential therapeutic strategy.
Multiple groups have identified distinct DIPG molecular subtypes, with H3.3K27M mutant tumors associated with shorter survival [22][23][24][25][26][27][28][29][30]. We found that H3.3K27M significantly affects the histone PTM landscape, with cell and tissue profiles clustering by mutation status. Further, a decrease in global H3K27me2/3 has been observed in H3.3K27M tumors, and our results corroborate these data, [5,7,8,10,12,13,31]. Of note, a new subgroup of DIPG lacking H3K27M mutation but still exhibiting genomic H3K27me3 loss was recently reported [32]: we did not observe this phenomenon in our cohort (Supplementary material 4). Our H3K27M mutant specimens also exhibit distinct modification states at H3K36 and H4K16. Previous work showed that PRC2 loss coincides with increased H3K36me2, and that an increase in H3K36 methylation impacts methylation at H3K27 [33][34][35]. Our results also confirm an increase in H3K36me2 in DIPG specimens on the H3.3K27M mutant peptide, as previously detected via mass spectrometry of H3.3K27M transgenic PDGF-driven murine glioblastoma [13]. Further, we detected high H4K16ac abundance in DIPG cells and tissue specimens. A recently published report characterizing high-throughput screening with nucleosome substrate methods demonstrates small-molecule inhibitors of the human histone lysine methyltransferase NSD2 [36], the writer of mono-and di-methylation of H3K36 [37]; another report identified DC_M01_7 as a novel inhibitor of MOF, the H4K16acspecific acetyltransferase [38,39]. Therefore, H3K36me2 and H4K16ac may represent potential therapeutic targets in DIPG, and further investigation of this strategy is warranted.
Using our method, we also characterized H3K27 and K36 modification states in cis, confirming me2/3 as the predominant modifications, and acetylation as the least abundant, on the H3.3K27MK36 peptide. We observed positive correlation between H3.3K27M and H3.1K27K36me2, and negative correlation between H3.3K27M and H3.1K27K36me3. Taken together, these results confirm that H3.3K27M affects PTMs on both mutant and WT H3.1 and H3.3 N-terminal tails. Importantly, our study is the first, to our best knowledge, to quantify the relative amount of mutant H3.3K27M peptide in DIPG tumor specimens. It is reported that relative abundance of mutant epidermal growth factor receptor (EGFR) can predict benefit from EGFR-tyrosine kinase inhibitors for advanced non-small-cell lung cancer [40]. Longitudinal quantification of the H3.3K27M peptide in patientderived specimens might similarly serve as a predictive biomarker for targeted therapy, and hence is worthy of consideration when designing future clinical trials.
Radiation therapy (RT) is standard treatment for DIPG, providing transient symptom relief but no change  in disease outcome [41]. Unfortunately, the mechanism of radiation resistance is unknown, and radio-sensitizers provide no benefit [42]. Patients with the H3.3K27M mutation also respond more poorly to RT [24]. Our study reveals RT can induce cell-specific and peptide-specific changes in histone acetylation, underscoring the heterogeneity of DIPG and warranting consideration when using histone deacetylase / acetyltransferase inhibitors in combination with RT. Importantly, H3K9me3 abundance increased in all four cell lines after RT, which kills rapidly proliferating cancer cells by inducing DNA damage beyond cellular capacity for repair, including DNA doublestrand breaks (DSBs) [43][44][45]. DSB repair is linked to rapid changes in epigenetic modifications, including increased H3K9 methylation, to recruit DNA repair proteins. Importantly, H3K9 methyltransferase G9a inhibition (BIX-01294) increases radio-sensitivity in glioma [46]. Indeed, on functional pathways analysis of cell gene expression profiles, we found DSB repair pathways implicated in RT treated cell lines relative to controls. These data suggest that combination therapy RT with G9a inhibition, aimed at reducing H3K9me3 levels and hence DNA DSB repair, may also be a useful therapeutic strategy in DIPG. We found distinct DIPG histone codes after treatment with JQ1, a BET/ bromodomain inhibitor. BET domain proteins are epigenetic readers that bind acetylated histones via their bromodomain region to regulate gene transcription [46][47][48]. BET inhibitors dislodge acetylated histone readers from chromatin, leading to oncogene repression [49][50][51]. JQ1 is the best known BET inhibitor, and induces a dose-dependent reduction in DIPG cell viability and genomic H3K27ac in vitro and in vivo [7,52,53]. Our results showed increased mono-and/or dimethylation of multiple residues with JQ1 treatment, consistent with previous reports that JQ1 globally reduces acetylation and increases methylation [7]. Since monomethylation of H3K27, H3K79, and H4K20 are all linked to gene activation, whereas tri-methylation of H3K27 and H3K79 are linked to repression [54], these findings may provide additional insight into the mechanism by which BRDi alters tumorigenic gene expression in DIPG. H4K20me, which is functional in DNA repair, represents a binding site for the 53BP1 protein; H3K9me3 and H4K20me3 represent epigenetic markers important to the function of 53BP1 in non-homologous end joining (NHEJ) repair [55]. In line with these known functions, pathways analysis of differential gene expression profiles implicated nucleotide excision repair and DNA mismatch repair via NHEJ as top activated pathways in JQ1 treated cells, compared to controls. As such, additional studies of the mechanism by which JQ1 may augment response to RT are necessary. Lastly, we also observed that reduction of histone acetylation in DIPG was greatest at 24 h, which may be important to consider when determining JQ1 dosing regimens.
Because fresh DIPG tissue is very rare, we chose to analyze archival post-mortem tissue specimens for this study, raising the possibility of treatment-induced effects. However, the patterns of peptide abundance observed in tumor tissues were similar to treatment-naïve cell lines. Analysis of treatment fresh tumor tissue specimens with this approach will be crucial to teasing out tumor biology from the effects of therapy. Optimization of our technical approach to minimize the amount of starting tumor tissue necessary for this purpose is currently underway.

Conclusions
In summary, targeted histone tail mass spectrometry of high grade glioma cell lines and tissues, including rare DIPG tumor tissue specimens, reveals distinct peptide modification states by H3.3K27M mutation status. Our study is the first to quantify the H3.3K27M peptide in DIPG, representing an exciting new opportunity to apply this technique on clinical specimens at diagnosis, and potentially throughout the course of treatment to guide clinical management. We also distinct epigenetic signatures of DIPG, including H3K36me2 and H4K16ac, which may represent novel therapeutic vulnerabilities for more effective treatment of this devastating disease.