Genetic background influences susceptibility to chemotherapy-induced hematotoxicity.

Hematotoxicity is a life-threatening side effect of many chemotherapy regimens. Although clinical factors influence patient responses, genetic factors may also play an important role. We sought to identify genomic loci that influence chemotherapy-induced hematotoxicity by dosing Diversity Outbred mice with one of three chemotherapy drugs; doxorubicin, cyclophosphamide or docetaxel. We observed that each drug had a distinct effect on both the changes in blood cell subpopulations and the underlying genetic architecture of hematotoxicity. For doxorubicin, we mapped the change in cell counts before and after dosing and found that alleles of ATP-binding cassette B1B (Abcb1b) on chromosome 5 influence all cell populations. For cyclophosphamide and docetaxel, we found that each cell population was influenced by distinct loci, none of which overlapped between drugs. These results suggest that susceptibility to chemotherapy-induced hematotoxicity is influenced by different genes for different chemotherapy drugs.

Introduction

Forty percent of people in the United States will develop cancer during their lifetimes[1]. The majority of these patients will be treated with repeated intravenous doses of systemic chemotherapy drugs. Many of these drugs cause adverse side effects that can negatively impact treatment. For example, regimens that include doxorubicin (DOX), cyclophosphamide (CYC) and docetaxel (TAX) cause febrile neutropenia (i.e. neutropenia with fever) in approximately 30% of patients[2]. Severe neutropenia can lead to treatment delays, increased infection rates, hospitalizations, and mortality[3, 4].

The occurrence and severity of neutropenia is influenced by non-genetic factors, such as age, disease comorbidities, cancer stage and type[3], as well as by genetic factors[5]. Genetic polymorphisms and the underlying genes that influence chemotherapy-induced toxicity are often interrogated via genome-wide association studies (GWAS), in which many patients receiving the same chemotherapy regimen are genotyped and an association is sought between genotype and the severity of toxicity. A GWAS of chemotherapy-induced neutropenia and leukopenia in a Japanese cohort implicated SNPs on different chromosomes for different drugs, implying that susceptibility to hematotoxicity may be due to distinct genes for different drugs[6]. Results in model organisms support these findings as well. A GWAS of anthracycline-induced cytotoxicity in mouse splenic cells also identified different loci for DOX versus idarubicin[7]. Another study of chemotherapy-induced toxicity in Drosophila melanogaster also found that the toxicity of each drug was modulated by different loci for each drug[8]. This implies that that the toxicity of different drugs may be influenced by different genes and that each drug, or drug class, may have distinct pharmacokinetic and pharmacodynamic mechanisms.

In this study, we sought to identify genomic loci that influence chemotherapy-induced hematotoxicity by mapping genes in Diversity Outbred (DO) mice. DO mice are an outbred stock derived from eight inbred founder strains (A/J, C57BL6/J, 129S1/SvImJ, NOD/ShiLtJ, NZO/HlLtJ, CAST/EiJ, PWK/PhJ and WSB/EiJ)[9] that accumulate recombination events with each outbreeding generation, resulting in fine mapping resolution. DO mice segregate at over 40 million genetic polymorphisms and have balanced allele frequencies[10], resulting in good power to detect the effects of genetic variants. These properties make the DO an excellent mapping population for the discovery of genes that influence xenobiotic-induced toxicities[11]. While we did not expect that the genetic polymorphisms in mice would be identical in humans, we expected to identify genes that provide mechanistic insight into the etiology of chemotherapy-induced hematotoxicity and may thereby provide translational relevance.

Materials and Methods

Mice Housing and Feeding

We obtained male and female DO mice, six to ten weeks of age, from outbreeding generations 3 through 5 (The Jackson Laboratory, Bar Harbor, ME) and allowed them to acclimate for at least one week. Mice were housed in specific pathogen free conditions in polycarbonate cages with Bed-o’Cob® corn cob bedding (The Andersons Inc., Maumee, OH). Females were housed five mice per cage and males were singly housed. Mice were kept on a 12 hour light/dark cycle in rooms maintained at approximately 22°C ± 4°C and 50% humidity (± 15%). Mice were fed LabDiet 5LL4 (St. Louis, MO) ad libitum and were provided with water acidified to a pH between 2.8 and 3.1. Mice were weighed on the first, third and last day of treatment. All animal procedures were reviewed and approved by The Jackson Laboratory’s Institutional Animal Care and Use Committee (Protocols JW10001 & 13005).

Blood Collection

We collected approximately 100 μl of blood from the mice via the retro-orbital vein before dosing. We euthanized the mice by CO2 asphyxiation[12] and collected post dose blood via cardiocentesis. We stored the blood in K2EDTA tubes and analyzed it within 24 hours on a Hemavet 950 FS hematology analyzer (Drew Scientific Group, Dallas, TX)[13].

Dosing

We aimed to model the first-cycle response to chemotherapy by dosing mice and measuring blood cell counts at an acute time point near the nadir of neutrophil counts. We performed preliminary dose-response studies in the founder strains and determined the dose and euthanasia time point for each drug (data not shown). Mice were randomly assigned to a single dosing group in batches of approximately 50 animals (25 females and 25 males). DOX, CYC and TAX were administered via the tail vein (i.v). We measured body weights daily and monitored mice for signs of clinical distress. Mice that experienced a weight loss of 20% or more were euthanized and were not included in this analysis.

Doxorubicin

A total of 396 DO mice (198 females and 198 males) were dosed i.v. with 20 mg/kg of a 2 mg/ml solution of DOX (Pfizer, New York, NY)[14]. We collected post-dose blood five days after dosing. Three mice died during dosing and 14 were removed from mapping analysis due to failed genotyping. We analyzed data for 379 mice (191 females and 188 males) and had power to detect peaks that explain 9.5% of the phenotypic variance.

Doxorubicin and CSF3

A total of 195 DO mice (100 females and 95 males) were dosed i.v. with 20 mg/kg of a 2 mg/ml solution of DOX. On days three through six after dosing, mice were injected subcutaneously twice per day with 150 μg/kg of a 30 μg/ml solution of Neupogen (Amgen Inc., Thousand Oaks, CA). Post-dose blood was collected six days after dosing. We removed 16 mice from the mapping analysis that failed genotyping. We analyzed data for 179 mice (84 females and 95 males) and had power to detect peaks that explain 19% of the phenotypic variance.

