Are Physiologically Based Pharmacokinetic Models Reporting the Right Cmax? Central Venous Versus Peripheral Sampling Site

Physiologically based pharmacokinetic (PBPK) models can over-predict maximum plasma concentrations (Cmax) following intravenous administration. A proposed explanation is that invariably PBPK models report the concentration in the central venous compartment, rather than the site where the samples are drawn. The purpose of this study was to identify and validate potential corrective models based on anatomy and physiology governing the blood supply at the site of sampling and incorporate them into a PBPK platform. Four models were developed and scrutinised for their corrective potential. All assumed the peripheral sampling site concentration could be described by contributions from surrounding tissues and utilised tissue-specific concentration-time profiles reported from the full-PBPK model within the Simcyp Simulator. Predicted concentrations for the peripheral site were compared to the observed Cmax. The models results were compared to clinical data for 15 studies over seven compounds (alprazolam, imipramine, metoprolol, midazolam, omeprazole, rosiglitazone and theophylline). The final model utilised tissue concentrations from adipose, skin, muscle and a contribution from artery. Predicted Cmax values considering the central venous compartment can over-predict the observed values up to 10-fold whereas the new sampling site predictions were within 2-fold of observed values. The model was particularly relevant for studies where traditional PBPK models over-predict early time point concentrations. A successful corrective model for Cmax prediction has been developed, subject to further validation. These models can be enrolled as built-up modules into PBPK platforms and potentially account for factors that may affect the initial mixing of the blood at the site of sampling. Electronic supplementary material The online version of this article (doi:10.1208/s12248-015-9796-7) contains supplementary material, which is available to authorized users.


INTRODUCTION
It has been observed that physiologically based pharmacokinetic (PBPK) models can often report over-predicted maximum plasma concentration (C max ) values compared to in vivo data. This is most obvious when administration is via the intravenous (i.v.) route over a short time period, for example, a fast infusion or bolus. It is suggested that this discrepancy, apart from practical limitations in drawing blood samples immediately after giving an i.v. dose, may be because PBPK models traditionally report the concentration in the central venous compartment, i.e. the Bpooled venous return^, whereas in clinical studies, sampling is usually taken from a peripheral vein in the arm due to ethical and logistical considerations (1). These two sites may not have the same concentrations, particularly at early time points where distribution throughout the body tissues is ongoing, and equilibration has not occurred.
With PBPK models increasingly used in assessments of pharmacokinetics, pharmacodynamics and toxicology, this is an issue that cannot be dismissed. PBPK approaches have been used in toxicology for the development of exposure limits (1) and to obtain predicted C max values to be utilised in the assessment of potential adverse effects (2). These models are also increasingly used in pharmaceutical regulatory submissions (3)(4)(5)(6). For modelling of pharmacodynamics, if arterial concentrations are considered to drive an observed effect, it is important to understand how this concentration relates to that observed from the clinical venous sampling. It is clearly critical to know what is being reported from PBPK models and how the underlying models relate to the in vivo situation.
The idea that over-prediction of C max , when reporting from the central venous compartment, is indicative of physiological differences between the central venous compartment and the peripheral sampling site is not new and not only theoretical; evidence to support the concept of local differences in distribution, mixing and concentration between potential sampling sites in vivo has been published. Concentration differences can be observed between arterial and venous sampling and can also be dependent on the location of the Electronic supplementary material The online version of this article (doi:10.1208/s12248-015-9796-7) contains supplementary material, which is available to authorized users. specific vessels within the body as reviewed by Chiou (7,8) as early as 1989. This review highlighted numerous examples where the choice of sampling site can lead to significant differences in the observed concentrations, both in human and in animal studies. The author also highlighted that utilising the Bwrong^observed concentration could have potentially serious implications.
Differences between the arterial and venous concentrations and potential models for inter-conversion of reported data have been described by a number of investigators (9)(10)(11)(12)(13)(14). Many of these publications focus on the implications of these concentration differences for modelling of pharmacodynamics, and these models are generally empirical or compartmental in structure and are not easily applied to PBPK models. In addition, some published models require prior knowledge of in vivo concentration-time or pharmacokinetic data, and/or parameter fitting, meaning that they could not be applied in a predictive scenario such as the estimation of first-in-man pharmacokinetics. A physiologically based recirculatory model for fentanyl was constructed and presented by Upton et al. (15). However, a number of the parameters were determined using prior knowledge of the compound distribution and hence are not applicable in a generic Bbottom-up^PBPK model. Levitt (16) provides the most pragmatic method to address the problem, by using an approach of creating a PBPK model for the arm to describe a peripheral sampling site. This model was developed to describe and predict the kinetics of organic solvents and may not be relevant for a wider variety of pharmaceutical compounds; although it does provide a useful starting point for model development.
The objective of this study was to develop and evaluate a potential corrective model to describe the C max at a peripheral sampling site, utilising the full PBPK model within the Simcyp population-based Simulator. A number of different options were assessed and compared to observed in vivo data. The model giving the best performance was selected for further validation and future development.

