Monocytic Subclones Confer Resistance to Venetoclax-Based Therapy in Patients with Acute Myeloid Leukemia

With a median age at diagnosis of 68, acute myeloid leukemia (AML) is predominantly a disease of the elderly (1). For decades, intensive induction chemotherapy has been the standard of care for patients with AML. However, most elderly patients are poor candidates for this type of therapy, given their higher treatment-related mortality and lower response rates (2, 3). The standard of care for these patients has therefore been to offer low-dose therapies such as hypomethylating agents (HMA). However, response to HMA therapy is limited, and long-term overall survival (OS) is negligible (4, 5). Recent clinical trials have reported that the addition of the highly specific BCL2 inhibitor venetoclax to the HMA backbone can greatly increase the response rates and potentially the OS for older, newly diagnosed patients with AML who are unfit for conventional chemotherapy (6, 7). These findings led to the recent FDA approval of this regimen for this population, and it is now considered to be the standard of care.

The combination of venetoclax and the HMA azacitidine results in a remission rate of approximately 70% (6, 7). However, a significant minority of patients do not achieve a remission and are refractory. In addition, the majority of patients who do achieve a remission ultimately relapse (6, 7). It is therefore critical to understand mechanisms of venetoclax resistance. Historically, adverse disease features that predict response to conventional therapies have been defined in the setting of chemotherapy-based regimens, and particular chromosomal and genomic abnormalities strongly predict poor outcomes in this context (8). Notably, though, with the exception of in vitro work which has demonstrated a link between p53 and venetoclax resistance (9), analyses of clinical response data have for the most part not identified traditional adverse risk features as poor-outcome predictors in patients with AML treated with venetoclax-based therapies (6, 7, 10).

One explanation for this observation is that venetoclax with azacitidine has a novel mechanism of action in AML (10, 11), necessitating a reassessment of biological features associated with prognostication. To this end, we considered the relative impact of the developmental stages of AML on BCL2-mediated metabolism. Previously, we identified BCL2 as differentially expressed in subpopulations of AML cells that are enriched for malignant stem/progenitor cells compared with more differentiated tumor cells (12). In addition, recent in vitro studies show reduced venetoclax sensitivity in primary AML cells with a monocytic phenotype (13). Thus, we hypothesized that clinical features of AML that are indicative of myeloid differentiation status may correlate with reduced BCL2 dependence in patients with AML. As reported herein, we discovered that more differentiated monocytic AML is much more likely to be refractory to venetoclax-based therapy. Further, we demonstrate that altered regulation of energy metabolism and survival is an inherent property of AML development and appears to underlie this resistance. Consequently, the selective pressure of venetoclax-based therapy mediates profound changes in the biology of leukemic cell populations that manifest in altered developmental and metabolic properties.

To test whether differentiation status may predict responsiveness to venetoclax + azacitidine (VEN+AZA) in the clinic, we retrospectively reviewed 100 consecutive, newly diagnosed, previously untreated patients with AML who received VEN+AZA at the University of Colorado between January 2015 and October 2019 (all baseline characteristics are listed in Supplementary Table S1). We analyzed several baseline factors to determine the ability of each to predict disease that was refractory to treatment as defined by the European Leukemia Network [ELN; lack of complete remission (CR), CR with incomplete recovery of peripheral blood counts (CRi), partial remission (PR), or morphologic leukemia free state (MLFS); ref. 8]. The median age of the cohort was 72 years; 20 patients (20%) had a documented antecedent hematologic disorder; 64 patients (64%) had adverse risk disease by ELN criteria (8).

To specifically examine features associated with myeloid differentiation, we initially employed the FAB (French, American, British) classification system. Although this system is no longer employed for clinical purposes, it provides a well-described and clinically associated means to segregate patients with AML by virtue of myeloid differentiation status. In our VEN+AZA-treated patient cohort, 13 patients (13%) were identified as the FAB-M5 subtype, which is defined as a more differentiated phenotype of monocytic AML, and 77 (77%) were FAB-M0 or M1, indicative of a less differentiated phenotype (Supplementary Fig. S1A). Univariate analysis revealed sex (P = 0.0495), presence of an RAS pathway mutation (P = 0.0039), and FAB-M5 maturation state (P < 0.0001) to be associated with disease that was refractory to VEN+AZA (Table 1). A multivariate analysis revealed only the FAB-M5 maturation state (P = 0.0066) to be predictive of refractory response (Table 1). Specifically, 62% of FAB-M5 patients were refractory to VEN+AZA, whereas 0% of FAB-M4 and only 8% of non–FAB-M5 patients were refractory (Supplementary Fig. S1B). In addition, the median OS in FAB-M5 patients was 89 days, compared with 518 days for non–FAB-M5 patients (P = 0.0039; Supplementary Fig. S1C). These findings indicate a strong correlation between myeloid differentiation status and resistance to venetoclax-based therapy.

To understand if the lack of response by monocytic AML to VEN+AZA is driven by intrinsic mechanisms, we sought to directly evaluate VEN+AZA sensitivity in vitro, where protection from extrinsic factors such as the microenvironment is minimal. Because the FAB system is no longer used in the AML field, we employed phenotypic markers that would serve as a surrogate for the FAB-M5 subtype. Previous studies have shown that FAB-M5 patients lose expression of the primitive marker CD117 and upregulate expression of the monocytic markers CD11b, CD68, and CD64 (14–18). This expression pattern was verified in our own analysis of CD117, CD11b, CD68, and CD64 gene-expression levels in all FAB subclasses in The Cancer Genome Atlas AML (TCGA-AML) dataset (ref. 19; Supplementary Fig. S1D). Therefore, we designed a multicolor flow cytometry panel including CD117, CD11b, CD68, and CD64 to distinguish patients with monocytic AML (FAB-M5) from patients with primitive AML (FAB-M0/M1/M2). As shown in Fig. 1, this approach readily distinguishes two predominant cell populations within patients with AML. For example, patient 51 (Pt-51; a typical FAB-M0/M1/M2) presented with a single dominant disease population that was phenotypically primitive as evidenced by CD45-medium/SSC-low/CD117⁺/CD11b⁻/CD68⁻ (Fig. 1A). This patient achieved CR with VEN+AZA treatment. In contrast, Pt-72 (a typical FAB-M5) was refractory to VEN+AZA and presented with dominant monocytic disease that was CD45-bright/SSC-high/CD117⁻/CD11b⁺/CD68⁺ (Fig. 1B). Analysis of an additional 12 primary AML specimens (Supplementary Table S2) confirmed the phenotypic profile for primitive versus monocytic specimens (Fig. 1C; Supplementary Fig. S1E and S1F). Hereafter, these AMLs are noted as “prim-AML” or “mono-AML,” respectively.