Cyclophosphamide

A total of 200 DO mice (102 females and 98 males) were dosed i.v. with 200 mg/kg of a 20 mg/ml solution of CYC (Baxter Healthcare, Deerfield, IL)[15]. Post-dose blood was collected five days after dosing. Three mice died during dosing and six were removed from mapping analysis due to failed genotyping. We analyzed data for 191 mice (97 females and 94 males) and had power to detect peaks that explain 18% of the phenotypic variance.

Docetaxel

A total of 181 DO mice (107 females and 74 males) were dosed i.v. with a total of 30 mg/kg by administering a 10 mg/ml solution of TAX (Sanofi Aventis, Bridgewater, NJ) for three consecutive days. Post-dose blood was collected seven days after the first dose. Four mice died during dosing and we removed 23 mice from mapping analysis that failed genotyping. We analyzed data for 154 DO mice (85 females and 69 males) and had power to detect peaks that explain 20% of the phenotypic variance.

Statistical Analyses

For pre-dose hematology values, we log-transformed all values and fit a mixed-effects model with sex and drug as fixed effects and dosing date as a random effect. We calculated the geometric mean and standard error on a log scale and transformed these values back to the natural scale. We tested for differences in pre-dose parameters between drugs using a mixed-effects model with sex and drug as fixed effects and experimental batch as a random effect. We tested for a change between the pre- and post-dose parameter values by fitting a mixed-effects model with sex and treatment as fixed effects and experimental batch as a random effect. We performed a likelihood ratio test between the null model (i.e. excluding treatment) and the full model. P-values were adjusted using a Bonferroni correction. The investigator was not blinded regarding the status of the mice.

Genotyping

We collected tail tips after euthanasia from each animal and isolated DNA. Geneseek (Lincoln, NE) genotyped each mouse at 7,854 markers using the Mouse Universal Genotyping Array (MUGA)[16]. Samples with allele call rates below 90% were removed due to low quality.

Quantitative Trait Locus Mapping

We reconstructed the genome of each DO mouse in terms of founder haplotypes using DOQTL, a software package that performs haplotype reconstruction and genetic mapping[17]. We transformed cell counts into normal scores and performed linkage mapping by fitting an additive haplotype model with a kinship adjustment. We regressed the post-dose cell counts on the pre-dose cell counts, sex, dosing batch and the genotype at each MUGA marker. We imputed all single high quality (i.e. FILTER = PASS) single nucleotide polymorphisms (SNP) from the Sanger Mouse Genome Project[18, 19] onto DO genomes[17] and fit an additive genotype model at each SNP using the same covariates as the additive haplotype model, but with the imputed Sanger genotypes rather than founder haplotypes. We calculated genome-wide p-values (pGW) via permutation of the phenotype values and selection of quantiles from the empirical null distribution of maximum –log10(p-values) in each permutation[20, 21].

Knockout Mice

We obtained three knockout mouse strains, each targeting a different gene or pair of genes: Abcb1a (FVB.129P2-Abcb1a N7, Taconic, Hudson, NY), Abcb1a/b (FVB.129P2-Abcb1a N12, Taconic, Hudson, NY) and Abcb4 (FVB.129P2-Abcb4/J, The Jackson Laboratory, Bar Harbor, ME). Male homozygous knockout mice were mated to two FVB/NJ female mice (FVB, The Jackson Laboratory, Bar Harbor, ME) to produce heterozygotes knockouts. We mated a pair of heterozygous knockouts to produce mice with wild type, heterozygous and homozygous knockout alleles. We dosed each of the cohorts with 20 mg/kg of DOX as described above and measured neutrophil counts before dosing and five days after dosing.

Repeat of Abcb4 Knockout

We repeated the Abcb4 knockout test using a different experimental design at the Saarland University Medical Center in Homberg, Germany. We obtained FVB.129P2-Abcb4/J and bred them to obtain 72 progeny (36 females and 36 males). We housed mice in individually ventilated cages (3 per cage) under standard conditions (12 hour light/dark cycle) and mice received water and a standard rodent diet (V1534, ssniff, Germany) ad libitum. Mice were divided by genotype into three groups of 24 mice each; wild type, heterozygous and homozygous knockout for Abcb4. We dosed 12 of the mice (6 males, 6 females) in each group with 0.9% saline solution and dosed the other 12 mice with 20 mg/kg of DOX; both were administered via tail vein injection. We collected blood from all mice five days after dosing and measured neutrophil counts on a Sysmex XE-5000 automated hematology analyzer. All animal experiments were approved by the respective government agency (Landesamt für Verbraucherschutz, Saarbrücken, Saarland;TV43/2012).

Comparison with Human GWAS

We obtained data from the Biobank Japan human GWAS of chemotherapy-induced neutropenia and leukopenia[22] which contained results for DOX, CYC and TAX. This study contained 13 122 patients treated with one or more drugs. There were 234, 758 and 523 patients dosed with DOX, CYC and TAX, respectively. We compared our candidate gene list for each drug by intersecting their gene symbols with the genes in our QTL support intervals.

Genome, Gene and Variant annotation

We used mouse genome build GRCm38 to report genome coordinates[23]. We obtained gene locations from Mouse Genome Database, (version MGI.20160711.gff3.gz[24]) and mouse genome variants from version 5 (REL-1505) of the Sanger Mouse Genomes Project[18, 19].

Data Availability

The phenotype and genotype data are available at ftp://ftp.jax.org/dgatti/TPJ201600234/.

Code Availability

The DOQTL software is freely available as an R package from the Bioconductor project at http://bioconductor.org/packages/release/bioc/html/DOQTL.html.

Results

