Bone Marrow Surveillance of Pediatric Cancer Survivors Identifies Clones that Predict Therapy-Related Leukemia

The adoption of dose-intensive chemotherapy in the treatment of childhood cancers has led to significant improvements in treatment response and overall survival (OS). However, therapy intensification has also led to an increase in therapy-related malignancies in pediatric cancer survivors. Therapy-related myelodysplastic syndrome and acute leukemias (t-MDS/AL) are common secondary malignancies in the early time period after treatment (1–3). They are associated with aggressive and chemoresistant diseases, resulting in very poor outcomes (4). With intensive therapy, including allogeneic hematopoietic stem cell transplant, OS reaches 40% to 50% for patients able to achieve remission (2, 5–7). However, outcomes remain dismal for those who do not, with OS as low as 5% at 5 years (6, 7). There is a pressing need for the development of surveillance protocols to detect the earliest evidence of t-MDS/AL in order to apply risk-adapted treatment strategies.

Neuroblastoma is the most common pediatric extracranial solid tumor (8). Approximately half of patients present with high-risk disease and require high-intensity multimodality therapy including chemotherapy, radiation, and immunotherapy (8–10). Despite intensive treatment, nearly half of high-risk patients with neuroblastoma relapse (11), necessitating further, and often prolonged, exposure to cytotoxic therapy. Consequently, rates of t-MDS/AL in these patients have been reported to be as high as 5% to 10% (1, 12, 13).

Historically, the pathophysiological basis for the development of t-MDS/AL has been linked to DNA damage from genotoxic therapy. Associations have been made between distinct types of oncologic therapies, specific genetic alterations, and patterns of transformation. For example, KMT2A rearrangements are associated with topoisomerase II inhibitors and early development of t-MDS/AL after oncologic therapy, whereas TP53 mutations have been associated with complex cytogenetics and later time to transformation (4, 14–17). In addition, some patients acquire transient cytogenetic abnormalities in bone marrow (BM) after chemotherapy without developing t-MDS/AL (18), whereas others develop t-MDS/AL without evidence of prior cytogenetic abnormalities. Together, these phenomena suggest diverse evolutionary paths leading to development of t-MDS/AL following cancer treatment. Recent studies in adult oncology patients have demonstrated that somatic mutations often precede cancer therapy (15, 16), likely associated with age-related clonal hematopoiesis (CH), and that patients with preexisting CH are at increased risk for developing t-MDS/AL (19–22). This increased risk may be due to a growth advantage of CH clones involving mutations in DNA-damage response pathway genes (TP53, PPM1D, and CHEK2) in patients undergoing oncologic therapy (15, 23). Limited data regarding the genetics of t-MDS/AL in pediatric patients suggest that chemotherapy-induced DNA damage may play a larger role in the oncogenesis in this patient population (24, 25). However, to our knowledge, there is no published systematic evaluation of the evolution of t-MDS/AL, and the role of CH as a predictor of t-MDS/AL in pediatric cancer cohorts has not been systematically investigated.

Here, we sought to study the genetic evolution of t-MDS/AL in patients with neuroblastoma using prospectively collected sequential BM samples archived over decades to evaluate the relationships between mutation acquisition, clonal dominance, and clinical presentation of t-MDS/AL. We use established clinical sequencing assays to evaluate their utility in a routine clinical surveillance setting. Last, we describe the prevalence of CH in a large pediatric oncology population.

Our primary study cohort consisted of a subset of patients with neuroblastoma treated at Memorial Sloan Kettering Cancer Center (MSKCC) diagnosed between 1988 and 2009 and enrolled on institutional protocol #00-109 (NCT 00588068). The study was conducted in accordance with U.S. Department of Health and Human Services approved assurance (FWA00004998) and with the approval of the Institutional Review Board (IRB). Informed written consent was obtained from each subject or the subject's guardian prior to study enrollment. All patients enrolled on this study had BM samples collected every 3 to 6 months for neuroblastoma and leukemia surveillance studies. Mononuclear cells were purified by Ficoll and cryopreserved in DMSO from these heparinized BM samples. Patients in the transformation cohort were chosen on the basis of development of t-MDS/AL. Given historically high rates of t-MDS/AL in this patient population, patients under the clinical care of the neuroblastoma team at MSKCC undergo annual BM surveillance for the presence of cytogenetic abnormalities to detect a preleukemic state. Patients who had cytogenetic abnormalities identified during routine surveillance but ultimately did not develop t-MDS/AL were studied as part of the transient cohort. Some patients in these cohorts had low-level cytogenetic abnormalities that did not meet strict criteria for clonality; however, they were included here due to such persistent abnormalities and/or multiple abnormalities on a single metaphase karyotype. We defined “clonal abnormalities” as having a burden by FISH ≥5% or at least two to three metaphase karyotypes as per ISCN standards (26) and noted “non-clonal” events wherever feasible. Last, a neuroblastoma-treated control cohort was selected from a group of patients also enrolled on one of two institutional treatment studies for patients with relapsed or refractory neuroblastoma and who did not have any history of cytogenetic or hematologic abnormalities to suggest t-MDS/AL.