Model Definitions
The base model in this study follows the assumption that the tissues surrounding a peripheral sampling site contribute, to varying degrees, to the observed venous concentration at that location. A similar approach was suggested by Levitt (16), where Bflow fractions^were used to describe the contribution of each tissue to the final overall concentration at the sampling site. These flow fractions were determined using a fitting approach where the values were adjusted within the PBPK model to recover the observed experimental data. The model described herein utilises the tissue concentration profiles predicted by the Simcyp Simulator full PBPK model (17), which allows simulation of tissue concentration data for all organs incorporated within the model. These tissue profiles, which are commonly reported with respect to blood, in combination with the tissue:plasma partition coefficient (Kps), can be used to calculate the emergent concentrations which leave the tissues. These can then be used in combination with drug data and the defined fractions to describe an overall concentration-time profile at a peripheral sampling site where a corrected C max value can be observed. The fractions determined by Levitt (16) provided the initial estimates that were incorporated in the model for validation.
The tissue concentrations selected for incorporation in the model were skin, adipose, muscle and a Bshunt^compartment. These tissues are logical selections since they can be found in the vicinity of the antecubital vein in the forearm, the most common clinical sampling site. The Bshunt^compartment describes the arterio-venous anastomoses in the skin of the hand, where there is a direct connection of arterioles to the venules (18). The model published by Levitt includes these tissues; however, they describe the shunt as a new tissue and state Bit is assumed the tissue supplied by the Shunt has the same pharmacokinetic properties as skin^, whereas for one of the models validated in this study a direct arterial contribution has been considered in combination with the skin contribution. In addition, the Levitt model also incorporates an Bother^compartment, representing the subcutaneous space and connective tissue, which was not available in the Simulator model. Assuming the summation of all fractions in the model is equal to a total of 1; different options were explored to assess the implications of assigning this fraction to components available in the current PBPK model. A schematic of the model is shown in Fig. 1.
Four different model options were evaluated, three of which were variations of the fractions provided in the publication by Levitt (16). In the following equations, adipose/muscle/skin/lung conc refers to the tissue concentration observed in the simulation results and (t) refers to each time point where a predicted concentration has been output for the simulation. As the arterial concentration is not a standard output from the simulator, the lung concentration was converted, by dividing by Kp Lung /BP and used as a surrogate for the purposes of this exercise (Eq. 7). Kp is the tissue:plasma partition coefficient for the specified tissue which can be predicted using either the corrected Poulin and Theil (19)(20)(21), or Rodgers and Rowland (22)(23)(24) methods or measured data can be used directly. BP is the blood to plasma ratio required for the conversion of tissue concentration data to the concentration of drug exiting each tissue and entering the venous plasma at the peripheral sampling site.