We dosed DO mice with one of three chemotherapy drugs that induce neutropenia in humans. Each drug has a different mode of action against tumor cells; DOX is a topoisomerase II inhibitor[25], CYC is a DNA alkylating agent[26], and TAX is a microtubule stabilizer[27, 28]. DOX is used to treat many tumor types, including breast, childhood solid tumors and lymphomas[29]. CYC is used to treat leukemias, lymphomas and solid tumors[30], and as a myeloablative agent[31]. TAX is used to treat breast cancer in both the adjuvant and metastatic setting[32]. These drugs are rarely administered alone and are commonly administered in combination with other chemotherapeutic agents. For example, the AC-T regimen in which DOX (A, Adriamycin®) and CYC (C) are followed by TAX (T, Taxotere®)[33]. Our goal in this study was not to mimic the exact details of clinical regimens. Rather, we aimed to determine whether genetic variants affect variation in hematotoxicity, whether genetic variants affect all hematopoietic cell sub-populations identically and to identify genomic loci that influence variation in hematotoxicity. We chose to model first-cycle responses because the first-cycle response in patients is often used to guide future dose modifications[34]. In addition to dosing three cohorts with either DOX, CYC or TAX, we dosed a cohort of mice with DOX followed by granulocyte colony stimulating factor (CSF3), which is commonly co-administered with chemotherapeutic agents to attenuate myelosuppression[35]. For each of these four dosing groups, we identified genomic loci associated with pre- and post-treatment blood cell counts.

Analysis of Pre-dose Parameters

We measured complete blood counts and body weight in each mouse before and after dosing with one of four dosing protocols: DOX, DOX with granulocyte colony stimulating factor (DOX+CSF3), CYC or TAX. We quantified the differences in mean pre-dose cell counts to determine if there were any baseline differences between drug groups (Table S1, Figure 1). We found no differences between pre-dose dosing groups for absolute neutrophil counts (NEUT, pBonf = 1.0), lymphocyte counts (LYMPH, pBonf= 0.44), red blood cell counts (RBC, pBonf = 0.22) and body weight (BW, pBonf = 1.0). We found that pre-dose values differed between drug groups for monocyte counts (MONO, pBonf = 0.041) and platelets (PLT, pBonf = 2.71×10−4). We categorized mice as neutropenic if NEUT fell below 120 cells/μl, which was the first percentile of the pre-dose NEUT distribution (i.e. 99% of pre-dose mice had NEUT > 120 cells/μl). This is similar to the neutropenia threshold of 110 cells/μl reported in ICR outbred mice[36]. We defined anemia as RBC < 7.88×106 cells/μl using the same 1% threshold in DO mice before dosing.

Figure 1

Changes in cell counts and body weight before and after dosing. The plots are organized with drugs in rows ((A) doxorubicin, (B) doxorubicin with CSF3, (C) cyclophosphamide and (D) docetaxel) and phenotypes in columns. Each panel plots the pre-dose (x-axis) value versus post-dose value (y-axis) for each mouse. Females are red and males are blue. Dashed lines plot the pre- and post-dose means for each sex. The black solid line along the diagonal shows where pre- and post-dose values would be equal. Points below the diagonal indicate a decrease in cell counts after dosing.

Effect of Chemotherapy Drugs on Blood Cell Counts

DOX had a general suppressive effect on the hematopoietic system, reducing mean LYMPH, NEUT, RBC and PLT (Figure 1A, Table S1). This effect is consistent with the neutropenic effects of DOX in the clinic[37]. The variance of post-dose BW was unchanged (pF > 0.05), but DOX increased the variance of all post-dose cell populations (pF < 10−10). There was no difference in neutropenia between sexes. The Pearson correlation between pre- and post-dose NEUT was 0.241, indicating that pre-dose values are not predictive of DOX response. DOX reduced mean BW by 11.7% in both sexes.

The addition of CSF3 to the DOX regimen increased mean NEUT by 54.6% (Figure 1B, Table S1), consistent with clinical experience[38]. The variance of NEUT increased after dosing with CSF3 (p = 1.80×10−17), indicating that genetic variation affects the response to CSF3 as well. The variance of BW and RBC was unchanged (pF > 0.05), but DOX increased the variance of post-dose LYMPH, MONO, NEUT and PLT (pF < 10−4). The Pearson correlation between pre- and post-dose NEUT was −0.039, again indicating that there was little predictive value in pre-dose counts. Mean LYMPH, MONO, PLT and RBC also decreased. Mean BW decreased by 11.9%.

CYC had a broad suppressive effect on the hematopoietic system (Figure 1C, Table S1), consistent with CYC’s role as a myeloablative agent[31]. Mean counts for all blood cell sub-populations decreased. The variance of BW and RBC was unchanged (pF > 0.05), but CYC increased the variance of post-dose LYMPH, MONO, NEUT and PLT (pF < 10−5). The Pearson correlation between pre- and post-dose NEUT was 0.125, again indicating little predictive value in pre-dose NEUT. Mean BW decreased by 3.9%. TAX produced a lesser myelosuppressive effect in DO mice when compared to the other drugs (Figure 1D, Table S1). There was a modest decrease mean NEUT and 2.8% of the mice were neutropenic. TAX suppressed mean RBC, causing 49.7% of the mice to become anemic. In contrast to the other drugs in this study, mean PLT increased in treated mice, which is consistent with limited clinical observations[39] and may be significant due to the relationship between PLT and metastasis[40]. The variance of BW and RBC was unchanged (pF > 0.05), but TAX increased the variance of post-dose LYMPH, MONO, NEUT and PLT (pF < 10−5). TAX produced a 9.2% decrease in mean BW.

Genetic Architecture of Chemotherapy-Induced Hematotoxicity

We performed genome-wide association (GWA) mapping on the change in cell counts to identify genomic loci that influence susceptibility to chemotherapy-induced hematotoxicity. We found that the hematotoxicity of each drug is influenced by a distinct set of genomic loci that have varying effects on the hematopoietic system. DOX toxicity is influenced by a pleiotropic locus on proximal chromosome 5 that affects cell counts for all of the blood cell populations as well as BW (Figure 2A). In the DOX+CSF3 cohort, the same locus on chromosome 5 had a major effect on white blood cell (WBC) populations and BW, but not on PLT or RBC (Figure 2B). For CYC, we found no pleiotropic loci; NEUT associated with a locus on chromosome 11, and PLT associated with two loci on the distal arm of chromosomes 5 and 19 (Figure 2C). Mice dosed with TAX showed peaks with less significant loci overall, which may be due to the modest change in mean NEUT at the dose we selected. There were loci that affect LYMPH on chromosome 3, a locus on chromosome 11 that affects RBC and loci on chromosomes 4 and 12 that affect NEUT (Figure 2D).

Figure 2