Multiple studies have suggested that leukemic stem cells (LSC) are an important target of AML therapies (20). Our previous studies have shown that a phenotype of low reactive oxygen species (ROS-low) enriches for functionally defined LSCs (12, 21). Therefore, to more directly assess drug responsiveness in the LSC subpopulation, we isolated ROS-low cells from prim-AML and mono-AML specimens. Because mono-AML has never been directly characterized by ROS level, we confirmed using colony-forming unit (CFU) assays that the ROS-low phenotype enriches for stem/progenitor potential in mono-AML (Supplementary Fig. S1G). These data indicate the ROS-low phenotype strongly enriches for stem/progenitor potential in mono-AML, similar to what was reported for prim-AML (12, 21). We then treated the ROS-low subpopulations from prim-AML or mono-AML with VEN+AZA in vitro. Our results show that ROS-low LSCs of the mono-AML specimens are significantly more resistant than those of the prim-AML specimens (Fig. 1D), suggesting the refractory responses seen in FAB-M5 patients can be at least partially attributed to intrinsic molecular mechanisms uniquely present in monocytic AML cells.

To identify intrinsic molecular mechanisms of mono-AML that may be responsible for resistance to VEN+AZA, we sorted ROS-low LSCs from mono-AML and prim-AML specimens and performed bulk RNA-sequencing (RNA-seq) analysis. After removal of low-expressing genes and normalization (Supplementary Fig. S2A and S2B), principal component analysis (PCA) showed that the mono-AMLs clustered separately from the prim-AMLs (Fig. 2A), highlighting distinct biological features. Indeed, the top 50 upregulated and downregulated genes in mono-AMLs include the monocytic markers MAFB, LYZ, and CD14 and the primitive marker CD34, respectively (refs. 18, 22–25; Fig. 2B). More broadly, using gene set enrichment analysis (GSEA), we observed that prim-AMLs are enriched for multiple LSC gene sets, whereas the mono-AMLs are enriched for monocytic differentiation and AML lysosome gene sets (18, 25–27), confirming their distinct transcriptome profiles (Fig. 2C; Supplementary Fig. S2C and S2D).

Interestingly, our GSEA also identified oxidative phosphorylation (OXPHOS) as the top upregulated gene set in mono-AML in the Kyoto Encyclopedia of Gene and Genomes (KEGG; Fig. 2D and E), and Seahorse functional assays confirmed that the basal respiration rate of OXPHOS is indeed significantly higher in ROS-low LSCs of mono-AML (Fig. 2F). Together, these data suggest mono-AMLs are transcriptionally distinct from prim-AMLs and exhibit elevated energy metabolism in the form of OXPHOS activity.

We next focused on expression of apoptosis family genes given that VEN is a BCL2-specific inhibitor, and several studies have shown that BCL2 expression strongly correlates with VEN sensitivity in vitro (28, 29). Among genes related to apoptosis regulation (Supplementary Fig. S2E), our analysis revealed significant and consistent loss of BCL2 in mono-AMLs (N = 5), compared with the prim-AMLs (N = 7; Fig. 2G). Analysis of the TCGA AML dataset also showed progressive loss of BCL2 gene expression through stages of AML morphologic maturation (FAB-M0 to FAB-M5; Supplementary Fig. S2F and S2G). As a result, significantly lower expression of BCL2 is observed in FAB-M5 relative to FAB-M0/M1/M2 in the TCGA AML dataset (Fig. 2H). Further, reduced expression of BCL2 in mono-AML was confirmed at the protein level (Fig. 2I).

Interestingly, loss of BCL2 also occurs during normal monocytic development (30, 31). We found consistent loss of BCL2 at the monocytic stage in both human and murine systems (Supplementary Fig. S2H–S2K). Together, these analyses indicate BCL2 loss is a conserved biological feature during both normal and malignant monocytic development. Further, the data suggest BCL2 loss in monocytic AML may drive resistance to venetoclax-based therapies.

We have previously reported that BCL2 mediates OXPHOS in ROS-low AML cells (12). Given that monocytic AML has loss of BCL2 expression yet elevated OXPHOS activity (Fig. 2D–I), we hypothesized that other members of the BH3 family may become more active as a means to compensate for the role of BCL2. Notably, our bulk RNA-seq data show that expression of MCL1 is comparable in seven prim-AML versus five mono-AML (Fig. 3A). We then analyzed the TCGA AML dataset and found that MCL1 expression is significantly higher in FAB-M5 patients than FAB-M0/M1/M2 patients (Fig. 3B; Supplementary Fig. S3A). These data suggest MCL1 may be relevant to the intrinsic biology of monocytic AML.

To evaluate our hypothesis, we first compared effects of the selective MCL1 inhibitor VU661013 (32) with VEN. Three primary mono-AML specimens were treated side by side with VU661013+AZA or VEN+AZA. We observed that VU661013 alone or in combination with AZA induced strong cell death in all monocytic specimens evaluated, and the effect was significantly greater than VEN alone or VEN+AZA (Fig. 3C; Supplementary Fig. S3B). Notably, VEN+AZA failed to inhibit OXPHOS in mono-AML, whereas VU661013+AZA did show significant inhibition of OXPHOS (Fig. 3D; Supplementary Fig. S3C and S3D). As expected, differential modulation of OXPHOS is reflected in production of ATP (Fig. 3E), demonstrating increased MCL1 dependency as a distinct feature of energy metabolism in monocytic AML.

