Impact of polymorphisms in drug pathway genes on disease-free survival in adults with acute myeloid leukemia.

Acute myeloid leukemia (AML) is a clinically heterogeneous disease, with a 5-year disease-free survival (DFS) ranging from under 10% to over 70% for distinct groups of patients. At our institution, cytarabine, etoposide and busulfan are used in first or second remission patients treated with a two-step approach to autologous stem cell transplantation (ASCT). In this study, we tested the hypothesis that polymorphisms in the pharmacokinetic and pharmacodynamic pathway genes of these drugs are associated with DFS in AML patients. A total of 1659 variants in 42 genes were analyzed for their association with DFS using a Cox-proportional hazards model. One hundred and fifty-four genetically European patients were used for the primary analysis. An intronic single nucleotide polymorphism (SNP) in ABCC3 (rs4148405) was associated with a significantly shorter DFS (hazard ratios (HR)=3.2, P=5.6 × 10(-6)) in our primary cohort. In addition, a SNP in the GSTM1-GSTM5 locus, rs3754446, was significantly associated with a shorter DFS in all patients (HR=1.8, P=0.001 for 154 European ancestry; HR=1.7, P=0.028 for 125 non-European patients). Thus, for the first time, genetic variants in drug pathway genes are shown to be associated with DFS in AML patients treated with chemotherapy-based autologous ASCT.

Introduction

Adult acute myeloid leukemia (AML) is a hematologic malignancy with widely heterogeneous clinical outcomes. New treatments for AML are increasingly being tested in clinical trials of patients with specific tumor cell mutations[1,2]. Although there have been substantial improvements in the number of patients who achieve complete remission, the choice of induction and post-remission therapy for adult AML is still based on the “one size fits all” principle. Most regimens incorporate antimetabolites (e.g. cytarabine, fludarabine), topoisomerase II inhibitors (e.g. etoposide, daunorubicin, idarubicin, mitoxantrone) and alkylating agents (e.g. busulfan, cyclophosphamide) for the treatment of AML. Prognostic factors for treatment response include age, prior exposure to chemotherapy, cytogenetic markers and expression profiles, and appearance of specific genetic mutations in tumor tissue, such as mutation and translocation of particular genes (e.g. FLT3, NPM1)[3,4]. However, these prognostic factors do not adequately capture the wide diversity of clinical outcomes in this disease. The percent of adults with AML who can survive 3 or more years and may be cured is approximately 5-70%[5,6].

One possible explanation for the difference in response to AML treatment is germline genetic variation. While the pharmacogenomics of AML drug response is an active area of research, there remain large challenges. These include the: (i) poor availability of uniform, well-collected and well-defined drug response phenotype information; (ii) lack of widely available germline DNA not contaminated with tumor cells (myeloblasts); (iii) limited availability of panels of genotype data in large patient cohorts; and (iv) inability to validate findings in replication studies. In this study, we overcome many of these challenges. In particular, to our knowledge, this is the first study that involved the analysis of large numbers of genetic polymorphisms in a cohort of AML patients treated with high-dose chemotherapy followed by autologous stem cell transplantation (ASCT). We present results from testing the association between germline variants in drug pathway and other genes with disease-free survival (DFS) in adult AML patients. We identified new associations between AML DFS and polymorphisms in several drug pathway genes for cytarabine, etoposide and busulfan, and also replicated SNPs previously reported to be associated with the AML response phenotypes.

Materials and Methods

Clinical Protocol, Study Criteria and Patient Cohorts

The population for the current study is AML patients who were enrolled in UCSF study protocols 9203 or 9303 between 1988 and 2010. The study protocols and patient selection criteria have been previously described[7-9]. This treatment protocol was used in patients with low and standard risk de novo AML, including APL (acute promyelocytic leukemia) in first or second complete remission. It was also used in a small number of patients with high risk AML (i.e. with secondary AML) if allogeneic SCT was not an option for the patient (e.g. unavailable donor)[7]. In step 1, patients were treated with consolidation chemotherapy including cytarabine 2000/mg/m[2] (i.v.) twice daily for 4 days concurrently with etoposide 40 mg/kg by i.v. infusion over the 4 days. During the recovery period from chemotherapy, peripheral blood stem cells were collected under granulocyte colony-stimulating factor stimulation. In step 2, patients underwent ASCT, which involved the preparative regimen of busulfan (total dose 16 mg/kg orally or 12.8 mg/kg intravenously, over 16 doses in 4 days) followed by etoposide 60 mg/kg (i.v. bolus) and reinfusion of blood or marrow stem cells. Patients had to be in complete remission for at least 30 days prior to step 2 (Figure 1). Complete remission was defined as normal bone marrow morphology with <5% blasts, resolution of previously abnormal cytogenetics and no evidence of extramedullary leukemia. In addition, patients must meet criteria for neutrophil and platelet counts, liver and kidney function[7-9]. Detailed procedures of patient enrollment, diagnosis, data collection and follow-up have been previously described[7-9]. Briefly, patients were actively followed up in the beginning within 6 months of diagnosis, with subsequent annual followed up by clinic visits. UCSF electronic medical records, the UCSF Blood and Bone Marrow Transplant Clinic database and patients' medical charts were abstracted to determine patients' remission status. The UCSF Committee on Human Research approved the research protocol (IRB no. 10-00649).