Heatmap of genetic mapping peaks shows that the hematotoxicity of each drug has a distinct genetic architecture. Each panel shows –log10(p-value) from genetic mapping of the change in cell counts or body weight for (A) doxorubicin, (B) doxorubicin + CSF3, (C) cyclophosphamide and (D) docetaxel. Brighter colors indicate greater significance. Values above –log10(p-value) = 6.2, which is equal to the pGW of 0.05, were set equal to 6.2 to reduce the effects of highly significant peaks.

We found a total of 21 loci met genome-wide significance (Table 2). We focused on genomic loci that affect NEUT because neutropenia is the most common dose-limiting toxicity of these drugs. NEUT in mice dosed with DOX or DOX+CSF3 showed a peak on proximal chromosome 5 (Figure 3A & B). We mapped the change in NEUT for mice dosed with CYC to a locus on chromosome 11 (Figure 3C) and for mice dosed with TAX to loci on chromosomes 4 and 12 (Figure 3D). We mapped NEUT using the pre-dose values from all 950 mice and did not find any peaks that overlapped with theses peaks, indicating that the peaks we found for each drug are related to treatment and not to constitutive NEUT (Figure 3E). The large peak on chromosome 1 contains chemokine receptor 4 (Cxcr4), which contains a missense SNP (rs8256191) and is involved in neutrophil trafficking from the bone marrow[9].

Table 2

List of significant QTL peaks. For each peak the drug regimen, cell type, chromosome, proximal and distal peak location, nominal p-value and adjusted p-values are listed.

Drug	Phenotype	Chromosome	Proximal Position	Distal Position	p-value	Adjusted p-value
Doxorubicin+GCSF	NEUT	2	148150899	148150944	6.00E-07	0.049
Doxorubicin+GCSF	MONO	2	148150899	148150944	8.70E-06	0.469
Docetaxel	LYMPH	3	88377902	88587928	3.49E-08	0.001
Docetaxel	MONO	3	90106020	90106020	4.27E-06	0.219
Docetaxel	NEUT	3	90106020	90106020	2.10E-05	0.598
Docetaxel	NEUT	4	142873065	142988627	7.27E-08	0.004
Doxorubicin+GCSF	EOS	5	5809501	6132452	4.53E-12	0.000
Doxorubicin	NEUT	5	6289554	6509240	1.51E-33	0.000
Doxorubicin+GCSF	LYMPH	5	6513094	7276358	5.06E-15	0.000
Doxorubicin+GCSF	MONO	5	6513094	7276358	2.80E-23	0.000
Doxorubicin+GCSF	NEUT	5	6513094	7276358	2.78E-17	0.000
Doxorubicin	LYMPH	5	8122454	8798106	1.73E-05	0.555
Doxorubicin	PLT	5	8125162	8797053	6.15E-10	0.000
Doxorubicin+GCSF	BW	5	8798831	9265645	4.88E-07	0.040
Doxorubicin	BW	5	9187736	9234163	7.17E-06	0.294
Doxorubicin	MONO	5	9285547	9302450	2.17E-18	0.000
Doxorubicin	RBC	5	10685282	10685282	2.25E-06	0.123
Cyclophosphamide	PLT	5	127467289	127475691	3.51E-07	0.029
Doxorubicin	RBC	6	27582992	27583007	6.98E-09	0.000
Doxorubicin+GCSF	MONO	6	92622580	92858576	2.19E-07	0.021
Doxorubicin+GCSF	EOS	6	92622580	92858576	4.09E-06	0.291
Doxorubicin+GCSF	MONO	7	143262081	143262081	1.16E-06	0.093
Docetaxel	RBC	11	44576416	44845071	1.19E-07	0.005
Cyclophosphamide	NEUT	11	69579943	69686580	2.83E-10	0.000
Docetaxel	NEUT	12	75408009	75837799	1.01E-06	0.062
Docetaxel	LYMPH	12	80024381	80024381	8.40E-07	0.054
Doxorubicin	NEUT	13	24678372	24872531	6.19E-07	0.041
Doxorubicin+GCSF	BW	13	62961059	63347403	5.26E-06	0.340
Doxorubicin+GCSF	LYMPH	13	63942114	63994071	5.13E-07	0.045
Doxorubicin+GCSF	LYMPH	16	39615002	39615002	2.81E-07	0.028
Doxorubicin	BW	17	84505914	84505914	2.79E-07	0.025
Cyclophosphamide	PLT	19	18875754	19340012	1.07E-07	0.006

Figure 3

Genome-wide association mapping plots of change in NEUT for (A) doxorubicin, (B) doxorubicin +CSF3, (C) cyclophosphamide, (D) docetaxel and (E) pre-dose NEUT. Each point plots the –log10(p-value) for the association between one SNP and the change in NEUT. Red lines show the pGW= 0.05 significance threshold. Orange triangles plot the locations of genes reported in the literature as being involved in the metabolism or transport of each drug[79–81]. The peaks on distal chromosome 4 in Figures 3D and 3E are 20 Mb apart and do not overlap.

Genetic Influences of Doxorubicin on Neutrophil Counts and Body Weight

We dosed 379 DO mice with DOX and mapped the change in NEUT before and after dosing. The additive heritability was 0.52 (p = 1.94 ×10−4), indicating that genetic background has a strong influence on the response to DOX. We mapped pre-dose NEUT and found that the highest –log10(p-value) in the interval from 3 to 30 Mb on chromosome 5 was 1.95, indicating that the chromosome 5 peak is not related to constitutive neutrophil levels (Figure 3E). We performed GWA mapping of the change in NEUT and found a peak on chromosome 5 at with an uncorrected p-value of 1.51 ×10−33 that explains 36.3% of the phenotypic variance (Figure 3A). This locus had pleiotropic effects on the hematopoietic system, influencing WBC, NEUT, MONO, PLT, RBC and body weight. We reconstructed the founder haplotypes of each DO mouse and estimated the founder allele effects on chromosome 5. DO mice carrying the NOD/ShiLtJ (NOD) or PWK/PhJ (PWK) alleles at the QTL on proximal chromosome 5 have a lower decrease in NEUT than mice carrying the other founder alleles (Figure 4A), while mice carrying the CAST/EiJ (CAST) allele have greatest decrease in NEUT. DO mice carrying one copy of the NOD allele had a 1-fold decrease in NEUT (i.e. NEUT decreased by half, Figure 4B) while mice carrying a CAST allele had a 2.7-fold decrease in NEUT. The NOD and PWK founder both carry a region of the M.m.musculus genome from 3 to 9.3 Mb[41] and we would expect that both alleles should have the same effect, but DO mice carrying the PWK allele had a 0.36-fold greater decrease in NEUT than mice carrying the NOD allele.

