DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML

Disruption of the epigenetic landscape is of particular interest in acute myeloid leukemia (AML) due to its relatively low mutational burden and frequent occurrence of mutations in epigenetic regulators. Here, we applied an information-theoretic analysis of methylation potential energy 5 landscapes, capturing changes in mean methylation level and methylation entropy, to comprehensively analyze DNA methylation stochasticity in subtypes of AML defined by mutually exclusive genetic mutations. We identified AML subtypes with CEBPA double mutation and those with IDH mutations as distinctly high-entropy subtypes, marked by methylation disruption over a convergent set of genes. We found a core program of epigenetic 10 landscape disruption across all AML subtypes, with discordant methylation stochasticity and transcriptional dysregulation converging on functionally important leukemic signatures, suggesting a genotype-independent role of stochastic disruption of the epigenetic landscape in mediating leukemogenesis. We further established a relationship between methylation entropy and gene expression variability, connecting the disruption of the epigenetic landscape to 15 transcription in AML. This approach identified a convergent program of epigenetic dysregulation in leukemia, clarifying the contribution of specific genetic mutations to stochastic disruption of the epigenetic and transcriptional landscapes of AML. 20 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 3 Main Text Acute myeloid leukemia (AML) is a hematopoietic malignancy driven by a combination of genetic lesions and epigenetic disruptions, leading to abnormal clonal expansion and blockade of myeloid differentiation (1, 2). While AML has a relatively low mutational burden, mutational 5 analysis has identified recurrent genetic lesions in epigenetic machinery like DNMT3A, IDH1 or IDH2, TET2, and MLL (KMT2A); in transcription factors like CEBPA and RUNX1; and in signaling molecules such as FLT3, KIT, and NRAS (1, 2). Interestingly, mutations in the epigenetic machinery often occur early in leukemogenesis, suggesting an important role of epigenetic regulation in leukemia (3, 4). Previous work has identified specific patterns of DNA 10 methylation linked to mutational profiles in AML, such as hypomethylation in DNMT3A-mutant AMLs and hypermethylation in IDH-mutant AMLs (5–8). Analysis of DNA methylation data has also revealed that some genetic subtypes of AML may be distinguished by methylation signatures, particularly subtypes with mutations in epigenetic machinery, and that this classification may hold prognostic relevance (9–13). Thus, the role of mutations in the epigenetic 15 machinery in mediating disruption of the epigenetic landscape of AML is a crucial topic of investigation. We and others have shown methylation stochasticity to be linked to plasticity, tumor evolution, and prognosis in leukemia and other cancers (14–22). To rigorously capture the higher-order 20 statistical properties of methylation, our group has previously developed information-theoretic

