The Role of Alleles in Susceptibility to Haloperidol-Induced Toxicity in Mice and Humans

Summary

Gary Peltz and colleagues examine the role of ABCB5 alleles in haloperidol-induced toxicity in a murine genetic model and humans treated with haloperidol.

Introduction

Relatively little is known about the mechanisms regulating the entry of many drugs into the brain, or about the genetic mechanisms affecting susceptibility to many drug-induced central nervous system (CNS) toxicities. Because of this lack of information, candidate drugs affecting the CNS have the lowest rate of survival during clinical development (<2%), and drug distribution across the blood-brain barrier has been identified as a contributing cause for drug candidate failure [1]. The blood-brain barrier is a very complex multicellular vascular structure, which separates the brain from the peripheral blood by actively controlling the entry and efflux of drugs [2]. Therefore, the mechanistic insight obtained by identification of genetic factors affecting a drug-induced CNS toxicity could aid in the development of new medicines and could improve our ability to optimally utilize existing treatments. For example, haloperidol is a potent dopamine receptor antagonist that is used to treat psychotic disorders [3]. However, haloperidol-induced alterations in the extrapyramidal motor system depress the ability to initiate voluntary movements. This effect triggers the characteristic extrapyramidal symptoms (EPS) that develop in 40%–76% of chronically treated human patients, which include tremors, Parkinsonian rigidity, and decreased spontaneous movement [4,5]. The same symptoms appear in some (but not all) of the 27 inbred strains that were treated with haloperidol [6]. In this murine genetic model, the time required (latency) for a mouse to move all four paws after being placed on an inclined wire-mesh screen was measured as an index of the haloperidol-induced Parkinsonian rigidity that is observed in haloperidol-treated human patients. Although all strains had plasma drug levels that were similar to those measured in human patients, the extent of the haloperidol-induced toxicity (HIT) varied in a highly (85%) heritable manner [6] among the inbred strains.

Because no pharmacogenetic factors for HIT have previously been identified in humans or mice, we wanted to analyze this murine genetic model of HIT. Of concern, mouse genetic studies have, overall, produced rather disappointing results; they have definitively identified <1% of the genes affecting the studied traits [7]. Beyond the well described problems associated with murine linkage studies [8], several analyses have predicted that the conventional methods used to analyze genome wide association studies (GWAS) using inbred mouse strains will produce even more disappointing results, owing to their low power and a high false positive rate [7,9,10]. However, the availability of next generation sequence (NGS) data, which can provide complete genome-wide variant information for an expanded number of inbred strains [11], could enable murine GWAS to identify causative genetic factors for many biomedical traits [12]. The detailed genetic variation information obtained using NGS has increased power for pinpointing underlying causal variants relative to studies using SNP array data [13]. Unlike the conventional murine GWAS methods, haplotype-based computational genetic mapping (HBCGM) [14] has provided an alternative method for analysis of murine GWAS data that has identified causative genetic factors for many biomedical traits (reviewed in [12]). HBCGM first organizes the allelic information into haplotype blocks, which facilitates the identification of the causative genetic factors through correlation of the phenotypic data with the haplotype-based pattern of genetic variation. Since our prior studies analyzed much less complex datasets, which had fewer SNPs and haplotype blocks, we did not know if HBCGM could analyze the far more complex datasets produced using NGS data. To determine if HBCGM could be successfully combined with an SNP database generated from whole genome sequence data, we analyzed a murine genetic model of HIT in order to identify causative genetic factors. A translational human genetic association study was then performed to determine if similar genetic factors affected human susceptibility to HIT.

Methods

Next-Generation Sequencing of the Genome of Nine Inbred Mouse Strains

