Genome-wide Interrogation of Germline Genetic Variation Associated With Treatment Response in Childhood Acute Lymphoblastic Leukemia FREE

Pediatric acute lymphoblastic leukemia (ALL) cure rates have increased from less than 10% in the 1960s to more than 80% today. Such advancement was partly derived from the identification of presenting clinical features (eg, molecular subtype, leukocyte count, age) associated with treatment outcome and subsequent implementation of risk-adapted therapy.^{1 - 2} The assessment of decreasing disease burden in response to therapy by sequential monitoring of minimal residual disease (MRD) status has now been integrated into risk stratification.^{3 - 5} Minimal residual disease assays provide a direct assessment of early treatment response and are associated with final treatment outcome.^{6 - 9}

Response to treatment varies during the 4- to 6-week phase of remission induction therapy, as exemplified by changes in early sequential MRD assays.^{4 ,8 - 9} Thus, some patients exhibit drastic depletion of leukemia (from 100% to less than 0.01% leukemia cells in bone marrow) within only 2 to 3 weeks of induction therapy, while others exhibit high levels of residual leukemia even after 4 to 6 weeks of therapy.

This interindividual variation in treatment response in cancer can arise from both tumor- and host-related factors; however, most prior studies focused on the former. Gene expression profiling of leukemic blasts has identified tumor genetic features associated with outcome^{10 - 11} and drug resistance.^{12 - 15} Much less is known about host genetic factors associated with cancer cure rates.^{16 - 19}

Taking a global approach to identify host genetic factors that may affect treatment response in ALL, we tested germline single-nucleotide polymorphisms (SNPs) for their association with MRD at the end of remission induction therapy in 2 independent cohorts of children treated for newly diagnosed ALL.

Two cohorts of patients were included (Table 1) with approval of their respective institutional review boards. Written informed consent from patients or their guardians (as appropriate) for genomic research was included as part of the treatment protocols at St Jude (Total Therapy protocols XIIIB and XV) and the treatment/biology protocol for the COG (Pediatric Oncology Group 9900 study).