Figure 1

Schematic of workflow applied to determine the association of genetic variations in 42 candidate genes with disease free survival (DFS) in AML patients treated with a two-step treatment protocol prior to autologous stem cell transplantation. The workflow begins with phenotype data collection and genotyping of DNA samples using a genomewide SNP array. Candidate genes were selected based on their roles in the: drugs' pharmacokinetic/pharmacodynamics pathway; DNA mismatch repair mechanism; association with the drug cytotoxicity in lymphoblastoic cell lines previously identified in literatures; and pharmacogenomics studies of AML drug response. After association of each SNPs with disease-free survival in the 154 AML patients of European ancestry, the SNPs with p<0.01 and their tag-SNPs were examined for their potential functional roles using databases to identify eQTL SNPs and predicted regulatory elements such as binding sites of transcription factors and for their replications of previous studies. The SNPs with p<0.01 were also examined for their associations with DFS in 125 non-European ancestry. Finally, imputation was performed to identify other SNPs with stronger associations with DFS in European ancestry.

DNA Isolation and Genotyping

DNA was isolated from peripheral blood stem cells, which were collected during the recovery from step 1 consolidation chemotherapy. As noted in the above section, patients were in complete remission prior to consolidation chemotherapy and hence the samples utilized in this step contained less than 5% leukemic cells. DNA was isolated at the UCSF DNA Banking and Extraction Services Lab. The lab followed standard DNA extraction protocol described in the Wizard® Genomic DNA Purification Kit (Promega). The DNA was then quantified using Picogreen and normalized to 50 ng/μL. For each sample, we genotyped 250 ng of DNA. The Illumina HumanOmniExpress v1.0 Beadchip was used, following the manufacturer's protocols, at the Center for Genomic Medicine, RIKEN, Yokohama, Japan. For quality control of the genotyping, we included one HapMap trio and three duplicates of the DNA samples from the AML patients. A total of 328 distinct DNA samples from patients were genotyped, along with 3 duplicates and 3 HapMap samples (trio).

Patients' Ancestral Origin

The genetic ancestral origin of patients was determined using principal components analysis implemented in Eigenstrat[10]. Genotype information on our 328 AML patients was analyzed in conjunction with SNP data from the HapMap project, which consist of Europeans (CEU and TSI), Asians (JPT and CHB), Africans (TSI and ASW) and Mexican (MXL). From these analyses we were able to distinguish 154 patients of European ancestry based on their close clustering with the European HapMap samples.

Pathway, Gene and SNP Selection

A total of 42 genes were selected for analysis based on the following criteria: (i) genes in the pharmacokinetic and/or pharmacodynamics pathway of the drugs administered (cytarabine, etoposide and busulfan)[11-15]; (ii) genes described in literature as having significant associations with the drug cytotoxicity in lymphoblastoid cell lines[16,17]; (iii) genes involved in DNA mismatch repair[18]; and (iv) genes found previously to be associated with the AML response phenotype[19] (see Supplementary Table 1, Figure 1). After filtering the SNPs with low call rates (< 90%) and SNPs with low minor allele frequencies (MAF) (MAF < 1%) in the 154 European ancestry patients, we selected the SNPs in the candidate genes and within 25000 bp upstream and downstream flanking regions.

Statistical Analysis of the Associations

The primary analysis was to estimate the association between SNPs in the selected candidate genes and disease free survival (DFS) in 154 patients of European ancestry. A Cox proportional hazard model was used to estimate the hazard ratios (HR) and 95% confidence limits for the effect of genotype on DFS. An additive coding of genotypes was used in all analyses. The genetic effect estimates were adjusted for levels of a clinical risk score (see Table 1). The SNP associations with p<0.01 were also tested in the 125 non-European patients using a Cox proportional hazard model. In light of the heterogeneous ethnicity of the samples, the genetic effect estimates were adjusted for the first 10 principal components calculated from the GWAS data, as well as the clinical risk score. The Cox proportional hazards function from the R-project (version 2.15.1) was used. We used 3 x 10-5 (= 0.05 / 1659) as the significance level after Bonferroni correction for multiple testing.