to model methylation potential energy landscapes (PELs), enabling the detection of stochastic disruption of the epigenetic landscape (23–26). These methods allow for direct comparisons of the underlying probability distributions of methylation, measuring differential methylation stochasticity by capturing changes in both mean methylation level and methylation 25 entropy, as well as other properties that are not detectable by mean-based analysis. We have previously applied these information-theoretic methods to investigate epigenetic landscape disruption in solid tumors (14, 23), in acute lymphoblastic leukemia (ALL) (16), and in MLL- rearranged AML (15), where we identified discordant methylation stochasticity over key regulators of the malignant phenotype, pointing to the role of stochastic epigenetic regulation in 30 oncogenesis. Li et al. subsequently carried out a comparative analysis of epigenetic heterogeneity in genetic subtypes of AML using empirical metrics of epiallele diversity, finding that specific somatic driver mutations were associated with epigenetic allele diversity over distinct loci (22). To further dissect the contribution of AML driver mutations to disruption of the epigenetic landscape, we applied the Correlated Potential Energy Landscape (CPEL) method (25, 26) to 35 characterize the methylation PELs of genetically defined AML subtypes, focusing on subtypes with mutually exclusive mutations in DNMT3A, IDH1/2, TET2, CEBPA double mutation (CEBPA-dm), KIT, or NRAS. We sought to carry out a comprehensive analysis of methylation stochasticity in genetically-defined AML subtypes, parsing out the contribution of individual mutational drivers to the disruption of the epigenetic landscape. By integrating epigenetic 40 landscapes with single-cell RNA sequencing, we investigated a functional relationship between DNA methylation stochasticity and gene expression. Information-theoretic methylation potential energy landscape analysis reveals discordant methylation stochasticity in AML 45 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 4 To investigate methylation stochasticity driven by genetic mutations in AML, we applied CPEL, an information-theoretic method for analysis of differences in methylation potential energy landscapes (26). CPEL detects differences in mean methylation level (MML), normalized methylation entropy (NME), and probability distributions of methylation (PDMs) based on methylation patterns within analysis regions (Figure 1A, S1). Differences in PDMs are 5 quantified by the uncertainty coefficient (UC), the geometric Jensen-Shannon divergence normalized by the cross-entropy between two groups, which captures both mean and entropy changes (Methods). We evaluated analysis regions containing an average of 6.5-9 CpGs per region (Data S1), requiring CpGs to have at least 10x coverage in every sample (Methods). 10 We employed a dataset of enhanced reduced representation bisulfite sequencing (ERRBS) data from a large cohort of primary AML samples (6). Samples were selected based on the presence of mutually exclusive mutations in DNMT3A (n=10), IDH1/2 (n=8), TET2 (n=6), CEBPA double mutation (CEBPA-dm, n=8), KIT (n=6), or NRAS (n=6), allowing for evaluation of epigenetic landscape disruption mediated by specific somatic mutations without the presence of 15 confounding mutations of other subtypes (Figure 1B, Data S2). ERRBS data from healthy CD34+ progenitors was used for normal reference (n=6) (27, 28). Strikingly, CEBPA-dm and IDH1/2-mut AMLs had distinct hypermethylation and increased methylation entropy over promoter regions (Figure 1C, S2-3). Principal component analysis (PCA) of MML and NME separated CEBPA-dm and IDH1/2-mut AMLs from other subtypes and normal samples (Figure 20 1D), indicating that these two subtypes have strong and convergent methylation disruption not found in other AML subtypes. DNMT3A-mut and TET2-mut AMLs had modest methylation stochasticity discordance (3,358 and 2,752 regions with significant UC, respectively) (Data S1), despite the canonical role of these enzymes as mediators of DNA methylation (29, 30). 25 We observed strong disruption of methylation over CpG islands (CGI), shores, and shelves (Figure S2, S3). Regions with significant differences in PDMs (UC DMRs) were enriched over CGI and promoter regions in CEBPA-dm, IDH1/2, TET2, KIT, and NRAS-mut AMLs (Figure S4). Interestingly, UC DMRs in DNMT3A-mut AMLs were depleted over these features and enriched over open seas and repetitive/CNV elements, indicating that DNMT3A mutations may 30 mediate changes in DNA methylation outside of CGI-associated features. UC DMRs were also enriched in heterochromatin regions, suggesting a role of increased methylation stochasticity in disruption of chromatin organization. HOMER motif analysis of DMRs revealed convergent enrichment of UC DMRs in CEBPA-dm 35 and IDH1/2-mut AMLs over transcription factor binding motifs for Homeobox family transcription factors such as HOXA1, HOXC6, CUX1, and OCT4 (POU5F1) (Figure S5, Data S4). CEBPA-dm and IDH1/2-mut DMRs also both showed enrichment over motifs for AP-2 (TFAP2A, TFAP2G), and LRF (ZBTB7A) transcription factors. These transcription factors play important roles in AML: disruption of Homeobox family members in AML is well-established 40 (31); AP-2 transcription factors have been implicated in disease progression through HOXA disruption (32); and ZBTB7A is necessary for myeloid differentiation and is frequently mutated in t(8;21) AML (33). In contrast, DMRs in DNMT3A-mut AMLs were enriched over motifs for ETS, RUNX, and bZIP family transcription factors, such as PU.1 (SPI1), a key regulator in hematopoiesis (34); RUNX2, implicated in differentiation and development of Cbfb-MYH11 45 AML (35); and AP-1 (JUNB), which has been linked to development of myeloproliferative (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 5 disease in mice (36) (Figure S5). Thus, genetic mutations, particularly mutations in epigenetic machinery, may contribute to leukemogenesis by mediating disruption of DNA methylation over binding motifs for key transcription factors. Together, these results indicate that genetic mutations are associated with highly stochastic 5 methylation in AML. Across subtypes, differential methylation and disruption of transcription factor binding motifs localizes to key regulatory elements involved in leukemia. Furthermore, we find that CEBPA-dm and IDH1/2-mut AMLs have distinct and highly entropic methylation profiles. Both of these mutations have previously been associated with hypermethylation in AML (2, 9). Li et al. identified characteristic distribution of epigenetic allele diversity associated 10 with CEBPA-dm AML, but observed less tightly linked profiles in IDH1/2-mut AMLs (22). The stochastic methylation landscapes of these subtypes suggests a crucial role of epigenetic regulation and disruption of methylation in AML. Discordant methylation stochasticity localizes to key regulators of the leukemic phenotype 15 across AML subtypes Next, we mapped regions with significant differences in PDM (UC DMRs), MML (dMML DMRs), or NME (dNME DMRs) to genes (Data S3). CEBPA-dm and IDH1/2-mut AMLs had strong convergence of genes containing DMRs (Figure 2A, S6A-B), indicating that the 20 increased methylation stochasticity observed in these AML subtypes occurs over similar targets. The overlap between UC DMR genes in CEBPA-dm and IDH1/2-mut AMLs consisted of 4,090 genes, including SETBP1, UHRF1, and MEIS2 (Data S5), all of which play roles in epigenetic regulation, hematopoiesis, or leukemia development (37–39). Our group has previously shown convergent methylation discordance over UHRF1 across cytogenetic subtypes of ALL (16). We 25 also identified a set of 24 genes with UC DMRs in all AML subtypes (Data S5). This gene set included IRX2, a member of the Iroquois homeobox gene family previously found to be predictive of outcome in infant ALL (40); DOK6, for which promoter methylation may serve as a prognostic biomarker in AML (41); and DUSP1, which is linked to therapy resistance in BCR- ABL1 chronic myeloid leukemia (42). 30 Surprisingly, DNMT3A-mutant AMLs had fewer genes associated with DMRs than TET2- mutants (837 versus 1,627 UC DMR genes; Figure 2A), despite having a greater number of DMRs in all (3,358 versus 2,752 UC DMRs; Data S1). This suggests that methylation discordance mediated by DNMT3A mutations in AML occurs mostly in non-promoter regions, a 35

that is supported by the unique enrichment of DNMT3A DMRs in open seas and repetitive elements (Figure S4). Previous reports have also found that methylation changes mediated by DNMT3A mutations target non-CGI features (43), corroborating these results. We observed that the set of genes associated with UC DMRs in CEBPA-dm and IDH1/2-mut 40 AMLs contained genes recurrently mutated in AML. Based on this observation, we explored whether other subtypes also exhibited methylation discordance over genes commonly altered in AML. We generated a list of genetic drivers of AML, encompassing genes frequently mutated involved in translocations in AML (1); 47 genetic drivers were covered in the CPEL analysis. Surprisingly, all AML subtypes demonstrated increased methylation stochasticity over key 45 drivers of AML, even in the absence of mutations in those genes (Figure 2B). In total, 38 of the (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 6 47 genetic drivers had increased methylation stochasticity in at least one AML subtype. We found significant enrichment of AML driver genes in the list of genes with UC DMRs for CEBPA-dm, IDH1/2-mut, KIT-mut, and TET2-mut AMLs (hypergeometric p < 0.01). This suggests that drivers of AML may be disrupted either genetically or epigenetically, thus contributing to leukemogenesis. This result mirrors previous work in ALL, which revealed 5 significant methylation discordance over chromosomal translocation genes across ALL cases, independent of cytogenetic status (16). Next, we performed over-representation analysis (ORA) and gene set enrichment analysis (GSEA) of genes with DMRs. ORA revealed enrichment of genes with UC DMRs in targets of 10 PRC2 complex subunits (EZH2, SUZ12) and other transcription factors such as REST, TRIM28, and FOXA1, which are involved in pluripotency and chromatin organization (44–48). Genes with UC DMRs in CEBPA-dm and IDH1/2-mut AMLs were also strongly enriched in NANOG targets, and all subtypes except DNMT3A- and NRAS-mut AMLs were enriched in SMAD4 targets (Figure S6C, Data S6). GSEA of the Molecular Signatures database (MSigDB) C2 15 collection revealed that genes with high UC were enriched in gene sets involving hematopoietic stem cell function, leukemia signatures, and H3K27me3 or H3K4me3 marks (Figure 2C, Data S7). We observed enrichment in MSigDB Hallmark gene sets relevant to oncogenesis such as epithelial-to-mesenchymal transition, hypoxia, IL2-STAT5 signaling, and KRAS signaling (Figure S6E, Data S8). Genes with dMML or dNME DMRs showed similar enrichment in 20 targets of the PRC2 complex and transcription factors regulating pluripotency, as well as in gene sets relevant to leukemogenesis and chromatin marks (Figure S6C-E, Data S6-8). These results demonstrate convergent disruption of methylation over crucial drivers of leukemia across all AML subtypes, regardless of genetic mutation. The enrichment of genes with high methylation stochasticity in gene sets relevant to leukemia biology suggests that stochastic 25 epigenetic landscape disruption occurs over a core set of leukemia regulators, potentially driving leukemogenesis across AMLs. CEBPA-dm and IDH1/2-mut AMLs had distinctive profiles of methylation stochasticity, indicating that these subtypes have additional factors mediating epigenetic disruption. Together, these findings highlight the role of methylation stochasticity in shaping the epigenetic landscape to mediate leukemic transformation. 30 Transcriptional dysregulation in AML occurs over epigenetically-disrupted regulators of leukemia To assess the influence of DNA methylation on gene expression in AML, we performed single- cell RNA sequencing on a cohort of primary AMLs obtained from Johns Hopkins Hospital and 35 Stanford University (Methods). This cohort included AMLs with mutations in DNMT3A (n=5), IDH1/2 (n=4), TET2 (n=7), CEBPA-dm (n=5), plus healthy CD34+ progenitors (CD34+ and GMP) for normal reference (n=3) (Data S2). Uniform manifold approximation and projection (UMAP) analysis of scRNA-seq data identified ten clusters. Samples were dominant in either Cluster 0 (progenitor-like) or Cluster 1 (monocyte-like) (Figure 3A, S7). All CEBPA-dm 40 samples were Cluster 0-dominant, suggesting that these AMLs have an undifferentiated, progenitor-like expression profile (Figure S7B-C). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 7 First, we performed differential expression analysis between each AML subtype and normal progenitors within UMAP clusters 0, 2, and 4 (Methods). These clusters are of particular interest due to their progenitor-like expression (Figure S7A). We observed the strongest transcriptional dysregulation in CEBPA-dm AMLs (1,081; 863; and 652 differentially expressed genes (DEGs) in clusters 0, 2, and 4 respectively) (Figure 3B, S8A-B; Data S9). Surprisingly, despite the high 5 methylation stochasticity observed in IDH-mutant AMLs, this subtype had minimal transcriptional disruption (162; 178; and 33 DEGs in clusters 0, 2, and 4 respectively), while DNMT3A-mutant AMLs, despite having only modest changes in methylation landscape, had stronger dysregulation of gene expression (697; 401; and 614 DEGs in clusters 0, 2, and 4 respectively) (Figure 3B, S8A-B; Data S9). 10 We identified five genes with differential expression across all AML subtypes, including ADGRE5 (CD97) upregulated in cluster 0 and BAHCC1 upregulated in cluster 4 (Figure 3C, S8C). Most of these genes had discordant methylation stochasticity in AMLs; for example, BAHCC1 was associated with UC DMRs in CEBPA-dm, TET2-mutant, and KIT-mutant AMLs 15 (Data S3). In fact, 4 of the 5 common differentially expressed genes had DMRs in CEBPA-dm AMLs. ADGRE5 has been identified as a crucial regulator of leukemic stem cells (LSCs) (49), and BAHCC1 was reported to sustain leukemogenesis through its function as a reader of the Polycomb mark H3K27me3 (50). Known genetic drivers of AML were also differentially expressed in CEBPA-dm and DNMT3A-mutant AMLs (Data S10). For example, CBLB, JAK2 20 and KDM6A were upregulated in both CEBPA-dm and DNMT3A-mutant AMLs in clusters 0, 2, and 4 (Figure S8D). Several of these genes had increased methylation stochasticity, such as JAK2 in CEBPA-dm and DNMT3A-mutant AMLs (Figure 2B). ORA revealed convergent enrichment of DEGs in targets of the PRC2 complex (SUZ12) and in 25 SMAD4 targets across clusters 0, 2, and 4 (Data S11). This parallels the enrichment of PRC2 and SMAD4 targets observed in DMR associated genes (Figure S6C). In agreement with the enrichment over genes with high methylation stochasticity (Figure 2C, S6D-E), GSEA revealed convergent enrichment of DEGs in all AML subtypes over leukemia signatures from the MSigDB C2 collection (Figure 3D, S8E; Data S12). We also observed a consistent upregulation 30 of gene sets relevant to oncogenesis from the MSigDB Hallmark collection, including PI3K/AKT, TGF-beta, and KRAS signaling (Figure 3E, S8F; Data S13). Interestingly, we MYC gene sets from the MSigDB Hallmark collection were downregulated across AML subtypes in all clusters analyzed, despite the known role of MYC deregulation in leukemogenesis (51). 35 Taken together, these results demonstrate the transcriptional disruption of key leukemia regulators across AML subtypes. The convergence of discordant methylation stochasticity and transcriptional dysregulation over similar gene targets and leukemia signatures highlights the important role of epigenetic landscape disruption on the transcriptional phenotype in AML, 40 suggesting that transcriptional changes linked to methylation stochasticity may be a defining feature of leukemogenesis. DNA methylation stochasticity is associated with transcriptional dysregulation in AML 45 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 8 To further investigate the connection between methylation discordance and transcriptional disruption in AML, we computed the mean-adjusted variability (MAV), a measure of stochastic gene expression, for each detected gene in each sample (Methods, Data S14). We observed a strong relationship between NME near the TSS and quartiles of genes ranked by MAV, where genes with higher MAV had increased NME near the TSS (Figure 4A, S9). Gene expression 5 mean was inversely related to mean methylation near the TSS, where genes with higher mean expression had lower MML near the TSS. This result links DNA methylation stochasticity to transcriptional variability, providing an explanation for the convergent enrichment of the two modes over key regulators of leukemia, and suggesting how these changes in DNA methylation may affect the transcriptional and therefore phenotypic state of AML. 10 Finally, we employed a multivariate information-theoretic method of gene regulatory network inference (52) to explore dependencies among the 38 genetic drivers of AML with significantly discordant methylation stochasticity (Figure 2B). We assessed CEBPA-dm, IDH1/2-mut, and TET2-mut AMLs, as these subtypes showed statistically significant overlaps between genetic 15 drivers of AML and genes with significant methylation stochasticity. This analysis revealed a complex regulatory network among 30 of these drivers of AML (Figure 4B). In this network, RAD21 was predicted to have the highest betweenness centrality, indicating that it exerts the greatest influence over the network. Interestingly, RAD21 was directly connected to IDH2, and IDH2 was directly connected to CEBPA, the genes defining AML subtypes with markedly high 20 methylation stochasticity. RAD21 is a subunit of the cohesin complex, with an essential role in mediating chromatin structure through formation of loops and topologically-associating domains (TADs) (53). The central position of this chromatin regulator in the gene regulatory network suggests that epigenetic and transcriptional dysregulation of key mediators of genome organization may influence leukemogenesis. 25 These results show that high methylation stochasticity is linked to transcriptional dysregulation in AML, with increased methylation entropy associated with higher gene expression variability. Additionally, we identify a gene regulatory network between non-mutated AML driver genes containing DMRs, suggesting that key drivers of AML may be genetically or epigenetically 30 disrupted in leukemia. Together, these findings provide evidence of how methylation discordance may directly contribute to transcriptional profiles driving leukemogenesis.

35 Methylation stochasticity has been shown to be associated with clinical outcome, chemoresistance, and disease progression in hematopoietic malignancies (16, 18–22). A recent study reported that somatic mutations are associated with epigenetic heterogeneity in AML, suggesting a relationship between increased epigenetic variability and inferior prognosis (22). Surprisingly, in this study IDH1/2 mutant AML had broad distributions of epialleles, without 40 characteristic clustering as seen for other genetic subtypes of AML (22). To investigate the contribution of specific, mutually-exclusive genetic mutations to stochastic disruption of the epigenetic landscape, we conducted an information-theoretic analysis of methylation potential energy landscapes in genetic subtypes of AML. Our analyses identified a core program of methylation stochasticity and transcriptional disruption involving key regulators of leukemia 45 biology across all AMLs. Furthermore, we link methylation entropy to gene expression (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 9 variability, suggesting that transcriptional dysregulation associated with methylation stochasticity may be a hallmark of leukemogenesis, regardless of genetic mutations. All subtypes of AML had significant methylation discordance over genetic drivers of AML, even in the absence of mutations in those genetic drivers. These genes formed a complex gene regulatory network, indicating that disruption of genetic drivers of AML by methylation stochasticity may 5 drive dysregulation of the transcriptional state of crucial regulators of leukemia. We additionally found that CEBPA-dm and IDH1/2-mutant AMLs have markedly increased DNA methylation entropy beyond this core leukemia program, suggesting that these subtypes have additional factors driving dysregulation of the epigenetic landscape. 10 Our observation of markedly increased methylation entropy in CEBPA-dm and IDH1/2-mutant AMLs is in agreement with findings by Li et al. (22), who reported a high degree of epigenetic heterogeneity in specific AML subtypes including CEBPA-dm and IDH1/2-mutated. Here, we describe a distinct and surprisingly convergent profile of methylation stochasticity in both CEBPA-dm and IDH1/2-mut AMLs. This complements previous findings of increased epiallele 15 heterogeneity in these AML subtypes by providing gene-level resolution, demonstrating the sensitivity of information-theoretic methods for methylation PEL analysis, which can capture higher-order statistical properties of methylation that may not be detected by empirical estimations of methylation entropy (23–25). In CEBPA-dm and IDH1/2-mut AMLs, we observed convergent enrichment of Homeobox, AP-2, and LRF transcription factor binding 20 motifs over regions with stochastic methylation disruption. These transcription factors are known to mediate leukemogenesis: overexpression or translocations involving Homeobox genes contribute to malignancy (31); AP-2 can negatively regulate CEBPA and may contribute to HOX dysregulation in AML (32, 54); and LRF (ZBTB7A) is a coordinator of hematopoietic differentiation (33). At the gene level, methylation landscape disruption in CEBPA-dm and 25 IDH1/2-mutant AMLs was also highly convergent, occurring over genes such as SETBP1, UHRF1, and MEIS2. These genes play crucial roles in epigenetic regulation, leukemia, or hematopoiesis: SETBP1 recruits chromatin modifiers to control gene expression, and is frequently mutated in leukemia (37); UHRF1 coordinates DNA methyltransferase DNMT1 to repressive chromatin marks, and was identified as a convergent target of stochastic DNA 30 methylation in ALL (16, 38); and MEIS2 is an important factor in myeloid differentiation (39). CEBPA- and IDH-mutant AMLs have previously been characterized as hypermethylated; IDH mutations in particular lead to a well-described hypermethylated profile (2, 7, 8). Mutant IDH enzymes produce the oncometabolite 2- hydroxyglutarate, which competitively inhibits alpha-35 ketoglutarate-dependent enzymes, including TET family members and histone lysine demethylases (55). The effect of IDH mutations on prognosis of AMLs is unclear, and may depend on co-occurring mutations, treatment, and other patient characteristics (30). Additionally, small-molecule inhibitors of IDH1 or IDH2 may be used in the treatment of IDH-mutant AMLs; these been reported to reduce methylation levels in human AML (55) and to lead to decreased 40 epigenetic allele diversity in a mouse model of AML (22). AMLs with biallelic or bZIP-domain mutations in CEBPA are a distinct subtype classified by the World Health Organization (56), and are associated with a favorable prognosis (57). Interestingly, mutations in CEBPA have been associated with resistance to IDH inhibitors (55). Aside from its canonical function as a transcription factor, CEBPA is also known to interact with epigenetic machinery, such as the 45 histone acetyltransferase p300 and the SWI/SNF complex (58–60). Recent work has also (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 10 implicated CEBPA in the control of DNA methylation through direct interaction with DNMT3A, where the interaction between wild-type CEBPA and DNMT3A blocks association with DNA; this interaction was shown to be lost upon CEBPA mutation, leading to aberrant hypermethylation (61). The disruption of interactions between CEBPA and epigenetic regulators upon CEBPA mutation may provide a basis for the observed increase in methylation 5 stochasticity, potentially contributing to the highly disordered epigenetic landscape of CEBPA- dm AMLs. A similar process may underlie the defined patterns of epigenetic alleles previously observed in AMLs with other canonical transcription factor mutations, such as t(8;21), inv(16), and t(15;17) (22). The disparate prognostic impact of CEBPA and IDH1/2 mutations indicates that increased methylation stochasticity is not always a marker of unfavorable outcome; rather, it 10 may reflect widespread deregulation of the epigenetic machinery, whether by oncometabolites (as in the case of IDH mutations) or by disruption of crucial protein interactions (as in the case of CEBPA mutations). DNMT3A-mut and TET2-mut AMLs had only modest discordance of methylation stochasticity, 15 despite the well-defined functions of these enzymes as direct modifiers of DNA methylation (29, 30). Nevertheless, methylation discordance in these subtypes localized to important regulators of malignancy. It is possible that mutations in these genes mediates changes in DNA methylation outside of regions covered by the ERRBS data analyzed here. Indeed, the enrichment of DMRs in DNMT3A-mut AMLs over open seas suggests that these mutations may predominantly affect 20 non-CGI regions. This observation aligns with previous studies showing differential methylation over open seas (43) and an absence of CpG island disruption (5) in DNMT3A-mutant AMLs. Likewise, TET2 mutations have been shown to affect primarily non-CGI methylation (62), particularly methylation over enhancer regions (63), which may not be captured by ERRBS. Further investigation into the effect of mutations in these enzymes on DNA methylation 25 stochasticity is warranted. We identified a surprising convergence of discordant methylation stochasticity over key regulators of leukemia and pluripotency across all AML subtypes. These included targets of the PRC2 complex, REST, and TRIM28, which are factors crucial for epigenetic regulation and 30 pluripotency (45–48). Genes with methylation discordance were also enriched in gene sets relevant to oncogenesis and leukemia across AML subtypes. Notably, we also observed significantly increased methylation stochasticity over genetic drivers of AML such as CREBBP, GATA2, and RUNX1T1 (1) in all AML subtypes, independent of mutations in those driver genes. In contrast, Li et al. found epigenetic heterogeneity over mediators of leukemogenesis only in 35 specific subtypes of AML with particularly high epigenetic allele diversity (22). Information- theoretic methods therefore may allow more sensitive detection of stochastic dysregulation of the epigenetic landscape. These results suggest disruption of the epigenetic landscape over critical drivers of leukemogenesis may be a core feature of leukemogenesis across AMLs. 40 Like methylation stochasticity, transcriptional dysregulation was convergently enriched over key regulators of leukemia across AML subtypes, suggesting a link between methylation disruption and transcriptional dysregulation. We further support this result by demonstrating a relationship between methylation stochasticity and gene expression, where increased methylation levels correlate with lower mean gene expression, and higher methylation entropy is associated with 45 greater gene expression variability. This corroborates the previous observation that higher (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 11 epiallele diversity corresponded with higher transcriptional variance across AML samples (22). The connection between methylation stochasticity and gene expression suggests that methylation discordance may set the stage for aberrant transcriptional programs in leukemia. We identified a complex gene regulatory network between known genetic drivers of AML. 5 These genes also had stochastic disruption of the methylation landscape, suggesting that epigenetic regulation of genetic drivers of AML, even in the absence of mutations in these genes, may contribute to leukemogenesis. This network was centered on RAD21, a subunit of the cohesin complex that is crucial for formation of chromatin loops and TADs (53), and is mutated in 3% of AMLs (1). In this network, RAD21 was closely connected to IDH2 and CEBPA, 10 suggesting a role of mutations in these genes in transcriptional disruption of RAD21. CEBPA stabilizes binding of RAD21 to maintain looping between enhancers and promoters in breast cancer cells (64). The central position of this chromatin regulator in the gene regulatory network suggests that chromatin organization may be disrupted by epigenetic and transcriptional dysregulation of key mediators of genome organization in AML. 15 In conclusion, we leverage high-coverage ERRBS and single-cell RNA-seq data to analyze DNA methylation stochasticity in distinct genetic subtypes of AML, employing an information- theoretic method to construct potential energy landscapes that encompass the higher-order statistical properties of methylation, including both mean methylation and methylation entropy. 20 This approach identified CEBPA-dm and IDH1/2-mutant AMLs as markedly high-entropy subtypes with distinctive and convergent profiles of epigenetic stochasticity. We further identified a core program of methylation stochasticity and gene expression dysregulation over key regulators of the leukemic phenotype across AML subtypes, along with a relationship between methylation entropy and gene expression variability, suggesting that leukemogenesis 25 may be mediated by stochastic disruption of the epigenetic landscape independent of genotype. Taken together, our results indicate that information-theoretic analysis of DNA methylation can elucidate the role of stochastic epigenetic regulation and provide insights into epigenetic drivers of genetically-defined subtypes of AML. 30 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 12

1. E. Papaemmanuil, M. Gerstung, L. Bullinger, V . I. Gaidzik, P. Paschka, N. D. Roberts, N. E. Potter, M. Heuser, F. Thol, N. Bolli, G. Gundem, P. Van Loo, I. Martincorena, P. Ganly, L. 5 Mudie, S. McLaren, S. O’Meara, K. Raine, D. R. Jones, J. W. Teague, A. P. Butler, M. F. Greaves, A. Ganser, K. Döhner, R. F. Schlenk, H. Döhner, P. J. Campbell, Genomic Classification and Prognosis in Acute Myeloid Leukemia. N. Engl. J. Med. 374, 2209–2221 (2016). 2. The Cancer Genome Atlas Research Network, Genomic and Epigenomic Landscapes of 10 Adult De Novo Acute Myeloid Leukemia. N. Engl. J. Med. 368, 2059–2074 (2013). 3. M. Schwede, K. Jahn, J. Kuipers, L. A. Miles, R. L. Bowman, T. Robinson, K. Furudate, H. Uryu, T. Tanaka, Y . Sasaki, A. Ediriwickrema, B. Benard, A. J. Gentles, R. Levine, N. Beerenwinkel, K. Takahashi, R. Majeti, Mutation order in acute myeloid leukemia identifies uncommon patterns of evolution and illuminates phenotypic heterogeneity. Leukemia 38, 15 1501–1510 (2024). 4. P. Desai, N. Mencia-Trinchant, O. Savenkov, M. S. Simon, G. Cheang, S. Lee, M. Samuel, E. K. Ritchie, M. L. Guzman, K. V . Ballman, G. J. Roboz, D. C. Hassane, Somatic mutations precede acute myeloid leukemia years before diagnosis. Nat. Med. 24, 1015–1023 (2018). 5. D. H. Spencer, D. A. Russler-Germain, S. Ketkar, N. M. Helton, T. L. Lamprecht, R. S. 20 Fulton, C. C. Fronick, M. O’Laughlin, S. E. Heath, M. Shinawi, P. Westervelt, J. E. Payton, L. D. Wartman, J. S. Welch, R. K. Wilson, M. J. Walter, D. C. Link, J. F. DiPersio, T. J. Ley, CpG Island Hypermethylation Mediated by DNMT3A Is a Consequence of AML Progression. Cell 168, 801-816.e13 (2017). 6. J. L. Glass, D. Hassane, B. J. Wouters, H. Kunimoto, R. Avellino, F. E. Garrett-Bakelman, O. 25 A. Guryanova, R. Bowman, S. Redlich, A. M. Intlekofer, C. Meydan, T. Qin, M. Fall, A. Alonso, M. L. Guzman, P. J. M. Valk, C. B. Thompson, R. Levine, O. Elemento, R. Delwel, A. Melnick, M. E. Figueroa, Epigenetic Identity in AML Depends on Disruption of Nonpromoter Regulatory Elements and Is Affected by Antagonistic Effects of Mutations in Epigenetic Modifiers. Cancer Discov. 7, 868–883 (2017). 30 7. E. R. Wilson, N. M. Helton, S. E. Heath, R. S. Fulton, J. E. Payton, J. S. Welch, M. J. Walter, P. Westervelt, J. F. DiPersio, D. C. Link, C. A. Miller, T. J. Ley, D. H. Spencer, Focal disruption of DNA methylation dynamics at enhancers in IDH-mutant AML cells. Leukemia 36, 935–945 (2022). 8. M. E. Figueroa, O. Abdel-Wahab, C. Lu, P. S. Ward, J. Patel, A. Shih, Y . Li, N. Bhagwat, A. 35 Vasanthakumar, H. F. Fernandez, M. S. Tallman, Z. Sun, K. Wolniak, J. K. Peeters, W. Liu, S. E. Choe, V . R. Fantin, E. Paietta, B. Löwenberg, J. D. Licht, L. A. Godley, R. Delwel, P. J. M. Valk, C. B. Thompson, R. L. Levine, A. Melnick, Leukemic IDH1 and IDH2 Mutations

in a Hypermethylation Phenotype, Disrupt TET2 Function, and Impair Hematopoietic Differentiation. Cancer Cell 18, 553–567 (2010). 40 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 13 9. M. E. Figueroa, S. Lugthart, Y . Li, C. Erpelinck-Verschueren, X. Deng, P. J. Christos, E. Schifano, J. Booth, W. van Putten, L. Skrabanek, F. Campagne, M. Mazumdar, J. M. Greally, P. J. M. Valk, B. Löwenberg, R. Delwel, A. Melnick, DNA Methylation Signatures Identify Biologically Distinct Subtypes in Acute Myeloid Leukemia. Cancer Cell 17, 13–27 (2010). 10. S. Deneberg, P. Guardiola, A. Lennartsson, Y . Qu, V . Gaidzik, O. Blanchet, M. Karimi, S. 5 Bengtzén, H. Nahi, B. Uggla, U. Tidefelt, M. Höglund, C. Paul, K. Ekwall, K. Döhner, S. Lehmann, Prognostic DNA methylation patterns in cytogenetically normal acute myeloid leukemia are predefined by stem cell chromatin marks. Blood 118, 5573–5582 (2011). 11. C. Gebhard, D. Glatz, L. Schwarzfischer, J. Wimmer, S. Stasik, M. Nuetzel, D. Heudobler, R. Andreesen, G. Ehninger, C. Thiede, M. Rehli, Profiling of aberrant DNA methylation in 10 acute myeloid leukemia reveals subclasses of CG-rich regions with epigenetic or genetic association. Leukemia 33, 26–36 (2019). 12. B. Giacopelli, M. Wang, A. Cleary, Y .-Z. Wu, A. R. Schultz, M. Schmutz, J. S. Blachly, A.-K. Eisfeld, B. Mundy-Bosse, S. V osberg, P. A. Greif, R. Claus, L. Bullinger, R. Garzon, K. R. Coombes, C. D. Bloomfield, B. J. Druker, J. W. Tyner, J. C. Byrd, C. C. Oakes, DNA 15 methylation epitypes highlight underlying developmental and disease pathways in acute myeloid leukemia. Genome Res. 31, 747–761 (2021). 13. G. Yamato, T. Kawai, N. Shiba, J. Ikeda, Y . Hara, K. Ohki, S.-I. Tsujimoto, T. Kaburagi, K. Yoshida, Y . Shiraishi, S. Miyano, N. Kiyokawa, D. Tomizawa, A. Shimada, M. Sotomatsu, H. Arakawa, S. Adachi, T. Taga, K. Horibe, S. Ogawa, K. Hata, Y . Hayashi, Genome-wide 20 DNA methylation analysis in pediatric acute myeloid leukemia. Blood Adv. 6, 3207–3219 (2022). 14. E. K. W. Lo, B. M. Mears, H. C. Maurer, A. Idrizi, K. D. Hansen, E. D. Thompson, R. H. Hruban, K. P. Olive, A. P. Feinberg, Comprehensive DNA Methylation Analysis Indicates That Pancreatic Intraepithelial Neoplasia Lesions Are Acinar-Derived and Epigenetically 25 Primed for Carcinogenesis. Cancer Res. 83, 1905–1916 (2023). 15. M. A. Koldobskiy, J. Abante, G. Jenkinson, E. Pujadas, A. Tetens, F. Zhao, R. Tryggvadottir, A. Idrizi, A. Reinisch, R. Majeti, J. Goutsias, A. P. Feinberg, A Dysregulated DNA Methylation Landscape Linked to Gene Expression in MLL-Rearranged AML. Epigenetics 15, 841–858 (2020). 30 16. M. A. Koldobskiy, G. Jenkinson, J. Abante, V . A. Rodriguez DiBlasi, W. Zhou, E. Pujadas, A. Idrizi, R. Tryggvadottir, C. Callahan, C. L. Bonifant, K. R. Rabin, P. A. Brown, H. Ji, J. Goutsias, A. P. Feinberg, Converging genetic and epigenetic drivers of paediatric acute lymphoblastic leukaemia identified by an information-theoretic analysis. Nat. Biomed. Eng. 5, 360–376 (2021). 35 17. A. P. Feinberg, A. Levchenko, Epigenetics as a mediator of plasticity in cancer. Science 379 (2023). 18. D. A. Landau, K. Clement, M. J. Ziller, P. Boyle, J. Fan, H. Gu, K. Stevenson, C. Sougnez, L. Wang, S. Li, D. Kotliar, W. Zhang, M. Ghandi, L. Garraway, S. M. Fernandes, K. J. Livak, (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 14 S. Gabriel, A. Gnirke, E. S. Lander, J. R. Brown, D. Neuberg, P. V . Kharchenko, N. Hacohen, G. Getz, A. Meissner, C. J. Wu, Locally Disordered Methylation Forms the Basis of Intratumor Methylome Variation in Chronic Lymphocytic Leukemia. Cancer Cell 26, 813– 825 (2014). 19. S. Li, F. E. Garrett-Bakelman, S. S. Chung, M. A. Sanders, T. Hricik, F. Rapaport, J. Patel, R. 5 Dillon, P. Vijay, A. L. Brown, A. E. Perl, J. Cannon, L. Bullinger, S. Luger, M. Becker, I. D. Lewis, L. B. To, R. Delwel, B. Löwenberg, H. Döhner, K. Döhner, M. L. Guzman, D. C. Hassane, G. J. Roboz, D. Grimwade, P. J. M. Valk, R. J. D’Andrea, M. Carroll, C. Y . Park, D. Neuberg, R. Levine, A. M. Melnick, C. E. Mason, Distinct evolution and dynamics of epigenetic and genetic heterogeneity in acute myeloid leukemia. Nat. Med. 22, 792–799 10 (2016). 20. P. Wahlberg, A. Lundmark, J. Nordlund, S. Busche, A. Raine, K. Tandre, L. Rönnblom, D. Sinnett, E. Forestier, T. Pastinen, G. Lönnerholm, A.-C. Syvänen, DNA Methylome Analysis of Acute Lymphoblastic Leukemia Cells Reveals Stochastic de novo DNA methylation in CpG islands. Epigenomics 8, 1367–1387 (2016). 15 21. A. Magi, G. Mattei, A. Mingrino, C. Caprioli, C. Ronchini, G. Frigè, R. Semeraro, D. Bolognini, A. Rambaldi, A. Candoni, E. Colombo, L. Mazzarella, P. G. Pelicci, High- resolution Nanopore methylome-maps reveal random hyper-methylation at CpG-poor regions as driver of chemoresistance in leukemias. Commun. Biol. 6, 1–11 (2023). 22. S. Li, X. Chen, J. Wang, C. Meydan, J. L. Glass, A. H. Shih, R. Delwel, R. L. Levine, C. E. 20 Mason, A. M. Melnick, Somatic Mutations Drive Specific, but Reversible, Epigenetic Heterogeneity States in AML. Cancer Discov. 10, 1934–1949 (2020). 23. G. Jenkinson, E. Pujadas, J. Goutsias, A. P. Feinberg, Potential energy landscapes identify the information-theoretic nature of the epigenome. Nat. Genet. 49, 719–729 (2017). 24. G. Jenkinson, J. Abante, A. P. Feinberg, J. Goutsias, An information-theoretic approach to the 25 modeling and analysis of whole-genome bisulfite sequencing data. BMC Bioinformatics 19, 87 (2018). 25. J. Abante, Y . Fang, A. P. Feinberg, J. Goutsias, Detection of haplotype-dependent allele- specific DNA methylation in WGBS data. Nat. Commun. 11, 5238 (2020). 26. J. Abante, J. Goutsias, CpelTdm.jl: a Julia package for targeted differential DNA methylation 30 analysis. bioRxiv [Preprint] (2020). https://doi.org/10.1101/2020.10.17.343020. 27. R. Rampal, A. Alkalin, J. Madzo, A. Vasanthakumar, E. Pronier, J. Patel, Y . Li, J. Ahn, O. Abdel-Wahab, A. Shih, C. Lu, P. S. Ward, J. J. Tsai, T. Hricik, V . Tosello, J. E. Tallman, X. Zhao, D. Daniels, Q. Dai, L. Ciminio, I. Aifantis, C. He, F. Fuks, M. S. Tallman, A. Ferrando, S. Nimer, E. Paietta, C. B. Thompson, J. D. Licht, C. E. Mason, L. A. Godley, A. Melnick, 35 M. E. Figueroa, R. L. Levine, DNA Hydroxymethylation Profiling Reveals that WT1 Mutations Result in Loss of TET2 Function in Acute Myeloid Leukemia. Cell Rep. 9, 1841– 1855 (2014). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 15 28. T. Qin, J. Sotzen, R. K. Rampal, F. T. Rapaport, R. L. Levine, V . Klimek, S. D. Nimer, M. E. Figueroa, Risk of disease progression in low-risk MDS is linked to distinct epigenetic subtypes. Leukemia 33, 2753–2757 (2019). 29. X. Wu, Y . Zhang, TET-mediated active DNA demethylation: mechanism, function and beyond. Nat. Rev. Genet. 18, 517–534 (2017). 5 30. A. P. Im, A. R. Sehgal, M. P. Carroll, B. D. Smith, A. Tefferi, D. E. Johnson, M. Boyiadzis, DNMT3A and IDH mutations in acute myeloid leukemia and other myeloid malignancies: associations with prognosis and potential treatment strategies. Leukemia 28, 1774–1783 (2014). 31. R. A. Alharbi, R. Pettengell, H. S. Pandha, R. Morgan, The role of HOX genes in normal 10 hematopoiesis and acute leukemia. Leukemia 27, 1000–1008 (2013). 32. X. Ding, Z. Yang, F. Zhou, F. Wang, X. Li, C. Chen, X. Li, X. Hu, S. Xiang, J. Zhang, Transcription factor AP-2α regulates acute myeloid leukemia cell proliferation by influencing Hoxa gene expression. Int. J. Biochem. Cell Biol. 45, 1647–1656 (2013). 33. E. Redondo Monte, A. Wilding, G. Leubolt, P. Kerbs, J. W. Bagnoli, L. Hartmann, W. 15 Hiddemann, L. Chen-Wichmann, S. Krebs, H. Blum, M. Cusan, B. Vick, I. Jeremias, W. Enard, S. Theurich, C. Wichmann, P. A. Greif, ZBTB7A prevents RUNX1-RUNX1T1- dependent clonal expansion of human hematopoietic stem and progenitor cells. Oncogene 39, 3195–3205 (2020). 34. R. C. Fisher, E. W. Scott, Role of PU.1 in Hematopoiesis. Stem Cells 16, 25–37 (1998). 20 35. Y.-H. Kuo, S. K. Zaidi, S. Gornostaeva, T. Komori, G. S. Stein, L. H. Castilla, Runx2 induces acute myeloid leukemia in cooperation with Cbfβ-SMMHC in mice. Blood 113, 3323–3332 (2009). 36. E. Passegué, E. F. Wagner, I. L. Weissman, JunB Deficiency Leads to a Myeloproliferative Disorder Arising from Hematopoietic Stem Cells. Cell 119, 431–443 (2004). 25 37. R. Piazza, V . Magistroni, S. Redaelli, M. Mauri, L. Massimino, A. Sessa, M. Peronaci, M. Lalowski, R. Soliymani, C. Mezzatesta, A. Pirola, F. Banfi, A. Rubio, D. Rea, F. Stagno, E. Usala, B. Martino, L. Campiotti, M. Merli, F. Passamonti, F. Onida, A. Morotti, F. Pavesi, M. Bregni, V . Broccoli, M. Baumann, C. Gambacorti-Passerini, SETBP1 induces transcription of a network of development genes by acting as an epigenetic hub. Nat. Commun. 9, 2192 30 (2018). 38. H. Sidhu, N. Capalash, UHRF1: The key regulator of epigenetics and molecular target for cancer therapeutics. Tumor Biol. 39, 1010428317692205 (2017). 39. T. Fujino, Y . Yamazaki, D. A. Largaespada, N. A. Jenkins, N. G. Copeland, K. Hirokawa, T. Nakamura, Inhibition of myeloid differentiation by Hoxa9, Hoxb8, and Meis homeobox 35 genes. Exp. Hematol. 29, 856–863 (2001). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 16 40. H. Kang, C. S. Wilson, R. C. Harvey, I.-M. Chen, M. H. Murphy, S. R. Atlas, E. J. Bedrick, M. Devidas, A. J. Carroll, B. W. Robinson, R. W. Stam, M. G. Valsecchi, R. Pieters, N. A. Heerema, J. M. Hilden, C. A. Felix, G. H. Reaman, B. Camitta, N. Winick, W. L. Carroll, Z. E. Dreyer, S. P. Hunger, C. L. Willman, Gene expression profiles predictive of outcome and age in infant acute lymphoblastic leukemia: a Children’s Oncology Group study. Blood 119, 5 1872–1881 (2012). 41. G.-K. Sun, L.-J. Tang, J.-D. Zhou, Z.-J. Xu, L. Yang, Q. Yuan, J.-C. Ma, X.-H. Liu, J. Lin, J. Qian, D.-M. Yao, DOK6 promoter methylation serves as a potential biomarker affecting prognosis in de novo acute myeloid leukemia. Cancer Med. 8, 6393–6402 (2019). 42. M. Kesarwani, Z. Kincaid, A. Gomaa, E. Huber, S. Rohrabaugh, Z. Siddiqui, M. F. Bouso, T. 10 Latif, M. Xu, K. Komurov, J. C. Mulloy, J. A. Cancelas, H. L. Grimes, M. Azam, Targeting c-FOS and DUSP1 abrogates intrinsic resistance to tyrosine-kinase inhibitor therapy in BCR- ABL-induced leukemia. Nat. Med. 23, 472–482 (2017). 43. Y . Qu, A. Lennartsson, V . I. Gaidzik, S. Deneberg, M. Karimi, S. Bengtzén, M. Höglund, L. Bullinger, K. Döhner, S. Lehmann, Differential methylation in CN-AML preferentially 15 targets non-CGI regions and is dictated by DNMT3A mutational status and associated with predominant hypomethylation of HOX genes. Epigenetics 9, 1108–1119 (2014). 44. L. A. Cirillo, F. R. Lin, I. Cuesta, D. Friedman, M. Jarnik, K. S. Zaret, Opening of Compacted Chromatin by Early Developmental Transcription Factors HNF3 (FoxA) and GATA-4. Mol. Cell 9, 279–289 (2002). 20 45. S. K. Singh, M. N. Kagalwala, J. Parker-Thornburg, H. Adams, S. Majumder, REST maintains self-renewal and pluripotency of embryonic stem cells. Nature 453, 223–227 (2008). 46. D. C. Miles, N. A. de Vries, S. Gisler, C. Lieftink, W. Akhtar, E. Gogola, I. Pawlitzky, D. Hulsman, E. Tanger, M. Koppens, R. L. Beijersbergen, M. van Lohuizen, TRIM28 Is an 25 Epigenetic Barrier to Induced Pluripotent Stem Cell Reprogramming. Stem Cells 35, 147– 157 (2017). 47. S. M. Jang, A. Kauzlaric, J.-P. Quivy, J. Pontis, B. Rauwel, A. Coluccio, S. Offner, J. Duc, P. Turelli, G. Almouzni, D. Trono, KAP1 facilitates reinstatement of heterochromatin after DNA replication. Nucleic Acids Res. 46, 8788–8802 (2018). 30 48. M. Sauvageau, G. Sauvageau, Polycomb Group Proteins: Multi-Faceted Regulators of Somatic Stem Cells and Cancer. Cell Stem Cell 7, 299–313 (2010). 49. G. H. Martin, N. Roy, S. Chakraborty, A. Desrichard, S. S. Chung, C. M. Woolthuis, W. Hu, I. Berezniuk, F. E. Garrett-Bakelman, J. Hamann, S. M. Devlin, T. A. Chan, C. Y . Park, CD97 is a critical regulator of acute myeloid leukemia stem cell function. J. Exp. Med. 216, 35 2362–2377 (2019). 50. H. Fan, J. Lu, Y . Guo, D. Li, Z.-M. Zhang, Y .-H. Tsai, W.-C. Pi, J. H. Ahn, W. Gong, Y . Xiang, D. F. Allison, H. Geng, S. He, Y . Diao, W.-Y . Chen, B. D. Strahl, L. Cai, J. Song, G. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 17 G. Wang, BAHCC1 binds H3K27me3 via a conserved BAH module to mediate gene silencing and oncogenesis. Nat. Genet. 52, 1384–1396 (2020). 51. B. Hoffman, A. Amanullah, M. Shafarenko, D. A. Liebermann, The proto-oncogene c-myc in hematopoietic development and leukemogenesis. Oncogene 21, 3414–3421 (2002). 52. T. E. Chan, M. P. H. Stumpf, A. C. Babtie, Gene Regulatory Network Inference from Single-5 Cell Data Using Multivariate Information Measures. Cell Syst. 5, 251-267.e3 (2017). 53. Y . Sun, X. Xu, W. Zhao, Y . Zhang, K. Chen, Y . Li, X. Wang, M. Zhang, B. Xue, W. Yu, Y . Hou, C. Wang, W. Xie, C. Li, D. Kong, S. Wang, Y . Sun, RAD21 is the core subunit of the cohesin complex involved in directing genome organization. Genome Biol. 24, 155 (2023). 54. K. Hilger-Eversheim, M. Moser, H. Schorle, R. Buettner, Regulatory roles of AP-2 10 transcription factors in vertebrate development, apoptosis and cell-cycle control. Gene 260, 1–12 (2000). 55. F. Wang, K. Morita, C. D. DiNardo, K. Furudate, T. Tanaka, Y . Yan, K. P. Patel, K. J. MacBeth, B. Wu, G. Liu, M. Frattini, J. A. Matthews, L. D. Little, C. Gumbs, X. Song, J. Zhang, E. J. Thompson, T. M. Kadia, G. Garcia-Manero, E. Jabbour, F. Ravandi, K. N. 15 Bhalla, M. Konopleva, H. M. Kantarjian, P. Andrew Futreal, K. Takahashi, Leukemia stemness and co-occurring mutations drive resistance to IDH inhibitors in acute myeloid leukemia. Nat. Commun. 12, 2607 (2021). 56. J. D. Khoury, E. Solary, O. Abla, Y . Akkari, R. Alaggio, J. F. Apperley, R. Bejar, E. Berti, L. Busque, J. K. C. Chan, W. Chen, X. Chen, W.-J. Chng, J. K. Choi, I. Colmenero, S. E. 20 Coupland, N. C. P. Cross, D. De Jong, M. T. Elghetany, E. Takahashi, J.-F. Emile, J. Ferry, L. Fogelstrand, M. Fontenay, U. Germing, S. Gujral, T. Haferlach, C. Harrison, J. C. Hodge, S. Hu, J. H. Jansen, R. Kanagal-Shamanna, H. M. Kantarjian, C. P. Kratz, X.-Q. Li, M. S. Lim, K. Loeb, S. Loghavi, A. Marcogliese, S. Meshinchi, P. Michaels, K. N. Naresh, Y . Natkunam, R. Nejati, G. Ott, E. Padron, K. P. Patel, N. Patkar, J. Picarsic, U. Platzbecker, I. 25 Roberts, A. Schuh, W. Sewell, R. Siebert, P. Tembhare, J. Tyner, S. Verstovsek, W. Wang, B. Wood, W. Xiao, C. Yeung, A. Hochhaus, The 5th edition of the World Health Organization Classification of Haematolymphoid Tumours: Myeloid and Histiocytic/Dendritic Neoplasms. Leukemia 36, 1703–1719 (2022). 57. J.-A. Georgi, S. Stasik, M. Kramer, M. Meggendorfer, C. Röllig, T. Haferlach, P. Valk, D. 30 Linch, T. Herold, N. Duployez, F. Taube, J. M. Middeke, U. Platzbecker, H. Serve, C. D. Baldus, C. Muller-Tidow, C. Haferlach, S. Koch, W. E. Berdel, B. J. Woermann, U. Krug, J. Braess, W. Hiddemann, K. Spiekermann, E. L. Boertjes, R. K. Hills, A. Burnett, G. Ehninger, K. Metzeler, M. Rothenberg-Thurley, A. Dufour, H. Dombret, C. Pautas, C. Preudhomme, L. Fenwarth, M. Bornhäuser, R. Gale, C. Thiede, Prognostic impact of CEBPA mutational 35 subgroups in adult AML. Leukemia 38, 281–290 (2024). 58. R. L. Erickson, N. Hemati, S. E. Ross, O. A. MacDougald, p300 Coactivates the Adipogenic Transcription Factor CCAAT/Enhancer-binding Protein α *. J. Biol. Chem. 276, 16348– 16355 (2001). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 18 59. T. Å. Pedersen, E. Kowenz-Leutz, A. Leutz, C. Nerlov, Cooperation between C/EBPα TBP/TFIIB and SWI/SNF recruiting domains is required for adipocyte differentiation. Genes Dev. 15, 3208–3216 (2001). 60. C. Müller, C. F. Calkhoven, X. Sha, A. Leutz, The CCAAT Enhancer-binding Protein α (C/EBPα) Requires a SWI/SNF Complex for Proliferation Arrest *. J. Biol. Chem. 279, 5 7353–7358 (2004). 61. X. Chen, W. Zhou, R.-H. Song, S. Liu, S. Wang, Y . Chen, C. Gao, C. He, J. Xiao, L. Zhang, T. Wang, P. Liu, K. Duan, Z. Cheng, C. Zhang, J. Zhang, Y . Sun, F. Jackson, F. Lan, Y . Liu, Y . Xu, J. J.-L. Wong, P. Wang, H. Yang, Y . Xiong, T. Chen, Y . Li, D. Ye, Tumor suppressor CEBPA interacts with and inhibits DNMT3A activity. Sci. Adv. 8, eabl5220 (2022). 10 62. J. Yamazaki, J. Jelinek, Y . Lu, M. Cesaroni, J. Madzo, F. Neumann, R. He, R. Taby, A. Vasanthakumar, T. Macrae, K. R. Ostler, H. M. Kantarjian, S. Liang, M. R. Estecio, L. A. Godley, J.-P. J. Issa, TET2 Mutations Affect Non-CpG Island DNA Methylation at Enhancers and Transcription Factor–Binding Sites in Chronic Myelomonocytic Leukemia. Cancer Res. 75, 2833–2843 (2015). 15 63. M. Tulstrup, M. Soerensen, J. W. Hansen, L. Gillberg, M. Needhamsen, K. Kaastrup, K. Helin, K. Christensen, J. Weischenfeldt, K. Grønbæk, TET2 mutations are associated with hypermethylation at key regulatory enhancers in normal and malignant hematopoiesis. Nat. Commun. 12, 6061 (2021). 64. A. S. Nacht, R. Ferrari, R. Zaurin, V . Scabia, J. Carbonell‐Caballero, F. Le Dily, J. Quilez, A. 20 Leopoldi, C. Brisken, M. Beato, G. P. Vicent, C/EBPα mediates the growth inhibitory effect of progestins on breast cancer cells. EMBO J. 38, e101426 (2019). 65. F. Krueger, S. R. Andrews, Bismark: a flexible aligner and methylation caller for Bisulfite- Seq applications. Bioinformatics 27, 1571–1572 (2011). 66. H. Li, B. Handsaker, A. Wysoker, T. Fennell, J. Ruan, N. Homer, G. Marth, G. Abecasis, R. 25 Durbin, 1000 Genome Project Data Processing Subgroup, The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078–2079 (2009). 67. H. Wu, B. Caffo, H. A. Jaffee, R. A. Irizarry, A. P. Feinberg, Redefining CpG islands using hidden Markov models. Biostatistics 11, 499–514 (2010). 68. Y . Hao, T. Stuart, M. H. Kowalski, S. Choudhary, P. Hoffman, A. Hartman, A. Srivastava, G. 30 Molla, S. Madad, C. Fernandez-Granda, R. Satija, Dictionary learning for integrative, multimodal and scalable single-cell analysis. Nat. Biotechnol. 42, 293–304 (2024). 69. I. Korsunsky, N. Millard, J. Fan, K. Slowikowski, F. Zhang, K. Wei, Y . Baglaenko, M. Brenner, P. Loh, S. Raychaudhuri, Fast, sensitive and accurate integration of single-cell data with Harmony. Nat. Methods 16, 1289–1296 (2019). 35 70. M. R. Corces, J. D. Buenrostro, B. Wu, P. G. Greenside, S. M. Chan, J. L. Koenig, M. P. Snyder, J. K. Pritchard, A. Kundaje, W. J. Greenleaf, R. Majeti, H. Y . Chang, Lineage- (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 19 specific and single-cell chromatin accessibility charts human hematopoiesis and leukemia evolution. Nat. Genet. 48, 1193–1203 (2016). 71. M. I. Love, W. Huber, S. Anders, Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550 (2014). 5 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 20

We thank all patients providing AML specimens for this study. We thank all members of the Majeti, Ji, and Feinberg labs for their helpful input and discussion. Funding: This work is supported by National Institutes of Health grant 5R01CA054358. 5 Author contributions: Conceptualization: EH, WZ, EMM, HJ, MK, APF Methodology: EH, WZ, EMM, LPG, RM, HJ, MK, APF Investigation: EH, WZ, EMM, AI, RT, LPG, RM Visualization: EH, WZ 10 Funding acquisition: APF Project administration: MK, APF Supervision: HJ, MK, APF Writing – original draft: EH, MK, WZ, APF Writing – review & editing: EH, MK, WZ, APF 15 Competing interests: A.P.F. is an inventor on patents to Johns Hopkins University 10,752,953, “Method of Detecting Cancer Through Generalized Loss of Stability of Epigenetic Domains, and Compositions Thereof,” and 16/310,176, “Potential Energy Landscapes Reveal the Information- Theoretic Nature of the Epigenome,” which are covered by a pathway licensing agreement with 20 Bristol-Myers Squibb. R.M. is on the Advisory Boards of Kodikaz Therapeutic Solutions, Orbital Therapeutics, Pheast Therapeutics, 858 Therapeutics, Prelude Therapeutics, Mubadala Capital, and Aculeus Therapeutics. R.M. is a co-founder and equity holder of Pheast Therapeutics, MyeloGene, and Orbital Therapeutics. No disclosures were reported by the other authors. 25 Data and materials availability: The ERRBS data used in this paper are publicly available (GSE98350). Single-cell sequencing data will be available in GEO upon publication. Processed data are included as suppl ementary data files , and o ther materials required for re -analysis are available upon request. 30 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 21 Supplementary Materials Figs. S1 to S9 Data S1 to S14 5 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 22 chr2:73,496,700−73,496,738 chr2:73,496,700−73,496,738 FBXO41 Normal Sample MML=0.22 Sample NME=0.49 D−2242 (CEBPA−dm) Sample MML=0.91 Sample NME=0.35 Differential Methylation dMML=0.8, dNME=−0.14, UC=0.62 chr19:49,623,173−49,623,192 chr19:49,623,173−49,623,192 PPFIA3 Normal Sample MML=0.03 Sample NME=0.18 D−2242 (CEBPA−dm) Sample MML=0.34 Sample NME=0.92 Differential Entropy dMML=0.25, dNME=0.63, UC=0.57 A t(15;17) SRSF2 PTPN11 PHF6 Del(5/7q) SETBP1 inv(16) t(8;21) RUNX1 BCOR TP53 NOTCH1 GATA2 ETV6 KDM6A SUZ12 EZH2 ASXL1 STAG2 SMC3 SMC1A RAD21 WT1 NPM1 FLT3.TKD FLT3.ITD NRAS KIT TET2 DNMT3A IDH2 IDH1 CEBPA−dm D−2240 D−2242 D−2253 D−2273 D−2545 D−2748 D−2753 D−7142 D−2268 D−5289 D−6887 D−7059 D−7180 D−7301 D−7418 D−7420 D−2185 D−2195 D−3331 D−5286 D−6239 D−6374 D−7131 D−7188 D−7313 D−7411 D−3330 D−6375 D−6456 D−6945 D−7115 D−7273 D−2208 D−2235 D−2549 D−6243 D−7143 D−7402 D−2190 D−2204 D−6882 D−7070 D−7150 D−7308 Sample Mutation B 0.2 0.4 0.6 −2000 −1000 0 1000 2000 Distance from TSS (bp) MML 0.1 0.2 0.3 −2000 −1000 0 1000 2000 Distance from TSS (bp) NME AML Subtype CEBPA−dm IDH1/2 DNMT3A TET2 KIT NRAS Normal C −20 0 20 −20 0 20 MML PC1: 24.2% MML PC2: 10.9% −20 0 20 40 −25 0 25 NME PC1: 15.9% NME PC2: 6.5% AML subtype CEBPA−dm IDH1/2 DNMT3A TET2 KIT NRAS Normal D Analysis region Analysis region Analysis region Differential Methylation Differential Entropy Analysis region Low methylation level High methylation level High methylation entropy Low methylation entropy Unmethylated Methylated (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 23 Figure 1: Potential energy landscape analysis identifies discordant methylation stochasticity in AML. (A) Illustration of differential mean methylation (left) and differential methylation entropy (right). Both cases have differences in probability distributions of methylation, but one is driven by differences in mean methylation level (left) and the other by differences in methylation entropy (right). (B) Mutational profiles (rows) of selected ERRBS 5 samples (columns) used in CPEL analysis. Mutations in genes used for definition of AML subtypes are indicated by color, and are mutually exclusive (CEBPA-dm, yellow; IDH1/2, green; DNMT3A, pink; TET2, dark blue; NRAS, purple; KIT, light blue). (C) Meta-region plots of CPEL MML (top) and NME (bottom) over promoter regions, colored by AML subtype. (D) Principal component analysis (PCA) of CPEL MML (top) and NME (bottom), colored by AML 10 subtype. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 24 Figure 2: Methylation discordance in AML localizes to key regulators of leukemia. (A) Overlap of genes associated with UC DMRs between AML subtypes. Genes with associated analysis regions are listed in Data S3. (B) Genetic drivers of AML have disruption of methylation (measured by CPEL UC) in AMLs, independent of mutations in that genetic driver. 5 Colors indicate the AML subtype (rows) with methylation disruption over the promoter of the given gene (columns). AML genetic drivers covered by CPEL analysis regions are included. *: significant enrichment of AML genetic driver genes in UC DMR genes for the given subtype (hypergeometric p-value < 0.01). (C) Gene set enrichment analysis of genes ranked by promoter UC over the MSigDB C2 gene set library (Methods). Select gene sets relevant to histone marks 10 and leukemic signatures are annotated. Positive NES = high UC (highly discordant methylation stochasticity); negative NES = low UC (lowly discordant methylation stochasticity). CEBPA−dm IDH1/2 KIT TET2 DNMT3A NRAS 21532034 966 793 47644037123617917517515615513711311110784 58 57 40 36 34 32 29 26 24 23 21 21 21 19 19 15 14 12 8 7 6 6 5 5 4 4 4 3 3 3 3 2 2 1 1 1 1 0 1000 2000 Intersection size 7181 5425 3161 1627 837 2800 2000 4000 6000 8000 # UC DMR genes A NRAS KIT* TET2* DNMT3A IDH1/2* CEBPA−dm* CREBBP DNMT3A CEBPA FBXW7 GATA2 IDH2 RB1 RUNX1T1 TP53 U2AF2 WT1 CUX1 ETV6 IKZF1 JAK2 MYH11 NF1 PRPF40B PTEN CBFB CBLB CDKN2A DEK EP300 EZH2 FLT3 GNAS IDH1 KIT KRAS MYC NPM1 PTPN11 RAD21 SF1 SH2B3 SRSF2 U2AF1 ASXL1 BRAF CBL KDM5A NUP214 PML RUNX1 SF3A1 TET2 AML genetic driver AML subtype B TET2 DNMT3A KIT IDH1/2 CEBPA−dm NRAS GENTLES_LEUKEMIC STEM_CELL_UP JAATINEN HEMATOPOIETIC_STEM CELL_UP EPPERT_HSC_R IVANOVA HEMATOPOIESIS_STEM CELL EPPERT_CE_HSC_LSC ALCALAY_AML_BY_NPM1 LOCALIZATION_UP IVANOVA HEMATOPOIESIS_STEM CELL_LONG_TERM WANG_MLL_TARGETS MIKKELSEN_NPC_HCP WITH_H3K27ME3 MEISSNER_BRAIN_HCP WITH_H3K27ME3 MIKKELSEN_NPC_HCP WITH_H3K4ME3_AND H3K27ME3 Normalized Enrichment Score −3 −2 −1 0 1 2 3 C (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 25 Figure 3: Single-cell RNA-seq reveals transcriptional dysregulation over leukemia signatures in AML. (A) UMAPs of scRNA-seq revealing 10 distinct clusters. Samples are pooled by genotype. (B) V olcano plots of differentially expressed genes in Cluster 0 for CEBPA-5 dm AMLs (left), IDH1/2-mutant AMLs (middle-left), DNMT3A-mutant AMLs (middle-right), and TET2-mutant AMLs (right) versus normal progenitors. (C) Boxplots of cluster-level pseudobulk expression profiles for ADGRE5 (in cluster 0) and BAHCC1 (in cluster 4), differentially expressed in all AML versus normal comparisons. *: Bonferroni adjusted p-value ≤ 0.05. (D) Gene set enrichment analysis over the MSigDB C2 gene set library of genes ranked by 10 average log2 fold change (versus normal samples) in cluster 0. Selected gene sets related to histone marks, leukemia, and hematopoietic signatures are annotated. (E) Gene set enrichment analysis over the MSigDB Curated Hallmark gene set library of genes ranked by average log2 A LINC00910 ARHGEF18 ATXN7 DUSP10 MDK PIK3R6 POLR2J3 TAPT1−AS1 MTX3VAMP5 SLC14A1 DDIT4 LINC01572 SORBS3 PFKMFHL1 ARL4C CCDC152 RBM38 TUBA4A CD96CKAP4 PRDX1 N_up: 862 N_down: 219 0 10 20 30 40 −5.0 −2.5 0.0 2.5 5.0 Log2FC −Log10 P CEBPA Cluster 0 total = 11173 variables DMXL2 BCAS4 MAP9 PLCG2RGS18 CRIP1LRRC61GIMAP7 REC8 PRICKLE1 PFKM C1QTNF4 ACSL1 IL1RAP ADGRE5 CA2 RASGRP3PIWIL4 N_up: 69 N_down: 93 0 10 20 30 −3 0 3 6 Log2FC −Log10 P IDH Cluster 0 total = 11173 variables PAWR MIR222HG CD96 CCDC26ADGRE5 NPDC1 IL1RAP MBNL1−AS1 SP140DUSP10FZD3 CLMNALDH2 PLCG2GAS7PTMSMPO DLEU2 C1QTNF4 N_up: 525 N_down: 172 0 10 20 30 40 −4 0 4 Log2FC −Log10 P DNMT3A Cluster 0 total = 11173 variables TNFRSF4 RAB11FIP1BEND7 F11R MPPED2 ATF3 ADGRE5 KLF2 TNFAIP2UHRF1 CRIP1 RASGRP3 TYROBPMCM2 CXCL8 CD36 S100A10 N_up: 57 N_down: 54 0 5 10 15 20 −6 −3 0 3 Log2FC −Log10 P TET2 Cluster 0 total = 11173 variables NS Log2FC p−adj p−adj and Log2FC B * ** * 0 1 2 3 4 CEBPA−dm IDH1/2 DNMT3A TET2 Normal Normalized counts Cluster 0 ADGRE5C * ** * 0 1 2 3 CEBPA−dm IDH1/2 DNMT3A TET2 Normal Normalized counts Cluster 4 BAHCC1 Cluster 0 IDH DNMT3A CEBPA TET2 WONG_EMBRYONIC_STEM_CELL_CORE RHEIN_ALL_GLUCOCORTICOID THERAPY_DN CROONQUIST_NRAS_SIGNALING_DN GRAHAM_CML_DIVIDING_VS_NORMAL QUIESCENT_UP BHATTACHARYA_EMBRYONIC_STEM CELL YAGI_AML_WITH_T_9_11 TRANSLOCATION HOEBEKE_L YMPHOID_STEM_CELL_UP REACTOME_CHROMATIN_MODIFYING ENZYMES EPPERT_HSC_R WP_HEMATOPOIETIC_STEM_CELL GENE_REGULATION_BY_GABP ALPHABETA_COMPLEX TAKEDA_TARGETS_OF_NUP98_HOXA9 FUSION_6HR_UP SENESE_HDAC1_TARGETS_UP KEGG_CHRONIC_MYELOID_LEUKEMIA SENESE_HDAC3_TARGETS_UP CASORELLI_ACUTE_PROMYELOCYTIC LEUKEMIA_UP Normalized Enrichment Score −4 −2 0 2 4 D *** *** *** ****** *** ***** * * ** * * *** ** ** ** ** ** **** ** ** * ** * ** * * * ** **** *** *** *** * * ** ******* *** ** ** * * **** *** **** *** *** **** *** *** ****** *** ****** ** *** *** Cluster 0 IDH TET2 DNMT3A CEBPA HALLMARK_MYC_TARGETS_V1 HALLMARK_OXIDATIVE_PHOSPHORYLATION HALLMARK_E2F_TARGETS HALLMARK_MYC_TARGETS_V2 HALLMARK_DNA_REPAIR HALLMARK_BILE_ACID_METABOLISM HALLMARK_XENOBIOTIC_METABOLISM HALLMARK_GL YCOL YSIS HALLMARK_MTORC1_SIGNALING HALLMARK_UNFOLDED_PROTEIN_RESPONSE HALLMARK_MITOTIC_SPINDLE HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION HALLMARK_HYPOXIA HALLMARK_KRAS_SIGNALING_UP HALLMARK_INFLAMMATORY_RESPONSE HALLMARK_COMPLEMENT HALLMARK_APOPTOSIS HALLMARK_IL6_JAK_STAT3_SIGNALING HALLMARK_TGF_BETA_SIGNALING HALLMARK_ANDROGEN_RESPONSE HALLMARK_PROTEIN_SECRETION HALLMARK_PI3K_AKT_MTOR_SIGNALING HALLMARK_TNFA_SIGNALING_VIA_NFKB HALLMARK_UV_RESPONSE_DN Normalized Enrichment Score −4 −2 0 2 4 E 0 2 4 5 6 1 3 7 8 9 0 2 4 5 6 1 3 7 8 9 0 2 4 5 6 1 3 7 8 9 0 2 4 5 6 1 3 7 8 9 0 2 4 5 6 1 3 7 8 9 −10 0 10 −10 0 10 −10 0 10 −10 0 10 CEBPA−dm IDH1/2 DNMT3A TET2 Normal Progenitor −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 0 10 UMAP 1 UMAP 2 0 1 2 3 4 5 6 7 8 9 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 26 fold change (versus normal samples) in cluster 0. *: Benjamini-Hochberg adjusted p-value ≤ 0.05; **: adjusted p-value ≤ 1e-3; ***: adjusted p-value ≤ 1e-4. Positive NES = upregulated in AML relative to normal samples; negative NES = downregulated in AML relative to normal. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 27 Figure 4: DNA methylation stochasticity is associated with gene expression in AML. (A) Relationship between MML or NME and gene expression mean (left), variance (middle), or variability (MAV; right) in IDH1/2-mut AMLs as a function of distance from the TSS. Lower methylation level is associated with higher expression level (top-left), while lower expression 5 levels are associated with lower levels of NME (bottom-left). Mean-adjusted expression variability is associated with higher NME near the TSS, connecting methylation stochasticity to expression variability. (B) Information-theoretic analysis of scRNA-seq data reveals a complex gene regulatory network between 30 genetic drivers of AML. All included genetic drivers show disrupted methylation landscapes by the CPEL model. The gene with the highest betweenness 10 centrality (RAD21) is marked. Quartile of expression rank product 0−25% 25−50% 50−75% 75−100% 0.0 0.2 0.4 0.6 0.8 −2000 −1000 0 1000 2000 Distance from TSS (bp) MML 0.1 0.2 0.3 −2000 −1000 0 1000 2000 Distance from TSS (bp) NME Quartile of variance rank product 0−25% 25−50% 50−75% 75−100% 0.0 0.2 0.4 0.6 0.8 −2000 −1000 0 1000 2000 Distance from TSS (bp) MML 0.1 0.2 0.3 −2000 −1000 0 1000 2000 Distance from TSS (bp) NME Quartile of MAV rank product 0−25% 25−50% 50−75% 75−100% 0.0 0.2 0.4 0.6 0.8 −2000 −1000 0 1000 2000 Distance from TSS (bp) MML 0.1 0.2 0.3 −2000 −1000 0 1000 2000 Distance from TSS (bp) NME Expression mean Expression variance Mean-adjusted variability A GNAS U2AF1 U2AF2 DEK IDH2 CUX1 ETV6 PTEN SRSF2 NPM1 RAD21 JAK2 RB1 SF1 FBXW7 KRAS EP300 CBLB FL T3 NF1 PTPN11 CREBBP CBFB KITEZH2 IKZF1 DNMT3A CEBP A SH2B3 IDH1 Betweenness centrality 0 100 200 300 400 Edge weight 1.75 1.80 1.85 1.90 1.95 B (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 28

ERRBS data download and processing FASTQ files containing enhanced reduced representation bisulfite sequencing (ERRBS) data from 44 primary patient AMLs were downloaded from GSE98350 (6). Patients were 5 annotated by 6 genetically defined AML subtypes, requiring mutually-exclusive mutations in each gene: DNMT3A (n=10), IDH1/2 (n=8), TET2 (n=6), CEBPA double mutation (CEBPA- dm, n=8), KIT (n=6), or NRAS (n=6) (Data S2). ERRBS data from 6 healthy CD34+ NBM samples was also downloaded for normal reference (GSE52945, GSE108247) (27, 28). Reads were trimmed and quality control was performed using TrimGalore v.0.6.6 10 (https://github.com/FelixKrueger/TrimGalore, RRID:SCR_011847) with default parameters. Alignment to the hg19 reference genome was performed using Bismark v.0.23.0 (RRID:SCR_005604) (65) with default parameters. The resulting BAM files were processed with samtools v.1.18 (RRID:SCR_002105) (66) for sorting, deduplication, and indexing. Finally, methylation calls were extracted and CpG_report files were generated using Bismark’s 15 methylation_extractor. Definition of analysis regions and potential energy landscape estimation Methylation analysis was performed using the correlated potential energy landscape (CPEL) model for targeted differential methylation using the Julia package CpelTdm.jl (26). 20 CpelTdm requires user-defined analysis regions as an input. These were generated by identifying autosomal regions containing at least 3 consecutive CpGs all passing a given coverage filter; coverage filters implemented were: requiring at least 10x coverage in every sample (all AML samples and normal); requiring at least 10x coverage in at least 4 samples per subtype (each AML subtype and normal); requiring at least 5x coverage in every sample (all AML samples and 25 normal); requiring at least 5x coverage in at least 4 samples per subtype (each AML subtype and normal). Regions over 200 bp were split into max(1, floor(region width/100)) approximately equally-sized subregions using the subdivideGRanges() function from the Bioconductor package exomeCopy (RRID:SCR_001276). The CPEL model was applied to ERRBS data for each AML subtype versus normal comparison using each coverage filter, run over all resulting regions or 30 regions subset to regions within promoters. The coverage filter requiring at least 10x coverage in every sample with all analysis regions as input was selected for downstream analysis, as this set of CPEL comparisons detected the most significant DMRs in TET2-mut AMLs, which had relatively few DMRs under most coverage filters (0-2752 UC DMRs at q ≤ 0.2) (Data S1). The CPEL model estimates the PDM within an analysis region R containing N CpG sites 35 as: π(𝑥) = 𝑃(𝑋 = 𝑥) (1) for a random vector 𝑿 = [𝑋! 𝑋# . . 𝑋$ ]% , where X = 0 represents an unmethylated CpG and X = 1 represents a methylated CpG. 40 After parameter estimation by maximum likelihood method, the mean methylation level (MML) can be computed as: µ(𝑋) = 𝐸 01 𝑁 2 𝑋& $ &'! 3 (2) (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 29 And the normalized methylation entropy (NME) can be computed as: ℎ(𝑥) = − 1 𝑁 2 π(𝑥) 𝑙𝑜𝑔#π(𝑥) ( (3) Both MML and NME are normalized to lie within [0, 1], with the highest methylation level taking MML=1 (fully methylated) and the lowest methylation taking MML = 0 (fully unmethylated). Similarly, an NME = 0 indicates there was only one methylation pattern found in 5 the analysis region (fully ordered), and NME = 1 indicates all methylation patterns are equally likely (fully disordered). Differential methylation, methylation entropy, and probability distribution of methylation For a group with M1 and another group with M2 values of CPEL statistics over an 10 analysis region R, the CPEL model calculates differential mean methylation level (dMML) with the test statistic: 𝑇))* = 1 𝑀! 2 µ=𝑋! (, )> ) ! , '! − 1 𝑀# 2 µ=𝑋# (, )> ) " , '! (4) Differential normalized methylation entropy (dNME) is calculated using the test statistic: 15 𝑇$). = 1 𝑀! 2 ℎ=𝑋! (, )> ) ! , '! − 1 𝑀# 2 ℎ=𝑋# (, )> ) " , '! (5) Differences between PDMs are quantified using the uncertainty coefficient (UC) with the test statistic: 𝑇/0) = 1 𝑀!𝑀# 2 2 𝑄 =𝑋! (, !), 𝑋# (, ") > ) " , "'! ) ! , !'! (6) 20 where 𝑋# (, !) and 𝑋# (, ") are the m-th methylation states in the first and second group over the analysis region, respectively. Here, the uncertainty coefficient is the geometric Jensen-Shannon divergence normalized by the cross-entropy between groups, given by: 𝑄(𝑿𝟏, 𝑿𝟐) = 𝐷(𝜋! ∥ 𝜋G) + 𝐷(𝜋# ∥ 𝜋G) 𝐶(𝑿𝟏, 𝑿J) + 𝐶( 𝑿𝟐, 𝑿J ) (7) 25 where 𝑿J is a random vector associated with the potential energy function 𝜋G, the average of the potential energy functions associated with 𝑿! and 𝑿#; the Kullback-Leibler divergence between two distributions 𝑓 and 𝑔 is 𝐷(𝑓 ∥ 𝑔) = 2 𝑓(𝒙) log𝟐 𝑓(𝒙) 𝑔(𝒙) 𝒙 (8) 30 and the cross-entropy between two random vectors 𝑿 and 𝒀 associated with 𝑓 and 𝑔 is 𝐶(𝑿, 𝒀) = − 2 𝑓(𝒙) log# 𝑔(𝒙) 𝒙 (9) Significance of differential statistics for a given analysis region is determined through permutation-based exact p-value computation or Monte Carlo-based permutation testing (26), (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 30 and p-values are corrected using the Benjamini-Hochberg (BH) procedure for controlling false discovery rate. For downstream analysis, differentially methylated regions (DMRs) were defined as regions with absolute dMML, dNME, or UC ≥ 0.1, and BH-adjusted p-value ≤ 0.2. Genomic annotations 5 Files and tracks have coordinates for hg19. Annotations for CpG islands (CGIs) were obtained from Wu et al. (67); CGI shores were defined as 2 kb regions flanking CGIs, CGI shelves were defined as 2 kb regions flanking shores, and open seas were defined as all other genomic regions. Annotations for promoters and gene bodies were obtained from the Bioconductor package TxDb.Hsapiens.UCSC.hg19.knownGene; promoters were defined as ± 2 10 kb flanking the transcriptional start site. Annotations for K562 chromatin states from ChromHMM (RRID:SCR_018141) and genic features were generated using the R package annotatr’s built-in hg19_K562-chromatin and hg19_basicgenes (DOI: 10.18129/B9.bioc.annotatr). 15 Enrichment analysis For enrichment analysis of DMRs over various genomic features, the odds ratio statistic and Fisher’s two-sided exact test were implemented on a 2 x 2 contingency table containing the number of regions overlapping and not overlapping the feature (such as DMR and overlapping promoter, not DMR and overlapping promoter, etc). 20 Motif analysis was performed using HOMER v.5.1 (RRID:SCR_010881) with size set to default and CpG normalization method. Input regions were CPEL dMML, dNME, or UC DMRs, and custom background region sets were defined as all CPEL analysis regions. Over-representation analysis (ORA) and gene set enrichment analysis (GSEA) were implemented using clusterProfiler v.4.12.0 (RRID:SCR_016884). To assess enrichment of genes 25 associated with DMRs or differentially expressed genes over transcription factor target sets, ORA was performed using the gene set library “ENCODE_and_ChEA_Consensus_TFs_from_ChIP-X” downloaded from Enrichr (https://maayanlab.cloud/Enrichr/, RRID:SCR_001575). For DMRs, a

gene set was defined as any gene containing a CPEL analysis region within its promoter; for DEGs, a background gene set was defined as all genes detected in the RNA-seq data. 30 GSEA was performed by ranking genes by the strength of the UC, dNME, or dMML of analysis regions within the gene promoter; or by log2FC reported by Seurat DEseq2 differential expression analysis (see Differential expression analysis ). When ranking by CPEL statistics (UC, dNME, dMML), most genes contained multiple analysis regions in the promoter ; in this case, the region with the strongest value was used for ranking. GSEA was performed over the MSigDB Hallmark 35 (H) and Curated (C2) gene set collections (Broad Institute, UCSD). Single-cell sequencing Primary patient AML samples for single-cell sequencing experiments were obtained from Johns Hopkins Hospital and Stanford University biobank. Two healthy CD34+ samples were 40 purchased from AllCells (Lot 3011884, 3011592) and one GMP sample was obtained from Stanford University for normal comparison. Samples were collected from peripheral blood or bone marrow; blast percentages ranged from 70-99%. Patients were annotated by 4 genetically defined AML subtypes with mutually exclusive mutations, corresponding to those used for methylation analysis: DNMT3A (n=5), IDH1/2 (n=4), TET2 (n=7), CEBPA-dm (n=5) (Data 45 S2). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 31 Cells were thawed in warmed RPMI + 10% FBS, then washed twice with 1X PBS + 0.04% BSA and strained using a 40 μm Flowmi Cell Strainer (SP Bel-Art). Cell concentration and viability were assessed via Trypan blue staining using a Luna-II cell counter (Logos Biosystems). All samples had greater than 85% viability. Samples were sequenced with either single-cell RNA-seq and single-cell ATAC-seq (separately); with single-cell multiome ATAC + 5 Gene Expression sequencing; or with all of: single-cell RNA-seq, single-cell ATAC-seq, and single-cell multiome ATAC + Gene Expression sequencing as biological replicates (Data S2B). We report here the single-cell RNA-seq and single-cell multiome Gene Expression portions only. For scRNA-seq, a single-cell suspension was used for construction of scRNA libraries, which were generated with Chromium NextGEM Single-cell 3’ Gene Expression reagents (v.3.1; 10 10x Genomics). For scATAC-seq, a single nuclei suspension was generated following the protocol Nuclei Isolation for Single Cell ATAC Sequencing (Demonstrated protocol CG000169, Rev D; 10x Genomics). Single-cell ATAC libraries were constructed with Chromium Next GEM Single Cell ATAC reagents (v.1.1; 10x Genomics). A 10x Chromium Controller instrument was used according to the manufacturer’s instructions for both scRNA-seq and scATAC-seq, 15 targeting a recovery of around 4,000 nuclei for each sample. For single-cell multiome sequencing, a single nuclei suspension was generated according to protocol Nuclei Isolation for Single Cell Multiome ATAC + Gene Expression Sequencing (Demonstrated protocol CG000365, Rev D; 10x Genomics). Single-cell multiome ATAC and gene expression libraries were generated using the Chromium Next GEM Single-cell Multiome 20 ATAC + Gene Expression reagents (10x Genomics) and a 10x Chromium Controller instrument according to the manufacturer’s instructions (User Guide G000338, Rev D), targeting a recovery of around 4,000 nuclei for each sample. Quality control on intermediate products and final libraries was performed using the Agilent Bioanalyzer High Sensitivity Kit (Agilent Technologies). The final libraries were 25 quantified via qPCR using KAPA Library Quantification Kits (KAPA Biosystems, #KK4824) and were sequenced on an Illumina NovaSeq6000 system or on an Illumina HiSeq 4000 system (Data S2) with a 1% spike-in of PhiX control library (Illumina). Single-cell alignment, normalization, and integration 30 Single-cell sequencing data were aligned to hg38. The single-omic scRNA-seq data were aligned using Cell Ranger (RRID:SCR_017344; 10x Genomics), and the multiomic data were aligned using Cell Ranger Arc (RRID:SCR_023897; 10x Genomics). The single-cell data were analyzed using the Seurat R package v.5.1.0 (RRID:SCR_016341) (68). For the single-omic scRNA-seq data, cells with number of non-zero 35 genes larger than 1,000, number of reads less than 60,000, and percentage of mitochondria reads less than 25% were retained. For the multi-omic data, only cells meet the following criteria were retained for further analysis: first, number of reads larger than 5,000 but smaller than 100,000 from the scATAC-seq data modality; second, number of reads larger than 1,000 but smaller than 100,000 from the scRNA-seq data modality; third, the percentage of mitochondria reads less than 40 30%. The gene expression (i.e., scRNA-seq) data were then normalized to the same library and log-transformed using the Seurat NormalizeData() function. Samples from different patients or technologies (single or multi-omics) were integrated using Harmony (RRID:SCR_022206) (69) based on the gene expression data. For one sample, we obtained replicate scRNA-seq data from Stanford University (Majeti Lab; sample SU758, 45 replicate SU758_stanford_single). For data integration, the top 2,000 most variable genes were (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 32 first obtained from each sample based on the mean-adjusted variance. Then the top 2,000 genes that were consistently ranked as the most highly variable genes were used for principal component analysis (PCA). The top 30 PCs were used for integration using Harmony. Low dimensional representation of the cells and cell clusters were obtained using UMAP and Louvain clustering based on the batch-adjusted PCs from Harmony. 5 To annotate the cell identity of the cell clusters in the integrated single-cell data, we obtained a bulk RNA-seq dataset that contains gene expression from 13 normal hematopoietic cell types from GEO (GSE74246) (70). Differential gene expression analyses were applied to each pair of cell types using DESeq2 (RRID:SCR_015687) (71) and the union of the top 100 differential genes from each comparison was obtained. Based on this differential gene set, we 10 calculated the Spearman’s correlation between the average gene expression from each cell cluster from the single-cell data and the quantile-normalized gene expression from the bulk data. The similarity between the cell cluster and the normal cell types were used to annotate their identities. 15 Differential expression analysis Differential expression analysis was performed using Seurat v.5.1.0 (68) between each AML subtype and normal samples (n=3, CD34-1, CD34-2, and GMP-1) within clusters 0, 2, and 4. Pseudobulk expression profiles were generated for each sample in each cluster using the Seurat function AggregateExpression(). Differential expression testing was performed using the 20 FindMarkers() function with test.use = “DESeq2”. Genes were called as differentially expressed if absolute average log2 fold change ≥ 0.5 and Bonferroni-adjusted p-value ≤ 0.1. Calculation of gene expression mean adjusted variability Let yij be the library-size normalized, imputed, and log2-transformed expression level for 25 gene i (i = 1, …, I) and cell j (j = 1, …, J). Let mi and si be the mean and standard deviation of the expression level for gene i across all cells respectively. A LOESS regression model was fitted across all genes where si is the response variable and mi is the independent variable. Let 𝑠4U be the fitted values of the standard deviations from the regression model. The gene expression mean- adjusted variability (MAV) of gene i, hi, is defined as the residual of the regression model, or 30 equivalently the difference between the observed and fitted standard deviation: ℎ5 = 𝑠5 − 𝑠4U⁠. MAV was calculated for each gene in each sample as described above. Genes were then ranked by MAV within each sample. A consensus ranking of genes for each AML subtype was generated by calculating the rank product of individual MAV ranks from each sample of the appropriate subtype. Similarly, a consensus ranking of expression level and variance for each 35 gene in each AML subtype was generated by taking the rank product of the ranks of the expression mean or variance within each sample for the appropriate subtype. We then divided genes into quartiles of gene expression rank product, variance rank product, and MAV rank product in each AML subtype, with a higher quartile indicating a high value of the associated statistic (ex. genes in the 75-100% quartile of MAV rank product have higher MAV than genes 40 in the 0-25% quartile). Finally, expression rank product, variance rank product, and MAV rank product were associated with values of CPEL MML and CPEL NME binned over promoter regions for visualization. Gene network inference 45 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint DNA methylation stochasticity is linked to transcriptional variability and identifies convergent epigenetic disruption across genetically-defined subtypes of AML 33 Gene network inference was performed from the single-cell gene expression data using a previously-described information-theoretic method (52) implemented via the Julia package NetworkInference.jl with default settings. Networks were visualized with the R package tidygraph v.1.3.1 (DOI: 10.32614/CRAN.package.tidygraph). NetworkInference.jl returns a list of every possible edge, and the confidence of each edge existing in the true network, so the top 5 10% of edges in the network were retained for plotting. For each gene in the network, the betweenness centrality was calculated using the function betweenness() from the R package igraph v.2.0.3 (RRID:SCR_021238). Betweenness centrality is computed as 2 𝑔567 𝑔57 ⁄ 586,587,786 (10) 10 where 𝑔57 is the number of shortest paths between vertices 𝑖 and 𝑗, while 𝑔567 is the number of shortest paths between 𝑖 and 𝑗 that also pass through vertex 𝑣. Betweenness centrality measures how the gene (node) lies on paths between other nodes in the network, so a gene with the highest betweenness centrality is predicted to exert the most influence over the network. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted October 29, 2024. ; https://doi.org/10.1101/2024.10.26.620422doi: bioRxiv preprint