A total of 52 patients with neuroblastoma with 199 samples (Table 1; Supplementary Fig. S1A and S1B) from three cohorts were included: (i) 17 patients with 88 samples in the transformation cohort; (ii) 14 patients with 74 samples in the transient cohort; and (iii) 21 patients with 37 samples in the neuroblastoma-treated control group. For the transformation cohort, t-MDS/AL diagnosis occurred at a median of 4.2 years (range, 0.9–18.2) from neuroblastoma diagnosis. In the transient group, median time to detection of a cytogenetic abnormality was 2.2 years from neuroblastoma diagnosis (range, 0.2–4.0), with no evidence of t-MDS/AL during follow up (median 7.4 years, range, 2.9–8.6). Neuroblastoma-control patients had a median leukemia-free follow-up of 8.1 years after diagnosis (range, 4.4–12.9).

All patients received similar induction therapy and similar doses of consolidative radiation to the primary site. Treatment for relapsed neuroblastoma consists of a fixed group of drugs, and consistent doses of radiation to resistant sites. Given that many patients received some treatment at outside institutions, making it more difficult to obtain exact drugs and dosages, and the nature of neuroblastoma therapy as described, cumulative treatment was determined based on total months of active therapy. All patients in the neuroblastoma-control cohort received total cumulative therapy at least as long as the transformation and transient cohorts as they were selected from a highly pretreated group of patients. Patient characteristics, including neuroblastoma treatment, are outlined in Table 1.

Patients in the transformation cohort commonly had high-risk cytogenetic abnormalities at diagnosis of t-MDS/AL (n = 7), whereas only 1 patient in the transient group had a clonal high-risk abnormality. High-risk abnormalities were defined as EVI1 rearrangements, abnormalities of chromosomes 5, 7, or 17p, or KMT2A rearrangements as per WHO classification and IPSS-R and ELN risk stratification schemas for myeloid malignancies and acute leukemias (27–29). Cytogenetic abnormalities in the transient group were predominantly not associated with high-risk hematologic malignancy.

Samples were collected throughout neuroblastoma treatment and during off-therapy follow-up. A subset of the originally collected samples was selected for study, targeted at beginning and ends of treatment intervals, including the time of leukemic transformation or cytogenetic abnormality whenever possible. Median time to first sample was 7.0 months (range, 0.7–48.3) and 4.6 (range, 0–28.2) months for the transformation and transient groups, respectively. Samples in the neuroblastoma-treated controls were selected more sparsely and directed at times later in each patient's treatment course. We aimed to select samples after a long therapy block at least 2 years after commencement of therapy and at the completion of all therapy. Median time to first sample in the control group was 38.3 months (range, 24.3–83.0).

Additional time points were identified to assess TP53 variants identified in transformation patients to better delineate the emergence of these clones relative to neuroblastoma diagnosis and treatment intervals. Every effort was made to obtain a sample from the earliest time point available, before the start of therapy whenever possible. One of the four mutations present in patient 118727 with a maximum variant allele frequency (VAF) of 3.8% on initial sequencing was not assessed by digital droplet PCR (ddPCR).

To evaluate the relevance of our findings across pediatric patients, we studied an additional cohort of pediatric and young adult patients with solid tumors. Eligible patients included those treated in the Department of Pediatrics at MSKCC who were ≤25 years old (median age, 11.1 years; range, 1.5 months–24.9 years), and who had their tumor profiled using our institutional sequencing platform (MSK-IMPACT: Integrated Mutation Profiling of Actionable Cancer Targets, see details below). Diagnoses included 144 patients with neuroblastoma and 513 patients representing diverse pediatric solid tumors (Supplementary Table S1). Most patients were enrolled on the IRB-approved institutional protocol #12-245 (NCT01775072). A subset of patients who underwent tumor-genomic profiling as standard of care did not directly consent, in which case an IRB waiver was obtained to allow for inclusion into this study. Patients included here were a subset of those previously described (23), and the samples included here were drawn from March 2014 to June 2018. MSK-IMPACT uses blood as a germline control to subtract private SNPs to identify tumor-specific variants in solid tumors. These blood controls were used to identify CH as described previously (23, 30, 31).

Targeted DNA sequencing was performed using HemePACT v3 and v4, which contain 585 and 576 genes, respectively, implicated in hematologic malignancies (full gene lists: Supplementary Tables S2 and S3). HemePACT data were used for the detection of small gene mutations, as well as focal- and arm-level copy-number abnormalities (CNA). RNA sequencing for recurrent gene fusions was performed using the ArcherDX FusionPlex PanHeme panel (full gene list: Supplementary Table S4), and high-confidence oncogenic fusions were confirmed through Archer's interface. Across all samples, median coverage for DNA capture of BM samples analyzed with HemePACT was 705x (mean, 717; range, 418–1,316; Supplementary Fig. S1C).