Figure 4

Association mapping of DOX-induced change in NEUT on chromosome 5 reveals two distinct loci with opposing effects. (A) X-axis shows the location in Mb on chromosome 5. Each colored line shows the effects for each founder allele as the log2-ratio of change in NEUT (see Methods). Founder allele effects on chromosome 5 for change in NEUT show that DO mice carrying the NOD (blue) or PWK (red) alleles have higher post-dose NEUT and mice carrying the CAST (green) allele have lower post-dose NEUT. (B) Estimated log2-ratio of the change in NEUT for each founder allele at 6.29 Mb on chromosome 5. Each bar and whisker shows the mean log2-ratio NEUT +/− standard error. Adding one NOD allele increases the log2-ratio NEUT by 1.3-fold and adding one CAST allele decreases log2-ratio NEUT by 0.61-fold. (C) Association mapping of the change in NEUT on proximal chromosome 5 reveals a set of significant SNPs for which NOD and PWK contribute the minor allele that is associated with higher post-dose NEUT. Each dot shows the –log10(p-value) for the association between change in NEUT and one SNP. SNPs above the dashed black line are highlighted in black and their minor allele is plotted in the top panel. The dotted line is the p = 0.05 significance threshold. The arrow denotes the most significant SNPs at which the effects are shown in panel D. (D) Log2-ratio of the change in NEUT versus the genotype at the most significant SNP in panel C (6.51 Mb on chromosome 5). Each dot represents one mouse. Red line show the least-squares fit and shaded region shows the 95% confidence interval. Dashed lines show the mean in each genotype group. Adding one NOD or PWK allele increases the log2-ratio NEUT by 1.3, which is the same as the increase in Figure 4B. (E) Association mapping of the change in NEUT after regressing out the first peak on proximal chromosome 5 reveals a set of significant SNPs for which CAST and PWK carry the minor allele that is associated with lower post-dose NEUT. Each dot shows the –log10(p-value) for the association between change in NEUT and one SNP. SNPs above the dashed black line are highlighted in black and their minor allele is plotted in the top panel. The dotted line is the p = 0.05 significance threshold. The arrow denotes the most significant SNPs at which the effects are shown in panel F. (F) Log2-ratio of the change in NEUT versus the genotype at the most significant SNP in panel E. (7.90 Mb on chromosome 5). Each dot represents one mouse. Red line show the least-squares fit and shaded region shows the 95% confidence interval. Adding one CAST or PWK allele decreases the log2-ratio NEUT by 0.44-fold, which is close to the values of 0.36 and 0.61 in Figure 4B.

We focused on the 15 266 SNPs with the most significant associations having –log10(p-values) ≥ 24 (Figure 4C). For all of these SNPs, the NOD and PWK strains contribute a resistant allele that is associated with a 1.3-fold increase in NEUT (Figure 4D), which is consistent with the increase of 1.3 that results from adding a NOD allele in Figure 4B. We added the genotype of the most significant SNP to the mapping model and rescanned the genome to determine if there were any other significant peaks. We identified another peak above the genome-wide significance level in the same location on chromosome 5. The most significant SNPs had a minor allele contributed by CAST and PWK (Figure 4E) and these alleles were associated with a NEUT decrease of 0.44-fold (Figure 4F). This combination of a NOD and PWK allele that increases NEUT, and a CAST and PWK allele that decreases NEUT explains the pattern of allele effects observed in Figures 4A and B. NOD and PWK contribute a resistant allele that increases the log2-ratio NEUT. CAST and PWK contribute a susceptible allele that decreases log2-ratio NEUT, lowering the PWK effects by 0.36 below NOD and lowering the CAST effect by 0.61 below the other strains (Figure 4B).

There are 178 annotated features in the region covered by the NOD/PWK allele SNPs in Figure 4C; 52 are protein coding genes, 56 are pseudogenes, 48 are non-coding RNAs and 22 are unclassified. We selected SNPs and indels that produce a coding or a splice site change in proteins and found 269 variants that intersect with 42 genes (Table S2). We accumulated annotation data about the genes under the NEUT association peak to select the most plausible candidates by combining evidence from missense and splice site variants, expression and protein quantitative trait loci (eQTL, pQTL) from separate liver data in DO mice[42] and genes annotated in the Comparative Toxicogenomics Database (CTD)[43]. Sorcin (Sri) contains one missense SNP (rs248635887) and is involved in calcium homeostasis. Sri is upregulated in DOX-resistant cell lines and Sri down-regulation in vitro induces apoptosis and renders cells less resistant to DOX[44-46]. The chromosomal region also contains three ATP-binding cassette genes; Abcb1a, Abcb1b and Abcb4. Abcb1a and Abcb1b are multidrug efflux transporters that are orthologs of the human multidrug-resistance gene 1 (MDR1 or P-glycoprotein). They are both upregulated in DOX-resistant tumors[47] and cell lines[48, 49] and are capable of exporting DOX from the cell. Abcb1a and Abcb1b carry either splice site or missense variants (rs31421794, rs8268091, rs6277917, rs16800179) for which NOD and PWK carry the alternate allele. Abcb4 is annotated as a phospholipid transporter in the liver[50], but it is also highly expressed in the megakaryocyte/erythroid progenitor[51, 52]. It is upregulated along with Abcb1a/b in cell lines exposed to DOX[53, 54] and the zebrafish ortholog has been shown to transport DOX from cells[55]. Abcb4 also contains both splice and missense SNPs (rs8279985, rs225705157) as well as an in-frame deletion that removes a lysine from the protein product (rs254540780).