Table 1

Demographic and characteristics of patients undergone autologous stem cells transplantation for acute myeloid leukemia from 1988 to 2010.

Variables	Patients of European ancestry (N=154)	Patients of non-European ancestry (N=125)
Sex, N (%)
Male	76 (49.3%)	67 (53.6%)
Female	78 (50.6%)	58 (46.4%)
Age (years)
Median (SD, range)
Age at Diagnosis (N=152)	47.0 (13.0, 18 - 72) (N=152)	40.5 (12.8, 17 - 68) (N=124)
Age at autologous bone marrow transplantation	47.0 (13.0, 19 - 72)	41.0 (12.8, 19 - 69)
Year of Transplantation
1988 to 1995	36 (23.4%)	38 (30.4%)
1996 to 2000	41 (26.6%)	31 (24.8%)
2001 to 2010	77 (50.0%)	56 (44.8%)
De novo/Secondary AML, N (%)
De novo	147 (95.5%)	110 (88.0%)
Secondary	7 (4.5%)	15 (12.0%)
Risk, N (%)
Acute promyeloid leukemia	15 (9.7%)	12 (9.6%)
Low	18 (11.7%)	18 (14.4%)
Standard	114 (74.0%)	80 (64.0%)
High	7 (4.5%)	15 (12.0%)
Disease Free Survival Time, mo
Median (SD)	21.4 (43.9)	16.7 (44.8)
Busulfan Area Under the Curve (AUC) (0-6 hr) (μmol*min)
Median (SD, range)	1164 (113, 973-1353) (N=12)	1173 (391, 777-1926) (N=9)

Fine Mapping of Associations via Imputation

