Mechanistic Modeling Reveals the Critical Knowledge Gaps in Bile Acid-Mediated DILI.

Bile salt export pump (BSEP) inhibition has been proposed to be an important mechanism for drug-induced liver injury (DILI). Modeling can prioritize knowledge gaps concerning bile acid (BA) homeostasis and thus help guide experimentation. A submodel of BA homeostasis in rats and humans was constructed within DILIsym, a mechanistic model of DILI. In vivo experiments in rats with glibenclamide were conducted, and data from these experiments were used to validate the model. The behavior of DILIsym was analyzed in the presence of a simulated theoretical BSEP inhibitor. BSEP inhibition in humans is predicted to increase liver concentrations of conjugated chenodeoxycholic acid (CDCA) and sulfate-conjugated lithocholic acid (LCA) while the concentration of other liver BAs remains constant or decreases. On the basis of a sensitivity analysis, the most important unknowns are the level of BSEP expression, the amount of intestinal synthesis of LCA, and the magnitude of farnesoid-X nuclear receptor (FXR)-mediated regulation.

Introduction

The capacity of drugs to inhibit the bile salt export pump (BSEP) has been correlated with clinical drug–induced liver injury (DILI).[1] Furthermore, known inhibitors of BSEP and other bile acid (BA) transporters have been implicated in the case studies of DILI. Glibenclamide has been associated with lobular liver injury in the clinic,[2,3] bosentan has been associated with many cases of hepatotoxicity and carries a black-box warning for hepatotoxicity,[4] and troglitazone was withdrawn from the market due to idiosyncratic DILI.[5,6] Currently marketed BSEP inhibitors, such as lapatinib, have been given black-box warnings due to issues with liver toxicity.[7] Interestingly, BSEP inhibitors often have additional proposed mechanisms for hepatotoxicity,[8,9,10,11] and this can often confound our understanding of how these drugs cause DILI. Understanding the effects of BSEP inhibitors on BA metabolism and processing will prove crucial in understanding the importance of BSEP inhibition in clinical DILI.

BAs are molecules that aid digestion of triglycerides and other fat-soluble nutrients. They are synthesized in the liver and undergo highly efficient enterohepatic recirculation.[12,13] While much is understood qualitatively about the basics of BA homeostasis, there are still many aspects of BA transport about which little quantitative information exists. For example, regulation of transporter activity by the farnesoid-X receptor (FXR) and other nuclear receptors generally responsible for BA synthesis has been demonstrated,[14,15,16,17] and chenodeoxycholic acid (CDCA) and its amide conjugates have been shown to decrease BA synthesis[18] and trigger FXR.[19] However, to our knowledge, FXR response and transporter activity have not been linked quantitatively to BA concentrations within hepatocytes. Secondary BAs, such as lithocholic acid (LCA) and deoxycholic acid (DCA), are synthesized in the gut by intestinal bacteria,[20] but to our knowledge, the rates of synthesis and the population variability of this synthesis have not been quantified in any species. Transporter kinetic and inhibition studies have used (or focused on) taurocholic acid, a metabolically stable model BA, but similar studies rarely have been performed on other BAs. Furthermore, while certain BA species, such as LCA and CDCA, are more cytotoxic than other BAs,[20,21,22,23,24] the effect that intrahepatic concentrations of these BAs have on hepatocytes remains to be elucidated. While BSEP transporter expression variability in humans has been quantified,[25] and underlying diseases have been shown to affect transporter expression,[26] the significance of this variability has not been explored fully.

The relative importance of these unknowns remains in question. The process of constructing a mathematical model of BA-induced DILI can be helpful in determining which data gaps should be filled first for maximum impact. While any mechanistic model of BA-induced DILI will be lacking in quantitative accuracy due to the significant data gaps, mechanistic modeling can be used to qualitatively describe elements of the system that would be more or less likely to contribute to BA-mediated DILI.

In this study, we use DILIsym (www.dilisym.com), a mechanistic mathematical model of DILI, to explore the nature of BA homeostasis and its disruption by BA transporter inhibitors. DILIsym has been used previously to assess the efficacy of different N-acetylcysteine treatment regimens in acetaminophen-overdosed patients[27] and to understand the species differences in methapyrilene toxicity using in vitro metabolism data.[28] For this study, we have constructed a module of BA homeostasis within DILIsym for rats and humans; information on the data used to construct this module is available in the Supplementary Materials online.