We obtained mice with one or more of the Abcb1a/b and Abcb4 genes knocked out and investigated whether these genes influence DOX-induced neutropenia. All knockouts were produced on the FVB/N (FVB) background, which is a strain with high fecundity[5], and is genetically similar to the A/J, C57BL/6J, 129S1/SvImJ and NZO/HlLtJ strains on proximal chromosome 5[41]. In each case, we dosed the mice using the same protocol used to map DOX-induced neutropenia in DO mice. We compared the response of homozygous and heterozygous knockouts of each of the Abcb* genes to wild type mice. We hypothesized that more copies of the causative gene would increase post-dose NEUT. The Abcb1a knockout mice did not show a change in NEUT with increasing allele count (Figure 5A). Both the Abcb1a/b and Abcb4 knockouts showed an increase in NEUT with increasing functional copies of the allele (Figure 5B,C). Each copy of the Abcb1a/b allele increased log2-ratio NEUT by 2.2-fold, which is larger than the effect size of 1.3 for the NOD allele in Figures 4B and D. This reflects the difference between complete knockout of Abcb1a/b versus the effect of an allele that increases the function or expression of Abcb1a/b. Each copy of the Abcb4 allele increased log2-ratio NEUT by 1.2-fold, which is similar to the effect size in Figure 4D. There is significant inter-individual variation in each allele group. Since these mice are inbred, with the exception of the knocked out gene, this indicates that non-genetic factors that affect NEUT as well. The human ortholog of Abcb1b (ABCB1, or P-glycoprotein) has been associated with DOX resistance[56, 57]. The CTD contains 138 associations between ABCB1 and DOX[58] and polymorphisms in ABCB1 have been associated with DOX-induced neutropenia[59]. This leads us to conclude that Abcb1b is a candidate gene for resistance to DOX-induced neutropenia.

Figure 5

Effect of knocking out ATP-binding cassette transporters on neutropenia. (A) The knockout alleles of Abcb1a are plotted on the x-axis versus the log2-ratio of post-dose over pre-dose NEUT. Each point represents one mouse. The red line is the least-squares fit and the pink shading is the 95% prediction interval for the fit. The p-value indicates whether the slope differs from zero. (B) Same as A, but for the Abcb1a and Abcb1b combined knockout. Each functional allele of Abcb1b produces a 2.2-fold increase in log2-ratio NEUT. (C) Same as A, but for the Abcb4 knockout. Each functional allele of Abcb4 produces a 1.2-fold increase in log2-ratio NEUT.

The Abcb1b allele from NOD and PWK mice contains several non-synonymous polymorphisms that may be causative variants. rs8268091 changes a glutamine to an asparagine in the cytoplasmic domain. rs6277917 changes an arginine to a glutamine in the cytoplasmic domain. rs16800179 changes a serine to an asparagine in the cytoplasmic domain. rs31409832 is contributed only by NOD and changes an arginine to a histidine in the transmembrane domain. An unnamed deletion at 8 825 108 bp lies in a splice site proximal to exon 14 in both the NOD and PWK alleles.

Abcb4 is a phospholipid transporter that is generally not associated with drug resistance[50]. However, the Abcb4 allele in both NOD and PWK contains a 3 base pair, in frame deletion (rs254540780) that removes a lysine from the translated protein. Abcb4 has been shown to transport several xenobiotics in zebra fish[60] and variants in ABCB4 have been associated with DOX-induced cardiotoxicity in children[61]. We repeated the Abcb4 knockout experiment in a separate facility in Germany and were unable to replicate the result (Figure S1). Therefore, our data are inconclusive regarding whether Abcb4 also plays a role in DOX-induced neutropenia.

There are many peaks that are below the significance threshold that may harbor genes with smaller effects. While looking under these peaks increases the risk of false positives, we found two genes previously linked to DOX’s mechanism of action. A peak on chromosome 13 contains tyrosyl-DNA phosphodiesterase 2 (Tdp2, Figure S2A), a gene involved in the release of topoisomerase II from the 5′ end of DNA[62], and a target for inhibition in cancers that are also treated with DOX[63]. Another peak on chromosome 14 contains topoisomerase (DNA) II beta (Top2b, Figure S2B), which is consistent with DOX’s role as a topoisomerase inhibitor[64].

Genetic Influences on Doxorubicin Myelosuppression with CSF3

We dosed 179 DO mice with DOX followed by CSF3. We mapped the change in NEUT before and after dosing and found a peak on chromosome 5 at 6.51 Mb with a p-value of 2.78 × 10−17 (pGW < 0.001) that explained 29.2% of the phenotypic variance (Figure 3B). This is the same interval that we found when mice were dosed with DOX alone, suggesting that Abcb1 may affect neutropenia even when CSF3 is co-administered with DOX. The pattern of founder allele effects on chromosome 5 was also the same as in the DOX dosing group (Figure 6A,B). LYMPH, MONO, EOS and BW also had a QTL in the same location on chromosome 5, suggesting that this locus has a pleiotropic effect on WBC and BW. The BW peak on chromosome 5 is in the same location as the NEUT peak (Figure 6C, D). However, the founder allele effects are reversed; the NOD allele associated with decreased weight the CAST allele associated with increased weight. There was no correlation between the change in NEUT and BW (ρ = −0.064, p = 0.439). The Bayesian credible interval for the BW peak covers almost 10 Mb (4.3 – 13.9 Mb) and there were no clear candidate genes.

Figure 6

DO founder allele effects for NEUT and BW on Chr 5. The x-axis shows the position in Mb on chromosome 5. The y-axis plots the centered founder allele effects in standard deviations for (A) NEUT in mice dosed with DOX, (B) NEUT in mice dosed with DOX+CSF3, (C) BW in mice dosed with DOX, (D) BW in mice dosed with DOX+CSF3. The high NOD allele effects (dark blue) in A and B indicate that the NOD allele is protective for neutropenia. The low NOD allele effects in C and D indicate that the NOD allele is deleterious for BW. The CAST allele (green) also shows reversed allele effects between NEUT and BW.

Genetic Influences on Cyclophosphamide Myelosuppression