Participants in the extended pediatric cohort had a tumor and blood sample (as a matched normal) sequenced using MSK-IMPACT, an FDA-authorized hybridization capture-based next-generation sequencing assay encompassing all protein-coding exons from the canonical transcript of 341, 410, or 468 cancer-associated genes (Supplementary Table S4). DNA was extracted from formalin-fixed paraffin-embedded tumor tissue and patient-matched blood samples and sheared, and DNA fragments were captured using custom probes. MSK-IMPACT contains most of the commonly reported CH genes with few exceptions. Earlier versions of the panel did not contain PPM1D or SRSF2. In addition, three genes commonly reported to be observed in patients with malignancies, SRCAP, BRCC3, and ZNF318, were not included, the first two belonging to the DNA damage response (DDR) pathway. All subsequent sequencing, variant calling, and analysis was performed as previously described (23).

Assays specific for the detection of C135S, V172F, R175H, R213P, R249S, L257Q, and G266R in TP53 were designed and ordered through Bio-Rad. Cycling conditions were tested to ensure optimal annealing/extension temperature as well as optimal separation of positive from empty droplets. Optimization was done with a known positive control.

After PicoGreen quantification and quality control (QC) by Agilent BioAnalyzer (RRID: SCR_018043), 56 to 800 ng of DNA were used to prepare libraries using the KAPA Hyper Prep Kit (Kapa Biosystems KK8504) with 7 to 8 cycles of PCR. One hundred to 900 ng of each barcoded library were captured by hybridization in equimolar pools of 11 to 19 samples using the HemePACT (Integrated Mutation Profiling of Actionable Cancer Targets related to Hematological Malignancies) assay (Nimblegen SeqCap), designed to capture all protein-coding exons and select introns of 585 or 576 commonly implicated oncogenes, tumor suppressor genes, and members of pathways deemed actionable by targeted therapies. Captured pools were sequenced on an Illumina HiSeq 4000 (RRID: SCR_016386) or HiSeq 2500 (RRID: SCR_020123) in High Output or Rapid mode in a PE100 or PE125 run using the HiSeq 3000/4000 SBS Kit, TruSeq SBS Kit v4, or HiSeq Rapid SBS Kit v2 (Illumina) producing an average of 1,009X coverage per tumor.

After RiboGreen quantification and QC by Agilent BioAnalyzer, 3 to 20 ng of total RNA were prepared using the MSKCC Pan Heme FusionPlex Kit for Illumina (ArcherDX) according to the manufacturer's instructions. Briefly, cDNA was synthesized using random priming for reverse transcription. Molecular barcode adapters were ligated to cDNA fragments and multiplex PCR with gene-specific primers was used to enrich for genes of interest. Barcoded samples were pooled equimolar and sequenced on a HiSeq 2500 in Rapid Mode or a NextSeq 500 (RRID: SCR_014983) in a PE125 (HiSeq) or PE150 (NextSeq) run, using the HiSeq Rapid SBS Kit v2 or NextSeq 500/550 High Output Kit v2.1 (300 Cycles; Illumina). Each sample yielded on average 7M reads, and FASTQ files were uploaded to the Archer Analysis bioinformatics suite for processing.