We have also collected novel in vivo data from the rats treated with the BSEP inhibitor glibenclamide and used these data to ensure the reasonability of the DILIsym model. We also compared the human representation in DILIsym to serum BA data in humans after acetaminophen overdose. This model was used to explore the changes in intrahepatic BA concentrations that are predicted to occur in humans when a BA transporter inhibitor is administered. Finally, a sample population was used to determine which system parameters have the most influence on the magnitude of the intrahepatic BA changes that occur after dosing with a BA transporter inhibitor and would thus be the system components that would be most useful to be illuminated by future experiments.

Results

The BSEP and Na+-taurocholate-cotransporting polypeptide inhibitor glibenclamide, when administered orally to rats, causes an increase in total plasma BAs; this increase in BA concentrations is associated with high systemic exposure to glibenclamide, but is not observed at lower exposure levels. shows the results from the single-dose studies. There was no relationship between dose and exposure. Glibenclamide is notoriously poorly absorbed[29]; as such, the systemic exposure for each individual rat is more important than the dose administered. At concentrations below a Cmax of 10,000 ng/ml, there was little correspondence between glibenclamide concentrations and serum BA fold change, suggesting that the variability in BA processing among rats exceeds the effects of transporter inhibition at these exposure levels. Rats with glibenclamide Cmax >10,000 ng/ml generally displayed a higher concentration of BAs, albeit over a wide range. As shown in , the simulated population (SimPops) covers the range of BA increases observed in the experimental rat data at each maximum simulated concentration tested, validating our model.

summarizes the results of the multiple-dose (300, 750, and 1,500 mg/kg) studies where individual BA species were measured. Serum data from day 1 are plotted in ; each data point represents an individual rat. Serum CDCA and LCA concentrations increased very soon after glibenclamide dosing with a return to the baseline predosing value by the end of 24 h; serum CDCA amide conjugates did not appear to change after glibenclamide dosing.

On day 7, plasma concentrations of glibenclamide did not vary widely in animals treated with the three different dose levels (data not shown); as a result, a single simulation was conducted on a rat SimPops at a simulated drug–exposure level similar to that observed in the rats in vivo. The results of this simulation were compared to the experimental data for LCA, CDCA, and CDCA–amide conjugates in . The simulated population showed similar dynamics to those observed in vivo; the simulated rats exhibited a mild increase in BAs between 1 and 2 h after dosing and returned to baseline values by 24 h. For unconjugated CDCA, the simulated rat population covered the range of observed responses; for LCA, the simulated responses covered the bulk of the experimental data but missed the highest BA concentrations, while for CDCA–amide, our simulated population covered a much wider range than the experimental data, likely due to a wider range of initial CDCA–amide concentrations in the rat SimPops. demonstrates that DILIsym provides a good approximation of the range of BA increases that were observed after administration of a transporter inhibitor.

As data detailing the effect of BA transporter inhibitors on BA concentrations in humans were scarce, the model was validated in humans by comparing the response of the system after acetaminophen dosing to serum BA data from the clinic.[30] Acetaminophen is not a BSEP inhibitor, so in this case, BA elevations are due to hepatocyte necrosis. The results of this validation are shown in . Our baseline human model had a similar response to that demonstrated by James et al., that is, little to no BA elevations in serum when injury was minimal, small but significant increases when there was moderate injury, and large elevations when the injury was severe. The small human SimPops also exhibited a similar range of BA concentration changes to the patients reported in the James study.

With the model structure validated using the rat and human experimental data, we proceeded to simulate a BSEP inhibitor in our baseline human model. The difference in liver BA area under the curves (AUCs) over the course of simulated drug administration in the presence and absence of the theoretical inhibitor is given in . In all zones of the liver, increases were observed in the AUC of bulk BAs, CDCA–amide conjugates, and LCA–sulfate conjugates. Meanwhile, the AUC of unconjugated CDCA and unconjugated LCA actually decreased, while the AUC of amide-conjugated LCA remained about the same.