In order to further clarify the association signals, we performed imputation on genes with p<0.005 in the primary analysis. For eight genes (ABCC3, DCK, GSTM1, GSTT1, MSH3, RRM1, SLC22A12, and SLC28A3) the genotypes at polymorphic sites known from the 1000 Genomes Project but not observed on the Illumina HumanOmniExpress v1.0 Beadchip were imputed using IMPUTE version 2 (version 2.3.0 for Mac OS X, http://http://mathgen.stats.ox.ac.uk/impute/impute_v2.html). The reference panel for the imputation was the 1000 Genomes Phase I integrated variant set referenced to NCBI build b37 (March 2012 relase, retrieved Jan. 20 2013 from http://http://mathgen.stats.ox.ac.uk/impute/impute_v2.html#reference). Imputed variants with imputation quality scores < 0.3 or MAF < 0.01 were excluded. The remaining imputed SNPs were each used in Cox proportional hazards models to predict disease free survival, as were the genotyped SNPs.

Functional Studies

The potential functional effects of the SNPs with p<0.01 associated with DFS in European population were examined using the following steps:

All tag-SNPs in linkage disequilibrium to the SNPs with p<0.01 in our primary analysis were identified using the Proxy Search in the Broad Institute SNAP (SNP Annotation and Proxy Search) (version 2.2), http://www.broadinstitute.org/mpg/snap/ldsearch.php. The search options used in this step were: SNP data set = 1000 Genomes Pilot 1 in CEU population panel, r2 threshold = 0.8 and distance limit = 500kbp.

Potential regulatory functions were identified by searching the following databases: Regulome Database[20] and eQTL browser (eqlt.uchicago.edu, http://eqtl.uchicago.edu/cgi-bin/gbrowse/eqtl/).

Literature searches were used to identify previous reports of associations with drug response phenotypes of any of the SNPs associated with AML response in this study.

Results

This study investigates the potential associations between variants in 42 candidate genes and disease-free survival (DFS) in adult acute myeloid leukemia (AML) patients treated with a 2-step approach to autologous stem cell transplantation (ASCT). High dose cytarabine, etoposide and busulfan were used in this treatment approach. Table 1 describes the demographic characteristics of this cohort, which consisted of 154 patients of European ancestry and 125 patients of non-European ancestry (African, Asian and Mexican). In the DFS analysis, 55 patients of European ancestry (35.7%) and 47 of non-European ancestry (37.6%) relapsed during the observation period, which extended from 1988 to 2010.

We applied quality control criteria to the SNP data, with a genotype call rate of 0.99 and minor allele frequency ≥ 0.01. After filtering, a total of 42 genes covered by 1659 SNPs were included in this association analysis. Furthermore, none of the DNA samples showed chromosomal abnormality by GenomeStudio (Illumina). Results in Figure 2 and Table 2 showed that among the 40 SNPs with p<0.01, the SNPs in the pharmacokinetic/ pharmacodynamic pathway genes have stronger associations compared to SNPs not in the drug pathway. After Bonferroni correction for multiple testing (1659 tests), a SNP in the first intron of ABCC3 (rs4148405) was significantly associated with DFS, with the minor allele (G) associated with shorter time to relapse (p-unadjusted=9.5 x 10-6, Figure 3a). Although other SNPs did not reach significance after Bonferroni-correction, overall there were 23 SNPs associated with DFS at p<0.005 (Table 2). These 23 SNPs are in or within 25000bp of 8 genes: SLC28A3, DCK, RRM1, GSTM1, ABCC3, MSH3, GSTT1, or SLC22A12. The majority of the minor alleles were associated with poor outcome (shorter DFS). Kaplan-Meier estimate plots of DFS are shown in Figure 3a-d for four of the top SNPs with MAF ≥ 3%. Other SNPs in the genes of the cytarabine pathway including NT5C2 and RRM2B were also associated with DFS but with weaker p-values (p<0.01, Table 2). Interestingly, SNPs in 3 out of 8 selected genes (SLC25A37, WNT5N and P2RX1), for which expression levels have previously been correlated with either etoposide or cytarabine IC50 values in lymphoblastoid cell lines (LCL), showed significant but weaker association (p<0.01) compared to genes in the drug pathways. Next, we examined the 40 SNPs in patients of non-European Ancestry. Only one SNP, in GSTM1-GSTM5 locus (rs3754446), was significantly associated with DFS in patients of non-European Ancestry (HR=1.7, p=0.028). Overall, in the entire cohort, the minor allele of the SNP (rs3754446) was significantly associated with shorter DFS (HR=1.7, p=0.00027).

Figure 2

Plot showing the significance (-log10 of the P-value) of associations of the 1659 SNPs with DFS in 154 AML patients of European ancestries. Only SNPs with minor allele frequencies of ≥1% in the selected 42 candidate genes are shown. Each dot represents a SNP. SNPs above the black dotted line are SNPs with p<0.01, and the SNP above the red dotted line reached a p-value,which was significant after correction for multiple testing (p<3x10-5).

Table 2

Significant SNPs (P<0.01) associated with disease free survival in 154 AML patients of European ancestry.

SNP	Chromosome	Gene	Feature	5′ Flanking Gene/3′ Flanking Gene	Risk allele	MAF	P (unadjusted)	HR	(95% CI) of the HR	Genotype Counts	Minor Allele	Major Allele	Reason for gene selection
rs4148405	17	ABCC3	intron	CACNA1G/ANKRD40	G	0.12	9.45E-06	3.12	(1.88-5.15)	2/34/117	G	T	Drug pathway
rs11140500	9	SLC28A3	intron	RMI1/NTRK2	T	0.01	0.00018	9.94	(2.99-33.09	0/3/151	T	C	Drug pathway
rs10805074	4	DCK	intron	MOBKL1A/LOC100128311	A	0.03	0.00039	4.59	(1.98-10.67	0/9/145	A	G	Drug pathway
rs7684954	4	DCK	intron	LOC100128311/LOC727995	A	0.03	0.00039	4.59	(1.98-10.67	0/9/145	A	G	Drug pathway
rs4593998	11	RRM1	intron	RRM1/LOC643244	A	0.14	0.00093	2.28	(1.40-3.71	4/34/116	A	G	Drug pathway
rs6842838	4	DCK	3′UTR of MOB1B	MOB1B/DCK	G	0.05	0.0010	3.21	(1.60-6.45)	0/15/139	G	T	Drug pathway
rs1385985	4	DCK	3′UTR of MOB1B	GRSF1/DCK	C	0.05	0.0010	3.21	(1.60-6.45)	0/15/139	C	T	Drug pathway
rs3754446	1	GSTM1	near-gene-5 [GSTM5]	GSTM1/GSTM5	G	0.38	0.0012	1.81	(1.26-2.59)	26/66/62	G	T	Drug pathway
rs7689093	4	DCK	near-gene-5 [DCK]	MOBKL1A/DCK	G	0.01	0.0016	5.34	1.89-15.13)	0/4/150	G	A	Drug pathway
rs1989983	17	ABCC3	near-gene-5 [ABCC3]	CACNA1G/ABCC3	A	0.11	0.0017	2.33	(1.37-3.96)	2/30/122	A	G	Drug pathway
rs2301835	17	ABCC3	synonymous variant in codi	CACNA1G/ABCC3	T	0.06	0.0029	2.60	(1.39-4.89)	1/16/137	T	C	Drug pathway
rs12515548	5	MSH3	intron	LOC100128458/RASGRF2	T	0.13	0.0029	2.07	(1.28-3.35)	6/27/121	T	C	DNA mismatch repair genes
rs2277624	17	ABCC3	synonymous variant	CACNA1G/ANKRD40	A	0.25	0.0040	1.75	(1.20-2.55)	12/54/87	A	G	Drug pathway
rs11090305	22	GSTT1	near-gene-5 [CABIN1]	GSTTP2/CABIN1	C	0.18	0.0044	1.98	(1.24-3.17)	3/50/101	C	T	Drug pathway
rs7130539	11	RRM1	intron	STIM1/OR55B1P	C	0.06	0.0047	2.56	(1.33-4.91)	1/16/137	C	T	Drug pathway
rs11031136	11	RRM1	intergenic	OR55B1P/LOC643244	G	0.06	0.0047	2.56	(1.33-4.91)	1/16/137	G	T	Drug pathway
rs528211	11	SLC22A12	intergenic	SLC22A11/SLC22A12	G	0.29	0.0048	0.51	(0.31-0.81)	13/62/79	A	G	Other AML study
rs2360872	11	SLC22A12	intergenic	SLC22A11/SLC22A12	C	0.29	0.0048	0.51	(0.31-0.81)	13/62/79	T	C	Other AML study
rs505802	11	SLC22A12	near-gene-5 [SLC22A12]	SLC22A11/SLC22A12	A	0.29	0.0048	0.51	(0.31-0.81)	13/62/79	G	A	Other AML study
rs524023	11	SLC22A12	near-gene-5 [SLC22A12]	SLC22A11/SLC22A12	G	0.29	0.0048	0.51	(0.31-0.81)	13/62/79	A	G	Other AML study
rs9734313	11	SLC22A12	5′UTR	SLC22A11/SLC22A12	T	0.29	0.0048	0.51	(0.31-0.81)	13/62/79	C	T	Other AML study
rs11231825	11	SLC22A12	synonymous variant	SLC22A11/SLC22A12	C	0.29	0.0048	0.51	(0.31-0.81)	13/62/79	T	C	Other AML study
rs2268166	11	RRM1	intron	STIM1/OR55B1P	G	0.06	0.0049	2.55	(0.20-0.75)	1/16/136	G	T	Drug pathway
rs11606370	11	SLC22A12	intron of NRXN2	SLC22A12/RASGRP2	A	0.30	0.0050	0.51	(0.32-0.82)	14/64/76	C	A	Other AML study
rs11191547	10	NT5C2	intergenic	CNNM2/NT5C2	T	0.31	0.0051	1.71	(1.17-2.49)	15/66/73	T	C	Drug pathway
rs11191549	10	NT5C2	near-gene-3 [NT5C2]	CNNM2/NT5C2	T	0.31	0.0051	1.71	(1.17-2.49)	15/66/73	T	C	Drug pathway
rs11191553	10	NT5C2	intron	CNNM2/LOC100128863	T	0.31	0.0051	1.71	(1.17-2.49)	15/66/73	T	G	Drug pathway
rs10883836	10	NT5C2	intron	LOC100128863/LOC729081	C	0.31	0.0051	1.71	(1.17-2.49)	15/66/73	C	T	Drug pathway
rs7095304	10	NT5C2	intergenic	NT5C2/LOC401648	A	0.31	0.0051	1.71	(1.17-2.49)	15/66/73	A	G	Drug pathway
rs6151816	5	MSH3	intron	LOC100128458/RASGRF2	T	0.12	0.0055	2.05	(1.24-3.41)	5/27/122	T	C	DNA mismatch repair genes
rs893006	11	SLC22A12	intron	SLC22A11/NRXN2	T	0.28	0.0055	0.51	(0.32-0.82)	13/61/80	G	T	Other AML study
rs7818607	8	SLC25A37	intergenic	SLC25A37/LOC646721	A	0.30	0.0057	1.75	(1.18-2.59)	16/59/79	A	C	Associated with drug cytotoxicy in LBL
rs2853229	8	RRM2B	intron	NCALD/UBR5	A	0.49	0.0062	1.68	(1.16-2.44)	38/74/42	A	C	Drug pathway
rs8534	8	SLC25A37	intergenic	SLC25A37/LOC646721	T	0.38	0.0067	1.68	(1.15-2.43)	24/69/61	T	C	Associated with drug cytotoxicy in LBL
rs8079740	17	ABCC3	intron of CACNA1G	SPATA20ANCC3	G	0.31	0.0078	0.56	(0.37-0.86)	20/56/78	A	G	Drug pathway
rs757420	17	ABCC3	intergenic	CACNA1G/ABCC3	T	0.30	0.0079	0.55	(0.36-0.86)	18/56/80	T	C	Drug pathway
rs2010851	12	WNT5B	near-gene-3 [WNT5B]	WNT5B/LOC100132548	A	0.30	0.0082	0.52	(0.32-0.84)	13/67/74	C	A	Associated with drug cytotoxicy in LBL
rs4995289	17	P2RX1	intergenic	P2RX1/ATP2A3	T	0.28	0.0088	0.51	(0.31-0.85)	14/59/81	C	T	Associated with drug cytotoxicy in LBL
rs1516801	17	P2RX1	intergenic	P2RX1/ATP2A3	G	0.28	0.0088	0.51	(0.31-0.85)	14/59/81	T	G	Associated with drug cytotoxicy in LBL
rs2607662	8	RRM2B	intron of UBR5	NCALD/UBR5	T	0.46	0.0095	1.64	(1.13-2.37)	31/79/44	T	C	Drug pathway

Note: The classifications near-gene-5 and near-gene-3, label SNPs that are outside transcribed regions, but within 2000 bp of a transcription region. Near-gene-5 includes upstream promoter region and untranslated 5′ mRNA.

LBL: Lymphoblastoid Cell Lines

Figure 3a – d

Kaplan-Meier estimate of disease-free survival (DFS) stratified by the top SNPs (a) rs4148405 ABCC3 (b) rs10805074 DCK (c) rs3754446 GSTM1 and (d) rs505802 SLC22A12 genotypes in patients of European ancestry.

Imputations of the eight candidate genes were performed to determine whether other SNPs in the regions have stronger association with DFS. The results showed that there are 234 imputed SNPs with p<0.01 (MAF ≥ 1%), and among these, there are 93 SNPs with improved p-values compare to the genotyped SNPs (Figure 4). Several SNPs in DCK have significant p-values <1.0 x 10-4 (Figure 4), with MAF 5-10%. Although none of these 93 SNPs are in exonic regions, on examination in the eQTL Browser, Regulome Database and GTEx eQTL Browser, we determined that SNPs in GSTM1 (rs929166, rs11101989) and MSH3 (6151896), are associated with their respective gene expression levels in liver[21] or lymphoblastoid cells[22] (data not shown). In the ABCC3 and SLC28A3 regions, imputation analysis did not identify other more significant SNPs in addition to the most significant genotyped SNPs, rs4148405 and rs11140500, respectively (Figure 4).

Figure 4

Plot showing the association of the genotyped and imputed SNPs with MAF ≥1% in the selected top 8 candidate genes (ABCC3, DCK, SLC28A3, SLC22AA12, MSH3, RRM1, GSTM1, GSTT1). Each colored diamond (not gray or black) represents a genotyped SNP and the grey/black dots represent imputed SNPs.

Using in silico analysis, we determined whether the 40 SNPs and their tag-SNPs were in known or predicted functional regions of the genome. Several of the SNPs were in DNA regions predicted to have binding sites for transcription factors. Some of these regions appear in the ENCODE Chip-Seq and DNase I peaks (Supplementary Table 2), suggesting that they could have regulatory functions. Interestingly, the GTEx (Genotype-Tissue Expression) eQTL Browser (http://www.ncbi.nlm.nih.gov/gtex/GTEX2/gtex.cgi) and eQTL Browser (http://eqtl.uchicago.edu/cgi-bin/gbrowse/eqtl/) showed that the SNP in GSTM1, rs3754446, is associated with GSTM1 and GSTM5 expression levels in the liver[21], brain and lymphoblastoid cell lines[22-24] (Supplementary Table 2).

Discussion

Previous pharmacogenomic studies of AML response to chemotherapy have been limited to a small number of candidate genes[2,25-27]. Though studies related to genes in drug pathways have been performed, many have not been replicated. In addition, some AML pharmacogenomics studies have been conducted using DNA from blast cells, which may have included somatic cell mutations in addition to germline polymorphisms[28-30].

To determine whether germline genetic variations are associated with AML response, we designed our own pharmacogenomic study in AML patients treated with a chemotherapeutic regimen consisting of cytarabine, etoposide and busulfan followed by autologous stem cell transplantation. Our study was focused primarily on 154 AML patients of European ancestry. The analysis was centered on 42 genes related to the pharmacokinetic and pharmacodynamic pathways of the chemotherapy. A few other genes that had previously been associated with drug response in AML were also included (Supplementary Table 1, Figure 1 and Figure 2). The goals of this study were to identify new associations with DFS in AML patients and to determine whether SNPs previously reported to be associated with AML response could be replicated.

Overall, the most significant SNP in our analysis was the intronic variant in ABCC3, rs4148405 (HR=3.1, p=9.5 x 10(-6)). ABCC3 is a multidrug resistance-associated protein (MRP3), which is known to transport the etoposide metabolite, etoposide glucuronide[31]. Several lines of evidence support a role for ABCC3 in DFS in AML patients. First in a previous study, a promoter variant in ABCC3, rs4793665 was associated with a shorter survival time in Israeli AML patients[29]. Though this variant was not associated with AML response in the current study, the data support a role of the transporter in response to chemotherapy in AML. Second, following etoposide administration to Abcc2 -/-, Abcc3 -/- mice, higher etoposide glucuronide levels were observed in the liver[32], consistent with a potentially important role of ABCC3 in etoposide pharmacokinetics. Third, cell lines transfected with ABCC3 show greater resistance to etoposide,[33] Finally, higher ABCC3 expression levels in leukemia cells are associated with poor outcome in children with leukemia[34,35].

Data from Regulome Database suggest that the ABCC3 variant, rs4148405, is in a functional location in the genome, as this region demonstrates direct evidence of binding through ChIP-seq studies. Transcription factors that play a role in hepatic gene regulation (e.g. CEBPB, USF1, FOXA1) have DNA response elements within this gene region (see http://regulome.stanford.edu/snp/chr17/48713567) [20,36]. Based on the results of our primary analysis with DFS, we speculate that the minor allele of rs4148405 is associated with higher expression levels of ABCC3 in the liver and/or leukemia cells, and thus reduced levels of etoposide in the tumor cells. Although the SNP rs4148405 is found at a MAF of greater than 10% in non-European populations, this SNP was not significantly associated with DFS in the AML patients of non-European ancestry. It is possible that different linkage disequilibrium patterns between rs4148405 and potential causative SNPs may have confounded the analysis.

In addition to a SNP in ABCC3, we also identified 23 SNPs in seven other genes (SLC28A3, DCK, RRM1, GSTM1, GSTT1, MSH3 and SLC22A12) that were associated with DFS in the European AML patients (with p< 0.005, Table 2). Interestingly, expression levels of these genes or other SNPs in these genes have been previously associated with response to chemotherapy in AML or other cancers[19,29,37-39]. Expression levels or SNPs in these genes have also been associated with IC50 values of various chemotherapy agents in cell lines[40,41]. Genetic polymorphisms in glutathione-S-transferases, such as GSTT1 and GSTM1, have been widely studied for their associations with drug toxicity, drug response and disease risk in leukemia patients[29,39,42]. Notably, GSTM1 and GSTT1 deletions have been implicated in various phenotypes associated with leukemia including drug response[43], busulfan pharmacokinetics[28] and disease risk[44]. Although we did not examine the effect of the GSTM1 deletion in this study, we identified several SNPs (genotyped and imputed) in the GSTM1-GSTM5 locus (rs3754446, rs929166 and rs11101989) associated with DFS in individuals of European ancestry that have not been previously reported (Table 2, Figure 2b). One of the SNPs in this locus, rs3754446, was also associated with AML response in the individuals of non-European ancestry (HR=1.7, p=0.028). Overall, in the entire cohort, this minor allele SNP, rs3754446, was significantly associated with shorter DFS (HR=1.7, p=0.00027). Perusal of eQTL databases (Regulome Database, GTEx Browser and eQTL Browser) suggests that the SNPs in the GSTM1-GSTM5 locus are associated with expression levels of GSTM1 and/or GSTM5 in liver[21], brain and lymphoblastoid cell lines[22-24]. Thus, the SNPs in these genes, which are involved in drug metabolism, could affect AML response by affecting the pharmacokinetics of the drugs used in the treatment of AML. We examined 21 AML patients from our overall cohort, where we have their first-dose busulfan AUC (area under the curve). Interestingly, we observed a significant association between rs3754446 and reduced busulfan AUC in these 21 AML patients (p = 0.03, Supplementary Figure 1). The minor allele, G, in rs3754446, was associated with lower busulfan plasma levels (AUCs), which was consistent with our observation that patients with the G allele had shorter DFS (Figure 3c). Because higher busulfan plasma levels have been associated with busulfan liver toxicity,[14,45] future studies are needed to determine whether the SNPs in GSTM1-GSTM5 are associated with liver toxicity. Collectively these data suggest that GSTM1 could play an important role in determining busulfan drug levels, drug response and/or drug toxicity.

Recently, using a drug-metabolizing enzyme/transporter genes (DMET) SNP array a synonymous variant in SLC22A12 (rs11231825) was found to be associated with response in 94 AML patients treated with a combination drug regimen of Gemtuzumab-Ozogamicin with Fludarabin-Cytarabine-Idarubicin[19]. Among the SNPs in the DMET genes that were found to be significantly associated with AML response, we were able to replicate the SNP (rs11231825) in the uric acid transporter, SLC22A12. Notably, this synonymous variant, which is in LD with a SNP upstream of SLC22A12, rs505802, has been found in various genomewide association studies to be associated with uric acid levels[46,47]. In these genomewide association studies [46,47], the minor allele T, has been associated with higher uric acid levels. Though speculative, our study and the previous study[19] which demonstrated that patients with the T allele have a better response to chemotherapy, suggest that higher uric acid levels may be associated with longer DFS time. In our study, the T allele, which is associated with higher uric acid levels[51,52], was associated with longer DFS time (see Figure 3d). Uric acid is a potent antioxidant, and it is possible that higher levels are beneficial for survival in AML patients.

Previous studies have shown that the nucleoside transporter, SLC28A3 (CNT3) plays a role in cytarabine cytotoxicity and resistance[40,48,49]. In this study, we observed several low allele frequency variants (MAF 1%) in SLC28A3 are associated with DFS in AML patients on cytarabine and other chemotherapy. Therefore, we hypothesized that CNT3 may transport cytarabine. Supplementary Figure 2a shows that radiolabeled cytarabine was taken up into CNT3 stably expressing cells and the uptake (over empty vector cells) was significantly enhanced in cell lines exposed to the equilibrative nucleoside transporter inhibitor S-(4-nitrobenzyl)-6-thioinosine (NBTI), which reduced background uptake of cytarabine in the cells. Notably, cytarabine uptake decreased significantly in CNT3 stable cells treated with the SLC28A3 inhibitor (Phloridzin) or with both inhibitors (Supplementary Figure 2a). Fludarabine, a known substrate of CNT3 was used as a positive control in this study (Supplementary Figure 2b)[50]. Because low allele frequency variant in CNT3 was found to be significantly associated with DFS (Table 2), we interpret the results with caution. Though our finding that SLC28A3 transports cytarabine supports the association, functional studies of variants in this region and/or a larger sample size are required to determine whether these uncommon variants are associated with cytarabine response.

Deoxycytidine kinase (DCK) plays an important role in activating cytarabine to its active metabolite, cytarabine triphosphate. Two tag-SNPs, rs4308342 and rs3775289, in DCK that were associated with DFS in our AML patients of European ancestry have been previously associated with the IC50 of another nucleoside analog, gemcitabine, in lymphoblastoid cell lines (LCL)[40]. Further, a previous study demonstrated that the level of cytarabine triphosphate in AML blast cells correlates with the ratio in expression levels of the cytarabine activating enzyme, DCK and the inactivating enzyme, 5′-nucleotidase, cytosolic II (NT5C2)[51]. In our study, in addition to the two tag SNPs in DCK associated with DFS in AML patients, five SNPs in the NT5C2 region were associated with DFS (Table 2). Though several of these SNPs are eQTLs (see Supplementary Table 2), further studies are required to determine whether these SNPs play important roles in determining the levels of cytarabine triphosphate in AML blast cells. Other genes, which play important roles in the cytarabine pharmacodynamics pathway are ribonucleotide reductase, RRM1 and RRM2, which are considered targets of nucleoside drugs such as cytarabine. The role of this enzyme is to regulate intracellular pools of ribonucleotides, such as deoxycytidine triphosphate (dCTP), which is important in building blocks for DNA replication. Studies have shown that AML blasts cells with high levels of dCTP are resistant to cytarabine, and that there is a significant correlation between RRM1 and RRM2 gene expression levels and dCTP levels after cytarabine treatment in AML blast cells[15]. Consistent with previous studies showing that SNPs in RRM1 are associated with response or toxicity to gemcitabine-based chemotherapy in lung and breast cancer patients[52,53], our findings suggest that the SNPs in RRM1 are associated with AML response to chemotherapy that include cytarabine.

Though our current association analysis supports the important roles of drug pathway genes, mainly transporters and enzymes, in AML response, we also selected 8 genes that have been associated with cytarabine or etoposide IC50 values in LCL[16,17]. A few SNPs in the 8 genes were significantly associated with response in AML patients suggesting that genes identified in in vitro assays in LCLs may also be important predictors of in vivo drug response in AML patients[16,17]

In summary, in this genetic association study of DFS in AML patients, we identified polymorphisms that have not been previously associated with AML response, including SNPs in ABCC3, DCK, GSTM1, MSH3, NT5C2, RRM1 and SLC28A3. A SNP in ABCC3, rs4148405, which remained significant after multiple testing, was associated with DFS, suggesting an important role of ABCC3 in determining etoposide levels in the liver and other tissues and hence AML response. Many of the significant SNPs or their tag-SNPs were eQTLs or located in functional regions in the genome. Finally, we determined for the first time that SLC28A3 (CNT3) transported cytarabine into cells suggesting an important role of this transporter in cytarabine cytotoxicity.