There were common and unique elements to the eligibility and treatment of the St Jude and COG cohorts (eFigure 1)^{20 - 21} (http://www.acor.org/ped-onc/diseases/ALLtrials/9906.html). Common elements included daily prednisone, weekly vincristine, weekly daunorubicin, thrice-weekly asparaginase, and intrathecal methotrexate therapy. After 28 days of therapy, St Jude patients received cytarabine plus etoposide (Total Therapy protocol XIIIB) or cytarabine plus cyclophosphamide and 6-mercaptopurine (Total Therapy protocol XV). Minimal residual disease was studied in bone marrow at days 19 and 46 by flow cytometry, with the latter time point corresponding to the end of induction treatment.⁸ In contrast, COG patients finished the induction phase after 28 days of therapy, and MRD status was assessed at day 8 (blood) and at end of induction at day 28 (bone marrow).⁷ For St Jude patients, MRD status was categorized as negative (<0.01%), positive (≥0.01% but <1%), and high-positive (≥1%). In COG patients, MRD classification was nearly identical: negative (≤0.01%), positive (>0.01%, but ≤1%), and high-positive (>1%).

Diagnostic immunophenotype and molecular subtype analyses were performed as described.^{7 - 8}

DNA (500 ng) was digested with restriction enzymes, amplified, labeled, and hybridized to the Affymetrix GeneChip Human Mapping 100K and 500K sets (Affymetrix, Santa Clara, California). Together, these 2 arrays interrogated genotypes at 588 920 SNPs.

Race was designated in mutually exclusive categories of white, black, and other by comparing patient SNP genotypes with reference populations in the HapMap resource (Supplemental Methods).

SNP genotypes were coded according to the number of B alleles in the genotype call as determined using BRLMM,²² with the AA, AB, and BB genotype calls coded as 0, 1, and 2, respectively. For genotypes that were not called by the BRLMM algorithm, we imputed the number of B alleles based on signal intensity and consistency with expected genotypes based on linkage disequilibrium.^{23 - 25}

SNPs with a minor allele frequency of less than 1% or call rates less than 95% (ie, the number of samples with definitive genotype call at this SNP is <95% of the total number of samples typed for this SNP) were excluded (Figure 1); patient samples that failed to achieve 95% call rates (ie, samples for which <95% of interrogated SNPs were successfully typed) were excluded (Figure 1 and Supplemental Methods).

Minimal residual disease was treated as an ordinal variable; ie, 1 for negative, 2 for positive, and 3 for high-positive, as defined above. To minimize confounding effects, patients with ALL subtypes that were strongly related to MRD status and that differed in frequency between the 2 cohorts (ie, E2A-PBX1, MLL rearrangements, BCR-ABL ALL) were excluded from the MRD analyses (Table 1, Figure 1, and Supplemental Methods). The final analysis included 476 796 SNPs among 318 St Jude and 169 COG patients (Table 1 and Figure 1).

SNPs associated with end-of-induction MRD were identified based on a bidirectional validation approach in both the St Jude and COG cohorts, comprising a 3-step analysis. Our goal was to find SNP genotypes that were associated with MRD in both cohorts—those that might be generalizable across induction treatment for ALL.

In step 1, we computed the statistical significance for each SNP genotype's association with MRD in each cohort separately. The Spearman rank correlation was used as the test statistic to account for the ordinal nature of MRD and the gene dosage effect of genotypes. P values were computed by a permutation-asymptotic hybrid method (Supplemental Methods). An additive model was assumed, although the trend test is also reasonably robust to moderate deviation from additivity.²⁶

In step 2, we determined the threshold for statistical significance by estimation of the false discovery rate (FDR) and an internal validation (Supplemental Methods) in each cohort. Using the P values obtained in step 1, in each cohort, FDR levels were estimated on a grid of per-test significance levels (P value cutoffs).²⁷ Based on the FDR estimates and the internal validation, a threshold (P ≤ .0125) was chosen for each cohort to declare a set of SNPs for further investigation.

In step 3, we used the COG MRD cohort to validate the top-ranked SNPs (P ≤ .0125) discovered in the St Jude MRD cohort and vice versa (bidirectional validation) using a rank-based inference procedure (Supplemental Methods). The 102 overlapping SNPs satisfying the significance threshold determined in step 2 (FDR estimation and internal validation) and step 3 (bidirectional validation) were prioritized for further investigation and association with additional relevant phenotypes. The absolute risk difference²⁸ was estimated for the 102 overlapping SNPs by comparing the frequency of MRD positivity in patients with 1 or 2 copies of the risk allele vs those homozygous for the alternative allele.

As a secondary approach, we used the St Jude cohort as a discovery set and the COG cohort as a test set (Supplemental Methods).

Operating characteristics of the Spearman rank correlation test were determined via a simulation study (Supplemental Methods and eFigure 4). The genotypes associated with MRD were also assessed by a pooled analysis that combined evidence across the 2 cohorts to provide a combined P value for each SNP (Supplemental Methods). The FDR, the false-positive report probability,^{29 - 33} and the population-attributable fraction (Supplemental Methods) for prioritized SNPs were estimated.

Statistical and computational analyses were performed using S-Plus software, version 7.0 (Insightful Corp, Seattle, Washington), R software, version 26.1 (http://www.r-project.org), and SAS software, version 9.1 (SAS Institute Inc, Cary, North Carolina). Data analyses were performed between 2006 and 2008.

Antileukemic Response. The relationship between the 102 overlapping MRD SNP genotypes and 2 additional leukemia response phenotypes was analyzed to prioritize SNPs and to minimize the risk of false discoveries.

For purposes of this analysis, patients were also retrospectively categorized into super-responders, responders, and poor responders based on consideration of MRD status at 2 time points.¹⁶ Minimal residual disease status was dichotomized as negative or positive as defined above. Super-responders were MRD-negative at both early (day 8 in COG, day 19 in St Jude) and later (day 28 in COG, day 46 in St Jude) time points; responders were MRD-positive at the early time point but became MRD-negative at the later time point; and poor responders had positive status at the later time point. The association between SNP genotypes and this MRD responsiveness phenotype was assessed by rank correlation in all evaluable patients in separate analyses of St Jude (n = 304) and COG (n = 154).

The cumulative incidence of hematologic relapse (including isolated and combined hematologic plus extramedullary relapses) as a function of SNP genotypes in the combined St Jude and COG cohorts was analyzed by the Gray test. Isolated extramedullary relapses, lineage switch, second malignancy, and death in remission were treated as competing events. Excluding individuals with E2A-PBX1, MLL rearrangements, or BCR-ABL ALL, 416 St Jude and 180 COG patients were included in this analysis, overlapping with but not identical to the MRD cohorts as defined in Figure 1 and Table 1. Of these, 33 in St Jude and 35 in COG experienced hematologic relapse. St Jude patients were divided into 4 strata according to their treatment protocol and risk classification, and COG patients formed the fifth stratum. The cumulative incidence hazard regression model of Fine and Gray³⁴ was used to confirm the directional association with relapse for SNPs that achieved P < .10 in the Gray test. Follow-up was confirmed as of March 2006 for St Jude patients and November 2006 for COG patients.

Pharmacokinetic Studies. Three pharmacokinetic phenotypes were available from a subset of St Jude patients for antileukemic agents used during remission induction. Patients in these 3 data sets overlapped with but were not identical to those studied in the primary St Jude cohort for MRD.

The first data set included plasma clearance of etoposide on day 29 of therapy in 97 patients enrolled in St Jude Total Therapy protocol XIIIB.³⁵ Although etoposide was a component of induction therapy for only a subset of the St Jude MRD cohort and none of the COG cohort, its elimination is mediated via cytochrome P450 3A4 (CYP3A4)³⁶ and P-glycoprotein,³⁷ a common mechanism of elimination that also affects prednisone,^{35 ,38} vincristine,^{39 - 40} and anthracyclines.^{41 - 42}

The second data set included intravenous methotrexate plasma clearance at day 1 in 319 patients treated in St Jude Total Therapy protocols XIIIB²⁰ and XV.²¹ The third data set included intracellular methotrexate polyglutamate accumulation in ALL blasts at 44 hours after receiving up-front methotrexate in 230 patients treated in St Jude trials.^{43 - 44} Although only a subset of the St Jude MRD cohort and none of the COG MRD cohort received intravenous methotrexate, all patients in both cohorts received intrathecal methotrexate, which is known to distribute from cerebrospinal fluid to blood systemically.^{45 - 47}

The relationship between SNP genotypes and pharmacokinetic variables was analyzed using linear regression.

From St Jude Total Therapy protocols XIIIB (accrual between 1994-1998) and XV (2000-2006), 371 children with newly diagnosed ALL had available germline DNA and evaluable MRD status at the end of induction therapy. Of the patients with ALL enrolled in the COG study P9906 (accrual between 2000 and 2003), 227 children had germline DNA and evaluable end-of-induction MRD status. The actual number of patients included in specific analyses is described below. We found no significant differences in the characteristics (age, initial leukocyte count, race, and sex) of the patients enrolled in the St Jude and COG trials who were not included in this analysis vs those who were. Partial overlap exists in the patients reported herein and those reported previously^{16 - 17 ,48 - 49} but not in genomic testing or analysis.

A total of 588 920 SNPs were genotyped in the germline DNA of 371 St Jude and 227 COG patients. After quality control filters were applied (Supplemental Methods), 476 796 SNPs were evaluated in 318 St Jude and 169 COG patients (Figure 1 and Table 1). We analyzed the association between germline SNP genotypes and MRD status independently in the St Jude and COG cohorts. A P value threshold of .0125 was established based on FDR estimates and an internal validation inference (Supplemental Methods and eFigure 2). Through a rank-based bidirectional validation, a significant impact of germline variation on MRD identified in the St Jude cohort was validated in the COG cohort (P = 2.2 × 10⁻⁶) and that identified in the COG cohort was validated in the St Jude cohort (P < 10⁻¹¹) (Supplemental Methods).

In total, 102 SNPs exhibited significant concordant association with end-of-induction MRD (P ≤ .0125) in both the St Jude and COG cohorts, with a median odds ratio of 2.18 and a median population-attributable fraction of 0.17 (Table 2 for top 25 SNPs by combined cohort P value; eTables 1 and 2 for full details). Among these 102 SNPs, 50 were annotated to genes. Because 45 SNPs were clustered at 15 genomic loci by linkage disequilibrium (pairwise r²>0.5), these 102 SNPs represented 72 unique genomic loci (eFigure 3). A SNP in the ST8SIA6 (NM_001004470.1) gene (for combined cohort, odds ratio, 3.91; absolute risk difference: 0.46; P = 9.6 × 10⁻⁸) had the strongest association with MRD but had no significant flanking SNPs and a relatively low minor allele frequency of 4% (Figure 2 [chromosome 10]). The next highest-ranked SNP (rs17007695) was in the interleukin 15 (IL15) (NM_172174.2) locus (Figure 2 [chromosome 4]; Table 2) and was notable for strong (for combined cohort, odds ratio, 2.67; absolute risk difference, 0.28;P = 8.8 × 10⁻⁷) and comparable association with MRD in both the St Jude (P = 4.4 × 10⁻⁴) and COG cohorts (P = 2.3 × 10⁻⁴). Moreover, this SNP was flanked by 4 IL15 SNPs (rs17015014, rs10519612, rs10519613, and rs35964658; Figure 3) that were also associated with MRD in both cohorts (Figure 3, Table 2, and eTable 1), and these 5 SNPs were in linkage disequilibrium with each other (pairwise r², 0.48-0.97).

Several of the highly ranked SNP genotypes have relatively small numbers of patients in the least common genotypic groups (eTable 2). Three (50%) of the 6 St Jude patients with the CC genotype, 35.6% of those with the CT genotype (n = 45), and only 15.8% of patients with the TT genotype (n = 267) at the IL15 SNP rs17007695 had detectable MRD at the end of induction therapy, with a similar finding observed in the COG cohort (Figure 4). The C allele at the IL15 germline SNP rs17007695 was weakly associated (P = .0701) with a higher IL15 expression in leukemic blasts, and overexpression of IL15 was associated with MRD in both cohorts (P = .0342 in St Jude and P = .0035 in COG; eTable 5).

All 102 SNPs remained significantly associated with MRD after adjustment for race, sex, leukocyte count at diagnosis, age, and ALL subtype (eTable 1). To further explore whether the genotypes had similar associations with MRD in each racial group, we examined the frequency of MRD positivity for each SNP genotype in each major racial group. An example is provided by the SNP rs13106616, illustrating that the GG genotype was similarly associated with a lower risk of MRD across 3 racial groups, although the allele frequency differed significantly by race (eTable 6). We also assessed the false-positive report probability for these 102 SNPs, and 82 (80.4%) exhibited a false-positive report probability of less than 0.5 (Table 1), a level associated with replicated associations in other contexts.^{26 ,29 - 33}

In addition to the bidirectional validation described herein, we also present a genome-wide analysis for SNPs associated with end-of-induction MRD by following the “discovery and validation” approach. In the discovery stage, we computed the statistical significance for each SNP genotype's association with MRD in the discovery cohort (St Jude), estimating permutation-asymptotic hybrid P values for association with MRD, as detailed in the Supplemental Methods. A P value threshold of 7 × 10⁻⁴ was arrived at by balancing the levels of false-negative and false-positive errors using the profile information criterion (eFigure 5)²⁷ ; 624 SNPs met this threshold. In the second stage, these SNPs were then tested in the validation cohort (COG). Of these, 39 exhibited concordant associations at P ≤ .05, more than what would be expected by chance (P = .021 by Fisher exact test), and those that overlap with the bidirectional approach are shown in Table 2 and eTable 1. At a P value threshold of .0125 for the discovery cohort, 8635 SNPs met this cutoff, 330 of which were validated in the COG cohort with P ≤ .05, exceeding what would be expected by chance (P = 1.8 × 10⁻⁹).

Although end-of-induction therapy MRD is highly associated with long-term treatment outcome, an earlier reduction of leukemic burden is also informative.⁵⁰ Thus, nearly all patients with negative MRD at early time points (day 19 in St Jude and day 8 in COG) remained leukemia-free. We examined which of the 102 overlapping SNPs could also distinguish patients who responded early (super-responders; n = 145 in St Jude and n = 26 in COG) vs those who remained MRD-positive at the end of induction therapy (poor responders; n = 59 in St Jude and n = 52 in COG) vs individuals who were MRD-positive at the early time point but MRD-negative later (responders; n = 100 in St Jude and n = 76 in COG). Of the 102 overlapping SNPs, 40 (40%) were also associated (P < .05) with early response in both cohorts (eTable 3).

Of the 102 SNPs, 21 were significantly associated with hematologic relapse by stratified Gray test and in a cumulative incidence hazard regression model (P < .05; eTable 3). An example is shown for rs1486649 (an intergenic SNP); there was a monotonic relationship between the number of copies of the C allele and the risk of hematologic relapse (Figure 5).

To understand mechanisms by which host genetic variation might affect treatment response, we tested whether the 102 overlapping SNP genotypes were related to antileukemic drug disposition in a set of St Jude patients evaluable for pharmacokinetics (Table 3; eTable 3). In total, 21 of the 102 MRD-related SNPs exhibited significant association with antileukemic agent pharmacokinetics, with 3 SNPs associated with more than 1 pharmacokinetic phenotype. Eight of 102 SNPs were associated with methotrexate clearance (at P < .05); all 8 genotypes were associated with positive MRD and greater drug clearance. Ten of the 102 SNPs were associated with the etoposide pharmacokinetics, with 7 of 10 associated with positive MRD and greater drug clearance. Similarly, 6 of the 102 SNPs were significantly associated with the leukemic cell accumulation of methotrexate polyglutamates, with 5 of 6 associated with positive MRD and lower methotrexate polyglutamates. Thus, of 24 significant associations, 20 were directly consistent with a pharmacokinetically intuitive association with MRD; ie, lower drug exposure translated into a higher level of MRD.

Eradication of cancer cells by chemotherapy is a composite phenotype that depends not only on the somatically acquired characteristics of the malignant cells but also on inherent patient characteristics. Childhood ALL has long served as a prototype for a malignancy that is curable with drugs. Early MRD measures are strongly associated with cure rates and are used to modify therapy.^{3 - 9 ,51} Eradication of MRD is affected by genetic characteristics of the blasts (eg, the Philadelphia chromosome) and by host characteristics such as age.^{7 - 8} Using a candidate gene approach, a few germline genetic variations have been shown to affect the level of MRD,^{16 ,52} but this has not been previously assessed on a genome-wide level. Herein, we used an agnostic genome-wide interrogation to identify 102 germline genetic variations associated with the level of residual leukemia in 2 independent cohorts and found that a high proportion (63 of 102 SNPs [61.7%]) also affected early response, relapse risk, or antileukemic drug disposition.

One of the strongest signals from the genome-wide scan came from 5 SNPs located in the IL15 locus, a proliferation-stimulatory cytokine.^{53 - 54} Interleukin 15 can protect hematologic tumors from glucocorticoid-induced apoptosis in vitro,⁵⁵ and IL15 expression in ALL blasts has been linked to risk and relapse of CNS leukemia.⁵⁶ Both higher IL15 expression (P = .0342 in St Jude and P = .003 in COG) and germline SNP genotypes were associated with an increased risk of positive MRD (eTable 5), and we found a trend (P = .0701) toward a relationship between IL15 SNP genotypes and IL15 expression in ALL leukemic blasts. Several of the IL15 SNPs that were associated with MRD have been linked to enhanced IL15 transcription/translation efficiency in vitro.⁵⁷ Thus, it is plausible that germline genetic variation in IL15 plays a role in treatment response in childhood ALL via affecting IL15's function or quantity in ALL blasts, and the fact that IL15 SNPs were prominent from unbiased genome scans in 2 independently treated cohorts points to its importance in ALL response, either as a prognostic marker or as a therapeutic target.

Because genome-wide interrogations for pharmacogenetics are still in their infancy, there are no published whole-genome data linking polymorphisms with anticancer drug response. We had the opportunity to couple the findings from our genome-wide SNP interrogation for MRD with 3 relevant host pharmacokinetic phenotypes: systemic clearance of 2 antileukemic agents (etoposide and methotrexate) and intracellular disposition of the latter. Although 4 to 8 different antileukemic agents were used in these 2 cohorts, remarkably, 21 of the 102 MRD-associated SNPs were also significantly associated with disposition of these 2 antileukemic agents. Although many additional genetic variations would be expected to be specific for antileukemic drugs other than methotrexate and etoposide and might therefore account for some of the remaining 81 MRD-associated SNPs, several of the pathways involved in methotrexate disposition and etoposide disposition (http://www.pharmgkb.org) are likely to be shared by other antileukemic agents. Particularly for etoposide, whose disposition involves CYP3A4 metabolism and P-glycoprotein excretion, it is likely that there is overlap in the genetic determinants of its disposition with those affecting anthracyclines, glucocorticoids, and vincristine.^{35 ,38 - 42} The majority (83.3%) of the associations between SNP genotypes and drug disposition were pharmacologically intuitive, with genotypes that were associated with increased drug exposure linked to lower levels of MRD. Together, these results suggest that more attention should be given to details of drug administration and risk factors for rapid drug clearance, in addition to the attention already placed on more granular risk classification of ALL.

There was also a high proportion (21/102) of SNPs that were associated with not only MRD but also with the risk of hematologic relapse in both cohorts. This high percentage is somewhat surprising in that the postremission therapy (which would ultimately be expected to have a significant effect on relapse risk) differed substantially in the COG and St Jude cohorts. This secondary analysis does lend credence to the hypothesis that we did identify true associations between SNP genotypes and poor response.

Like all risk features, genotypes that are informative for pharmacogenetic phenotypes are likely to be highly dependent on therapy. For this reason, we purposefully chose 2 cohorts (St Jude and COG) that had received somewhat different remission induction regimens, with slightly different time points for the primary phenotype (MRD), to identify polymorphisms more likely to have prognostic significance across multiple therapeutic regimens. Compared with the traditional “discovery and validation” approach, this bidirectional approach minimizes bias against the discovery cohort (more stringent P value cutoff). The disadvantage of this approach is that we might have missed SNPs more specific to the few elements of therapy that differed between the cohorts, and may have a higher FDR.

It is important to consider race, both from the standpoint of its possible effects on antileukemic drug efficacy^{58 - 60} and from its influence on germline SNP allele frequency.⁶¹ We found good agreement between self-declared race and that determined using ancestry-informative SNPs, and the 102 MRD-associated SNPs remained significant after adjusting for ancestry (eTable 1). The fact that SNP genotypes maintained significance after adjusting for race, despite substantial differences in some allele frequencies by race, suggests that inherent differences in ALL prognosis among racial groups are partly influenced by differences in allele frequencies among racial groups, which could in the future lead to “race-neutral” (but genomically based) individualization of therapy.

We acknowledge that despite the fact that these SNP genotypes were associated with MRD in 2 independent cohorts, there is a danger of false-negative and false-positive findings, especially when sample size is relatively small. However, phenotypes of interest in pharmacogenetic studies (eg, CYP2C9/VKORC for warfarin^{62 - 63} and TPMT for thiopurine^{52 ,64} ) may have effect sizes that exceed those likely to be observed for multigenic common diseases,²⁴ and, therefore, smaller sample sizes may suffice in the former. By identifying 102 SNPs based on association in 2 independent cohorts and by further validation of 62% of these SNPs (eTable 3) to be associated with the related phenotypes of relapse, “super response” at days 8 or 19, and antileukemic drug pharmacokinetics, we have further decreased the chance of false discoveries. The SNPs we identified may be in linkage disequilibrium⁶⁵ with the truly causative genetic variants that have not yet been directly genotyped (eTable 4). Importantly, few of the 102 polymorphisms we identified have previously been suggested as candidates for affecting anticancer drug efficacy, and approximately half of the genomic variants are not annotated to genes at all, illustrating the need to further explore mechanisms by which germline genomic variation affects interindividual variability in antileukemic drug response.

Although the acquired genetic characteristics of tumor cells play a critical role in drug responsiveness, our results show that inherited genetic variation of the patient also affects effectiveness of anticancer therapy, and that genome-wide approaches can identify novel and yet plausible pharmacogenetic variation. Such variation may be factored into treatment decisions in the future by placing additional emphasis on optimizing drug delivery to overcome host genetic variation, in addition to the current emphasis on tumor genetic variation.