A. Modified Levitt arm (LAC) model
In this model, the fractions were based on those reported in the Levitt publication with the 10% assigned to 'other' being split equally between the skin, adipose, muscle and shunt: Peripheral Site Concentration t ð Þ with respect to plasma The Levitt model assigns the same distribution parameters and volume/blood flow ratio to the shunt as the skin and therefore the skin concentration was initially used as a surrogate for the shunt: B. Modified Levitt Bskin^(LAS) model In this instance, the 10% of the arm concentration assigned to 'other' by Levitt was assigned completely to the shunt (as skin concentration): Peripheral Site Concentration t ð Þ with respect to plasma The shunt equation is the same in this model as described in the LAC model. C. Modified Levitt arm Barterial^(LAA) model The lack of arterial contribution in the shunt for the LAC and LAS model was not considered to be physiologically representative, so in the third model, the fraction for Bother^was assigned to the shunt but used arterial concentration, as shown in Eq. 6. This 10% was removed from the skin fraction that was representing the shunt in the LAS model. Physiological arm concentration (PAC) model Literature searches were used to identify relevant physiological data that could describe the tissue concentration fractions required for this model. These fractions were derived from reported relative blood flows in the forearm (25)(26)(27)(28)(29). Unfortunately, these data are limited and no description of the contribution of the anastomoses to the antecubital vein blood flow could be obtained. The described fractions for this model therefore only refer to the adipose, muscle and skin concentrations: Since the sampling volume is usually small in relation to the total blood volume, it was assumed to be negligible in this case, and therefore mass balance is maintained in the underlying PBPK model. Calculations of the peripheral site profile were performed in Microsoft Excel 2010 (Microsoft Corporation, Redmond, WA, USA) using the mean concentration-time profiles for each tissue. Following the initial validation steps of these models, the arterial contributing fraction in the LAA was varied from 0.02 to 0.625, with skin fraction in the shunt being increased/decreased to maintain a total shunt fraction of 0.625, to assess the impact on C max prediction and potentially further refine the model.