Viably frozen cells were thawed at 37°C, pelleted, and resuspended in RLT (cell pellets were suspended directly), and nucleic acids were extracted using the AllPrep DNA/RNA Mini Kit (QIAGEN; catalog # 80204) according to the manufacturer's instructions. RNA was eluted in nuclease-free water and DNA in 0.5X Buffer EB.

Some samples were lysed and cryopreserved in TRIzol Reagent (Thermo Fisher; catalog # 15596018). Phase separation for these samples was induced with chloroform. RNA was precipitated with isopropanol and linear acrylamide and washed with 75% ethanol. The samples were resuspended in RNase-free water.

Viably frozen cells were thawed and pelleted and incubated for at least 30 minutes in 360 μL Buffer ATL + 40 μL proteinase K at 55°C, and DNA was isolated with the DNeasy Blood & Tissue Kit (QIAGEN; catalog # 69504) according to the manufacturer's protocol. DNA was eluted in 0.5X Buffer AE.

After PicoGreen quantification, 10 ng DNA were combined with locus-specific primers, FAM- and HEX-labeled probes, restriction enzyme (Mse I, Hind III, or Hae III), and digital PCR Supermix for probes (no dUTP). All reactions were performed on a QX200 ddPCR system (Bio-Rad; catalog # 1864001), and each sample was evaluated in two to six wells. Reactions were partitioned into a median of approximately 14,000 droplets per well using the QX200 droplet generator. Emulsified PCRs were run on a 96-well thermal cycler using cycling conditions identified during the optimization step (95°C 10′; 40 cycles of 94°C 30′ and 52–55°C 1′; 98°C 10′; 4°C hold). Plates were read and analyzed with the QuantaSoft software to assess the number of droplets positive for mutant DNA, wild-type DNA, both, or neither.

Raw sequence data were aligned to the human genome (NCBI build 37) using BWA (http://bio-bwa.sourceforge.net/; RRID: SCR_010910) version 0.7.17. PCR duplicate reads were marked with Picard tools (https://broadinstitute.github.io/picard/; RRID: SCR_006525) version 2.18.2. For alignment, we used the pcap-core dockerized pipeline version 4.2.1 available at https://github.com/cancerit/PCAP-core/wiki/Scripts-Reference-implementations.

QC of the fastq data and bam data was performed with FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/; RRID: SCR_014583) version 0.11.5 and Picard tools, respectively. In addition, sequential samples across patients were evaluated with CONPAIR (https://github.com/nygenome/Conpair). Mismatched samples were reviewed and resolved or excluded from further analysis.

Variants were called using a combination of variant callers. For single-nucleotide variants (SNV), we used CaVEMan (http://cancerit.github.io/CaVEMan/; RRID: SCR_017089) version 1.7.4, Mutect (https://software.broadinstitute.org/cancer/cga/mutect; RRID: SCR_000559) version 4.0.1.2, and Strelka (https://github.com/Illumina/strelka; RRID: SCR_005109) version 2.9.1. For small insertions and deletions, we used Pindel (http://gmt.genome.wustl.edu/packages/pindel/; RRID: SCR_000560) version 1.5.4, Mutect version 4.0.1.2, and Strelka version 2.9.1.

VAF was uniformly reported across all called variants using a realignment procedure (https://github.com/cancerit/vafCorrect). Likely artifact variants were filtered out based on: (i) The number of callers calling a given variant and the combination of filters from the triple callers. (ii) Variants with VAF <2% or less than 3 mutant supporting reads were excluded. In some cases, where specifically mentioned in text, VAFs <2% were reported if there was adequate coverage and the variant had been previously identified in a different sample from the same patient.

Samples sequenced using the ArcherDX Pan Heme panel as described above were uploaded to Archer's bioinformatics suite for analysis. Variants and fusions were visualized, manually reviewed, and filtered for artifacts. Details regarding the included fusions were downloaded from the suite and can be found in Supplementary Table S5.

Copy-number analysis based on targeted DNA sequencing was performed using methods previously described for calling mosaic copy-number CH using targeted sequencing data (32). Copy-number changes identified in a single sample across consecutive samples for any given patient were excluded from representation in any main or supplemental summary figures. Large arm-level events were filtered for using a size threshold of log(nbases)>16, and single events spanning both chromosome arms were reported separately, whereas multiple events on a single arm were merged.

Analysis was performed using R v3.5.0 (RRID: SCR_001905). Differences in cohort characteristics were evaluated using ANOVA, and outcomes between groups were evaluated using t tests (pairwise where appropriate), unless otherwise noted. Relative risk was calculated using the riskratio function from the fmsb R package (https://cran.r-project.org/web/packages/fmsb/fmsb.pdf) version 0.70, comparing the transformation with the transient cohorts.

Growth rates reported were calculated based on the change in VAF over change in time. Mathematically, this was calculated as log(VAF t₁/VAF t₀)/(t₁-t₀).

Code to generate all figures is available at https://github.com/spitzerb/tMN_CCR_2022.

We detected at least one leukemia-associated alteration in all samples collected at t-MDS/AL diagnosis, including a total of 39 putative oncogenic events (median, 2 per patient; range, 1–4), and 32 presumed somatic variants of uncertain significance. A summary of all oncogenic events is depicted in Fig. 1A and B and Supplementary Fig. S2. Complete details are outlined in Supplementary Tables S5, S6, S7, and S8.

The most frequent driver events at transformation were KMT2A rearrangements (n = 6 patients) and TP53 mutations (n = 9 variants, n = 5 patients). The remaining patients had mutations in common myeloid malignancy genes including NPM1, IDH1, PTPN11, NRAS, CUX1, STAG2, WT1, and PPM1D, at a median VAF 0.21 (range, 0.026–0.75), and 1 patient had only structural abnormalities. KMT2A rearrangements were mutually exclusive with TP53 mutations (with the exception of 1 patient with a clonal KMT2A translocation and subclonal TP53 mutation). At the time of transformation, 13 (76.5%) of t-MDS/AL cases had abnormalities identified on conventional cytogenetics, whereas all patients (100%) had events identifiable by our integrated capture-based sequencing strategy.

KMT2A rearrangements were identified in 6 patients with t-MDS/AL (Supplementary Table S5). All KMT2A rearrangements occurred within the canonical breakpoint region and included common fusion partners (ELL and AFF1), less common partners (FRYL and GAS7), and two novel partner genes (PRDM10 and DDX6). The fusion with PRDM10 was identified by both capture-based approaches and conventional cytogenetics at t-AML diagnosis, whereas the previously undescribed rearrangement involving DDX6 was only captured through RNA capture-based techniques. Both novel fusion partners are involved in histone methyltransferase activity (33, 34).

Recent studies in adult patients with cancer (15, 21, 22) have shown that early preleukemic mutations, in the form of CH, are often established before their primary cancer diagnosis (15). However, CH is remarkably less prevalent in younger patients (23, 30, 35, 36), and thus, it is not yet clear whether the same processes occur in pediatric cancer patients. In our cohort, 5 patients harbored TP53 mutations at t-MDS/AL diagnosis. We used ddPCR assays, which survey 10⁴ molecules to assess the timing of emergence of the TP53 mutations relative to neuroblastoma diagnosis and therapy (Fig. 1C; Supplementary Fig. S3A). A total of 28 samples were tested, including samples from additional time points not included in the initial sequencing cohort (Fig. 1C; Supplementary Fig. S3A; Supplementary Table S9). In all 5 patients with TP53 mutations, the earliest time point of mutation detection occurred after the first treatment window. For 4 of 5 patients, the first sample available was obtained within 7 months (range, 0.2–6.5 months) of neuroblastoma diagnosis, and we did not find evidence of TP53 mutations in these samples. This suggests that these mutations were acquired on therapy or were below the limit of ddPCR detection before systemic therapy for neuroblastoma. For the final patient, the first available BM sample was collected >4 years after neuroblastoma diagnosis and was positive for the two TP53 mutations identified at transformation. In this group, the burden of TP53 mutation increased both on and off therapy; however, there was a greater increase in VAF during periods of chemotherapy and/or radiation compared with no therapy or immunotherapy alone (Supplementary Fig. S3B).

We analyzed serial BM samples from individual patients in the transformation cohort to study the evolution toward leukemia. De novo variant calling was performed for all samples and identified mutations were also evaluated across all related samples. Disease-defining gene mutations were identified in 11 of 14 patients with any gene mutation present at diagnosis at a median of 11.8 months (range, 1.3–32.4) before t-MDS/AL presentation at a clinically detectable VAF of 0.02, corresponding to a median of 3 years (36.8 months; range, 10.8–78.1) following neuroblastoma diagnosis. Median VAF at initial detection was 0.08 (range, 0.03–0.33) and generally increased over time leading to transformation (Fig. 2A–C; Supplementary Fig. S4A and S4B). Although all but one KMT2A rearrangement was identified clinically at t-MDS/AL diagnosis, routine cytogenetic surveillance captured only two of six rearrangements prior to disease transformation (at 2 and 28.5 months before overt disease). However, our RNA-capture approach identified KMT2A rearrangements in 5 of 6 patients prior to transformation (median 28.5 months before transformation, range, 2.3–61.8, and 21.9 months after neuroblastoma diagnosis, range, 10.8–28.4; Fig. 2B; Supplementary Fig. S4B). Altogether, detection of mutations and KMT2A rearrangements identified aberrations associated with a subsequent t-MDS/AL in 13 of 17 patients with neuroblastoma at a median of 15.0 months before transformation (range, 1.3–32.4) and 26.9 months after initial neuroblastoma diagnosis (range, 10.8–78.1).

We observed diverse patterns of disease evolution (Fig. 2A–C; Supplementary Figs. S4 and S5). In 9 patients, acquisition of genomic alterations occurred sequentially leading up to transformation. These events included secondary mutations as well as oncogenic CNAs. For example, patient 118725 had a clonal hotspot NPM1 mutation 36 months before transformation, followed by acquisition of IDH1^R132H 12 months later, ultimately leading to leukemic transformation (Fig. 3A). In patients with mutations in TP53 (e.g., 107649, 118732), allelic imbalances were detected at the time of transformation (Fig. 3B; Supplementary Fig. S5), consistent with loss of heterozygosity at this locus. In addition, we observed parallel and convergent evolution toward the same oncogenic pathway. For example, several patients had multiple TP53 mutations (118727 and 118731), or mutations in TP53 and PPM1D (118730 and 118731), likely representing distinct clones (Supplementary Fig. S5; Fig. 3C). One patient (118728) had several RAS-pathway mutations and another (118725) had independent clones with IDH1 and TET2 mutations (Supplementary Fig. S5; Fig. 3C).

We also noted molecular abnormalities emerge and regress over time, with such oncogenic events occurring in 5 patients preceding t-MDS/AL. For example, patient 118731 had a clonal oncogenic BCOR mutation 2 years before transformation, which fell below the limit of clinical detection at the time of t-MDS (Fig. 3C). Patient 118729 had a transient IDH2^R140Q mutation 1 year prior to transformation in the context of a longstanding (total >2 years) preleukemic KMT2A rearrangement (Supplementary Fig. S5). Other genes involved in transient molecular abnormalities in the transformation cohort included TET2, BCOR, RUNX1, and RB1. This reveals a more complex molecular and clonal landscape with distinct clonal dynamics in emerging neoplasms.

All patients' treatment courses, clinical cytogenetic, and current molecular analyses are depicted in detail in Supplementary Fig. S5.

There was a notable difference in conventional cytogenetic abnormalities between the transformation and transient groups. More transformation patients (n = 7) than transient (n = 1) had high-risk clonal events (27), P < 0.05 (Fig. 4A), and partial gains of 1q were prevalent in the transient group (5/14) but not in the transformation cohort (0/17, P = 0.01). By definition, the neuroblastoma-treated control group did not have any documented cytogenetic abnormalities, nor did they have any abnormalities detected by our capture-based CNA analysis.

In contrast to the transformation cohort, patients in the transient and neuroblastoma-control groups had a paucity of acquired gene mutations compared with patients who developed t-MDS/AL (Fig. 1A and B). Only 4 of 14 patients from the transient cohort and 1 of 21 from the neuroblastoma-control cohort had somatic gene mutations at any time point evaluated. The mutations in the transient cohort included known oncogenic variants in TP53, KRAS, PPM1D, and DNMT3A (non-R882; Supplementary Table S7). The VAF for mutations detected in the transient cohort was five-fold lower than that of t-MDS/AL diagnostic samples (mean VAF: 0.042, range: 0.025–0.085 vs. 0.22, range 0.026–0.75; P < 1e-5) and 2.5-fold lower than those of premalignant samples in t-MDS/AL patients (mean VAF 0.10, range: 0.026–0.33, P < 1e-4; Supplementary Fig. S6). VAFs in all transient patients resolved over time or decreased to <0.05 at last evaluation with long follow-up (n = 2, nearly 3 years and >12 years). One control case had a mutation in BCOR; this patient died of neuroblastoma 7 years after initial diagnosis though only 18 months after the last evaluated BM sample (Supplementary Fig. S5). Although patients in the transformation cohort had a significantly higher number of oncogenic mutations (P < 0.001; Fig. 4B and C), they had received less chemotherapy and radiation than either the transient or neuroblastoma-control cohorts at the median time of transformation (P = 0.04 and 0.0076, respectively). Relative risk for the development of t-MDS/AL in the presence of an oncogenic mutation was 8.8 (95% confidence interval, 1.3–57.8, P < 0.001).

Taken together, our findings show that detection of CH in patients with neuroblastoma identifies those at higher risk of developing t-MDS/AL. To evaluate the frequency of potentially preleukemic events across pediatric oncology patients, we evaluated 657 pediatric patients (≤25 years-old and seen in the Department of Pediatrics at MSKCC) whose tumor was sequenced on our institutional platform, MSK-IMPACT, as described above (cohort characteristics in Fig. 5A and B and Supplementary Table S1). Median time to collection for all patients in this cohort was 155 days following primary malignancy diagnosis (Supplementary Fig. S7A). Twenty-six (3.9%) patients were identified as having CH, defined as at least one mutation with VAF ≥0.02, as previously described (refs. 30, 37; Fig. 5D and E). Median VAF for all mutations was 0.033 (range, 0.020–0.176). Of these, 7 patients (1.05% of the total cohort) each had one CH variant classified as having oncogenic potential in a gene implicated in myeloid disease (CH-PD), with a median VAF of 0.027 (range, 0.023–0.081; Fig. 5C). In contrast, the frequency of CH-PD was 17.6% in patients >25 years old (Supplementary Fig. S7B), further demonstrating the rarity of CH in pediatric oncology patients. Average age at time of collection and time to sample collection from diagnosis was not significantly different between patients with CH-PD compared with all patients (Mann-Whitney-Wilcoxon test). Of the 7 patients with CH-PD, 3 subsequently died from their primary malignancy and 4 are still alive a median of 3.7 years (range, 2.5–4.2) from CH analysis with no evidence of hematologic malignancy. Four additional patients in this cohort developed t-MDS/AL, and 1 patient an unrelated hematologic malignancy. None of these patients had evidence of CH-PD with VAF ≥0.02 at the time of their MSK-IMPACT analysis, at a median of 17.5 months (range, 6.1–33.0) prior to diagnosis of hematologic malignancy. Patient characteristics from both with CH and secondary hematologic malignancy from the MSK-IMPACT cohort are outlined in Supplementary Table S10.

Historically, t-MDS/AL was considered to be a direct consequence of cytotoxic damage induced by cancer therapy. However, in recent years, molecular studies have provided a cascade of data demonstrating an increased risk of t-MDS/AL in patients with pre-existing CH (21, 22, 38, 39) and poorer OS in cancer patients with CH (30, 40). Here we study a unique cohort of similarly treated patients with high-risk neuroblastoma longitudinally, as well as 657 pediatric cancer patients with samples at a single time point. We characterized disease-defining alterations at t-MDS/AL diagnosis in patients with neuroblastoma and evaluated the emergence of these mutations retrospectively in longitudinal samples.

We demonstrate that t-MDS/AL-associated mutations are detectable up to 3 years before transformation (range, 1.3–32.4 months), with a median detection time of approximately 1 year (11.8 months) before diagnosis. Importantly, capture-based DNA and RNA assays enable the detection of acquired gene mutations, copy-number alterations, and fusion events not detected by conventional cytogenetics. Although cytogenetic abnormalities can be transiently acquired, the presence of mutations in myeloid genes is strongly enriched in patients who subsequently developed t-MDS/AL. Acquisition of secondary oncogenic events or increase in burden of abnormalities often heralded the onset of overt hematologic disease. This is consistent with recent studies in healthy adult populations demonstrating that both mutation number and clone size are significant predictors of t-MDS/AL (19–22). We also observed clones both developing and expanding during therapy. It is possible that these clones escaped our level of detection in the earliest samples, even with the use of ddPCR in the case of TP53 mutations, and were thus present but undetectable prior to therapy. Recent work in adults with myeloproliferative neoplasms in fact suggests that JAK2 mutations responsible for these diseases may derive from prenatal clones in some cases (41, 42). However, the absence of observed mutations early in treatment among transformation patients in our cohort raises a possible difference in etiology of t-MDS/AL in pediatric and adult patients. Mutations in pediatric patients may be more likely to result from a consequence of genotoxic therapy, whereas such mutations in adults have been demonstrated to preexist before the start of therapy (15, 21, 22, 38, 39) and subsequently expand in the presence of therapeutic selective pressure (37).

The presence of hematopoietic clones with myeloid gene mutations as isolated CH outside the setting of hematologic malignancy is rare but does occur in children with advanced cancers receiving systemic therapy (1% of pediatric MSK-IMPACT cohort and 4/14 transient cohort). Clones in these patients usually involved single oncogenic mutations at low VAF, and these patients have thus far not progressed to t-MDS/AL. This confirms that although most CH clones represent the premalignant stages of t-MDS/AL, some may represent a state of indolent CH. In addition, 4 patients in this cohort did develop subsequent t-MDS/AL at a median of 18.9 months (and all but one at ≥15 months) following MSK-IMPACT evaluation. None of these patients had evidence of high-risk CH at the time of analysis. Therefore, distinct from adult populations, where CH likely pre-exists at the time of malignancy diagnosis and evaluation at a single time point may identify high-risk patients, optimal monitoring of pediatric patients should include longitudinal surveillance and close follow-up for patients undergoing high-risk treatment (Supplementary Fig. S8). Moreover, longer follow-up of pediatric CH patients may reveal specific genetic or clinical parameters that predict for subsequent hematologic malignancies. Future multicenter studies are warranted to provide important insights into leukemic predisposition in the pediatric context.

This study has some limitations. Given that our center has a large referral base, our primary study cohort may represent a group of higher-than-average risk neuroblastoma patients. However, this potential bias would have affected all cohorts of neuroblastoma patients alike, and any between-group comparisons should remain unaffected. As noted above, our transformation cohort actually received less therapy than the transient or neuroblastoma control cohorts at the average time to transformation (Fig. 4A). In addition, our sample size for evaluating longitudinal samples is relatively small, and we have limited power to establish a fully predictive model relating specific VAF or change in frequency to t-MDS/AL progression. Thus, further prospective work should focus on validating the findings in our current study and attempting to establish more specific parameters for high-risk patients. Last, our study primarily consisted of BM samples, which are not as readily available as peripheral blood (PB). However, it should be noted that BM and PB have been demonstrated to be comparable in multiple settings, including MDS (43), acute leukemia (44), and CH (45), and thus, we recommend the use of either source for evaluation of CH in high-risk patients.

Last, it remains to be seen how emerging sequencing technologies, such as error-corrected sequencing and high-sensitivity SNP arrays, which improve our ability to identify low-level variants, will alter our predictions regarding patients at high-risk for developing t-MDS/AL. Nonetheless, this study focuses on the utility of established clinical assays and VAFs within the range of current clinical reporting guidelines (≥0.02; refs. 30, 35, 36). Interestingly, among the 3 patients in the neuroblastoma transformation cohort with t-MDS/AL-associated mutations whose mutation was not identified at time points prior to t-MDS/AL diagnosis at the clinical level of detection used in this study (VAF ≥0.02), all 3 had evidence of these mutations at lower allele burdens present in prior samples (VAF 0.0012–0.006 by targeted panel sequencing and frequency 9.2e-4 to 0.012 by ddPCR). These findings should be interpreted with caution, but should be further investigated in future work.

In conclusion, our data demonstrate the clinical utility of routine surveillance of pediatric cancer patients and survivors to identify patients at risk of developing t-MDS/AL. We provide strong support for incorporating deep molecular profiling approaches to comprehensively capture genetic alterations that are implicated in t-MDS/AL pathogenesis, as well as integrating longitudinal sampling into surveillance protocols. Larger studies are warranted to validate clinical guidelines for patients with identifiable CH at time of diagnosis or during therapy. However, in our cohort, the median time to detection of a clonal event was >2 years after initiation of therapy. In addition, 3 of 4 patients in our larger pediatric oncology cohort who underwent profiling at single time points and developed subsequent t-MDS/AL, did so >1 year from analysis. This suggests that pediatric survivors without CH at first surveillance would also benefit from annual follow-up (Supplementary Fig. S8). Last, implementation of routine surveillance practices in clinical practice will also enable studies into mechanisms of leukemic transformation.

G. Gundem reports personal fees from Isabel Inc. outside the submitted work. J.S. Medina-Martínez reports other support from Isabel Inc. outside the submitted work. K.L. Bolton reports grants from Bristol Myers Squib outside the submitted work. M.E. Arcilla reports personal fees from Invivoscribe, Biocartis, AstraZeneca, Bristol Myers Squibb, Clinical Care Options, Janssen Global Services, LLC, Physicians' Education Resource, LLC, and Roche outside the submitted work. S. Modak reports personal fees from YMAbs Therapeutics and Illumina RP outside the submitted work. A.L. Kung reports other support from Isabel Technologies, personal fees and other support from Imago Biosciences and Emendo Biotherapeutics, and personal fees from DarwinHealth outside the submitted work. R.L. Levine reports being on the supervisory board of QIAGEN and is a scientific advisor to Imago, Mission Bio, Zentalis, Ajax, Auron, Prelude, C4 Therapeutics, and Isoplexis; receiving research support from and consulting for Celgene and Roche and has consulted for Incyte, Janssen, Astellas, Morphosys, and Novartis; and receiving honoraria from AstraZeneca, Roche, Lilly, and Amgen for invited lectures and from Gilead for grant reviews. S.A. Armstrong reports personal fees and nonfinancial support from Neomorph, Imago Biosciences, Cyteir Therapeutics, Accent Therapeutics, Mana Therapeutics, and C4 Therapeutics, as well as grants from Novartis, Syndax, and Janssen outside the submitted work. N.K.V. Cheung reports personal fees from YMAbs Therapeutics, Abpro Lab, Eureka Therapeutics, and Biotec Pharmacon, as well as other support from YMAbs Therapeutics outside the submitted work. E. Papaemmanuil reports other support from Isabel Inc. outside the submitted work. No disclosures were reported by the other authors.

The data generated in this study from the neuroblastoma cohort are available within the article and its Supplementary Data Files.

The data utilized in this study from the IMPACT cohort has been previously published. The patient subset reported in this manuscript is available on cbioportal at: http://www.cbioportal.org/study/summary?id=msk_ch_ped_2021.

B. Spitzer: Conceptualization, data curation, formal analysis, investigation, visualization, methodology, writing–original draft, writing–review and editing. K.D. Rutherford: Formal analysis, investigation, visualization, methodology, writing–original draft, writing–review and editing. G. Gundem: Data curation, formal analysis, methodology, project administration. E.M. McGovern: Data curation, project administration. N.E. Millard: Conceptualization, data curation. J.E. Arango Ossa: Software, formal analysis. I.Y. Cheung: Resources, project administration. T. Gao: Formal analysis, methodology. M.F. Levine: Software, investigation. Y. Zhang: Data curation, central cytogenetic review. J.S. Medina-Martínez: Software, methodology. Y. Feng: Resources. R.N. Ptashkin: Formal analysis. K.L. Bolton: Formal analysis. N. Farnoud: Investigation. Y. Zhou: Methodology. M.A. Patel: Project administration. G. Asimomitis: Investigation. C.C. Cobbs: Methodology. N. Mohibullah: Methodology. K.H. Huberman: Methodology. M.E. Arcilla: Methodology. B.H. Kushner: Resources, data curation, writing–review and editing. S. Modak: Resources, data curation, writing–review and editing. A.L. Kung: Resources, funding acquisition. A. Zehir: Formal analysis. R.L. Levine: Resources, supervision. S.A. Armstrong: Conceptualization, supervision, funding acquisition. N.K.V. Cheung: Conceptualization, resources, supervision, funding acquisition, writing–review and editing. E. Papaemmanuil: Conceptualization, resources, supervision, funding acquisition, methodology, writing–review and editing.

The authors acknowledge the use of the Integrated Genomics Operation Core, funded by the NCI Cancer Center Support Grant (CCSG, P30 CA08748), Cycle for Survival, and the Marie-Josée and Henry R. Kravis Center for Molecular Oncology.

Research funding for this work included awards from Geoffrey Beene Cancer Research Center (E. Papaemmanuil and N.K.V. Cheung), Gabrielle's Angel Foundation (E. Papaemmanuil), V Foundation (E. Papaemmanuil), and the Damon Runyon Cancer Research Foundation (Damon Runyon-Rachleff Innovation Award, E. Papaemmanuil).

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.