We dosed 191 DO mice with CYC and mapped the change in NEUT before and after dosing. We found a peak on chromosome 11 at 70.82 Mb with a p-value of 2.78 × 10−17 (pGW < 0.001) that explained 15.7% of the variance (Figure 3C). We found 232 SNPs with p-values below 1 × 10−8 (Figure 7A). Mice carrying the A/J, CAST, NZO/HlLtJ (NZO) or PWK alleles at these SNPs had lower post-dose NEUT. These SNPs intersect with the exons or untranslated regions of 7 genes (Table S3). Two of these genes, dynein axonemal heavy chain 2 (Dnah1) and SRY-box 15 (Sox15), contain missense SNPs. While there are no high significance SNPs that intersect with the exons of transformation related protein 53 (Trp53), this genes is within the confidence interval and is a promising functional candidate because it is a DNA check point gene that induces cell cycle arrest and apoptosis in the presence of DNA damage[65]. Trp53 may be influenced by nearby regulatory SNPs. Trp53 protein levels are upregulated in a rat model of CYC-induced liver injury[66] and CYC increases the levels of activated TP53 in HepG2 cells[67] and thus may be associated with the survival of replicating cells.

Figure 7

Association mapping plots for (A) CYC on chromosome 11 and (B) TAX on chromosome 4. The top panel in each plot shows the minor allele for the SNPs that are highlighted in red in the middle panel. SNPs in the middle panel are highlighted in red if they are above the pGW = 0.05 significance threshold. Trp53 is highlighted in panel A as a candidate gene.

Among the cell count traits that we measured, only NEUT had a significant peak on chromosome 11, suggesting that the gene that protects some mice from CYC-induced toxicity is acting primarily in NEUT. The action of Trp53 as a DNA-repair gene would be consistent with this hypothesis. We regressed out the most significant SNP and mapped the change in NEUT. The lowest p-value in the interval between 60 and 80 Mb was 4.74 × 10−3, which suggests that there are no other loci in this interval.

Genetic Influences on Docetaxel Myelosuppression

We dosed 154 DO mice with TAX and mapped the change in NEUT before and after dosing. We mapped a peak on chromosome 4 at 143.3 Mb with a p-value of 5.12 × 10−7 (pGW < 0.03) that explained 11.4% of the variance (Figure 3D). The p-value for pre-dose NEUT at this locus was not significant (p = 0.011). There was one SNP at 142 870 495 bp on chromosome 4 with a p-value below the significance threshold (rs225386879) and it did not overlap with the exons or untranslated regions of any genes. This SNP is not located in a conserved region but is in a CpG island and thus may have a regulatory role on one or more genes nearby. The closest gene is PR domain containing 2, with ZNF domain (Prdm2), which is a nuclear protein methyltransferase that has tumor-suppressing functions[68] that may impact cell survival in the context of chemotherapy. There are several genes with plausible functional roles in the QTL interval. Kazrin (Kazn) is a protein that co-localizes with acetylated microtubules[69]. Caspase 9 (Caps9) is part of the apoptosis cascade and has been associated with increased DNA damage in hematopoietic cells after exposure to an alkylating agent[70]. TAX also activates Casp9 and increases apoptosis in prostate cancer cell lines[71]. Another functional candidate is kinesin family member 1B (Kif1b), a molecular motor that is involved in vesicle transport and axon growth[72]. The peak on chromosome 4 also weakly affected MONO (Figure 2D) and we hypothesize that the causative gene(s) acts in the direct precursor of NEUT and MONO. We regressed out the most significant SNP and did not find any other significant SNPs.

Overlap with Human Neutropenia Loci

We obtained data from a Japan-based GWAS of chemotherapy-induced neutropenia that dosed patients with DOX, CYC and TAX[6]. We looked for overlap between the genes that they identified for each drug and found none consistent with our study. Variants in ABCB1 were associated with severity of DOX-induced bone marrow toxicity in a Danish pediatric cohort[73]. In another study of breast cancer patients treated with CYC and DOX along with filgrastim, variants in ABCB1 were not associated with hematototoxicity[74].

Discussion

In this study, we dosed DO mice with three chemotherapy drugs, each with distinct modes of action, as well as CSF3, a stimulator of NEUT production. As seen in humans, we observed differences in the neutropenic response between subjects. For both CYC and DOX, we observed neutropenia levels on par with Grade 4 neutropenia in humans[3]. While mean NEUT did not decrease significantly in the TAX group, it did decrease in a subset of mice, recapitulating the variation in response seen in human patients[75]. The combination of DOX with CSF3 increased neutrophil counts in many mice, but some still experienced neutropenia. This variation in response to CSF3 may have a genetic basis and further research in this area is warranted.

The genetic architecture of hematotoxicity was quite different between drugs. DOX-induced hematotoxicity was influenced by a pleiotropic locus on chromosome 5 that affected all cell populations as well as BW. This suggests that the gene or genes in this region are acting in a manner that protects the entire hematopoietic system and one of them may be a therapeutic target to protect patients from hematotoxicity. There also appear to be at least two genes in this region with different allelic effects, as shown in Figure 4. We showed that Abcb1b influenced DOX-induced hematotoxicity and Abcb4 or Sri may play a role as well. One possible explanation for the lack of replication of the Abcb4 knockout experiment may be due to diet and the influence of the microbiome on xenobiotic metabolizing enzymes. We have included a comparison of the two diets in Table S4.

When we dosed DO mice with DOX, we identified a pleiotropic locus on chromosome 5 that influenced most hematopoietic cell populations. The DOX-associated locus contains Abcb1a and Abcb1b, which are the murine orthologs of ATP-binding cassette, sub-family B, member 1 (ABCB1), an efflux transporter that has been associated with DOX resistance in human cancers[76]. However, there is currently no evidence from epidemiological or pharmacogenomics studies that ABCB1 variants affect susceptibility to myelosuppression. We showed that the presence of functional copies of Abcb1b is positively correlated with resistance to neutropenia after treatment. The protective effect of Abcb1b on the hematopoietic system suggests that it is acting in the liver or in one of the early hematopoietic progenitor cells. If these genes are functioning in the hematopoietic system, they may be affecting the survival of uncommitted progenitor cells in the hematopoietic hierarchy. Abcb1b is constitutively expressed in hematopoietic progenitor cells and may be acting at that level[52]. Abcb1b may also be induced in liver and may be acting to protect the entire body by excreting DOX metabolites directly into the bile.