Compound Data and Simulations
Literature searches for drugs with an available Simcyp Simulator compound dataset were conducted to identify clinical studies where the dose was administered by i.v. infusion, and in vivo concentration-time data were available. Data extracted from the publications included the reported steady-state volume of distribution (Vss), in addition to the in vivo concentration-time data and study details. Where concentration-time data were only presented in figures, the data were extracted using the GetData graph digitiser version 2.22 (Get Data Graph Digitizer, 2012, http://getdata-graph-digitizer.com/) for comparative purposes. Where Vss was not reported for a specific study, this was obtained by non-compartmental analysis using Phoenix version 1.3 (Pharsight), utilising the extracted concentration-time profile data. Physicochemical and elimination data already present in the Simcyp compound file were utilised.
It was assumed that reasonable prediction of Vss can be used as an indication of the correctness of the predicted Kp values used in the PBPK model. However, clearly this is necessary but may not guarantee that all Kp values are correct. Since Vss is a description of the distribution under steady-state conditions, it should theoretically be unaffected by early distribution differences. For the purposes of these simulations, Kp values and the overall Vss were initially predicted using BMethod 2^within the simulator. In some cases, the predicted Vss did not match the observed values; therefore, all tissue Kp values were scaled up or down to recover the observed Vss prior to the simulations being performed. Suitable values for the universal Kp scalar were determined using the inbuilt sensitivity analysis or parameter estimation tools within the platform.
Simulations were performed using the Simcyp Simulator (Version 13 release 1) (Simcyp, Sheffield, UK), with the provided healthy volunteer population and concentration time points equally divided over 2000 time points for the length of the simulation period. Study design and dosing regimen were matched to the published trial design for all studies.

Analysis and Comparisons
Visual comparisons of reported in vivo and predicted concentration-time profiles were made. Predicted C max values were compared with observed values. As the time point of the highest concentration (T max ) in simulations can differ from the reported T max , predicted and observed concentrations were compared at the in vivo T max . This is due to the fact that the simulation reports 2000 time points whereas in clinical studies, a finite number of samples can be taken. The fold difference between these observed and predicted values for C max was calculated as a measure of prediction accuracy, with the number of results between 0.8-and 1.25-fold or within ±2-fold recorded. Predicted C max values within 2-fold of reported in vivo values were considered to be reasonable predictions. The prediction accuracy of all the models described above was compared for selection of the most suitable peripheral site model.

Compound Data
A total of 15 studies with i.v. administration via infusion were identified from the published literature, covering seven compounds, where the relevant compound data were available in the Simcyp Simulator. The compounds investigated were alprazolam, imipramine, metoprolol, midazolam, omeprazole, rosiglitazone and theophylline (n=4, 3, 1, 2, 2, 1 and 2 sets of observed data, respectively). The studies investigated are summarised in Table I.

Comparison of Models
A visual comparison was made of the early-stage concentration-time profiles from published studies with predicted profiles using the different models described above. Selected key profiles are shown in Fig. 2, with results for all studies available in the supplementary material. Clear improvements in predictions can be seen for a number of studies, metoprolol and omeprazole in particular, with the LAA simulated C max generally closest to the observed C max . However, it can be noted from the profiles (see supplementary material) that in a number of cases, the simulated profiles for all of the models have converged at an earlier time point than an observation has been made in the in vivo studies. This is true even when a difference was predicted between the initial phase of the peripheral sampling site and central venous profiles and makes validation of the models difficult. The fold difference between the observed and predicted concentrations at the observed T max for all models is shown in Table II, with summary statistics in Table III and a graphical summary of the prediction accuracy for the LAA is given in Fig. 3. The LAA is considered to be the most accurate model, with all of the predicted concentrations within 2-fold of the observed values. This comparison confirms the utility of the corrective model, particularly for those compounds identified previously where early sampling time points are available. The final model selected based on these results was the LAA, which incorporates the arterial fraction in the shunt.

Systematic Alteration of the Arterial Contribution
Changing the fractional contribution of the arterial concentration to the shunt resulted in an improvement for some compounds, bringing the C max predictions within a 0.8-to 1.25-fold range, while other compounds moved outside the 2-fold range (Table IV). A C max prediction within 2-fold of the in vivo value was observed for all studies when using the LAA with a fraction of 0.1 for the arterial contribution. Figure 4 gives a graphical representation of the consequences of changing this fraction, though due to extremely poor results, the arterial fraction of 0.625 is excluded from this figure. These observations suggest that the original fraction of 0.1 for arterial contribution was preferable and should be maintained in the final model. However, further validation with a larger observed dataset is required.
The final selected fractions for use in the LAA are shown in Table V, and the final equations are shown below.
Peripheral Site Concentration t ð Þ with respect to plasma

DISCUSSION
A systematic validation and comparison of a number of methods to improve the prediction of C max from PBPK models has been performed and reported. Models that were evaluated were based on the full PBPK model in the Simcyp Simulator (Version 13 Release 1) and include standard tissue concentration outputs to describe the concentration-time profile for venous plasma at a peripheral sampling site, i.e. the antecubital vein. While the consideration was made at the inception of this model that construction of a PBPK model or compartment for arm might be a possible solution, the model selected allows the use of the existing physiologically based model as implemented in the Simulator and the associated outputs. This model provides an elegant solution where mass balance is maintained in the original model, and comparisons can be easily drawn between the concentrations in the central venous compartment and at the new peripheral sampling site. It is important to note that all of the validation is performed based on peripheral sampling data. The perfect scenario to validate this model would be the use of both peripheral and central venous or arterial sampling data for comparative purposes; however, such data are experimentally difficult and ethically questionable to obtain in humans.
Based on the current validation, the LAA model was considered to give the best performance of all of the models investigated, with the predicted C max closest to the observed values, and within the specified 2-fold range for all studies. This 2-fold range is a commonly used criterion for assessing the accuracy of predictions in in vitro-in vivo extrapolation and has recently been subject to some debate, with additional models for success criteria being proposed, though there is no consensus on the most appropriate method for validation (30). The LAA model is also considered to be the most physiologically accurate, where the shunt considers a contribution from the arterial blood in line with the role of the arterio-venous anastomoses in the skin of the hand, this is consistent with descriptions of the circulatory physiology of the skin (18). This contributing factor is absent from all other models investigated in this study (LAC, LAS and PAC). The lack of physiological data, particularly relating to the shunt, that could be identified in the literature when building the alternative PAC model is disappointing. This factor potentially contributes to the poor performance of this model in comparison to the LAA. The final variations of the arterial fractions made some improvement to the selected model for specific compounds, but led to a worsening prediction for others. Overall, it was considered that a model giving a reasonable prediction for all studies was preferred to improving the prediction of one or two compounds at the expense of accuracy for others, potentially skewing the model. Usually, over-fitting of a model for one or two cases can reduce its generalisation ability. In addition, the final value selected is in line with observed data cited by Chiou (7) where arterio-venous shunting in dog leg is observed to be approximately 4% of blood flow, which gives confidence that the final value of 10% is physiologically reasonable in the absence of any observed human data.
An assumption of the current approach is that accurate prediction of Vss means the Kp values are correct, so a drawback to the model is that it relies on the accurate prediction of Vss values within the underlying PBPK model. On this basis and using the in vivo Vss values, a fitting strategy was used to ensure the distribution volume was reasonably described (by an accurate Vss prediction) before proceeding to predict the peripheral sampling site concentration. Scaling of Kp values was required for all compounds apart from metoprolol. However, if in vivo data are lacking to validate the Vss or Kp values or if a bottomup approach is intended, the confidence in the prediction of distribution cannot be validated. If animal data are available, these can be utilised within the model to ensure the distribution is reasonably described or the physicochemical properties of the compound can be used to guide the selection of a predictive model ( (31)(32)(33), personal communication from Iain Gardner).
In the initial stages of this project, the use of different fractions in the models for compounds with different Vss values or physicochemical parameters was considered. This was ultimately dismissed, as the reason for the differences in the central venous and the peripheral sampling site concentrations depends on both the compound and physiological parameters, and therefore, the model structure should not be compound dependent. However, compound-specific differences may be evident in the concentration-time profiles, as lipophilic drugs will be expected to distribute more extensively into tissues leading to bigger site-specific concentration differences at early times after drug administration (7). These compound-specific differences are already incorporated in the underlying PBPK model, with lipophilicity and protein-binding considered as contributing factors for the prediction of Kp values. The resultant differences can be observed in the profiles provided in Fig. 2 and in the supplementary material, e.g. theophylline with low logP and Vss values has minimal differences between the concentration profiles for the central venous compartment and the peripheral sampling site model, whereas  imipramine with both high Vss and logP shows significant differences. For other compounds, the contributions of these factors are less clear. The lack of reported early time points in the observed data for these compounds limits the conclusions that can be drawn as the convergence of predicted concentration profiles prior to the observed data hinders the determination of which model fits best in these circumstances. It is worth noting that omeprazole does not have a particularly high logP or Vss, and significant differences are noted in the profiles-suggesting other factors may also be contributing to this phenomenon (e.g. charge/ compound type). With such a limited clinical dataset, it is challenging to explore the potential causes further and draw any strong conclusions. It should also be noted that no consideration was given to the accuracy of clearance values in this study, and the elimination phase of each simulated profile was not compared to the in vivo situation. The assumption is made that for low and moderate clearance compounds, there will be little impact on the C max value, meaning the simulated results remain valid for all of the compounds that were investigated in this study. The selection of the predicted concentration at the time point of the highest observed concentration for the calculation of fold differences is a key factor in understanding the results. Comparisons can be drawn between the highest observed and predicted concentrations, but these could be potentially meaningless as the T max reported from the simulation can differ substantially from that observed. Over the course of each simulation in this study, 2000 data points are reported whereas the clinical/observed data is usually limited to around 10-20 samples over several hours. This can be particularly problematic for comparisons when the first sampling time point is much later than the end of infusion and/or the T max observed in the simulation, as there is no way of knowing what is happening in vivo where there are no reported observations. While all compounds showed an improvement in terms of fold error in C max prediction when using the LAA model, it is obvious from the comparison of the observed and simulated concentration-time profiles that the utility of the model can only be fully tested when concentration data for early sampling time points are available. This is in agreement with Chiou (7) and Shankaran (1) who highlighted the differences occurring in the initial phases of the concentration-time profiles. These time points are expected to be within the initial distribution phase (e.g. 5-10 min after administration), and samples taken at >30 min after administration are too late to observe the   (Fig. 2), it seems that only three studies, for metoprolol and omeprazole, had early enough time points to fully illustrate the benefit of using the corrective model and in most of these cases, only one point  is relevant. As can be seen in Fig. 2 and also in the supplementary material, most studies have clearly reached an equilibration point before the first in vivo sampling time and the predicted profiles have converged. The differences between the models are mainly limited to the initial distribution phase, which is in line with the observations made by others in this arena (1,7). A few exceptions arise in this dataset; in study 4 for alprazolam and study 7 for imipramine, the first sampling point is not until 30 min after infusion, yet a significant correction using the peripheral sampling site model was still observed at this time, suggesting that in some cases, distribution differences can persist for longer than the initial 5-10 min.
Another consideration is variability, as the model has been compared only to the highest concentration point extracted digitally from the mean concentration-time profile for each dataset. C max may not be at the same time for all individuals, hence biasing the C max of the mean concentration profile. Out of the 15 datasets studied, only two (34, 35) had tabulated C max and T max data for the i.v. study with a measure of variability (range or standard deviation) available. Variability in the reported C max values was <30%, and the reported mean C max (taken as the mean of individual C max values) was similar to the C max value determined from the mean concentration profile (12.6 vs 12.6 ng/mL and 146 vs 115 ng/mL for alprazolam and rosiglitazone, respectively). In addition, the numbers of subjects considered in the studies are relatively low (n values range from 1 to 14). Both of these issues could lead to uncertainty in the validation of the model. This is due to a paucity of data; in an ideal scenario, individual data for a large number of subjects would be available, ensuring that variability within the wider population is fully captured with bestand worst-case scenarios fully represented. This is a general problem with validating models and does not only apply to this study, but it becomes more pertinent when dealing with such a limited number of datasets. The Simulator reports concentration-time points for each simulated individual within a virtual trial for plasma and individual tissues, so if clinical data were available, further validation in this area would be feasible.
Expansion of the underlying models could also be considered in the future to include the effect of certain drugs on the cardiovascular system or the impact of the application of heat on the distribution at the sampling site. An example of this is metoprolol, which alters heart rate, cardiac output and consequently tissue blood flows (36,37). If metoprolol distribution is perfusion-limited then it could be affected by the action of the drug on the blood flows. Similarly, the arterial contribution could be increased to mimic the effect of a temperature increase in the vicinity of the sampling site. While in the current study, it is assumed, in the absence of further information, that the conditions of a standard temperature and a sedentary subject are met, the use of a heat-box to induce an Barterialised concentration^at a venous site have been documented previously (38,39). However, the development of such models is outside the scope of the current study.
The initial constraints of this investigation include the availability of both a library compound within the Simulator and an in vivo study utilising an infusion administration, the validation was therefore limited by the sampling points in those studies. Further i.v. infusion studies with relevant sampling time points have been identified for compounds not currently available in the Simcyp Simulator. This extended validation requires the development of additional compound files which was outside the scope of this initial study. As proof of concept, the improvement shown in this report seems to be adequate and allows a basis for further validation as well as investigation into potential applications within a number of different arenas.

CONCLUSION
A peripheral sampling site model that improves the prediction of C max values has been developed and validated in this study. The tissue fractions that constitute the sampling site concentration were tested for seven compounds and provide a good concordance between the observed and predicted C max values. Further validation is required with more compounds and clinical studies, particularly those with early sampling time points, to confirm the utility of the model. However, the potential for the use of such a model has been adequately illustrated herein. It is envisaged these models can be enrolled as built-up modules within PBPK platforms to address potential affections and variability of the initial mixing of the blood at the site of sampling (such as effect of heat, variation in adipose and muscle content of the body), and hence, they may give a more realistic simulated drug concentrations in blood or plasma at early time points.

ACKNOWLEDGMENTS
The authors would like to thank Eleanor Savill and Jessica Waite for their help in preparing the manuscript.