The changes in BA concentrations with inhibitor as a function of the zonality of the liver were also examined using DILIsym (). In the periportal region, where most BAs are concentrated, the results were similar to those of the overall liver. However, in the centrilobular region, all BA species except for unconjugated LCA increased significantly. For CDCA–amide and LCA–sulfate conjugates, the increases in AUC were far greater than those predicted in the periportal region. However, the absolute centrilobular AUC remained lower than the absolute periportal AUC for all BA species.

The large model human population was used to understand the sensitivity of the model parameters to BA accumulation. The most important variables for the accumulation of each BA species in simulated human liver are shown in . Since the baseline human model simulations suggested that CDCA–amide and LCA–sulfate were the most likely BA species to accumulate in the liver, we assessed the correlation between changes in system parameters and increases in CDCA–amide and LCA–sulfate concentrations. The sensitivity analysis suggested that BSEP expression, as represented by CDCA–amide canalicular transport Vmax, was the most important variable for CDCA–amide accumulation, while the magnitude of LCA synthesis in the gut was the most important variable for LCA–sulfate accumulation.

The difference in hepatic CDCA–amide accumulation when CDCA–amide canalicular efflux Vmax changes compared to when the canalicular transport regulation scaling factor changes (which roughly corresponds to the magnitude of the response of BA transporter expression to FXR activation), as seen in , is particularly instructive. When CDCA–amide canalicular efflux Vmax was varied, the baseline CDCA–amide concentration of the individual changed; however, the percentage increase in CDCA–amide concentration after administration of the inhibitor was approximately constant across the variable range. In contrast, when the canalicular regulation scaling factor was varied, the baseline concentration of CDCA–amide remained the same (although meal fluctuations were more pronounced with less regulation), while the percentage increase in CDCA–amide concentration was much higher when canalicular transport regulation was lower.

Discussion