To more directly interrogate the functional role of MCL1, we also employed siRNA-mediated knockdown of MCL1 in three primary mono-AML specimens. As shown, two independent siMCL1 sequences (siMCL-#B and siMCL-#C) successfully reduced MCL1 protein (Fig. 3F; Supplementary Fig. S3E) and impaired OXPHOS activity (Fig. 3G; Supplementary Fig. S3F and S3G). In addition, production of the tricarboxylic acid cycle intermediates citrate, α-ketoglutarate, and malate was significantly lowered in a subset of monocytic AML specimens (Supplementary Fig. S3H). Notably, in using siRNAs to reduce MCL1 expression, significant cell death was observed with or without the addition of AZA in all three mono-AML specimens tested (Fig. 3H; Supplementary Fig. S3I). These data suggest that MCL1 is a key mediator of OXPHOS and viability for the mono-AML cell type.

Given that OXPHOS is critical for LSCs, as well as other types of cancer stem cells (12, 33–35), we next examined whether preferential dependence on MCL1 is also reflected in the stem/progenitor function of monocytic AML. In two independent mono-AML specimens, we observed that treatment with VU661013+AZA was significantly better than VEN+AZA in reducing CFU of mono-AML (Fig. 3I). We repeated the same treatments on three normal specimens isolated from umbilical cord blood (CBMC-1, -2, and -3) and observed that VU661013+AZA had only a modest impact on CFU potential, similar to VEN+AZA (Supplementary Fig. S3J), implying preferential reliance on MCL1 in malignant cells.

To more directly assess the increased reliance of monocytic LSCs on MCL1, we performed xenograft studies using prim-AMLs and mono-AMLs treated with VEN+AZA or VU661013+AZA followed by transplantation into immune-deficient NSG-S mice (21). Our results showed VEN+AZA effectively impaired the LSC engraftment potential in primitive AML as expected (10), but less so in monocytic AML (Fig. 3J, left). Conversely, the VU661013+AZA treatment was significantly better in reducing the LSC engraftment potential of monocytic AMLs (Fig. 3J, right). Together, these data demonstrate that monocytic AML displays greater reliance on MCL1 than BCL2 for energy production, stem/progenitor potential, and survival.

Based on the above findings, we next investigated the extent to which monocytic disease is evident in patients who initially responded but then relapsed on VEN+AZA therapy. In analyzing patients with AML prior to VEN+AZA treatment, we noted that the majority of patients actually present with tumors showing a mixture of the monocytic and primitive phenotype, which we term “MMP-AML” (for mixed monocytic/primitive-AML). We analyzed the characteristics of two patients with MMP-AML (Pt-12, Pt-65) during the course of treatment (Fig. 4A and B). Upon relapse after an initial CR, both patients showed almost complete loss of the primitive subpopulation and emergence of a dominant monocytic phenotype (CD45-bright/SSC-high/CD117⁻/CD11b⁺/CD68⁺). Thus, VEN+AZA treatment appeared to induce striking in vivo selection for the monocytic subpopulation in each patient (Fig. 4A and B).

Of note, this monocytic selection phenotype seems to be a unique clinical characteristic of VEN+AZA therapy. Indeed, previous analyses of patients treated with conventional chemotherapy have shown consistent enrichment of more primitive LSC phenotypes (36). To further corroborate this finding, we analyzed RNA-seq data of 11 pairs of diagnostic and relapsed specimens after conventional chemotherapy from a separate study by Shlush and colleagues (37). In this setting, we observed a gain of the LSC gene-expression signature, and loss of monocytic markers (CD11b and CD68) and a monocytic gene-expression signature at relapse, suggesting suppression of the myeloid phenotype following chemotherapy (Supplementary Fig. S4A and S4B). Lastly, we compared paired diagnosis versus relapse specimens from 6 patients with AML treated with conventional chemotherapy at our institution. In no case was a monocytic phenotype evident at relapse. In fact, for 2 patients with monocytic characteristics at diagnosis, conversion to a more primitive phenotype at relapse was observed (Supplementary Fig. S4C). Together, these data suggest that relapse following conventional chemotherapy strongly favors a primitive phenotype, and that selection of a monocytic phenotype at relapse appears to be a distinct characteristic of VEN+AZA therapy.

To investigate the origin of relapsed monocytic subpopulations in Pt-12 and Pt-65, we performed whole-exome sequencing (WES) analysis on sorted primitive and monocytic subpopulations from paired diagnosis and relapse specimens (sorting strategy outlined in Supplementary Fig. S4D). For each patient, three subpopulations of cells were isolated: diagnosis primitive (“Dx-prim”), diagnosis monocytic (“Dx-mono”), and relapse monocytic (“Rl-mono”). In Pt-12, 400× WES detected five unique nonsynonymous cancer-related mutations in 49 commonly mutated genes in AML (Supplementary Table S3). Although the variant allele frequencies (VAF) for three of five mutations are similar among all three sorted subpopulations, two mutations, SMC1A^R807H and NRAS^Q61K, presented unique patterns (Supplementary Fig. S4E). The SMC1A^R807H mutation was detected only in Dx-prim, suggesting that the Dx-prim subpopulation represents a genetically distinct subclone that was effectively eradicated by VEN+AZA treatment (Fig. 4C, clone 1 in teal). In contrast, the NRAS^Q61K mutation was detected only in Dx-mono and Rl-mono, suggesting the monocytic subpopulation that emerges at relapse represents a genetically distinct subclone that preexisted at diagnosis (Fig. 4C, clone 2 in pink).

In Pt-65, our analysis detected three unique nonsynonymous cancer-related mutations (Supplementary Fig. S4F). Interestingly, the mutation profile between the Dx-prim and Dx-mono subpopulations was identical, except for a low-level (∼1%) RAS^G13D mutation in the Dx-mono subpopulation. These data suggest the majority of cells in the Dx-mono subpopulation arose from the same genetic subclone from which Dx-prim was derived, despite their phenotypic differences. In particular, an EZH2^D185H mutation presented at an identical VAF of 46% in Dx-prim and Dx-mono but 0% in Rl-mono, suggesting the Dx-prim and Dx-mono cells containing the EZH2^D185H mutation arose from the same genetic subclone and were eradicated by VEN+AZA treatment (Fig. 4D, clone 1 in teal and brown, the gradual transition of color from teal to brown illustrates phenotypic progression from primitive to monocytic). Intriguingly, two KRAS mutations, G13D and G12V, emerged in Rl-mono at a VAF of 21% and 31%, respectively. The KRAS^G13D mutation was confirmed at 1% VAF in Dx-mono, demonstrating that this subclone preexisted at a low level in the monocytic subpopulation at diagnosis (Fig. 4D, clone 2 in pink). Although not confirmed due to limitation of sequencing depth, we infer that the KRAS^G12V mutation likely also preexisted as a phenotypically monocytic subclone at diagnosis (Fig. 4D, clone 3 in red). Therefore, overall, our WES analyses of Pt-12 and Pt-65 suggest dominant monocytic disease at relapse preexisted within the monocytic subpopulation of cells at diagnosis.

For patient 12 (Pt-12), the WES analysis shows the same monocytic clone at diagnosis and relapse (Fig. 4C, clone 2). However, based on the clinical course of the patient (i.e., initial CR followed by relapse), we inferred that biological properties of the monocytic subclone had at least partially evolved between diagnosis and relapse. To investigate further, we performed Cellular Indexing of Transcriptomes and Epitopes by Sequencing (CITE-seq; ref. 38) to simultaneously profile surface antigen expression and transcriptomics at single-cell resolution. Analysis of the CITE-seq data revealed three major clusters at diagnosis and two major clusters at relapse for Pt-12 (Fig. 5A; Supplementary Fig. S5A). Shown by surface antigen expression, the three clusters at diagnosis consisted of a lymphocytic cluster (CD3⁺/CD33⁻), a primitive myeloid cluster (CD34⁺/CD33⁺), and a monocytic myeloid cluster (CD11b⁺/CD33⁺; Supplementary Fig. S5A, left, and S5B). In contrast, the relapse specimen was comprised of a lymphocytic cluster (CD3⁺/CD33⁻) and a monocytic myeloid cluster (CD11b⁺/CD33⁺), but no primitive cluster, consistent with our flow cytometric analysis (Supplementary Fig. S5A, right, and S5B). Importantly, the transcriptomes of lymphoid clusters observed at diagnosis and relapse were quite similar, as depicted by their close proximity in the UMAP projection (Fig. 5A). In contrast, the mono-AML populations at diagnosis versus relapse were completely distinct (Fig. 5A), suggesting a significant change in the transcriptional profile of relapsed disease.

To characterize the transcriptional changes at relapse, we first analyzed the CITE-seq data using “clustifyr” (bioRxiv doi: 10.1101/855064), a program that assigns each cluster to its closest normal hematopoietic lineage counterpart according to transcriptomic similarity. The clustifyr analysis showed that the primitive cluster at diagnosis is indeed transcriptionally similar to hematopoietic stem cells (HSC), and the monocytic cluster at diagnosis is most similar to mature monocytes (Fig. 5B, blue and pink). Interestingly, within the relapse monocytic cluster, clustifyr identified two subpopulations including one that is most similar to mature monocytes and one that is closer to monocytic progenitors (CFU-monocytes; Fig. 5B, pink and red). These data suggest that a significant subpopulation of monocytic cells at relapse have acquired unique biology that clearly differs from the monocytic cells at diagnosis.

Next, using “gprofiler” pathway analysis, we observed significant upregulation of an LSC signature in the Rl-mono population (Fig. 5C and D). Intriguingly, the signature is derived specifically from MLL-rearranged leukemia (MLL⁺). This is notable because unlike most other previously reported LSC signatures, MLL-specific LSCs have a myeloid phenotype (39, 40). Thus, in comparison with Dx-mono, the Rl-mono population appears to have dedifferentiated to a more stem-like phenotype, while retaining myeloid characteristics. Furthermore, the Rl-mono subpopulation shows increased HOXA9 and OXPHOS signatures, characteristics of MLL⁺ leukemia and VEN+AZA-resistant leukemias, respectively (Fig. 5D). Similar results were observed when the gprofiler analysis was used to compare the “CFU-monocytes” subpopulation with the “Monocytes” subpopulation within the parent Rl-mono cluster (Supplementary Fig. S5C). Of note, Dx-prim uniquely enriches for HSC (normal stem cells) and non-MLL LSC gene-expression signatures including the GOODELL_HSC, EPPERT_HSC, and NG_LSC gene sets relative to Dx-mono, consistent with their expected non-MLL LSC nature (Supplementary Fig. S5D). Further, we plotted the enrichment score of a curated MLL-specific LSC signature obtained by Somervaille and colleagues and two other non-MLL LSC signatures obtained by Eppert and colleagues and Ng and colleagues (Supplementary Table S4; refs. 25, 26, 39). This analysis showed that although the Rl-mono cluster loses non-MLL LSC gene-expression signatures (Supplementary Fig. S5E), the “CFU-monocytes” subpopulation highly expresses the MLL-specific LSC signature (Fig. 5E). Consistently, expressions of the HOXA9 and MEIS1 genes were up in the Rl-mono cluster relative to Dx-mono and Dx-prim (Fig. 5F and G), further suggesting activation of an MLL-specific LSC program in VEN+AZA relapsed monocytic AML.

Finally, we asked if relapsed monocytic AML retains MCL1 dependence. Despite loss of BCL2 in both Dx-mono and Rl-mono clusters, MCL1 expression is sustained in both and presents a trend of increase at relapse (Fig. 5H). This result was corroborated by analysis of the single-cell RNA-seq data from van Galen and colleagues showing simultaneous loss of BCL2 and gain of MCL1 expression in monocytic AML subpopulations (Supplementary Fig. S5F; ref. 41). Importantly, viability assays showed that VU661013+AZA performed significantly better than VEN+AZA in eradicating relapsed monocytic cells from Pt-12 and Pt-65 (Supplementary Fig. S5G). Further, CFU assays showed that the VU661013+AZA regimen is significantly better than VEN+AZA in reducing the CFU capacity of relapsed monocytic AML from Pt-69 (Supplementary Fig S5H and S5I). Together, these findings indicate that relapsed monocytic AML retains dependence on MCL1.

Our analysis of patients treated with venetoclax-based therapy has revealed several previously unrecognized characteristics of AML. First, genetic subclones can manifest in distinct developmental stages, where the intrinsic properties of monocytic differentiation mediate venetoclax resistance. Second, leukemic cells at different developmental stages rely on distinct mechanisms to mediate energy metabolism. In particular, regulation of OXPHOS through BCL2 is more prevalent in primitive AML (12, 42–44), whereas reliance on MCL1 appears to be more important in monocytic AML. Third, intrapatient heterogeneity among LSCs underlies therapeutic resistance. Specifically, primitive or monocytic LSC phenotypes can coexist in the same patient with AML (i.e., MMP-AML) and demonstrate differential response trajectories to venetoclax-based therapy. Together, these data clearly indicate that to improve outcomes, AML therapies must be designed so as to target unique biological properties found in distinct AML subclones, developmental stages, and LSC subpopulations.

A notable aspect of these findings is the observation that AML cells can apparently switch from BCL2 to MCL1 dependence to drive energy metabolism as cells acquire a more differentiated developmental state. This finding suggests developmental heterogeneity is a previously underappreciated factor in determining therapeutic response. Based on our findings in AML, we note that similar considerations may be relevant in other forms of cancer. As we have shown, loss of BCL2 and sustained MCL1 dependency is a conserved property of monocytic development in both normal and malignant hematopoiesis. Others have shown similar switches can be found in the context of neuronal development and hematopoietic T/B-cell development (45, 46). Therefore, it is possible that in cancer types where differing developmental stages are evident, the relative reliance on BCL2 or MCL1 may vary in a manner analogous to what we have observed in AML.

Although the correlation between myeloid differentiation and venetoclax resistance seems clear, the underlying genetic events that may drive this process are as yet not well characterized. Intriguingly, the 2 patients profiled in this article had relapsed disease arising from RAS mutations, and we report here that RAS pathway mutations are a univariate, but not multivariate, predictor of refractory response to VEN+AZA (unlike FAB-M5 being highly significant in both univariate and multivariate analyses). These data suggest a possible link between the RAS pathway and myeloid development. To this point, analysis of data from the BEAT-AML study (47) shows coenrichment between FAB-M5 and RAS mutations (data not shown). In addition, mutations in RAS pathway members such as PTPN11 have been strongly associated with FAB-M5 in previous studies (48). Further assessment of the RAS pathway and other genetic factors is ongoing, but evidence thus far suggests that genetic factors such as RAS mutations may contribute to venetoclax resistance either through epistatic mechanisms and/or by driving monocytic differentiation.

An intriguing finding revealed by the WES analysis is the possibility that genetically distinct subclones with independent LSC populations may co-occur in some patients. This phenomenon was observed for both patients in our study (Fig. 4), albeit with differing presentations. For one patient (Pt-12), the diagnostic specimen had two predominant genetically independent subclones, each with its own LSC population (clones 1 and 2, Fig. 4C). Upon treatment with VEN+AZA, the monocytic clone was strongly selected. However, it is important to note that the patient did initially respond to therapy, implying that the relapse population must have evolved to some degree. In fact, CITE-seq analysis in this patient confirmed that the transcriptional profile of monocytic cells at relapse had markedly changed, with the selection and/or acquisition of a more stem-like profile. In contrast, at presentation the second patient (Pt-65) had primitive and monocytic populations with identical genetic characteristics (clone 1, Fig. 4D), implying a parent–progeny relationship. Notably, VEN+AZA therapy appears to have been highly effective in eradicating clone 1 but resulted in selection for rare genetically distinct monocytic subclones 2 and 3 that were resistant. Thus, in both patients, the data clearly indicate at least two genetically distinct LSC populations, one of which was ultimately resistant to venetoclax-based therapy. The prevalence of this type of dual LSC pathogenesis remains to be determined, and of course additional resistance mechanisms may emerge as relapse to VEN+AZA continues to be characterized in the future.

Our data suggest that prospective identification of VEN+AZA-resistant subpopulations will be of clinical value. Analogous to the previously employed FAB classification system, our approach thus far has been to use flow cytometry to define stages of differentiation, and we utilized SSC, CD45, CD117, CD11b, CD64, and CD68 to distinguish monocytic AML versus primitive AML. This preliminary panel provides a starting point toward developing a rapid tool to make treatment decisions related to the use of VEN+AZA. A more refined method will likely include additional features such as genetic mutations and/or metabolic properties.

Lastly, our studies have additional implications for the design of improved AML therapies. For elderly or otherwise unfit patients, VEN+AZA is well tolerated and provides deep and durable remissions for the majority of patients. For patients who respond well, we believe that identifying and targeting the VEN+AZA-resistant subpopulations is key to improved remission duration. To this end, well-tolerated consolidation regimens that selectively target VEN+AZA resistance are an important future objective. Based on the findings presented here, MCL1 inhibition may be one such strategy. Similarly, we have recently reported that inhibition of fatty-acid oxidation may also be a strategy to enhance the efficacy of VEN+AZA (11). Regardless of the specific strategy, treating patients prior to the emergence of resistant disease is likely to be highly beneficial. For the problem of disease that is refractory to VEN+AZA, our findings suggest it should be possible to prospectively identify such individuals and treat them with alternative therapies that enhance, or perhaps replace, venetoclax.

One hundred consecutive patients with newly diagnosed AML and no prior therapy who were prescribed VEN+AZA at our institution and had taken at least one dose of either therapy were identified. Patients were treated over a period from January 2015 to October 2019. Six patients were not included in the analysis because they had no documented follow-up. Of the 100 evaluable patients, all received the same regimen of VEN+AZA, as previously reported (6); some were treated in the context of a clinical trial, whereas others were treated outside of a clinical trial, with off-label use (prior to approval of venetoclax in November 2018) or on-label use after approval. Baseline details of all 100 patients are summarized in Supplementary Table S1. The University of Colorado Institutional Review Board approved a request to retrospectively analyze these patients (#19-0115). A single experienced hematopathologist (J. Schowinsky) categorized all patients into appropriate FAB groups by retrospectively reviewing the microscopic descriptions of all baseline pathology reports, including all of the flow cytometry histograms and the IHC studies, if performed. Responses were assessed in accordance with the European Leukemia Network (8). Nonresponding patients were defined as those who failed to achieve a CR, CRi, PR, or MLFS as described previously (6).

Logistic regression was used to determine the effects of age, sex, prior hematologic disorder, complex cytogenetics, ELN risk group, RAS pathway, TP53, IDH1/IDH2, NPM1, FLT3-ITD, and ASXL1 mutation, and FAB-M5 classification on refractory status. Univariate logistic regression was used to assess the effect of each predictor on its own, whereas multivariate logistic regression was used to assess their effects in aggregate. All univariate predictors with a P ≤ 0.10 were included in a multivariate model.

Primary human AML specimens were obtained from apheresis product, peripheral blood, or bone marrow of patients with AML who gave written informed consent for sample procurement on the University of Colorado tissue procurement protocols (Colorado Multiple Institutional Review Board Protocol #12-0173 & #06-0720). See Supplementary Table S2 for clinical information of primary AML specimens.

Primary human AML cells were stored in freezing media composed of 50% FBS (Corning), 10% DMSO (Sigma), and 40% Iscove's Modified Dulbecco's Medium (IMDM; Gibco) and then cryopreserved in liquid nitrogen. Freshly thawed or sorted cells were cultured in cytokine-added Serum Free Media (cytokine⁺ SFM) at 37°C, in 5% CO₂ incubator. SFM is composed of IMDM (Gibco), 20% BIT-9500 (Stemcell Technologies), 10 μg/mL LDL (Low Density Lipoprotein, Millipore), 55 μmol/L 2-Mercaptoethanol (GIBCO), and 1% Pen/Strep (Gibco). Cytokine⁺ SFM was made by supplementing the SFM with FLT3, IL3, and SCF cytokines (PeproTech), each at 10 ng/mL. All reagents are listed in Supplementary Table S5.

0.5 × 10⁶ freshly thawed primary AML cells were stained with our phenotyping panel containing antibodies against human CD45, CD34, CD117, CD11b, and CD64 at 4°C for 15 minutes, washed with ice-cold FACS buffer, fix/perm for 20 minutes using the BD Cytofix/Cytoperm kit (BD), washed with perm/wash, stained with intracellular antibody against human CD68 in perm/wash buffer for 30 minutes, washed with FACS buffer, resuspended in FACS buffer, and analyzed on BD FACSCelesta flow cytometer (BD). FCS files were analyzed on Flowjo 10.5.3 (Flowjo).

Primary AML cells were stained with anti-CD45, CD3, and CD19 in FACS buffer (PBS + 0.5% FBS) at 4°C for 15 minutes. The cells were then washed in FACS buffer and stained with 5 μmol/L CellRox-DeepRed dye (Life Technologies) in FACS buffer at 37°C in CO₂ incubator for 30 minutes. After CellRox staining, the cells were washed twice and resuspended in FACS buffer containing 10 μmol/L DAPI. The BD FACSAria II Cell Sorter was used to analyze ROS profile and sort. For sorting ROS-low LSCs from primitive AML cells, we first gated on live (DAPI⁻), nonlymphocytic (CD3⁻, CD19⁻), and primitive subpopulation (CD45-medium, SSC-low), and then gated on ROS-low cells (lower 20% of CellRox histogram). For sorting ROS-low cells from monocytic AML cells, we first gated on live (DAPI⁻), nonlymphocytic (CD3⁻, CD19⁻), and monocytic subpopulation (CD45-bright, SSC-high), and then gated on ROS-low cells (lower 20% of CellRox histogram). ROS-high cells were sorted similarly except higher 20% of CellRox was gated in the second step. The primitive and monocytic phenotypes were validated using our phenotyping panel on BD FACSCelesta (BD).

Freshly thawed and/or sorted AML cells were preincubated in cytokine⁺ SFM for 1 hour before treatment with various regimens prior to functional or viability assays. All drugs in the combo regimens were added simultaneously. Stock concentration of venetoclax and VU661013 is 10 mmol/L in DMSO; azacitidine is 40 mmol/L in saline.

All viability assays were performed after 24 hours of drug treatment unless otherwise specified. Cells were pelleted and stained in 1 x Annexin V staining buffer containing fluorophore conjugated Annexin V for 15 minutes in 4°C. Stained cells were then resuspended in 1 x Annexin V buffer containing DAPI and 0.5% FBS and analyzed on BD FACSCelesta (BD). Viability was determined by percentage of Annexin V-, DAPI- cells within the parent population.

The two independent siMCL1 sequences were purchased directly from Dharmacon's ON-TARGETplus siRNA Reagents collection: siMCL1-#B (J-004501-16) and siMCL1-#C (J-004501-15). The lyophilized siRNA products were resuspended in RNAse/DNAse-free water at 5 μmol/L and used as stock solution. Electroporation of siRNA was performed using the Neon Electroporation Transfection System (Thermo) according to the manufacturer's protocol. Specifically, for every 2 × 10⁶ primary AML cells, the cells were spun and resuspended in 80 μL Buffer T with the addition of 20 μL siRNA stock solution in an Eppendorf tube. Electroporation was then performed using microtips (10 μL capacity) in a cuvette containing 3 mL of Buffer E with the following electroporation settings: 1,600 V, 10 ms, 3 pulses. To complete the 100 μL volume containing 2 × 10⁶ cells, 10 continuous repeats of 10 μL each were performed, the cells from all 10 runs were sequentially combined into one well containing media, and the electroporated cells were then cultured for subsequent analysis.

Extracellular flux assay XF96 kit was used to measure oxygen consumption rate (OCR). Drug-treated or siRNA-exposed primary AML cells were plated in Cell-Tak (Corning)–coated XF96 cell culture microplates at 200 K cells/well in five replicates and measured according to the manufacturer's protocol. OCR was measured at basal and after injection of 5 μg/mL oligomycin (Sigma), 2 μmol/L FCCP (Sigma), 5 μmol/L Antimycin A (Sigma), and 5 μmol/L Rotenone (Sigma).

Cells were harvested and lysed in 1x Laemmli Sample Buffer (Bio-Rad) at 5 × 10⁶/mL. About 50K to 200K cells per lane were loaded and resolved on 6% to 12% SDS-PAGE gels, transferred to PVDF membranes, and blocked in 5% milk in TBS with 0.1% Triton X-100 (Sigma). After incubation, the PVDF membranes were incubated with primary antibodies at 4°C overnight, washed, incubated with secondary antibodies at room temperature for 2 hours, and subjected for imaging on the ChemiDoc Imaging System (Bio-Rad).

Primary AML or normal CBMC cells were seeded in human methylcellulose complete media (R&D Systems) at 50K to 250K cells/mLor 10K cells/mL, respectively, and grown in 37°C, 5% CO₂ incubator for 2 to 4 weeks before counting on stereomicroscope. For drug treatment, primary AML cells were treated in culture dishes overnight with various inhibitors in cytokine⁺ SFM.

One day prior to transplant, freshly thawed primary AML cells were treated in culture dishes overnight with various inhibitors in cytokine⁺ SFM. NSG-S mice were conditioned with 25 mg/kg busulfan via i.p. injection. Second day at injection, overnight-treated primary AML cells were washed with FACS buffer and resuspended in saline at 50 × 10⁶/mL concentration. Anti-human CD3 antibody (BioXCell) was added at a final concentration of 1 μg/10⁶ cells to avoid potential graft-versus-host disease. Per mouse, 5 × 10⁶ cells in 0.1 mL saline were tail-vein injected; there were 8 to 10 mice per experiment group. Mice engrafted with primary AML cells were sacrificed after 5 to 8 weeks. All animal experiments were approved by the University of Colorado Anschutz (Institutional Animal Care and Use Committee) under protocol number 0308.

mRNA isolation, quality check, library construction, and sequencing were performed according to same protocol used previously (21). Single-end reads of 100 nt were generated for each sample on the Illumina HiSeq2500 platform. For analysis, raw reads were demultiplexed using bcl2fastq version 2.19.0. Quality filtering and adapter removal were performed using Trimmomatic-0.36 with the following parameters: “TRAILING:13 LEADING:13 ILLUMINACLIP:adapters.fasta:2:30:10 SLIDINGWINDOW:4:20 MINLEN:35.” Processed/cleaned reads were then mapped to the Homo sapiens reference genome sequence (GRCh38, hg38) with STAR-2.6.0c given the following parameters: “–twopassMode Basic –runMode alignReads –genomeDir $ –readFilesIn ${SAMPLE} –outSAMtype BAM SortedByCoordinate –outSAMstrandField intronMotif –outFilterIntronMotifs RemoveNoncanonical.” The subread-1.6.1 package (featureCounts) was used to derive gene counts given the following parameters: “-s 2 -t exon -g gene_name.” Differential expression analysis and data normalization were performed using DESeq2-1.16.1 with an adjusted P value threshold of 0.05 within an R-3.4.1 environment. A batch factor was given to the differential expression model in order to control for batch differences.

Design of oligo-antibody panel is detailed in Supplementary Table S6. For constructing CITE-seq libraries, the following modifications were added to the standard 10x Genomics V3 protocol. During cDNA amplification, an additive primer was added at 0.2 μmol/L to selectively amplify the antibody tags. In the following purification, the supernatant is removed in the first cleanup step to be further purified, whereas the beads are cleaned to generate the gene-expression libraries. The supernatant is cleaned with 1.4x Ampure Select beads (Agencourt) and eluted in 40 μL of EB buffer. CITE-seq libraries are amplified with Illumina compatible primers for 14 cycles and purified with 1.2x Ampure Select beads. All libraries were normalized by Qubit (Invitrogen) and Tapestation (Agilent) analysis. Libraries were sequenced on a NovaSeq 6000 (Illumina) for a depth of 5,000 reads per cell for CITE-seq libraries. For analysis, CITE-seq data were first processed through the 10x Genomics Cell Ranger pipeline (v3.1) and then analyzed in R (v3.6) with the Seurat (v3.1.1) package 1. Cell barcodes with 25% or more mitochondrial RNA content were excluded from analysis. RNA counts were normalized using scran to estimate normalization size factors, and ADT counts were normalized using centered log-ratio normalization in Seurat. Cell type similarity and gene signature analyses were conducted using R package “clustifyr” (https://github.com/rnabioco/clustifyR). Cell clusters from CITE-seq data were assigned to cell-type identities based on Spearman correlation against a reference dataset generated from microarray data of purified populations of human hematopoietic cells (GSE24759), comparing the expression of variable genes determined by the VST method in Seurat. Gene signature scoring for both positive and negative markers was calculated by clustifyr, first on a per-cell basis by Spearman correlation, and then averaged over each cluster. Comparisons between curated gene sets were performed using hypergeometric tests or the Jaccard index. Gene set and pathway enrichment analysis was performed in R via the “gprofiler2” package, an interface to the gProfiler web suite, using previously curated AML gene sets.

WES libraries were generated using the Agilent Sure Select XT Exome prep kit with 200 ng of input as per protocol (Agilent). The probe used was SureSelect XT Human All Exon V7 (Agilent). Libraries were normalized by Qubit (Invitrogen) and Tapestation (Agilent) and sequenced on a Novaseq 6000(Illumina) to obtain 400x coverage. For analysis, we adopted the IMPACT (Integrating Molecular Profiles with Actionable Therapeutics) pipeline that we previously developed (49). IMPACT links variants detected from WES to actionable therapeutics. Briefly, IMPACT takes sequence data as input and outputs a VCF file containing predicted deleterious mutations. The sequencing reads were mapped to the human hg19 reference exome using the Burrows–Wheeler Aligner. SAMTools and BCFtools (v1.1) were utilized to detect variants from the BAM file and output into a VCF file. In the IMPACT pipeline, we used ANNOVAR (v2014-11-12) to annotate the variants. Synonymous and intronic variants were removed. Variants were further analyzed by deleterious prediction tools such as SIFT and PolyPhen2. We also focused on 49 genes commonly mutated in AML (http://raindancetech.com/thunderbolts-myeloid-panel/) and cancer-related mutations reported in the Catalogue of Somatic Mutations in Cancer (COSMIC; https://cancer.sanger.ac.uk/cosmic) to infer clonal dynamics in paired diagnosis and relapse specimens for VEN+AZA trial patients.

GSEA was performed using GSEA version 3.0 (Broad). Normalized count matrix produced from DESeq2 analysis was formatted into GCT files containing expression values. CLS files were built to label biological states. When performing GSEA, several gene set databases were used including: c2.cp.kegg.v6.0.symbols.gmt, c2.cp.reactome.v6.0.symbols.gmt, and an in-house database containing a collection of 791 gene sets of our interest. Following parameters were used: Number of permutations = 1,000, permutation type = gene_set. Otherparameters were left at default values.

After 16 hours of knockdown with siSCR control or siMCL1, primary AML cells were washed in ice-cold PBS and snap-frozen for metabolomic analysis. Technical replicates of four were used per condition, 100K cells per replicate. Metabolomic methods are detailed previously (11, 50). Metabolite levels were normalized to protein quantification.

Methods used for statistical analyses were detailed in figure legends. In all figures, ns indicates a not significant P value of >0.05; *, **, ***, and **** indicate P < 0.05, P < 0.01, P < 0.001, and P < 0.0001, respectively.

The data that support the findings of this study are available from the corresponding author upon reasonable request. Patient-related clinical data are included as Supplementary tables. Bulk RNA-seq (GSE132511), CITE-seq (GSE143363), and WES (PRJNA600769) data have been uploaded to public databases.

D.A. Pollyea is consultant at AbbVie, Genentech, and Servier, and reports receiving a commercial research grant from AbbVie. M.R. Savona reports receiving commercial research grants from Astex Corporation, Incyte Corporation, and Takeda Pharmaceutical Company; has ownership interest (including patents) in Karyopharm Therapeutics; has a consultant/advisory board relationship with Bristol-Myers Squibb, Celgene Therapeutics, TG Therapeutics, AbbVie, and Ryvu; and has received other remuneration from Boehringer Ingelheim. C.A. Smith is a consultant for AbbVie. No potential conflicts of interest were disclosed by the other authors.

Conception and design: S. Pei, D.A. Pollyea, A. Gustafson, M.R. Savona, C.A. Smith, C.T. Jordan

Development of methodology: S. Pei, D.A. Pollyea, M.L. Amaya, A. D'Alessandro, H.E. Ramsey, C.T. Jordan

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): S. Pei, D.A. Pollyea, A. Gustafson, B.M. Stevens, M. Minhajuddin, A. Inguva, A. Winters, H. Ye, A. Krug, C.L. Jones, N. Khan, J. Ponder, J.M. Ashton, T. Nemkov, A. D'Alessandro, J.A. Gutman, M.R. Savona, C.A. Smith

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): S. Pei, D.A. Pollyea, A. Gustafson, M. Minhajuddin, R. Fu, K.A. Riemondy, A.E. Gillen, R.M. Sheridan, J. Kim, J.C. Costello, J. Schowinsky, D. Abbott, A. Hammes, J.R. Myers, J.M. Ashton, T. Nemkov, A. D'Alessandro, M.R. Savona, C.T. Jordan

Writing, review, and/or revision of the manuscript: S. Pei, D.A. Pollyea, J. Kim, J.C. Costello, M.L. Amaya, B. Adane, J. Schowinsky, D. Abbott, A. Hammes, M.R. Savona, C.T. Jordan

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): M. Minhajuddin, A. Krug, J. Ponder, T. Nemkov, C.T. Jordan

Study supervision: D.A. Pollyea, C.A. Smith, C.T. Jordan

We thank all of the patients and their families, as well as the nurses, pharmacists, and advanced practice practitioners who were involved in their care. C.T. Jordan is generously supported by the Nancy Carroll Allen Chair in Hematology Research and the NIH (R01CA200707). D.A. Pollyea is supported by the University of Colorado Department of Medicine Outstanding Early Career Scholar Program, the Robert H. Allen MD Chair in Hematology Research, and the Leukemia and Lymphoma Society Scholar in Clinical Research Award. M.R. Savona is a Leukemia and Lymphoma Society Clinical Scholar and is supported by an E.P. Evans Foundation Discovery Research Grant and the Biff Ruttenberg Foundation. A. Inguva is supported by NIH MSTP T32 (GM008497). A. D'Alessandro is supported by funds from the Boettcher Foundation (Webb-Waring Early Career 2017). C.L. Jones is supported by a Leukemia and Lymphoma Society Special Fellow Award and the American Cancer Society (25A5072). We appreciate bioinformatic support provided by the Informatics Fellows of the RNA Bioscience Initiative at the University of Colorado School of Medicine. We also appreciate funding support for the Biostatistics and Bioinformatics Shared Resource and the Genomics Shared Resource of the University of Colorado Cancer Center (P30CA046934).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.