Abstract
The routine assessment of xenobiotic in vivo kinetic behavior is currently dependent upon data obtained through animal experimentation, although in vitro surrogates for determining key absorption, distribution, metabolism, and elimination properties are available. Here we present a unique, generic, physiologically based pharmacokinetic (PBPK) model and demonstrate its application to the estimation of rat plasma pharmacokinetics, following intravenous dosing, from in vitro data alone. The model was parameterized through an optimization process, using a training set of in vivo data taken from the literature and validated using a separate test set of in vivo discovery compound data. On average, the vertical divergence of the predicted plasma concentrations from the observed data, on a semilog concentration-time plot, was approximately 0.5 log unit. Around 70% of all the predicted values of a standardized measure of area under the concentration-time curve (AUC) were within 3-fold of the observed values, as were over 90% of the training set t1/2 predictions and 60% of those for the test set; however, there was a tendency to overpredict t1/2 for the test set compounds. The capability of the model to rank compounds according to a given criterion was also assessed: of the 25% of the test set compounds ranked by the model as having the largest values for AUC, 61% were correctly identified. These validation results lead us to conclude that the generic PBPK model is potentially a powerful and cost-effective tool for predicting the mammalian pharmacokinetics of a wide range of organic compounds, from readily available in vitro inputs only.
Physiologically based pharmacokinetic (PBPK) models are mathematical descriptions of the flow of blood throughout the body, developed for the simulation of xenobiotic absorption, distribution, and elimination. The essential concepts were outlined over 60 years ago in a farsighted paper (Teorell, 1937) that presented many of the mathematical relationships required to simulate blood flow and tissue distribution.
Simulation modeling ideas were developed further by Mapleson (1973), to explain the effect of anesthetics, and early attempts to apply the approach to drugs were published in the 1960s by Bellman et al. (1961). Probably the most important contributions in that period were made by Bischoff and Dedrick (1968), who demonstrated that PBPK models could be used for the a priori prediction of the pharmacokinetics of thiopental. During the following decades, developments were made by academics such as Rowland (Rowland, 1986), Sugiyama (Sugiyama and Ito, 1998), and Amidon (Yu and Amidon, 1999), as well as scientists working in the environmental health field, in particular Anderson and Clewell (Andersen et al., 2002). Recent reviews (Grass and Sinko, 2002; Leahy, 2003) have discussed the application of these approaches to the prediction of pharmacokinetics in drug discovery.
It is of interest to us to apply the PBPK approach to the estimation of plasma levels in animals from in vitro data alone. If this can be achieved, with sufficient confidence in the outcome, a generic model for each species, including humans, could be used to address a wide range of problems for which the estimation of pharmacological effects, safety margins, and exposure limits is required. In particular, we believe that the exploitation of data from the growing number of validated in vitro toxicology assays requires a sound approach to the estimation of plasma and tissue levels within the whole animal. Since it is also our goal that predictive methods should be applicable to a wide variety of chemicals prior to any animal experiment, we are interested in establishing generic PBPK models, which are parameterized for the physiology of the animal, independently of any specific xenobiotic. The properties of the xenobiotic that determine its overall kinetic properties could then be determined through the use of separate in vitro surrogates of the important absorption, distribution, metabolism, and elimination (ADME) properties. Coupled with in vitro efficacy or toxicity data, simulation of the overall biological effect in vivo could become routine.
This paper describes the work we have done to parameterize a generic PBPK model for the rat, and to assess the reliability of the model in estimating plasma levels of xenobiotics, where these data are available from experimentation. To date, no other truly generic PBPK model has apparently been published, although there are numerous examples of compound-specific PBPK models for the rat, such as those for simulating the pharmacokinetics of tolbutamide (Sugita et al., 1982), diazepam (Igari et al., 1983), β-lactam antibiotics (Tsuji et al., 1983), cyclosporin (Bernareggi and Rowland, 1991), and a homologous series of barbiturates (Blakey et al., 1997). These and many other such models rely upon data derived from in vivo studies to estimate the extent of tissue distribution. Although computational models for predicting in vivo tissue/plasma partition coefficients from physicochemical and biochemical properties have been developed (Yokogawa et al., 1990, 2002; Poulin and Theil, 2000; Rodgers, 2003), these are not generally applicable, being appropriate only for a particular class of basic compounds (Rodgers, 2003), for bases in general (Yokogawa et al., 1990, 2002), or requiring modification of the fundamental model, depending upon distribution characteristics that cannot be established in advance with any certainty for novel compounds (Poulin and Theil, 2000). Furthermore, many of the published PBPK models for the rat also utilize estimates of clearance determined in vivo (Tsuji et al., 1983; Bernareggi and Rowland, 1991; Blakey et al., 1997). The PBPK model presented herein predicts distribution and elimination kinetics in the rat from in vitro data alone, in a manner that is independent of ionization status and does not require a priori knowledge of in vivo pharmacokinetic properties.
In this paper, we have concentrated on predicting the in vivo pharmacokinetics of compounds for which plasma levels have been determined following an intravenous dose. Work that we have done to extend the model to predict plasma levels following an oral dose will be reported separately.
Materials and Methods
Model Inputs. A generic PBPK model, which enables prediction of the pharmacokinetic behavior of any given compound dosed intravenously in a specified rat population, without recourse to data generated through in vivo studies, is described herein. The only compound-dependent inputs required by this model are: molecular weight; the octanol/water partition coefficient (logP); the octanol/water distribution coefficient at pH 7.4 (logD7.4); all pKa values that affect the ionization status at pH 7.4; the hepatic intrinsic metabolic clearance (CLint); and the fraction unbound in plasma (fup).
Model Description. The PBPK model is based upon that published by Bernareggi and Rowland (1991), as shown in their Fig. 1, but with substantial modification of the tissue distribution (Fig. 1) and elimination (Fig. 2) components, and comprises a series of compartments representing 14 major organs and tissues in the body, interconnected by further compartments representing arterial and venous blood pools, according to the principles developed by Bischoff and others (Bischoff, 1975).
Administration of xenobiotics is via intravenous infusion, and their transport between the compartments that represent the organs and tissues occurs exclusively via blood flow. Distribution into the organ and tissue compartments is “diffusion limited” (Fig. 1). Thus, the cellular membrane represents a diffusion barrier, and movement between blood and tissue is modeled dynamically, rather than assuming that compound in the effluent blood of the tissue is in equilibrium with that within the tissue (Lutz et al., 1980), as in the “flow-limited” PBPK model of Bernareggi and Rowland (1991). The tissue is viewed as consisting of three, well stirred subcompartments (capillary bed, interstitial fluid, and intracellular space), so that restricted diffusion across the capillary wall (for example, for modeling the blood-brain barrier) can also be represented. This arrangement is the minimum required to model adequately all possible restrictions to passive transport between blood and tissue. The possibility that movement across the plasma membrane is restricted to the un-ionized free compound is also accommodated in the model.
Elimination occurs from the compartment representing the liver, as in the model of Bernareggi and Rowland (1991), but also from that representing the kidneys. Hepatic metabolism is modeled as a first-order process, the rate of which is determined by CLint and the unbound compound concentration at the metabolic site. Renal excretion is represented by a physiologically based model (Fig. 2) that has been developed from the work of Komiya and coworkers (Komiya, 1986, 1987) and Katayama et al. (1990). According to this model, free compound within the plasma is filtered, within the glomerulus, into the lumen of the renal tubule; plus there may be additional active secretion of compound into the renal tubule lumen. Active renal secretion is characterized by an equation of the Michaelis-Menten form, where the maximum rate of tubular secretion and the Michaelis constant govern the rate of secretion of unbound compound into the tubular fluid. Both ionized and un-ionized species within the lumen of the renal tubule may be reabsorbed through solvent drag; however, only the un-ionized fraction is subject to reabsorption through passive diffusion (Komiya, 1986). Inherent assumptions are that the permeability of the renal tubule to both water and solutes remains constant along its length (Komiya, 1986), and the renal tubule lumen and kidney tissue/capillary bed behave as well stirred compartments, so that any compound reabsorbed from more distal portions of the tubule is made available for active secretion into the proximal tubule (consistent with the morphology of the nephron and its vasculature) (Katayama et al., 1990). Any compound within the renal tubule that escapes reabsorption is excreted unto the urine.
Model Parameters. The physiological parameters used in the model were obtained from the literature and are given in the Appendix (which is available online as supplemental data); these were scaled according to the actual body weights of the animals used in the in vivo studies being simulated. Tissue and organ volumes were derived from two comprehensive compilations of physiological data for use in pharmacokinetic models (ILSI, 1994; Brown et al., 1997) and represent the extravascular (combined interstitial fluid and intracellular space subcompartments) volumes only. Blood flow rates for stomach, gut, pancreas, spleen, and hepatic artery (and hence total liver blood flow) were obtained from a number of sources (Sasaki and Wagner, 1971; Malik et al., 1976; Nishiyama et al., 1976; Brown et al., 1997); all other blood flows were from Bernareggi and Rowland (1991). The glomerular filtration rate and urine flow rate were taken from Davies and Morris (1993), and the renal tubular lumen volume was derived from a textbook of physiology (Pitts, 1974). A hematocrit of 0.503 (Altman and Dittmer, 1971) was assumed.
Parameterization of the distribution and elimination components of the generic PBPK model for rat required the development of a number of correlation models. Two such models, for the prediction of parameters corresponding to the effective in vivo lipophilicity and plasma protein binding, were derived through a process of optimization of the performance of the PBPK model, in terms of the estimation of experimental plasma levels. A comprehensive training set of in vivo data was used for this purpose. Both the optimization process and the training set data are described in greater detail below.
Instantaneous equilibration of the compounds between erythrocytes and plasma and between plasma and the interstitial fluid was assumed; i.e., there was no effective barrier to the transfer of compounds between the capillary bed and interstitial fluid model subcompartments, and PvSv (Fig. 1) was set to a nonlimiting value. The parameter PmSm (Fig. 1) was optimized across the training set and fixed at the same value for all compounds.
Intracellular space/interstitial fluid (unbound) partition coefficients (KpIC:IFu) were calculated for each compound using a model developed by Yokogawa et al. (1990, 2002), which predicts tissue partition coefficients from a measure of lipophilicity and appropriate estimates of tissue lipid content. In this work, an estimate of the effective in vivo lipophilicity was used in the calculation, rather than a measure of actual lipophilicity (i.e., an experimental or calculated logP). Similarly, an estimate of the effective in vivo plasma protein binding was used as an input to the PBPK model, in place of a measured value for fup. The effective lipophilicity and fup for each training set compound were simultaneously determined through a grid search for the pair of values that generated the best estimate of experimental plasma levels. Regression analysis between the effective lipophilicity values for the training set compounds and the available input variables (listed above) showed that effective lipophilicities are best described as a simple linear combination of logP and the logarithm (to base 10) of the fraction un-ionized at pH 7.4. The best predictive model for the effective protein binding was found by linear regression, using a linearized protein binding parameter (fupLin) and the logarithm of the fraction un-ionized at pH 7 to predict the logarithm of the effective fup values. The parameter fupLin is calculated from the in vitro fup as follows:
An estimate of the blood cell/plasma partition coefficient was derived from a correlation between a set of measured blood cell/plasma partition coefficients and the variables logP and fupLin. This estimate was utilized along with the hematocrit and in vitro fup to calculate the blood/plasma concentration ratio (R).
The in vitro microsomal CLint, scaled to ml/min/g liver, formed a direct input to the PBPK model as an estimate of the in vivo CLint. Binding to albumin and lipoproteins also occurs in the interstitial fluid, and the fraction unbound in this subcompartment (fut) was calculated from the effective in vivo fup, according to a methodology proposed by Poulin and Theil (2000). To simplify the model, specific binding to intracellular components was not considered.
The parameters describing passive reabsorption from the renal tubule (the permeability-surface area product) and reabsorption through solvent drag (the reflection coefficient) were calculated for each compound from logP and molecular weight, respectively. Although the PBPK model makes provision for representing active renal secretion, this capability was not used.
For the results presented here, stochastic simulations were performed to incorporate known variability in the animal body weights or imprecision in the values of the physicochemical parameters, microsomal CLint and in vitro fup. For every set of input data (corresponding to dosing of a single compound in a specified animal population), multiple iterations of the simulation were performed, with the values of the input parameters at each iteration being sampled randomly from the assumed distributions, these being either uniform or normal. Thus, variability was generated from the ranges or the means and standard deviations of the input data.
Training Dataset. The set of in vivo data used in training the model comprised 214 instances of data (where an instance corresponds to a single plasma concentration-time profile) for 82 different compounds, and was derived from numerous published studies of intravenous dosing in rat. The compounds in this training set were drawn from many therapeutic areas and were selected for diversity in compound parameters. No attempt was made to avoid or eliminate compounds with nonlinear pharmacokinetics or those subject to processes not explicitly modeled, such as active transport.
For the majority of the training set compounds, the required model inputs were determined at Cyprotex. Values for logP, logD7.4, and pKa were obtained primarily using a Sirius GLpKa apparatus (Sirius Analytical Instruments Ltd., East Sussex, UK), but in a small number of cases were generated using predictive software (ACD/PhysChem Batch, version 7.10; Advanced Chemistry Development, Inc., Toronto, ON, Canada). Microsomal CLint and fup were determined in vitro through incubation with hepatic microsomes and equilibrium dialysis, respectively. These data were supplemented by values from the literature where these were available.
Test Dataset. To objectively evaluate the performance of the model, an independent set of in vivo test data was constructed. This consisted of 194 instances of plasma concentration-time data for 134 discovery compounds dosed intravenously in rat. These data were supplied by several separate pharmaceutical/biotechnology companies, having been generated through their internal pharmacokinetic studies, and the results presented here were derived in the absence of any prior knowledge of the in vivo pharmacokinetics or the chemical structures.
The test set compounds were varied in terms of physicochemical properties and represented diverse therapeutic areas. Many of the requisite model inputs for these compounds were determined at Cyprotex, although the companies supplying the in vivo data also provided certain of the in vitro data, which were obtained through a variety of methods.
Calculation of the Plasma Concentration Weighted Mean Log -Fold Error (wMLFE). For each pair of in vivo and simulated plasma concentration-time profiles, the log -fold prediction error was determined at each simulated time point for which there were corresponding in vivo data, and the mean of these errors over all time points was calculated, to give an overall mean prediction ratio for each instance of simulated data. The wMLFE represents the weighted mean of these individual means. The weights used in the calculation arise from there being multiple instances and/or sources of in vivo data for several compounds, and hence the contribution of each individual log -fold prediction error to the overall mean is weighted accordingly; i.e., so that each compound contributes equally, whatever the number of instances of in vivo data for that compound.
Principal Components Analysis. The values of seven variables for the training set compounds were transformed by subtracting the mean of the value, and dividing by the standard deviation, so that each transformed variable had a mean of 0 and a standard deviation of 1. The selected variables were: fractional charge at pH 7.4; fraction un-ionized at pH 7.4; logD7.4; logP; fupLin; log10(CLint); and log10(mol. wt.). Principal components analysis (Mardia et al., 1979) was performed on the transformed data, and the scores of the training set compounds were recorded. The first three principal components explained 81% of the variance of the dataset. The variable values for the test set were normalized in the same way, using the means and standard deviations of the values for the training set compounds. The scores of the test set compounds in the coordinates of the principal components of the training set were recorded.
Results
For any given set of input data, output from the PBPK model is in the form of a predicted plasma concentration-time profile. When stochastic simulations are performed for a single set of input data, each iteration generates a predicted profile, and hence the total output consists of a population of profiles that reflect the inherent uncertainty in the input data. Examples of typical simulation results, plotted on the same axes as the analogous in vivo data, are given in Fig. 3 for selected training set compounds and in Fig. 4 for a similar selection of test set compounds.
The simulated profiles in Fig. 3, A and B, and Fig. 4, A and B, illustrate accurate estimation of plasma concentrations over time, for selected training set and test set compounds, respectively. There is little variation within the population of profiles generated for either clozapine (Fig. 3A) or two example compounds from the test set (Fig. 4, A and B), and the fit of the model output to the single set of observed data is almost exact in each case. In contrast, Fig. 3B demonstrates some observed in vivo variability, in this instance for erythromycin, and a combination of predicted variability and uncertainty, generated from both known variability in the animal weights and multiple estimates of in vivo fup and CLint.
Some other simulation results are shown in Fig. 3, C and D, and Fig. 4, C and D. The median predicted profile for pentazocine deviates from the observed profile, but the range of predicted profiles encompasses the in vivo data (Fig. 3C); hence, the model results can still be considered acceptable. Figure 3D indicates somewhat inaccurate estimation of the in vivo tissue distribution of phenytoin, resulting in a tendency to underpredict plasma levels, although there is clearly a high degree of variability in the observed data for this compound, which is reflected in the model output. Figure 4D shows a similar, but more pronounced overestimation of tissue distribution for an example test set compound. However, the elimination half-life has been accurately predicted for the test set compound represented in Fig. 4C.
For the remainder of the results presented here, the median of the population of predicted profiles generated from each set of input data was used as an individual estimate of the plasma concentration-time course. To assess the overall performance of the model in terms of successfully predicting in vivo plasma levels, the plasma concentration wMLFE was determined for both the training and test sets; this statistic corresponds to the mean vertical deviation (in log units) of a simulated data point from a corresponding observed data point on a semilog plot of plasma concentration versus time. The plasma concentration wMLFE values calculated for the training and test sets were 0.46 and 0.53, respectively, and hence the predicted plasma concentrations deviate from the observed values on average by around 0.5 log unit. By way of illustration, the median predicted profiles shown in Fig. 3C and Fig. 4C both have an associated MLFE of approximately 0.5.
The frequency distributions of the actual mean -fold errors in plasma concentration prediction for both sets of compounds are shown in Fig. 5, A and B. A high proportion (72%) of plasma concentration predictions for the training set compounds are on average within a factor of 5 above or below the observed data points, and just a few (9%) are more than 20-fold in error overall (Fig. 5A). Whereas the majority (65%) of the predictions for the test set compounds are again within 5-fold of the observed data on average (Fig. 5B), the mean prediction error is less than 2-fold for a lower percentage of the test set compounds, and more than 20-fold for a greater proportion of these, compared to the training set.
The nature of the test set compounds for which plasma concentrations were poorly predicted by the PBPK model was explored in more detail by means of a descriptive classification model, developed using the RPART method of the R statistical software system (http://cran.us.r-project.org/doc/packages/rpart.pdf), that provides some indication of biochemical and physicochemical commonalities between them (Fig. 6). Compounds demonstrating the characteristics of those that might be expected with a reasonably high degree of confidence to be poorly predicted, according to the classification model, were removed from the test set, including that for which the observed and predicted profiles are plotted in Fig. 4D. Thus, compounds with both logD < 1.95 and fup ≤ 0.065 were eliminated; these criteria define a particular area of property space within which the likelihood of obtaining a poor plasma concentration prediction (mean -fold error greater than 20) is approximately 60%, containing just 13% of the test set compounds overall, but 65% of those for which the predictions were in error, on average, by a factor of more than 20. Outside this area, the likelihood of a plasma concentration prediction being similarly poor is just 5%.
The frequency distribution of the plasma concentration mean -fold prediction errors for the reduced test set (comprising 117 compounds) is shown in Fig. 5C. The error distribution pattern for this reduced set is somewhat different from that of the full set, with 73% of the predictions being within 5-fold of the observed data, and a lower proportion being more than 10-fold in error. The plasma concentration wMLFE for the reduced test set was 0.45 and, therefore, approximately the same as for the training set.
Although the primary outputs from the PBPK model are the predicted in vivo plasma concentration-time profiles that are generated, these can be utilized for estimation of standard pharmacokinetic (PK) parameters of interest, including area under the concentration-time curve (AUC) and elimination half-life (t1/2), allowing direct comparisons to be made with analogous in vivo data. Since different methods of extrapolating AUC from zero time to the first time point and from the last time point to infinity can vary in the estimates they yield, the simulation results were compared to observed data in terms of a standardized parameter, the dose-normalized AUC from the first to the last recorded time points (AUCt1-tlast-DN), as well as t1/2.
The capability of the PBPK model to accurately predict the selected PK parameters in rat has been evaluated in terms of the median values and interquartile (IQ) ranges of the predicted/observed ratios, for both the training and test sets, as given in Table 1. The respective summary data indicate that the prediction of both AUCt1-tlast-DN and t1/2 for the training set compounds is generally successful. The median predicted/observed ratio is close to 1.0 for both parameters, and half of the predictions are, on average, within a range of approximately 0.5 to 1.5 times the observed values, although the statistics do indicate a slight tendency to underpredict AUCt1-tlast-DN. The selected parameters are rather less accurately predicted for the test set compounds, as demonstrated by the greater deviation of the median predicted/observed ratios from a value of 1.0. Furthermore, the skew in the distribution of AUCt1-tlast-DN predictions is more apparent, while a notable overprediction of t1/2 for many of the compounds is also indicated.
The frequency distributions of the predicted/observed ratios of AUCt1-tlast-DN and t1/2, for both the training and test sets, are more clearly demonstrated by the histograms shown in Figs. 7 and 8. The predicted/observed distributions for AUCt1-tlast-DN are similar for both sets of compounds, in both qualitative and quantitative terms (Fig. 7). Around half of all the predictions are within a factor of 2 above or below the observed values, and almost 70% are within a factor of 3. The predicted/observed ratios conform closely to a normal distribution, although a slight tendency to underpredict AUCt1-tlast-DN is evident. Conversely, there is an obvious discrepancy between the training set and test set predicted/observed distribution patterns for t1/2 (Fig. 8). Estimation of this parameter for the training set is highly successful, with over 70% of the predictions being within 2-fold of the observed values, and more than 90% within 3-fold. Furthermore, the predicted/observed ratios are normally distributed. However, only 44% of the predicted test set t½ values are within a factor of 2 of the observed values, although 60% are within a factor of 3. There is also a clear bias toward overprediction of this parameter, with a significant number of predictions being between 2 and 5 times the observed values, and a few even higher.
The performance of a predictive method, according to a given criterion, can also be represented graphically by means of a lift chart (Witten and Frank, 2000), such as that shown in Fig. 9. This plot illustrates the success rate of the PBPK model at selecting the 25% of the test set compounds with the largest values for AUCt1-tlast-DN (and, therefore, with the lowest plasma clearances). The horizontal axis shows the sample size as a percentage of the total number of compounds comprising the test set, while the vertical axis indicates the percentage of compounds correctly selected according to the stated criterion. The diagonal line represents the expected success rate if the selections were made at random, whereas the plot that would be generated by a 100% success rate is indicated by the line on the far left, which reaches a maximum (all selections correct) at a sample size equal to 25% of the total number of compounds. The actual performance of the PBPK model at ranking compounds according to AUCt1-tlast-DN is shown by the central line. The plot demonstrates that of the 25% of the test set compounds ranked by the PBPK model as having the largest values for AUCt1-tlast-DN, 61% have been correctly identified. This is clearly substantially greater than the level of success to be expected (25%) if the compounds were selected at random.
To visualize the multivariate distributions of the training and test set compounds in the input space of the PBPK model, the scores for the first versus second and third versus second principal components were plotted for the training and test sets (Fig. 10, A and B, respectively). It can be seen from Fig. 10B that the training set compounds fall into three reasonably well separated groups. Furthermore, distinct physicochemical attributes can be ascribed to each of these groups. The compounds having their second principal component (PC) score less than zero (PC2 < 0) are strong bases, all having fractional charge at pH 7.4 of 0.56 or greater. The compounds with PC2 > 0 are acids and weaker bases, having fractional charge at pH 7.4 of 0.355 or less (see Fig. 10, A and B). This set of compounds can be subdivided on the value of the third PC. The group of 12 compounds with PC3 > 1 are all strong acids. The other compounds, with PC3 < 1, are bases or weaker acids, with a fractional charge at pH 7.4 greater than –0.545 (Fig. 10B). Thus, the three groups correspond to 1) strong bases, 2) strong acids, and 3) weak/moderate acids and bases. The first PC is dominated by lipophilicity (logP, logD7.4), plasma protein binding, and intrinsic clearance, and displays a continuum of values for the training set (Fig. 10A).
The test compounds show some similarities to, but also significant differences from, the training set scores. The first PC is again mostly a continuum, but with three outliers having PC1 > 4, compared to one for the training set. The distributions of the second and third PC values are, however, clearly different from those of the training set. The distributions for the test set are more continuous, so the separation into the three groups observed in the training set is not so clearly seen in the test set. Furthermore, although there are nascent clusters that correspond to those in the training set, the test set centers, particularly for groups 2 and 3, are displaced relative to those of the training set.
Discussion
The PBPK model is potentially the most powerful tool currently available for predicting in vivo PK, investigating the physiological and chemical interactions that give rise to organism- and compound-dependent PK, and linking with pharmacodynamics to provide an integrated PK/pharmacodynamics view of therapeutic and/or toxic effects of xenobiotics. However, the lack of a generic model that can be applied to a wide range of compounds and/or can utilize inexpensively and reliably determined inputs for any given compound has, thus far, prevented widespread adoption of the PBPK modeling approach for PK prediction. We have developed a generic PBPK model for the rat that overcomes this obstacle and thereby facilitates the routine and cost-effective prediction of PK. In principle, a similar model can be derived for any species for which the relevant physiological data are available; we have described an analogous model for the human in the accompanying paper (Brightman et al., 2006).
Recent developments in both the pharmaceutical industry and those sectors of the manufacturing industry that produce environmental chemicals have increased the need for facile prediction of PK in humans and other species. Within the pharmaceutical industry, the realization that consideration of pharmacokinetic and toxicity criteria during drug discovery will reduce the likelihood of costly failure during development has driven a greater requirement for PK data. Within the manufacturing industry, increasingly stringent international legislation concerning the safety of new and existing compounds necessitates quantification of the putative extent and duration of exposure in a variety of possible scenarios. At variance with the need to generate more PK data is the desire, both for cost and ethical reasons, to reduce animal experimentation where possible. The integration of in vitro ADME data and in silico methods can address this dilemma, with PBPK modeling offering several advantages over alternative in silico or computational methods for the prediction of PK.
The first such advantage is that a PBPK model is likely to be more robust than alternative predictive methods. A large proportion of the variation in the pharmacokinetic behavior of compounds can be explained by the physiological and biochemical processes involved in xenobiotic distribution, metabolism, and elimination that are, at least partially, independent of the properties of the compounds. These processes are represented by the parameters and differential equations of the PBPK model. Thus, for example, the approach of in vivo plasma clearance to a maximal value with increasing hepatic elimination is determined by the limiting effect of hepatic blood flow. Similarly, although the partitioning of compounds from plasma into tissues is compound-dependent, the relative extent of the distribution of lipophilic and hydrophilic compounds is dependent on the volumes of adipose and muscle tissue, respectively, whereas differences in their distribution rates are determined, in part, by the specific perfusion rates of these tissues. Many alternative approaches to predicting pharmacokinetic behavior rely upon explaining all of the variation between compounds through a statistical model, with consequent dependence on the training data composition. Application of such methods can be expected to result in inaccurate predictions when applied to areas of chemistry outside of, or poorly represented within, the training data. As we have illustrated (Fig. 10), there are distinct differences in physicochemical and ADME-related properties between marketed drugs and compounds typical of those in drug discovery. Consequently, a statistical model that has been trained using data for existing drugs will most likely have a limited ability to generalize to typical discovery compounds. Within the generic PBPK model, the use of such models has been restricted to calculation of a small number of internal variables. In addition, these internal models have been built using data from different, though overlapping, sets of compounds. These measures reduce the risk of poor prediction when the PBPK model is applied to novel areas of chemistry. The model validation results described herein demonstrate that the pharmacokinetic properties of the outlying test set compounds (see Fig. 10A) are well predicted, despite these compounds being outside the parameter space of the training set.
A second advantage is that the PBPK model predicts plasma concentration-time profiles, from which all required PK parameters can be calculated using standard formulae. By contrast, other approaches are usually restricted, in practice, to the prediction of one parameter. Consequently, if a number of parameters are to be considered in decision-making, corresponding models must be obtained for each. Assuming reliable models can be derived for the parameter space being studied, the individual parameters predicted by each model (for example, clearance, half-life, and volume of distribution) must be consistent with each other to be of value in decision-making. A PBPK model predicts a consistent set of standard parameters for a given compound, calculated from a single predicted plasma profile. In addition, the simulated plasma or tissue profiles can be used as a basis for determining additional, nonstandard parameters, such as the time following dosing for which the plasma or tissue concentration is greater than a specified value, which can be related to therapeutic and/or toxic thresholds and therefore used to estimate relevant exposure.
One possible role of the PBPK model in drug discovery is to enable drug metabolism/pharmacokinetics scientists to reliably identify those compounds that are likely to have suitable in vivo PK. The model provides two sets of data to this end: the predicted plasma and tissue concentration-time profiles, and PK parameters derived from the plasma profile. Either, or both, can be used in compound selection. If, for instance, the in vivo potency is known or can be estimated, then the tissue or plasma profiles can be used to select compounds that have the longest durations at concentrations greater than that required for therapeutic effect. Alternatively, total exposure, as measured by dose-normalized AUC, may be used as a selection criterion.
The effectiveness of PK parameter prediction can be determined by relative accuracy and/or the ability to rank compounds successfully. Anecdotally, we have observed that drug discovery drug metabolism/pharmacokinetics scientists would value tools enabling predictions to be made within a 3- or sometimes 2-fold margin of error. We have shown that the PBPK model predicts t1/2 and AUCt1-tlast-DN within 3-fold of the observed values for, respectively, 60% and almost 70% of a test set composed of discovery compounds. Assessment of the ranking capability of a predictive method can be performed by means of a lift chart, which can be used to quantitatively predict the likely consequences of making a particular operational decision. We have illustrated (see Fig. 9) an arbitrary scenario in which the identification of compounds with AUCt1-tlast-DN in the top quartile of the test set is required. Using the generic rat PBPK model to select the 25% of the test set with the highest AUCt1-tlast-DN results in 61% of the required compounds being selected, corresponding to an enrichment of more than 2.4-fold compared to random selection. In contrast, making the selection based on the 25% of compounds with the lowest intrinsic clearance results in 39% of the required compounds being selected, an enrichment of less than 1.6-fold (not shown). Specifying different selection criteria enables a series of curves to be generated, from which the potential benefits of the different strategies, in terms of the enrichment of compound selection, can be compared.
We have described our findings concerning the current version of a generic rat PBPK model, but model improvement is ongoing. As for statistical model development, additional compounds can be added to the optimization training set and additional descriptors assessed for their ability to predict the optimization parameters. However, a further advantage of the PBPK model is that its predictive capability can be improved by adding new features that represent additional physiological or biochemical processes, or by improving extant features of the model. The motivation can be to understand specific determinants of observed PK and/or to address shortcomings in the predictions. In the latter case, prioritization can be guided by identified weaknesses in the model. Thus, the observation that hydrophilicity coupled with high plasma protein binding is a strong indicator of poor plasma concentration prediction has guided us to improve the modeling of renal clearance. Additionally, we are developing a model of the blood-brain barrier, initially as a passive permeability barrier. These developments should lead to both significantly more accurate prediction of plasma PK and more realistic prediction of the time course of compound concentration in the brain. Other factors such as nonlinearity, specific transport processes, and population effects, which are not explicitly described in the current model, could also be explored using the PBPK approach.
Developing the PBPK model for the rat has enabled us to verify that it is possible for such a model to generalize reliably outside the property space of the training set. However, the greater value to industry lies in being able to predict PK for human subjects exposed to xenobiotics. The work we have carried out in developing and validating a similar generic PBPK model for humans is described in the accompanying paper (Brightman et al., 2006).
Acknowledgments
We thank Dr. Darwin Cheney for helpful comments on the manuscript.
Footnotes
-
Article, publication date, and citation information can be found at http://dmd.aspetjournals.org.
-
doi:10.1124/dmd.105.004804.
-
ABBREVIATIONS: PBPK, physiologically based pharmacokinetics; ADME, absorption, distribution, metabolism, and elimination; AUC, area under the concentration-time curve; AUCt1-tlast-DN, dose-normalized AUC from the first to the last recorded time points; CLint, hepatic intrinsic metabolic clearance; fup, fraction unbound in plasma; fupLin, linearized plasma protein binding parameter; fut, fraction unbound in the interstitial fluid; IQ, interquartile; KpIC:IFu, intracellular space/interstitial fluid (unbound) partition coefficient; MLFE, mean log -fold error; PC, principal component; PK, pharmacokinetic(s); PmSm, permeability-surface area product for movement across the cell membrane; PvSv, permeability-surface area product for movement across the capillary wall; wMLFE, weighted MLFE.
-
↵ The online version of this article (available at http://dmd.aspetjournals.org) contains supplemental material.
- Received March 23, 2005.
- Accepted October 11, 2005.
- The American Society for Pharmacology and Experimental Therapeutics