The novel in vivo data included herein demonstrate the importance of considering population variability and individual characteristics in any model of BA homeostasis and its disruption. The variability between rats shrouded the detection of changes in BA due to low concentrations of glibenclamide (Cmax: 500–10,000 ng/ml) and the resultant transporter inhibition. The lack of an apparent exposure/response effect in the individual BA data is due to the population variability. Even with higher doses of glibenclamide, when an exposure/response relationship was observed, we found that there was a wide range of responses, as measured by a fold increase in BA concentrations, among the rats with glibenclamide Cmax higher than 10,000 ng/ml. Few rats reached the higher drug-exposure levels; further experiments, perhaps using a more readily absorbed formulation of glibenclamide, could help illuminate this effect. The variability in the results also makes model validation difficult, as our ability to reproduce the experimental conditions (including individual rats' initial plasma and liver BA concentrations) within the context of a simulation is challenging; while our model is reasonably validated, further experimentation in this area will continue to improve our representation.

Our results suggest that due to the population variability inherent in the BA system, a baseline model is of limited utility for predicting the response of rats (or patients) to BA transporter inhibitors. A baseline model is useful for gleaning interesting information about the basic behavior of the system, but quantitative predictions of BA increases and decreases should not be based on the behavior of a baseline BA model. Constructing an appropriate population for the human model is, therefore, of the utmost importance for the prediction of BA-mediated DILI. To that end, we used the BA model to identify the most important parameters that are likely to vary among individuals. These parameters should be included in the creation of a reasonable population for the prediction of BA toxicity.

The model suggests that the end-product BA conjugates are the BA species that are most likely to accumulate in the liver. This is likely due to the fact that the conjugation of BAs is fairly rapid in the liver, and these pathways do not easily saturate. Previous research has placed the Km for BA amidation in the 1 mmol/l range,[31] which is several orders of magnitude higher than typical in vivo unconjugated BA concentrations. As a result of the predicted accumulation of BA conjugates in the liver, BA amide and sulfate conjugates are potential liver toxicophores; the toxicity of these molecules should be explored and quantified with future experimentation.

The model also suggests that an interesting zonal effect exists in BA accumulation. Centrilobular cells, which generally have a lower BA concentration than their periportal counterparts,[32] have a much higher predicted relative accumulation of all BA species. In fact, the AUC of unconjugated BAs was predicted to increase in the centrilobular zone, even though the same BAs' AUC decreased in the periportal zone. Despite the differential BA changes in the two zones, the absolute AUC of BAs in the centrilobular zone was predicted to remain lower than the periportal zone. However, if centrilobular cells are more sensitive to BA toxicity than periportal cells, an increase in BA concentrations could lead to zonal centrilobular toxicity. Indeed, some case reports of liver injury due to glibenclamide demonstrate centrilobular injury.[2,3] The investigation of the differential effects of BAs on centrilobular and periportal liver cell toxicity would be another interesting and worthwhile experiment.

Simulations using a model population provided further insight into the most valuable experiments that could be performed to improve our understanding of BA-mediated DILI. Among the many system unknowns varied within the model population, the most important variables included LCA synthesis in the gut and the regulatory response to BA accumulation. Therefore, quantifying LCA synthesis in the gut in addition to changes in BSEP expression and activity with increasing hepatocyte BA concentrations would be important and impactful experiments. The importance of BSEP expression levels to BA accumulation in hepatocytes, for example, could be quantified using existing animal or in vitro models.

The modeling also revealed the importance of transporter expression levels on the behavior of the system; however, it is important to note a key difference between the behavior of the system when these expression levels are varied and its behavior when the regulatory response is varied. When transporter expression levels—especially canalicular transporter levels—are decreased, the baseline BA concentration in the liver increases; the percentage increase in BAs after dosing the inhibitor remains somewhat constant. Varying the magnitude of the regulatory response (FXR activation), by contrast, does not have a large effect on baseline BA concentrations in the liver; however, the percentage increase in BAs after inhibitor dosing is highly variable. Both transporter expression levels and the magnitude of the regulatory response are expected to vary within a population; however, only the variability in transporter expression levels may be apparent in untreated individuals. This result suggests both the promise and the limitations of prescreening individuals for certain risk factors before prescribing known BSEP inhibitors.

Methods

The BA submodel in DILIsym comprises several components. It includes (i) the synthesis and metabolism of BAs in hepatocytes, (ii) the basolateral and canalicular active transport of BAs, (iii) the release of BAs from the gallbladder in humans, (iv) the synthesis of secondary BAs and deconjugation of BAs in the gut, (v) the recirculation of BAs from the gut and subsequent active uptake by the liver, and (vi) the regulatory effects of BAs on transporter expression and BA synthesis. The model contains representations of LCA and CDCA and its conjugates, the BA species most frequently linked to toxicity in in vitro experiments,[20,21,33,34] and a “bulk” BA representation that contains the other BAs. A more detailed description of the BA homeostasis model and its parameterization is presented in the Supplementary Materials online.

To represent BA dynamics in humans, the model was optimized to the known profile of BAs in serum published by Trottier et al.,[35] as well as to overall concentrations of BAs in the liver measured by Setchell[36] and García-Cañaveras et al.[37] There were 47 system variables that were either unknown or expected to vary within the human population that were fit to the BA profile; a list of the variables used in the optimization is provided in the Supplementary Materials online. The optimization of the model based on the published BA profile was performed in a manner similar to the SimPops method outlined in previous publications, wherein a genetic algorithm was used to generate values of variable parameters that lead to model outcomes (in this case, BA profiles) that are within the range of experimental data.[27,28] Using this optimization method allowed the selection of a population of 2,400 individuals with reasonable baseline BA concentrations for the large human sample population. A smaller, 10-parameter, 331-individual SimPops was also constructed for the purpose of model validation. While the large population was given artificially wide parameter ranges for the purpose of the sensitivity analysis, the small human sample population was constructed with more constrained parameter ranges for the purpose of approximating a plausible population of humans. For example, the transporter Vmax ranges included four orders of magnitude in the large population, but were constrained to ranges suggested by transporter expression profiles from Meier et al.[25] in the small population.

A smaller SimPops was constructed in rats for the validation against the experimental data. A population of 191 rats using a more limited 11 variable set (Supplementary Table S4 online) was generated. This population was intended to represent rats with both plausible serum BA concentrations and plausible values for the parameters that were varied. The baseline rat model was optimized to the BA profile from the control rats in the present experiment, and the SimPops was constructed around this baseline.

Details on the simulations performed for the model validation and for the model exploration can be found in the Supplementary Materials online. Multivariate analysis on the population sample was performed using JMP 9 from SAS (Cary, NC).

Multiple-dose glibenclamide study in rats

For the BA-profiling experiments, 16 male 8- to 9-week-old CD-1 rats from Charles River Laboratories (Raleigh, NC) weighing between 200 and 300 g were randomized into four groups of four animals each. These groups were administered daily doses of either the vehicle control (0.5% hydroxypropylmethylcellulose /0.1% Tween 80 in water) or glibenclamide (Sigma-Aldrich, St Louis, MO) in vehicle via oral gavage for 7 days at three dose levels (300, 750, and 1,500 mg/kg/day). Details on the treatment of the rats are available in the Supplementary materials online.

Rats underwent a viability check twice per day, and detailed clinical observations were taken at least twice over the course of the study. Blood was drawn for serum BA profiles at 1, 3, 6, and 24 h after dosing on day 1 and day 7, and glibenclamide concentrations were measured for toxicokinetic analysis using the same blood samples from day 7. The animals were euthanized under isoflourane anesthesia on day 8 and underwent necropsy. Clinical and anatomic pathology data were collected from the animals, and these results are reported in the Supplementary Materials online. The serum BA concentrations from day 1 were used for comparison to the simulation results.

BA profiling in serum was performed using liquid chromatography–tandem mass spectrometry analysis at GlaxoSmithKline (Ware, UK). BAs in liver tissue 24 h after the final dose also were profiled; results from this analysis, as well as the liquid chromatography–tandem mass spectrometry analytical method, are presented in the Supplemental Materials online. This study was conducted in accordance with the GlaxoSmithKline Policy on the Care, Welfare and Treatment of Laboratory Animals and was reviewed by the Institutional Animal Care and Use Committee.

For this experiment, doses of glibenclamide were administered every 24 h for 7 days. Systemic exposure data were compared to simulation results for day 7; day 1 BA data were compared to simulated BA concentrations on day 1.

Short-term glibenclamide studies in rats

For the short-term studies, male Han Wistar rats (substrain AlpkHsdBrlHan:WIST; AstraZeneca, Macclesfield, UK) of 10–12-week age (300–400 g) were used. Details on the treatment of the rats used in this study are available in the Supplemental Materials online. All animals were treated in accordance with approved UK Home Office license requirements.

Two experiments were performed in which total plasma BA concentrations and glibenclamide plasma concentrations were determined. Glibenclamide (Sigma-Aldrich) was formulated as a solution or suspension in hydroxypropyl-β-cyclodextrin (Acros Organics, distributed by Fisher Scientific, Loughborough, UK) in aqueous 0.2 mol/l Na2CO3/NaHCO3 buffer (pH 10). In the first experiment (29 animals), four groups of five animals received a single dose via oral gavage of 50, 250, and 500 mg/kg glibenclamide or 10% (w/v) hydroxypropyl-β-cyclodextrin vehicle alone, and blood samples were taken at 1, 6, and 24 h after dosing; an additional three animals per glibenclamide-treated group received a single dose of 50, 250, and 500 mg/kg glibenclamide, and blood samples were taken at 1, 3, 6, 12, and 24 h after dosing. In the second experiment (20 animals), two groups of five animals each received two oral doses of 250 and 500 mg/kg glibenclamide in 20% (w/v) hydroxypropyl-β-cyclodextrin vehicle at 0 and 4 h, and blood samples were taken predose and at 0.5, 1, 2, and 4 h after dosing; blood samples from an additional five animals per glibenclamide-treated group were taken predose and 1 h after dosing to increase the dataset for this time point. Details on the blood-sampling procedure are located in the Supplementary Materials online. Animals were dosed 2 h into the light cycle.

Author Contributions

J.L.W. wrote the manuscript. J.L.W., K.Y., S.Q.S., S.H.S., J.L.A., K.L.R.B., P.B.W., and B.A.H. designed the research. J.L.W. and K.Y. performed the research. J.L.W., K.Y., S.Q.S., S.H.S., J.L.A., D.B., and B.A.H. analyzed the data. D.B. contributed new analytical tools.

Conflict of Interest

The authors declared no conflict of interest.

Study Highlights

Table 1

Bile acid changes in a model of human bile acid inhibition in the periportal and centrilobular zones of the liver

Table 2

Correlation coefficients for the 10 most correlated system parameters to hepatic CDCA–amide accumulation and LCA–sulfate accumulation