Genomic DNA was obtained from the Jackson Laboratory (http://www.jax.org). A total of 3 μg of DNA from each of nine inbred mouse strains was sonicated using a Covaris instrument with the following settings: duty cycle, 10%; intensity, 4; cycles per burst, 200; and time, 60 seconds. After end repair of the DNA fragments and addition of “A” base to the 3′ ends, standard Illumina genomic DNA sequencing adaptors were ligated. The samples were run on a 3% agarose gel, gel slices between 300–500 bp were excised, and the DNA was purified with a Qiagen column. The DNA was divided into three aliquots of equal quantity, and amplified separately using one of the following three DNA polymerases: Phusion HF (Finnzymes), Kapa HiFi (Kapa Biosystems), and Herculase (Agilent). For analysis of the SJL strain, Phusion GC (Finnzymes) was also used to amplify a separate aliquot. The PCR conditions for amplification were as follows: 98°C for 3 minutes, followed by nine cycles of 98°C for 20 seconds, 62°C for 30 seconds, and 72°C for 45 seconds, and then a final 72°C for 5 minutes. The amplified products were separated on an agarose gel, slices between 400 and 600 bp were excised, and the DNA was eluted from the gel slices and sequenced on an Illumina HiSeq 2000 instrument.

Sequence Analysis and Variant Calling

The quality statistics for the allele calls at each position were evaluated using the FASTX-Toolkit (http://hannonlab.cshl.edu/fastx_toolkit/). A rapid decay in sequence quality was observed toward the end of each read, usually within the last 10 (or sometimes 20) bps. To avoid artifacts caused by low quality sequence data, the last 10 or 20 bps of each read were trimmed, based upon a visual inspection of the quality statistics plot for the reads in each lane. The extent of trimming was designed to ensure that the quality scores should be ≥30 in the 1st quartile and ≥20 in the 5th percentile. All sequences were trimmed by 10 bps; with the exceptions being one of the three lanes for the MA/MyJ and MRL, and two of the three SM/J lanes. Since a much larger number of sequences were obtained for SJL, SJL sequences were trimmed by 20 bps to ensure the overall quality of the base calls. The read sequences were then aligned to the reference C57BL/6J mouse genome sequence (MGSCv37 assembly) using BWA [15]. Optical and PCR duplicates were then marked using the Picard tool (http://picard.sourceforge.net/). The sequence variants (SNPs and INDELs) for each strain were identified using Samtools [16], and the default settings were applied. The raw sequence data for the 17 inbred mouse strains described in [11] were obtained from the European Nucleotide Archive (http://www.ebi.ac.uk/ena/), and the same data processing procedures were applied, but these sequences were not trimmed.

To ensure high quality variant calls, the variants were further filtered using two additional criteria. (i) To ensure a high level of confidence in the existence of the variant, we require the QUAL score, which represents the (Phred-scaled) probability that all observed samples are homozygous reference alleles, should be ≥50. (ii) To ensure a high level of confidence that the variant allele is indeed homozygous, the genotype call (the “GT” value) must be a homozygous alternative call, and that the Phred-scaled genotype likelihood (the “PL” value) for the homozygous alternative call is at least 20 units above that of any other PL value. If the variant allele call for a strain does not meet these criteria, the call is converted to an “N” (undetermined) instead of being disregarded, which ensures that we don’t mistakenly assume that it is the reference allele. Information about the fold-coverage and other aspects of the sequencing methods are provided in S1 Text.

Haplotype Block Construction and Genetic Mapping in Mice

HBCGM utilizes bi-allelic SNPs that are polymorphic among the strains. Only bi-allelic SNPs with at least one definitive homozygous alternative call among the 23 strains were included, while all other variants (including INDELs) were marked and excluded from HBCGM analysis. The locations and potential codon-changes caused by the SNPs are then annotated relative to the encoded genes using predictive gene models from Ensembl version 65. Only SNPs meeting the following criteria were used for haplotype block construction: (i) polymorphic among the strains with input trait data; and (ii) there were at least eight strains (or one-half of the number of input strains) with unambiguous allele calls. Haplotype blocks with 2, 3, 4, or haplotypes were then dynamically produced as described in [17], and the correlation between the input phenotypic data and the haplotype pattern within each identified block was evaluated as described [14]. The genes are then sorted based upon the ANOVA p-value (in increasing order) for numeric data or by the F statistic (in decreasing order) for categorical data. For a gene that is covered by multiple blocks, the smallest p-value (or the largest F statistic) obtained for all blocks within that gene was used.

The genetic effect size is calculated using our previously described method [18]: Genetic effect size = (SSB − (k − 1) * MSE)/(SStotal + MSE), where SSB is the between-group sum-of-squares of the ANOVA model, SStotal is the total sum-of-squares, k is the number of groups, and MSE is within-group sum-of-squares divided by (n − k), where n is the total number of objects in the model. Given the very large number of haplotype blocks produced by the HBCGM method, all results were automatically filtered according to pre-determined criteria using software that selected the correlated genes that were expressed in a designated target organ and had missense SNPs.

Robustness Assessment

To assess the robustness of the HBCGM output obtained using the data obtained from different sets of input strains, we developed software that would automatically perform a re-sampling analysis. For each iteration, a number (1, 2, or 3) was randomly selected to indicate the number of strains whose data would be randomly excluded from the 17 strains with available latency data. Next, the selected number of strains (and each of the individual strains to be deleted were also randomly selected) was excluded for this iteration. Then, HBCGM were run on this re-sampled data and the result was recorded (using p = 0.01 as the cutoff). This process was repeated 100 times. If a re-sampled strain panel was identical to one that had already been re-sampled, this iteration was disregarded, and another subset was randomly chosen using the same procedure. To measure the overall association of a gene with the phenotype among the 100 random re-sampled sets of output, we defined a robustness score based on the p-values obtained, which was calculated as: Score= ∑ i=1 100 -log 10 ( p i ) , where p_i was the resulting p-value for that gene in the ith re-sampled experiment (if the p-value was bigger than the cutoff and therefore censored in the result, p_i was set to 0.01, i.e., the cutoff). Genes with a larger robustness index have stronger overall correlation among the different subsets analyzed. Their calculated scores then were used to rank the genes, and those with the highest scores were output.

Statistical Analysis of Haloperidol Latency Data in Mice

The haloperidol-induced latency data measured on days 0, 3, 7, 30, 60, and 120 (MPD datasets: 39407, 39408, 39409, 39410, 39411, and 39445, respectively) and the plasma haloperidol concentrations measured on days 30, 60, 90, and 120 (MPD numbers: 39403–39406, respectively) were obtained from the Mouse PHENOME Database (http://phenome.jax.org/). The correlation analyses were performed using the phenotypic data for all 27 strains that were evaluated. The non-parametric Spearman’s correlation coefficient (rho) was calculated [19], and the statistical significance of the correlation between the tested datasets was evaluated using this result. This non-parametric test was used to avoid producing spurious results that could be caused by the non-normal distribution of the data and by the presence of extreme outlying values. The p-values were then adjusted for multiple testing using the Benjamin-Hochberg method [20]. To minimize the effect produced by the strain with an extreme outlying value and to better visualize the relationship of the variables, the haloperidol-induced latency data from day 7 and day 30 were log-transformed (10-based). (It should be noted that the log-transformation has no effect on the Spearman’s rho, nor did it affect the corresponding test result.)

Haloperidol Toxicity Measurements

All animal experiments were performed according to protocols approved by the Stanford Institutional Animal Care and Use Committee. All mice were obtained from Jackson Labs, and were used at 10–12 weeks of age, and the results are reported according to the ARRIVE guidelines (S1 Text) [21]. Haloperidol was administered to the mice by one of two methods: (i) a haloperidol pellet (Innovative Research of America) was subcutaneously implanted, which continuously released a dosage equivalent to 3 mg/kg/day, using published methods [22]; or (ii) with 10 mg/kg haloperidol (Sigma) IP qd for measurement of drug and metabolite concentrations. The haloperidol-induced latency was measured as the time required for a mouse to move all four paws after being placed on a vertical metal mesh screen as previously described [22].

Haloperidol Quantitation

Haloperidol was purchased from Sigma-Aldrich; haloperidol-d4 was purchased from Toronto Research Chemicals Inc.; HPLC grade acetonitrile and water were from Honeywell Burdick and Jackson International Inc.; and formic acid was purchased from Pierce Thermo Scientific.

Mouse brain tissues of known weight were homogenized in methanol using a Precelly 24 (Bertin Technologies) homogenizer. 50 μl of homogenate was aliquoted, and haloperidol-d4 was added as internal standard to a final concentration of 100 ng/ml. After the mixture was incubated with three volumes of cold acetonitrile at −20°, it was centrifuged at 14,000 rpm for 10 minutes. The supernatant was dried in a speed vacuum and re-suspended in the original volume. A standard curve was prepared using brain homogenate from an untreated mouse in a linear range from 0.5 ng/ml to 250 ng/ml. Quantitative analysis of haloperidol and its metabolite (HPP+) was performed on an Agilent accurate mass QTOF 6520A coupled with an UHPLC Infinity 1290. The analytes were separated on a Phenomenex Kinetex XB-C18 column 2.1 × 100 mm. The flow rate was 0.5 ml/min with a gradient of solvent B (acetonitrile with 0.1% formic acid; solvent A is 0.1% formic acid) from 5% to 35% in 10 minutes, then 35% to 95% in 5 minutes. Positive electrospray full scan in the range of m/z 110–1,000 spectra were collected. Data analysis was done using Agilent quantitative analysis software.

Desorption Electrospray Ionization Mass Spectrometroscopic Imaging

High mass resolution/high mass accuracy mass spectrometry tissue imaging was performed using a lab-built desorption electrospray ionization mass spectrometry imaging (DESI-MSI) source coupled to an LTQ-Orbitrap XL mass spectrometer (Thermo Scientific). DESI-MSI was performed in the positive ion mode from m/z 200–500, using the orbitrap as the mass analyzer at a resolving power of 60,000. The spatial resolution of the imaging experiments was of 150 μm. A histologically compatible solvent system dimethylformamide:acetonitrile (DMF:ACN) 1:1 (v/v) was used for analysis [23] at a flow rate of 0.8 μl/min. The N₂ pressure was set to 175 psi. After DESI-MSI, the same tissue section was subjected to hematoxylin–eosin (HE) staining. The software ImgGenerator (freeware, http://www.msimaging.net/) was used for converting raw files into 2-D images. Spatially accurate ion images were assembled using BioMap software (freeware, http://www.maldi-msi.org/). The protonated form of haloperidol was detected in the tissue sections at m/z 376.14679 (mass error of −1.7 ppm). The isotopic distribution of the ion at m/z 376.14679 was also found to agree with the molecular formula of haloperidol. DESI-MS ion images of haloperidol were compared to optical images of the same tissue in HE-stained tissue sections. The total abundance of haloperidol in each pixel from a DESI-MS ion image was extracted and averaged to obtain the average abundance of haloperidol per tissue section. This procedure was performed for each tissue section analyzed for each strain of mice. Based on the average total abundance of haloperidol for each strain, the fold change was then calculated.

RNA In Situ Hybridization

Abcb5 and Drd2 probes were designed to hybridize to their corresponding mRNAs. The 741 bp Abcb5 probe was PCR amplified from RIKEN_cDNA_clone B020023O04 using the primer set: AGGCCAGAAACAGAGGATTG; and CCTGGTAGAGCATGGCTTTG. The 1,126 bp Drd2 probe was PCR amplified from total cDNA from a C57BL/6J mouse using the primer set: CCGTGAACCCCATCATCTAT and GGTTTGGTGCATGTATGGTG. PCR reactions were performed using Phusion High-fidelity DNA polymerase (Thermo Scientific) and the products were cloned into pCR-BluntII-TOPO vector (Life Technology). Clones with correct insertions were linearized using either BamHI or NotI (New England Biolabs) for in vitro transcription. Anti-sense and sense RNA probes were then in vitro transcribed using Riboprobe System-SP6/T7 (Promega) and DIG RNA labeling mix (Roche), DNase treated, precipitated in lithium chloride (Ambion), and stored at −80°C in aliquots. Freshly dissected mouse brain was immediately bisected sagittally, and fixed in 4% paraformaldehyde in DEPC treated phosphate buffered saline at 4°C overnight. Fixed tissues were then saturated in DEPC-treated 20% sucrose at 4°C overnight, then flash frozen in OCT compound (Tissue-Tek) and stored at −80°C. The brain sections were cut at 20 μm thickness, and serial sections were collected on Superfrost Plus glass slides (Fisher Scientific). Sets of serial sections across the whole brain were then hybridized anti-sense probes, and adjacent sections were hybridized with the sense probe (as negative controls), at 63°C in a moist chamber overnight. After washing and blocking, tissue sections were then incubated with AP-conjugated anti-DIG antibody (Roche) at 4°C overnight. Tissue sections were then washed and chromogenic reactions were then performed by incubation with the AP substrate BM Purple (Roche) at room temperature. The sections were then examined using a light microscope, and a purple/blue color indicates the localization of the targeted mRNA within the tissue.

Statistical Analysis of Data Obtained from Inbred and Chromosome Substitution Strains

The latency measurements from chromosome substitution strain (CSS) mice were log-transformed, and a two-sample t-test was applied to compare the latency of CSS12 to that of C57BL/6. The differences between the latency on the log-scale and the corresponding 95% confidence intervals were evaluated, and then transformed back to the original scale. Similarly, the measured haloperidol and HPP⁺ concentrations were log-transformed and a two-sample t-test was applied to compare the concentration in A/J to that in C57BL/6. The differences between the latency on the log-scale and the corresponding 95% confidence intervals were evaluated and then transformed back to the original scale. For comparison of the HPP⁺ levels in A/J, CSS12, and C57BL/6 strains, the concentrations were log-transformed and an F-test was applied. The differences between the latency on the log-scale and the corresponding 95% confidence intervals were evaluated using SAS (version 9.3), and were then transformed back to the original scale.

Analysis of Abcb5 Knockout Mice

Mice with a homozygous Abcb5 gene deletion were prepared and bred onto the C57BL/6 by nine backcrosses, and these mice did not express Abcb5 mRNA (JPG, MMG, and colleagues, personal communication). The latency in control wild-type littermate and Abcb5 knockout mice were analyzed after 0 to 7 days of haloperidol 10 mg/kg/day IP administration. The statistical significance of the difference in the latency measurements between wild-type and Abcb5 knockout mice was determined using a two-sample two-sided t-test on log-transformed data for each day as described in the preceding paragraph. Since gender could affect the haloperidol response, we performed a two-way ANOVA to assess the effect of gender. However, the gender effect was only significant on day 4 (p = 0.046) and was not significant for all of the other days. Therefore, an adjustment for the gender effect was not performed.

Human Genetic Association Analysis: Human ABCB5 Alleles and Haloperidol Toxicity

The study was approved by the ethics committee at the University of Munich, and was carried out in accordance with the ethical standards laid down in the 1964 Declaration of Helsinki and later revisions (S2 Text). Additional description of these participants and their characteristics can be found elsewhere [24]. Written informed consent was obtained at the time of recruitment for each participant. In brief, this human cohort consists of 101 patients of European descent that were treated with haloperidol as a monotherapy for a period of one month during the acute phase of their psychotic illness. The patients were enrolled over a 3-year period (1999–2002) at the Department of Psychiatry, Ludwig-Maximilians-University of Munich, Germany. Patients were excluded if they had a known contra-indication for treatment with haloperidol, had tardive dyskinesia previously, or severe neurological or medical disorders, organic brain diseases, pregnancy, or acute suicidal tendencies. Furthermore, the patients were excluded if they received a co-medication, such as beta blockers, or antidepressants that could possibly influence the antipsychotic treatment or its side effects. The patients were treated with haloperidol during the first phase of the illness’s episode, and then shifted to another antipsychotic treatment if there was a lack of response or if severe side effects developed. Patients were treated with haloperidol according to clinical need without any dose limitation during the acute phase of the illness. The mean dose of haloperidol was 8.76 mg/d (SD = 4.4) for day 1, 10.54 mg/d (SD = 6.0) for day 3, 10.88 mg/d (SD = 5.9) for day 7, 10.43 mg/d (SD = 6.2) for day 14, 10.95 mg/d (SD = 7.3) for day 21. As this was a naturalistic study doses were adapted to clinical requirements. Diagnosis was obtained through the Structured Clinical Interview for DSM-IV (SCID-I) interview, and psychopathological measurements were administered by two psychiatrists with highly reliable inter-rater evaluation results (k > 0.80). Haloperidol plasma levels, positive and negative syndrome scale (PANSS), ESRS, and UKU scores were assessed at baseline and on days 3, 7, 12, 21, and 28 as described [24]. Since genotyping was performed after the clinical part of the study was closed, it was not possible for any of the evaluators to be aware of the genetic results.

The alleles at 38 SNPs in ABCB5 were genotyped on a MassARRAY platform (Sequenom) for 85 patients, which represents all of the available DNA samples obtained from the individuals in this cohort. Briefly, after a multiplex PCR and shrimp alkaline phosphatase treatment, a single base extension reaction was performed. Following the de-salting of extension products with SpectroCLEAN resin (Sequenom), samples were spotted on SpecroCHIPs GenII (Sequenom), and analyzed with a MassARRAY MALDI-TOF mass spectrometer. Allele specific extension products and the resulting genotypes were identified using Typer 3.4 Software (Sequenom). For genotyping quality assurance, the CEU HapMap Trios (Coriell Institute for Medical Research) were included and compared with the HapMap database (www.hapmap.org). Nine of the SNPs were not polymorphic in these 85 patients, and were excluded from further analysis. None of the remaining 29 SNPs violated the HWE condition, since the smallest p-value obtained using a traditional Chi-squared test [25] was 0.11. Of note, 26 SNPs exhibited all three possible genotypes (i.e., homozygous major allele, heterozygous and homozygous minor allele), while three SNPs had only two genotypes. For the 26 SNPs with all three possible genotypes, we used an additive model to assess the association using ANOVA model. For the three SNPs with only two observed genotypes (homozygous major allele and heterozygous), only one grouping was tested using two-sample t-test. For each patient, the degree of parkinsonoid toxicity was rated as 0 (absent), 1 (mild), 2 (moderate), or 3 (severe) for each time point. A patient was classified as exhibiting parkinsonoid toxicity if the score was ≥2, as not exhibiting this toxicity if the score was ≤1 at each time point. Patients were similarly rated and classified for dyskinesia and akathesia, the other two toxicities. A combined measure of HIT at each time point was defined as follows: a patient was designated as exhibiting toxicity for (i.e., 1) if any of the three toxicity measurements was ≥2, and as non-toxic (i.e., 0) if the three toxicity measurements were all ≤1. Finally, as a cumulative index of the toxic effect of haloperidol, the area-under-curve (AUC) of the combined toxicity measures for days 0, 1, 3, and 7 was calculated. This value represents the cumulative toxicity induced by haloperidol over the indicated period, and it assumes one of the following values: 0, 1, 2, 3, or 4.

Results

Whole Genome-Genetic Association Mapping of Susceptibility to HIT

As described in S1 Text, we analyzed NGS data to produce a 16 million SNP database with alleles covering 26 strains, including the 12 strains sequenced here, 13 strains analyzed in [11], and the reference C57BL/6 strain (S1 Fig.; S1 Table). Since our prior studies analyzed much less complex datasets with fewer SNPs and blocks, we did not know if HBCGM could analyze this far more complex (~900-fold more blocks) whole genome NGS-based dataset. However, we found that whole genome HBCGM could identify the known genetic factors for four biomedical traits (S2 Fig.; S2 Table; S3 Text). Moreover, we also found that its performance was superior to that of another genetic association method, which analyzed the correlation with a selected set of individual SNPs, especially when the inbred strains exhibited three or more distinct phenotypic responses (S3 Table).

We then examined available data on HIT (measured after 3, 7, 30, 60, or 120 days of haloperidol treatment) for the 17 strains with available genomic sequence (S2B Table). As described in the S2B Table, our analysis of these data indicated that the prolonged latency appearing in some of the inbred strains (i) was due to haloperidol treatment; (ii) was not correlated with plasma drug levels; and (iii) had its relative magnitude in the strains maintained across all treatment periods (ranging from 3 to 120 days) analyzed. The lack of correlation between latency and plasma haloperidol level was expected, since the effect of CNS-acting drugs is dependent upon the brain (and not the plasma) drug concentration. HBCGM analysis indicated that genetic variation within Abcb5 had the highest correlation (p = 8.3 × 10⁻¹¹) with the day 30 latency values (Fig. 1A). HBCGM analysis also indicated that Abcb5 was highly correlated with the latencies measured after 60 (rank number 3, p-value 2.9 × 10⁻¹³) or 120 (rank number 2, p-value 1.7 × 10⁻¹⁰) days of haloperidol treatment. Several other factors also indicated that Abcb5 was a likely candidate gene: (i) When the robustness of the mapping results were assessed by randomly deleting sets of strains and iteratively repeating the analysis, Abcb5 was among the top genes identified after the latency data were analyzed for treatment days 30, 60, or 120 (S4 Table). (ii) Among the 440 SNPs within or near Abcb5, eight missense SNPs divided the 17 strains into three haplotypic groupings with different latencies: group 1 (range: 3–44 seconds; 14 strains); group 2 (144–171 seconds; NZW, NZO); and group 3 (234 seconds, A/J). The three strains with extreme latency values share alleles at three missense SNPs that are within a nucleotide-binding region of the protein (Fig. 1B). Also, another codon changing SNP (cSNP) (rs30781912) converted a serine, which is conserved in all vertebrates, to an alanine in the strain (A/J) with the highest latency. (iii) We investigated whether Abcb5 alleles could predict the measured latency for two inbred strains with extreme latency values (NON/ShiLtJ, 2.8 seconds; and NZL/LtJ, 211 seconds) whose genomic sequence was not yet characterized. After characterizing their Abcb5 alleles, we found that their measured latencies were well within the range predicted by their Abcb5 haplotype (Table 1). Lastly, genetic variation within Abcb5 could affect susceptibility to HIT, since it is a member of the ABC-transporter gene family. Although little is known about murine Abcb5, its human homologue (ABCB5) is an energy-dependent efflux transporter that mediates the resistance to multiple drugs and chemotherapeutic agents [26–30]. Drug efflux pumps located at the luminal membrane of capillary endothelial cells in the brain, which export a drug back into the circulating blood after it has entered the brain, are major determinants of drug levels in the brain [31].

Analysis of Chromosome Substitution Strains and Brain Haloperidol Levels

A series of C57BL/6J-A/J CSS have been produced where each strain is homosomic for a single specified A/J chromosome on an otherwise C57BL/6 genetic background [32]. Since A/J and C57BL/6 exhibited extreme latency values, we tested the haloperidol-induced latency in four CSS. Of these, only CSS12 (which have A/J alleles for every gene on Chromosome 12) had a significantly increased latency after 5 (6.7-fold, t value = 12.98, p < 0.0001, 95% CI 4.9–9.3) and 10 days (8.5-fold, t = 8.41, p < 0.0001, 95% CI 4.8–14.8) of haloperidol treatment (Fig. 2A). Among the 34 genes on Chromosome 12 that were correlated with susceptibility to HIT, Abcb5 had the highest level of correlation and was the most plausible candidate based on its biological function (S5 Table).

Haloperidol is oxidatively metabolized to a pyridinium species (HPP⁺) that is present in brain tissue obtained from haloperidol-treated humans [33] or rodents [34]. HPP⁺ is a structural analog of the oxidative metabolite (MPP⁺) of a Parkinson-inducing neurotoxin N-methyl-4-phenyl-1,2,3,6-tetrahydropyridine (MPTP), and HPP⁺ has been shown to induce toxicity in dopaminergic neurons [33,35,36]. To determine if differences in brain drug (or metabolite) levels were associated with susceptibility, we measured the concentration of haloperidol and HPP⁺ in brain tissue obtained from C57BL/6J and A/J mice after treatment with haloperidol 10 mg/kg/day for 4 days. This higher dose of haloperidol was used to ensure that we could measure the brain concentrations of both haloperidol and HPP⁺, since the metabolite is present at much lower (~100-fold) concentrations than the parent drug. C57BL/6 mice treated with this haloperidol dose do not exhibit an increased latency (see below). Consistent with HIT occurring in A/J but not C57BL/6J mice, the brain levels of haloperidol (2.2-fold and 2.3-fold, t values 4.88 and 4.99, p-values 0.008 and 0.002, 95% CI 1.4–3.5 and 1.5–3.5, respectively) and HPP⁺ (1.4-fold and 1.7-fold, t values 4.75 and 2.39, p-values 0.009 and 0.05, 95% CI 1.1–1.7 and 1.0–2.6, respectively) were both higher in A/J than in C57BL/6J mice in two independently performed experiments (Fig. 2B).

We also examined HPP⁺ levels in brains obtained from A/J, C57BL/6J, and CSS12 mice after a prolonged period (10 days) of haloperidol administration (3 mg/kg/day). The brain HPP⁺ level in CSS12 mice was significantly above (1.3-fold, t = 3.22, p = 0.004, 95% CI 1.1–1.5) that in C57BL/6 mice, but was below that of A/J mice, which had 1.9-fold (t = 7.65, p = 3 × 10⁻⁵, 95% CI 1.6–2.2) increase in HPP⁺ relative to C57BL/6 (Fig. 2C). These results suggest that the differential strain susceptibility to HIT could result from differences in Abcb5-regulated drug transport, and that an A/J Chromosome 12 on an otherwise C57BL/6 genetic background increases brain HPP⁺ levels. The correlation between the relative brain abundance of HPP⁺ in A/J, C57BL/6, and CSS12 mice and their haloperidol-induced latencies (Fig. 2A) is consistent with the other studies [33,35,36] indicating that HPP⁺ is causative of HIT.

Desorption Electrospray Ionization-Mass Spectrometry Imaging of Haloperidol Abundance in Brain Sections

In DESI-MSI, thin tissue sections are bombarded with solvent microdroplets to dissolve hundreds of chemicals from the sample surface; the secondary microdroplets formed enter a mass spectrometer, and are analyzed for their chemical content. A movable stage enables multiple pixels (in a row by column matrix) encompassing the entire thin section to be sampled, which generates a detailed 2-D chemical map of the distribution of molecules within the sample surface. The relative amount of haloperidol in each region was analyzed by this method to produce a 2-D ion image of its relative abundance throughout the entire brain section. We analyzed independently prepared sagittal brain sections (two per mouse) from A/J and C57BL/6 mice (three per strain) after treatment with haloperidol (10 mg/kg/d IP) for 5 days. This analysis revealed that haloperidol was diffusely distributed throughout all brain regions, but the amount of haloperidol in A/J brain tissue was consistently greater (~6-fold) than in C57BL/6 brain tissue (Fig. 3). Because of its low abundance, we could not detect the oxidative metabolite by this method. The inter-strain differences indicated by DESI-MSI were greater than those identified by quantitative MS analysis of extracted whole brain. This difference could be attributed to the fact that DESI-MSI was performed on tissue sections in a specific sagittal plane, which does not account for the abundance of haloperidol in other brain regions that were sampled by quantitative MS analysis of extracted whole brain. Nevertheless, direct analysis of brain tissue and DESI-MS imaging both clearly show that haloperidol was more abundant in A/J mouse brain.

Abcb5 mRNA Is Expressed near Brain Capillaries and by Purkinje Cells

Since we did not know if Abcb5 is expressed in brain, in situ hybridization studies were performed on tissue sections from mouse brain. The hybridizations were performed using an anti-sense probe that hybridizes to Abcb5 mRNA, and with an Abcb5 sense probe that serves as a negative control. Abcb5 mRNA expression was observed in the perivascular region of brain capillaries (Fig. 4A). Since this is the anatomic location of the blood-brain barrier function, this expression pattern is consistent with Abcb5 having a role in regulating brain haloperidol levels. Of importance, Purkinje cells in the cerebellum also expressed Abcb5 mRNA; and these cells also expressed dopamine receptor D2 (Drd2) mRNA, which encodes the receptor that binds haloperidol (Fig. 4B). Interestingly, Abcb5 mRNA was not detected in the substantia nigra (S3 Fig.). A low level of Abcb5 mRNA expression was also detected in regions of the hippocampus and frontal cortex (S4 Fig.), but not by the Drd2-positive cells in the substantia nigra (S3 Fig.).

Abcb5 Knockout Mice Have a Prolonged Haloperidol-Induced Latency and Increased Brain Haloperidol Levels

To more precisely determine whether Abcb5 plays a role in HIT, we measured the latencies in homozygous Abcb5 knockout mice (on a C57BL/6 background) before and after treatment with haloperidol (10 mg/kg/day IP). Before drug treatment, there was no difference in the latency of Abcb5 knockout mice relative to wild-type littermates (t = 1.54, p = 0.185) (Fig. 5A). Although there were only marginal differences on treatment day 4 (t = 2.15, p = 0.05), the latencies in Abcb5 knockout mice were significantly prolonged after 5 (10.1-fold, t = 3.09, p-value = 0.03, 95% CI 1.5–69.2), 6 (13.6-fold, t = 4.79, p-value = 0.005, 95% CI 3.5–55.1), and 7 (50.3-fold, t = 7.91, p-value = 0.0005, 95% CI 14.1–179.6) days of haloperidol treatment relative to wild-type mice (Fig. 5). In a second independently performed experiment, haloperidol treatment induced a markedly prolonged latency in Abcb5 knockout mice (Fig. 5B). In contrast, the latencies measured in wild-type littermates did not increase with haloperidol treatment and were comparable to those of C57BL/6 mice. Furthermore, DESI-MSI revealed that the amount of haloperidol in brain tissue obtained from Abcb5 knockout mice (n = 3, 6,282 ± 2,279) was 6.8-fold higher than in wild-type littermates (n = 3, 935 ± 200, p = 0.02) (Fig. 5C). The prolonged latency and increased brain haloperidol levels observed in Abcb5 knockout mice confirms that Abcb5 affects susceptibility to HIT by increasing brain haloperidol levels.

Human ABCB5 Alleles Associate with Haloperidol Toxicity

A human cohort with 85 patients of European descent was treated with haloperidol as a monotherapy for one month during the acute phase of their psychotic illness, and their response to therapy and the appearance of HIT was assessed at baseline and after 1, 3, 7, 14, and 21 days of treatment as described [24]. Their clinical and demographic characteristics are described in Table 2. We investigated whether ABCB5 SNP alleles were associated with the appearance of the three types of symptoms of HIT (Parkinsonoid, Dyskinesia, and Akathisia scores) that were measured in this cohort. Of the 29 SNPs analyzed: six were codon changing (cSNPs), three were in the 5′ UTR, and the others were intronic (Table 3). We first examined the effect that treatment duration had on the incidence of each of the three types of HIT toxicities throughout the course of treatment (S6 Table). The incidence of akathesia and parkinsonoid features increased through day 14, while dyskinesia decreased after day 3 (Fig. 6A), which may be due to fact that the onset of the other toxicities masked the dyskinesia. Therefore, we analyzed a composite index of all three symptoms of HIT, since our mouse genetic analysis identified a genetic factor that alters the brain haloperidol concentration, which will affect susceptibility to all of these symptoms. In examining the allelic associations, we found a clear time-dependent allelic effect on HIT that was maximal on days 1, 3, and 7 (Figs. 6B and S5). Therefore, it is likely that the protective allelic effect would dissipate after prolonged dosing owing to haloperidol (or metabolite) accumulation, but individuals with a protective ABCB5 allele would exhibit less toxicity during the early treatment period (i.e., days 3–7).

Since each of these three different toxicities could all be caused by an increased drug (or drug metabolite) concentration in the brain, we examined whether ABCB5 alleles were associated with a combined toxicity index, which indicated whether a patient manifested any of the three measured toxicities on a given day. To reduce the number of hypotheses that were tested simultaneously, the area under the curve (AUC) for the cumulative toxicity measurements on days 0 to 7 was calculated, and its association with each of the 29 SNPs was tested. (A truncated dataset was used because many of the values for day 14 and 21 were not available for many patients.) Remarkably, a SNP (rs17143212) that was a missense SNP (Thr131Ile) had a particularly strong association (raw p = 0.002) with the toxicity index, and a total of ten (out of 29) SNPs showed a plausible association (with raw p < 0.05) with the toxicity index (Table 2). Among these SNPs, four were consecutive SNPs within the second identified region of linkage disequilibrium in the ABCB5 gene (Fig. 6C). The probability that ten SNPs would randomly have these p-values (< 0.05) is <8 × 10⁻⁷, indicating that it is highly unlikely that all of these associations were false positives resulting from a random match. Since the AUC values didn’t follow a Gaussian distribution and the sample size may not be large enough, the resulting p-values from the ANOVA models could be biased. Therefore, we performed a permutation test, which accounts for multiple testing and makes no assumption about the data distribution, to rigorously assess the statistical significance of each association. To do this, the combined toxicity value for the 85 patients was permutated, while their genotypes remained unchanged; this enabled the linkage disequilibrium among the SNPs to be retained. The distribution of the minimum p-value for the 29 SNPs was evaluated by performing this permutation experiment 1,000 times. The observed raw p-value for each SNP was then compared to this distribution to obtain the true p-value for that SNP under the null hypothesis that this region was not associated with HIT. Strikingly, rs17143212, which caused a significant (131Thr → ILE) amino acid change in the predicted protein, still exhibited a significant association (permutation test p-value = 0.034) after the permutation test. The four consecutive SNPs that had significant raw p-values didn’t have significant permutation test p-values. However, if additional patients were available for analysis, we would be able to better assess the significance of their associations. Nevertheless, at least one human ABCB5 SNP shows a significant association (after correction for multiple testing) with a combined measure of HIT.

Discussion

This study identifies murine Abcb5 alleles as pharmacogenetic factors affecting susceptibility to HIT. Because Abcb5 is expressed in brain capillaries, it can function as part of the blood-brain barrier to regulate the level of drugs in the brain [37]. This function is consistent with the measured differences in the brain levels of haloperidol and its metabolite (HPP⁺) in inbred and chromosome substitution strains. These findings suggest a new model for the genetic mechanisms underlying inter-individual and cell type-specific differences in susceptibility to haloperidol-induced CNS toxicity (Fig. 7).

After haloperidol enters the brain (because of its lipophilicity or other intrinsic properties), the brain concentrations of haloperidol and its toxic metabolite will be affected by Abcb5 function at the brain-blood barrier. Thus, a strain with an Abcb5 allele with reduced transporter function is more susceptible to HIT because of the increased haloperidol concentration in the brain. Like the causative toxins (MPTP, rotenone) used in the other animal models of toxin-induced Parkinson’s disease [38–41], the oxidative metabolite of haloperidol (HPP⁺) is a potent inhibitor of NADH-linked mitochondrial respiration [36]. Brain enzymes convert the Parkinsonian neurotoxin MPTP to its toxic metabolite [42–44], which is concentrated in cells with a high affinity dopamine uptake system [41,45], which causes the Parkinsonian-like symptoms to develop. HPP+ is a more potent inhibitor of mitochondrial electron transport than the oxidative metabolite of MPTP [36]. The toxicity of paraquat, which is another neurotoxin that selectively damages dopaminergic neurons, is also dependent upon dopamine transporter activity [46]. Thus, we propose that haloperidol in Drd2-expressing dopaminergic neurons is converted to its oxidative metabolite HPP⁺ by the drug metabolizing enzymes present in these neurons [47,48]; and the mitochondrial toxicity induced by this metabolite produces the Parkinsonian-like symptoms.

Of interest, we found that Abcb5 and Drd2 were also expressed in cerebellar Purkinje cells. Given their large size, high metabolic requirements, and Drd2 expression, Purkinje cells should also be affected by HIT, since they are particularly susceptible to toxin-induced injury [49]. However, the symptoms of HIT observed in mice and humans are predominantly extrapyramidal and do not have a substantial cerebellar component. Why does HIT produce Parkinsonian symptoms, but not cerebellar toxicity? We believe that Abcb5 expression by Purkinje cells protects them from HIT by their ability to transport haloperidol out of the cells. In contrast, haloperidol-targeted neurons in the substantia nigra, which lack Abcb5 expression, will be susceptible to HIT (Fig. 7).

Although a relatively small cohort of patients who were treated with haloperidol was evaluated here, a human ABCB5 SNP, which causes a significant amino acid change in the protein, was also associated with susceptibility to HIT. However, the human allelic affect is exerted during the early period of drug treatment. The time-dependent effect observed in humans may be caused by the fact that continued haloperidol treatment may overwhelm the ability of the protective allele to clear the toxic metabolite. Despite some evidence for familial clustering [50] and scattered reports (from studies that evaluated a small number of patients) of candidate gene associations [51,52], we did not have any solid information about human pharmacogenetic factors for HIT. A recently completed, large human GWAS did not reveal any significant associations with antipsychotic drug-induced toxicities [53]. However, this study did not perform the detailed analysis of HIT that was performed here. Our results suggest that an examination of the time dependence of individual and the aggregation of cumulative toxicity measurements may be required to identify pharmacogenetic factors for antipsychotic agents. Like the murine gene, human ABCB5 is also highly polymorphic; alleles at ten codon changing SNPs (cSNPs) are predicted to affect human ABCB5 transporter function, including two (Iso648Thre, Glu679Lys) within the nucleotide-binding region [54]. So, it is certainly possible that analysis of a larger cohort could identify other human ABCB5 alleles that also affect ABCB5 function, and possibly susceptibility to drug toxicity. It is also possible that analysis of a larger human cohort could reveal SNP alleles in other genes that also affect susceptibility to HIT. Although haloperidol is a first-generation antipsychotic agent, the more commonly used second-generation antipsychotics (aripiprazole, clozapine, olanzapine, risperidone) are also associated with a significant incidence of EPS-related toxicities, which is higher than was originally appreciated [55]. It would be of interest to determine whether genetic variation in ABCB5, or other efflux transporters, could affect susceptibility to their toxicity. Simulation studies have convincingly demonstrated that using pharmacogenetic information for drug selection [56] or dose adjustment [57] could reduce drug toxicity. Characterization of genetic factors affecting the toxicity caused by haloperidol (or second-generation antipsychotics) could lead to pharmacogenetic strategies that would improve the outcome for patients taking these drugs.

Lastly, analysis of chromosomal substitution strains indicated that the presence of A/J alleles on Chromosome 12 could partially induce susceptibility to HIT. However, the magnitude of the haloperidol-induced latency and the brain HPP⁺ levels in CSS12 mice were significantly below those of A/J mice. Thus, other genetic factors could also affect susceptibility to HIT in mice, which raises the possibility that there are other human pharmacogenetic factors to be found.

Conclusions

These findings indicate that Abcb5 is a component of the blood-brain barrier in mice and suggest that genetic variants in this gene underlie, at least in part, the differences in susceptibility to haloperidol-induced toxicity among inbred mice strains. Moreover, the human genetic association study indicates that a specific ABCB5 allele also affects the susceptibility of people to haloperidol-induced toxicity.

Supporting Information

Attachment 1

Attachment 2

Attachment 3

Attachment 4

Attachment 5

Attachment 6

Attachment 7

Attachment 8

Attachment 9

Attachment 10

Attachment 11

Attachment 12

Attachment 13