When we dosed mice with both DOX and CSF3, we found an association for neutropenia in the same region of chromosome 5 as DOX alone. Although CSF3 attenuated neutropenia in many mice, the chromosome 5 locus still affected susceptibility to neutropenia. It is interesting to note that the loss in BW did not correlate with the degree of neutropenia. This suggests that different mechanisms of toxicity may affect different organ systems. For example, the NOD allele of Abcb1b may confer protection from DOX-induced neutropenia, but it does not appear to protect mice from weight loss. Conversely, mice carrying the CAST allele on proximal chromosome 5 had the most severe neutropenia, but the least weight loss. If Abcb1b is acting in the liver, then it should export DOX metabolites into the bile and out of the body, precluding any other toxicity. Alternately, if Abcb1b is protecting cells in the hematopoietic system, then it is exporting DOX metabolites back into the circulation, where they may have a deleterious effect on other organ systems.

In contrast to DOX, different loci affected each cell population when we dosed mice with CYC or TAX. CYC is a myeloablative drug that reduced all cell counts and BW, yet changes in each cell population were associated with distinct loci. The same pattern of differential association of cell counts with distinct loci was true for TAX. This suggests that different chemotherapy drugs may affect different hematopoietic progenitor cells and their committed cell populations in a distinct manner. There may be lineage-specific genes that are expressed in different progenitor populations that confer resistance to toxicity and a greater understanding of the genetics underlying these alleles may help to improve patient outcomes.

When we mapped the change in NEUT after CYC administration, we found a broad peak on chromosome 11 containing many genes. We narrowed the candidate gene list to Trp53, a plausible candidate gene involved in DNA repair and cell cycle regulation. CYC is a DNA alkylating agent that damages DNA and also leads to cell cycle arrest[77], consistent with the action of Trp53. The protective allele of Trp53 may be acting in committed neutrophil progenitor cells, since the QTL appears only for NEUT and not for other myeloid derived cells. When we regressed out the effect of the chromosome 11 locus, we did not find significant peaks over any phase I or II metabolizing genes or DNA repair genes that have been shown to influence differential metabolism in humans[77, 78]. It is important to note that while there may be more than one gene at this locus, we may be unable to dissect it due to the width of the linkage disequilibrium blocks in early generations of the DO.

For each drug that we tested, we identified a different locus with chemotherapy-induced neutropenia. This is consistent with findings in both humans[6] and model organisms[7, 8]. This result highlights the importance of developing tractable translational models to understand the role of genetic variation in chemotherapy-induced toxicities. Inbred mouse strains are valuable tools, but these results would not have been discovered using a single inbred strain. The broad genetic diversity in DO mice allowed us to sample alleles from a large number of genes to interrogate their effects. DO mice also have the advantage that they are a whole body system that includes complex interactions between organs that may influence toxicity.

We were initially concerned that our results provide no overlap with the results of a clinical GWAS for chemotherapy-induced neutropenia[22]. However, we agree with the authors of the human GWAS that their study had low statistical power and was further weakened by uncontrolled confounding variables. Each of the three drug groups had a small number of cases, ranging from 168 for doxorubicin to 423 for cyclophosphamide. These are very low sample sizes for an unbiased human GWAS. Further, the authors did not use measured neutrophil counts; rather, they used neutropenia grades binned into two groups (Grades 1 & 2 or 3 & 4). This artificial discretization further reduces power. In our view, this highlights the strength of using whole animal dosing to perform GWAS in murine models. We can afford adequate sample sizes, we can control confounding covariates and we can impute full genome sequences to find candidate genes. Further, we can control the minor allele frequency in DO mice to ensure that we do not have rare variants (i.e. with a MAF < 5%), which is impossible in human GWAS studies. The high MAF in the DO (an average of 1/8) allowed us to identify large QTL peaks above the α = 0.05 genome-wide significance level for all three drugs and identified plausible and translationally relevant genes for two. This stands in contrast to the clinical GWAS in which the authors “failed to identify genetic variants that surpassed the genome-wide significance level”[22].

One of the limitations of this study was the relatively small numbers of mice in each cohort. We dosed between 200 and 400 mice per drug and, while we identified peaks for each drug, the width and resolution of the peaks was limited by early DO outbreeding generation with low numbers of recombinations. Many chemotherapy drugs have strong effects and this produces correspondingly strong heritabilities and peaks. The mice in this study were also from generations 3 through 5 of outbreeding and hence they had large recombination blocks (median = 6.47 Mb), which limits the mapping resolution of the experiment. DO mice are currently at generation 23 of outbreeding and we estimate that their median recombination block size is now 2.9 Mb. Investigators who wish to estimate the number of mice required to power a study using DO mice are referred to the power simulations and associated figure in Gatti et al.[9].

These findings leave us with several questions for future research. In which organ or cell population is Abcb1b acting? What is the other gene on proximal chromosome 5 that makes DO mice susceptible to neutropenia? What are the genes that regulate the response to CSF3? How can we translate these findings to the clinic? While we do not expect genetic polymorphisms in mice to be identical to polymorphisms in the human population, the pathways that affect the variation in susceptibility to myelosuppression may translate between species. An improved mechanistic understanding of the pathways that affect susceptibility to chemotherapy-induced myelosuppression may lead to new treatment options in the clinic.

Table 1

Additive heritability estimates for circulating cell populations and body weights. The p-values were adjusted using a Bonferroni correction.

	Doxorubicin			Doxorubicin+Neupogen			Cyclophosphamide			Docetaxel
Phenotype	Heritability	p-value	Adj. p-value	Heritability	p-value	Adj. p-value	Heritability	p-value	Adj. p-value	Heritability	p-value	Adj. p-value
NEUT	0.516	8.07E-06	1.94E-04	0.395	0.044	1.000	0.272	0.052	1.000	0.352	0.028	0.681
LYMPH	0.369	4.94E-04	1.19E-02	0.486	0.019	0.459	0.143	0.175	1.000	0.240	0.084	1.000
MONO	0.141	0.067	1.000	0.175	0.209	1.000	0.117	0.217	1.000	0.249	0.078	1.000
PLT	0.424	1.13E-04	2.72E-03	0.263	0.068	1.000	0.255	0.061	1.000	0.213	0.106	1.000
RBC	0.240	0.010	0.230	0.137	0.201	1.000	0.341	0.025	0.591	0.401	0.017	0.407
BW	0.104	0.126	1.000	0.302	0.043	1.000	0.000	0.500	1.000	0.256	0.072	1.000