Predictive Performance of Next Generation Human Physiologically Based Kinetic (PBK) Model Predictions Based on In Vitro and In Silico Input Data

The goal of the present study was to assess the predictive performance of a generic human physiologically based kinetic (PBK) model based on in vitro and in silico input data and the effect of using different input approaches for chemical parametrization on those predictions. To this purpose, a dataset was created of 38,772 Cmax predictions for 44 compounds by applying different combinations of in vitro and in silico approaches for chemical parameterization and these predicted Cmax values were compared to reported in vivo data. Best results were achieved when the hepatic clearance was parameterized based on in vitro (i.e. hepatocytes or liver S9) measured intrinsic clearance values, the method of Rodgers and Rowland for calculating tissue:plasma partition coefficients, and the method of Lobell and Sivarajah for calculating the fraction unbound in plasma. With these parameters, the median Cmax values of 34 out of the 44 compounds were predicted within 5-fold of the observed Cmax and the Cmax of 19 compounds were predicted within 2-fold. The median Cmax values of 10 compounds were more than 5-fold overestimated. Underestimations (>5-fold) did not occur. A comparison of the current generic PBK model structure with chemical-specific PBK models available in literature was made to identify possible kinetic processes not included in the generic PBK model that might explain the overestimations. Overall, the results provide crucial insights in the predictive performance of PBK models based on in vitro and in silico input and the influence of different input approaches on the model predictions.


Introduction
Physiologically based kinetic (PBK) modelling has a crucial role in next-generation (non-animal) risk evaluations to predict internal human plasma (or tissue) concentrations, and to relate these concentrations to in vitro biological effect concentrations (Blaauboer, 2010;Fabian et al., 2019;Louisse et al., 2017;Punt et al., 2019Punt et al., , 2021Wetmore et al., 2015;Yoon et al., 2015). The ultimate applicability of PBK models in next generation (non-animal) risk evaluations will, however, depend on the possibility to make predictions without the support of animal in vivo kinetic data (e.g. plasma or tissue concentrations) (Paini et al., 2019;Peters and Dolgos, 2019). When available, human kinetic data may be applied for PBK model development as well as data from chemical analogues (OECD, 2021), but such data are generally scarce or unavailable for nonpharmaceuticals. Even though PBK models are increasingly parameterized based on in vitro and in silico input data, in vivo data are currently still needed to assess the validity of the model predictions for a given chemical. Moreover, when certain kinetic processes cannot be parameterized based on in vitro or in silico experiments, they are usually obtained by fitting model predictions to in vivo data (Peters and Dolgos, 2019;Tsamandouras et al., 2015). Given that human kinetic data are difficult to obtain, other strategies to evaluate the adequacy of in vitro-and in silico-based PBK models to estimate in vivo kinetics are needed. Initial estimates of plasma and tissue concentrations of orally consumed compounds can effectively be made with minimal generic PBK models that are defined based on 1) a first order intestinal absorption rate (ka) and a fraction absorbed (fa), 2) intrinsic hepatic clearance (CLint), 3) tissue:plasma partition coefficients, 4) the fraction unbound in plasma (fup), and 5) passive renal excretion (defined as the glomerular filtration rate times the unbound concentration of the compound in plasma (GFR x cPlasma x fup) (Jones and Rowland-Yeo, 2013). Gaining confidence in PBK model predictions requires confidence in the quality of the input data that that are used to parameterize these kinetic processes (Gouliarmou et al., 2018;Louisse et al., 2020;Utsey et al., 2020). In addition, it is important to determine for which compounds a minimal PBK model as described above provides a sufficient level of prediction, and for which compounds additional kinetic processes need to be added to the model (for example, extrahepatic metabolism and/or transporter-mediated kinetics).
Recently, we evaluated the predictions made for Cmax in rats (upon a single exposure) for 44 compounds with a minimal rat PBK model, applying different methods to obtain chemical-specific input parameters (Punt et al., 2022). This analysis revealed that Cmax predictions ranged up to six orders of magnitude for some of the compounds as a result of the input parameters that were used. Cmax predictions were typically found to be within 5-fold of the observed Cmax for those compounds for which the unbound intrinsic clearance was low (<1 L/h) or relatively high (>20 L/h). However, one must be cautious with generalization of these findings to other compounds with different physicochemical properties, and findings obtained with the generic rat PBK model may not be applicable to a generic human PBK model. Since within next generation risk assessment predominantly human PBK model predictions will be required, it is crucial to particularly find means to evaluate human PBK model predictions in the absence of human in vivo data.
The goal of the present study was to assess the predictive performance of a minimal generic human PBK model, as described above, based on in vitro and in silico input data to predict plasma Cmax. The Cmax was selected as most critical kinetic parameter that is most frequently used in PBK modelling-based reverse dosimetry approaches to predict in vivo toxicity from in vitro data (Louisse et al., 2017), particularly when the biological effect has a threshold-related mechanism of action (Groothuis et al., 2015).To evaluate the predictive performance of the minimal human generic PBK model to predict Cmax values of a chemical, a literature search was performed to gather human in vivo plasma concentrations for a range of compounds. Cmax predictions were subsequently made upon a single oral dose based on a variety of input approaches for estimating the chemical-specific parameters (CLint, fup, partition coefficients, ka and fa) and compared to the observed Cmax values as reported in the collected literature. Based on the data obtained, we characterized the contribution of different input approaches to the variation in Cmax outcomes for individual compounds. In addition, for the compounds that were not predicted within 5-fold of the observed Cmax, a literature study was performed to find existing chemical-specific PBK models that have been developed for these compounds, allowing to evaluate which additional kinetic processes may need to be added to the generic PBK models to better describe their in vivo kinetics and improve plasma Cmax estimations.

Chemical selection and collection of human in vivo Cmax values
A dataset of 44 model compounds was created to evaluate the performance of the Cmax predictions by the PBK model based on different in vitro and/or in silico input approaches. The starting point for the selection of the model compounds, was the availability of in vitro human intrinsic hepatic metabolic clearance data as key input to PBK model development. The R httk package (version 2.0.4.) by Pearce et al. (2017) was used as primary source for these in vitro metabolic clearance data, based on which the clearance data for 40 compounds could be derived. The dataset was extended with rosuvastatin and fluvastatin, being structurally related compounds for which transporter-mediated processes in liver and kidney play a main role in the kinetics (Chan, 2019). In addition, ochratoxin A and coumarin were added to extend the dataset to also include nonpharmaceuticals. For the latter 4 compounds, in vivo human kinetic studies are available, and the intrinsic clearance was measured in the present study in incubations with human liver S9.
For the selected compounds, a literature search was subsequently performed in Scopus TM1 to identify human kinetic studies that report peak plasma concentrations of the compounds after a single oral dose or within the first 24 hour of a repeated oral dose study. For this literature study the following key words were used: ( ( TITLE ( "compound name" ) AND ALL ( bioavailability OR pharmacokinetics OR kinetics ) ) AND ( ( human OR man OR volunteer OR subject ) ) AND (Cmax OR "c max" OR "maximal concentration" OR "maximum concentration" OR "peak concentration" ) ). The studies that were identified for each compound were subsequently filtered to exclude 1) results obtained for specific patient groups like people with renal impairment or gastric by-pass, 2) studies with children, and 3) studies using slow or extended-release formulations. Description of kinetics in those situations require specific model adjustments, which was beyond the scope of the present study that focusses on the applicability domain of a generic PBK model for an average healthy adult. The final list of model compounds and related in vivo Cmax data (and related oral doses) from 421 different literature studies that were obtained are provided in the supplementary Excel file 2 (SM1).

2.2
Generic PBK model code and input parameters PBK model predictions were performed based on a published generic human PBK model code by Jones and Rowland-Yeo (2013) that was implemented as R script 3 by Punt et al. (2021) with minor modifications as summarized in the model code that is provided on GitHub 4 . The model consists of 13 compartments, corresponding to the major organs in the body and an arterial and venous blood compartment. The model requires chemical-specific parameters for intestinal uptake, distribution (i.e., partition coefficients, blood:plasma ratio (assumed to be 1 in the present study for all compounds), fraction unbound in plasma), hepatic clearance, and renal clearance (assumed to be the glomerular filtration rate times the free plasma concentration). Table  1 provides an overview on how these different input parameters were parameterized using a range of in vitro and/or in silico methods. Further details on these input approaches are given in the text below. The differential equations of the model are solved with the deSolve package (Soetaert et al., 2010).

Fraction unbound plasma
In vitro with rapid equilibrium dialysis Data derived from the httk package (Pearce et al., 2017) In vitro 44 In silico predicted based on log P, log D and pKa Lobell and Sivarajah, 2003 Lobell and Sivarajah 44 In silico predicted based on the SMILES

ADMET predictor ADMET 44
Absorption from the gastrointestinal tract was described in the model by a first-order uptake process from the intestine to the liver compartment and requires an absorption rate constant (ka) and fraction absorbed (fa) as input (Jones and Rowland-Yeo, 2013). For the parameterization of these input constants an in silico approach based on a QSAR from Hou et al. (2004) was applied that predicts the Caco-2 apparent permeability (Papp) based on the topological polar surface area (TPSA) of the compounds (Equation 1). For 30 of the 44 compounds, the QSAR-based approach was compared with in vitro measured Papp coefficients in Caco-2 transwell absorption experiments. These in vitro measured Caco-2 absorption data were partly obtained from Punt et al. (2022) and partly newly generated in the present study. Details of the QSAR calculations and Caco-2 experiments are provided in the supplementary information (SM2 5 ). Both the QSAR-derived Papp values and the in vitro measured values were scaled to an uptake rate constant (ka) and fraction absorbed (fa) based on the following equations.
Log Papp (cm/s) = -4.36 -0.01*TPSA (1) Log Peff,human (10 -4 cm/s) = 0.4926*log Papp (10 -6 cm/s) -0.1454 (2) ka (/h) = Peff*2 (cm/s) /R (cm) * 3600 (s/h) (3) fa = 1−(1+(2*Peff*<Tsi>)/(7*R)) −7 (4), in which equation 2 scales the Caco-2 apparent permeability to an human effective permeability based on Sun et al. (2002) and equations 3 and 4 describe how the effective permeabilities are converted to ka and fa as described by Yu and Amidon (1999). For the calculation of ka and fa with equations 2 and 3, an intestinal radius (R) of 1 cm was used and a small intestinal transit time <Tsi> of 3.32 hr (Grandoni et al., 2019). Physicochemical data (log P, log D and pKa values, TPSA), which are used as input to calculate the fup and tissue:plasma partition coefficients and intestinal uptake were derived with ADMET Predictor software 6 (v9.0, Simulation Plus, Lancaster, CA, USA) and with ChemAxon 7 (ChemAxon Ltd., Budapest, Hungary). Given that slight differences occur between the results of these two software packages with respect to the log P and pKa estimates, the influence of these differences on the PBK model predictions was evaluated. The log P, log D and pKa(s) that were obtained for the 44 compounds with each of the two software packages are provided in the Github repository 4 .
For the parameterization of fup values two in silico approaches were compared with in vitro measured values. One in silico approach for the calculations of fraction unbound in plasma, is a method of Lobell and Sivarajah (2003). Log P and information on the pKa(s) are required input parameters for this calculation. The R code for this calculation is provided in the in the Github repository 4 . As second in silico approach the fup calculations were obtained with the ADMET Predictor software 5 . The in vitro-derived human-specific fup values for 39 compounds were taken from the httk package with the original data  Wetmore et al. (2015) and Shibata et al. (2000). In case of the calculation of partition coefficients, three approaches were compared including the in silico approaches of (i) Rodgers and Rowland (2006), (ii) Berezhkovskiy (2004), which corresponds to the corrected method of Poulin and Theil (2002), and (iii) the in silico approach of Schmitt (2008). Log P and information on the pKa(s), and fup are required input parameters for these calculations. The R codes for these different calculation were obtained from Utsey et al. (2020) and were adjusted to fit the pipeline of the PBK model calculations of the current study. The codes can be found in the Github repository 4 . In case of the method of Rodgers and Rowland, the effect of including lysosomal trapping in the calculation method, as described by ( Schmitt et al., 2021), was tested. However, as this adjustment of the Rodgers and Rowland method did not lead to substantial differences in predicted partition coefficients (SM1), only the original approach of Rodgers and Rowland was included in the final dataset. In case of the method of Schmitt (2008), the membrane affinity (log MA) was calculated from the log P based on a QSAR from (Yun and Edginton, 2013) as provided in the code. The harmonized tissue composition data from Utsey et al. (2020) were used as input.
Three different approaches to obtain model parameter values for hepatic intrinsic clearance were compared. These include i) an in silico approach using ADMET Predictor 5 , ii) an in vitro approach based on clearance studies with primary human hepatocytes, and iii) an in vitro approach based on clearance studies with human liver S9. For the in silico based approach, the cytochrome P450-dependent human hepatic clearance rates of the 44 compounds (CYP_HLM_CLint) were predicted with the ADMET Predictor 5 . These predicted clearance rates, expressed as µL/min/mg microsomal protein, were scaled in the model to the whole liver based on a microsomal protein yield of 40 mg microsomes/gram liver (Barter et al., 2007). The primary human hepatocyte clearance data were derived from the database of the R httk package, containing the values that were originally measured by Wambaugh et al, 2019, Sohlenius-Sternbeck et al, 2012, Obach, 1999, Lombardo et al, 2018, Ito et al, 2004, Wetmore et al, 2015and Shibata et al, 2000. For 18 out of the 44 compounds, the intrinsic hepatic clearance was also measured in incubations with human liver S9 in the present study. The protocol for these incubations is provided in the Supplementary file (SM2 5 ). The primary hepatocyte intrinsic clearance data were scaled to the whole liver based on a hepatocellularity (number of hepatocytes per gram liver) of 117.5 x 10 6 ( Barter et al., 2007), whereas the S9-derived clearance data were scaled based on an S9 protein yield of 121 mg/g liver . Corrections for unspecific binding of the compounds to the hepatocytes or S9 in the in vitro incubations was applied based on a calculation method of Kilford et al. (2008) for primary hepatocytes and a method of Hallifax and Houston (2006) for the S9. Although the latter calculation method was developed to predict the unbound intrinsic clearance in microsomal incubations, it was assumed to also be applicable to S9 incubations.

2.3
Human PBK model predictions and data analysis By combining different input approaches, a total of 38,772 Cmax predictions were made for the different model compounds at the same exposure conditions as applied the in vivo studies from which the reported Cmax values were obtained. For each chemical the predicted Cmax was divided by the observed Cmax as marker of the quality of the of PBK model prediction for that compound. As a result of the different input combinations a range in predicted:observed ratios are obtained for each chemical of which the median is calculated. Median predicted Cmax values that are more than 5-fold higher than the observed Cmax values were considered as overestimated, and median Cmax values that are more than 5-fold lower than the observed Cmax values were considered as underestimated, though the latter did not occur in the present data set (see Results section). We also assessed the number of chemicals that are predicted within 2-fold, which is often used as a requirement when a PBK model is used within a regulatory context to demonstrate that the proposed model is fit for purpose (Shebley et al., 2018). The effect of different input approaches on the Cmax predictions was determined, by comparing for each input approach and compound the median Cmax and predicted:observed ratios and determining the differences between the input approaches in predicted median Cmax values.
A sensitivity analysis was performed for the predictions by changing the input value of a parameter by 1% and determining the relative change in Cmax, expressed as the normalized sensitivity coefficient (NSC) according to equation 5.
where C is the initial value of the model output, C′ is the modified value of the model output resulting from an increase in parameter value, P is the initial parameter value, and P′ is the modified parameter value. The sensitivity analysis was performed at an equal oral dose of 1 mg/kg bw for all compounds and input combinations. The R codes for the above analyses can be found in the Github repository 4 . To better understand the potential causes of why the Cmax of certain compounds could not be predicted within 5-fold of the observed Cmax, a literature study was performed to explore critical differences between the current minimal PBK model and existing PBK models for the different compounds. For this literature study the following key words were used: ( ( ALL ( "compound name") AND ALL (PBPK OR PBK OR PBBK OR PBTK) ). The obtained results were manually screened for the presence of PBK model codes for the given substance.

3.1
Evaluation of the collected data on human plasma concentrations For 41 of the 44 compounds, two or more in vivo studies were available from which the peak plasma concentrations could be derived. For these compounds the variation in reported plasma concentrations between studies was assessed by calculating the normalized Cmax (Cmax /dose) for each study and assessing the distribution of these normalized Cmax values. For 22 of the 44 ALTEX, accepted manuscript published January 19, 2022 doi:10.14573/altex.2108301 5 compounds, the difference in normalized Cmax values between the different studies ranged between 1.3-and 3-fold and for 10 compounds between 3 and 5-fold (see SM1). For the remaining 12 compounds larger differences between studies are observed, with the highest difference occurring for metoprolol (16-fold difference normalized Cmax values between studies) and dextromethorphan (13-fold difference normalized Cmax values between studies). These compounds are both CYP2D6 substrates (Frank et al., 2007) and are therefore prone to large interindividual human variation. The PBK model predictions (based on combinations of input approaches) were compared to each of the in vivo reported Cmax from all the different studies of the 44 compounds.

3.2
Performance of the generic PBK model based on different input approaches In Fig. 1 the ratios between PBK model-predicted Cmax values and in vivo observed Cmax values are shown. A large variation (1-4 orders of magnitude) in predicted:observed ratios can be observed, which is predominantly the result of the different input approaches that were applied and to a lesser extent reflects the experimental variation between the in vivo studies to which the PBK-model predictions are compared (as described above). The highest range in predicted:observed ratios is observed for curcumin. Depending on the input combinations, a 12-to 16,654-fold overprediction of the reported Cmax values is obtained. Despite the variation in predicted:observed ratios for the different model compounds in the dataset, the median Cmax of the majority of the compounds (32 out of the 44) are within 5-fold of the observed Cmax values, and the Cmax of 19 out of 44 compounds could be predicted within 2-fold. The Cmax values 12 compounds were predominantly overestimated (i.e. having a median predicted Cmax that is >5-fold higher than the observed Cmax). Although some of the input combinations led to significant underpredictions for various compounds, only the median Cmax of ochratoxin A was more than 5-fold underpredicted.

3.3
Sensitivity analysis Improving the predictive value of the minimal PBK model requires insights in the critical input parameters that affect the model predictions. To this end, a sensitivity analysis was performed for the predictions by changing the input value of a parameter by 1% and determining the relative change in Cmax, expressed as the normalized sensitivity coefficient (NSC) at a dose of 1 mg/kg bw. Fig. 2A shows the results for the most influential input parameters that affect the Cmax predictions (max NSCs>0.5 in absolute value). The NSCs of remaining input parameters can be found in Figure S1 5 . Figure 3 reveals that the scaled unbound intrinsic clearance (CLint,u), the fraction unbound (fup) and parameters that determine the blood flow to and the volume of the liver (FVli, FQh, Fqgu) are among the most critical input parameters. All these parameters determine the availability of the compounds for metabolic clearance. Other important input parameters that affect the Cmax predictions relate the oral absorption (ka and fa) and the blood flow (QC). The B:P ratio also showed to be a sensitive parameter. This parameter was set to a default value of 1 for all compounds in the present study, since measured data on B:P ratios are generally lacking and no in silico tools are available to estimate the B:P ratio. Altogether, these most influential parameters can be summarized as factors that affect the bioavailability of the compound, such as the fraction absorbed and the liver clearance, and factors that that influence the volume of distribution, such as the fraction unbound in plasma, B:P and to a lesser extent the tissue:plasma partition coefficients (Fig. S1 5 ). Figure 2A, reveals that not all compounds are equally sensitive to the different input parameters. The observed differences between compounds in sensitivity were found to relate to the extent of metabolic clearance: the sensitivity increases with increasing CLint,u until a maximum sensitivity is reached for compounds with a high CLint,u (Fig. 2B). A similar relation is observed between the CLint,u of the compounds and the sensitivity towards the other input parameters like ka, fa, and QC ( Fig. S2 5 ).

Fig. 1: Ratios between PBK model-predicted Cmax values and in vivo observed Cmax values observed for 44 reference compounds in human
Per chemical, different predicted Cmax values are obtained by running simulations with the different input approaches (in vitro or in silico approaches to parameterise input parameters) as presented in Table 1. Each predicted Cmax is then compared with the in vivo Cmax values for the chemical in the dataset. The median of these predicted:observed ratios is depicted along the individual datapoints. Datapoints between the dashed horizontal lines are within 5-fold of the observed Cmax and the datapoints between the dotted horizontal lines are within 2-fold of the observed Cmax. Compounds for which the median predicted Cmax is more than 5-fold overestimated are depicted in light grey.

Fig. 2: A) Normalized sensitivity coefficients (NSCs) of the Cmax predictions to different input parameters for the different compounds. B) NSCs for CLint,u plotted against CLint,u values used as input.
The datapoints in the figures correspond to the NSCs for a random selection of 12 Cmax simulations based on different input approaches per chemical. BP = blood plasma ratio, CLint,u = unbound intrinsic liver clearance , fa = fraction absorbed, FQgu = blood flow fraction to the gut, FQh = blood flow fraction to the liver, fup = fraction unbound plasma, FVli = volume fraction liver, ka = intestinal uptake rate, QC = cardiac output.

Effect of the input approaches on Cmax predictions
Given the sensitivity of the Cmax predictions to chemical-specific input parameters like CLint, fup, fa, and ka, and to a lesser extent to the partition coefficients, the quality of the input approaches for these parameters can have a substantial effect on the Cmax predictions. It is therefore important to understand if certain input approaches perform better than others. To this end, we determined whether differences occur in median Cmax predictions between the different input approaches. In Figure 3 the Cmax predictions with more than 3-fold differences between the applied input approaches are highlighted. These results reveal that differences in Cmax predictions (and corresponding predicted:observed ratios) most frequently occur as a result of differences in calculation methods for the partition coefficients ( Figure 3A), with the method of Rodgers and Rowland performing best. The median Cmax of 19 compounds was predicted within 2-fold of the observed Cmax values with the Rodgers and Rowland method. With the methods of Berezhkovskiy and Schmitt 12 and 13 compounds were predicted within 2-fold of the observed Cmax, respectively. Particularly for acidic compounds (pKa < S6), like bosentan, diclofenac, fluvastatin, ibuprofen, ochratoxin A, naproxen, S-warfarin, and tolbutamide, the method of Berezhkovskiy resulted in relatively high partition coefficients and low Cmax predictions compared with the other input approaches ( Figure 3A). The method of Schmitt resulted in overpredictions of the Cmax of particularly lipophilic compounds, like bisphenol A, bosentan, clozapine, and midazolam.
In case of the of the parameterization of the intrinsic hepatic clearance, the in silico calculated clearance values were found to only provide rough estimates of the intrinsic clearance values as no direct correlation was observed between the in vitro (hepatocyte) measured clearance data and in silico (ADMET microsomal) predicted ones (R 2 = 0.06, Figure S3A). Yet, substantial differences (>3-fold) in median Cmax predictions between the in silico and in vitro input values for intrinsic clearance only occurred for chlorpromazine, coumarin, curcumin, and naloxone ( Figure 3B). For those compounds the in silico calculated clearance values led to (higher) overestimations of the observed Cmax values compared with the in vitro measured clearance values. For the compounds for which S9 data were available, the clearance values strongly correlated with the hepatocyte clearance data (R 2 = 0.75, Figure 3C). Significant differences (>3-fold) in Cmax predictions between the latter two approaches were only observed for naloxone ( Figure 3B). A relatively lower number of compounds could be predicted within 2-fold when parameterized based on S9 clearance data. This was, however, expected as the S9 clearance data are not equally distributed over the dataset, and contain a relatively high number of chemicals that are >5-fold overpredicted (irrespective of the input approach). Figure 3C reveals that the highest number of compounds were predicted within 2-fold when the in silico predicted fup values based on the method of Lobell and Sivarajah (2003) or the in vitro measured fup values were used. Although these two approaches appeared to perform equally well, the use of in vitro measured unbound fractions particularly did lead to overestimations of the Cmax of buspirone, carvedilol, coumarin, prazosin, and resveratrol. The Cmax predictions based on the ADMET predicted fup frequently led to underpredictions. Regarding the sources for the physicochemical characteristics (ADMET versus ChemAxon) and intestinal uptake data (in vitro Caco-2 data versus in silico predictions) no substantial (>3fold) differences in Cmax predictions occurred between the different input approaches. Figure 4 depicts the results of the dataset in which the most significant outliers as described above are removed. This includes removal of the simulations based on 1) the methods of Berezhkovskiy and Schmitt for the partition coefficients, 2) the in silico intrinsic hepatic clearance data, and 3) the ADMET predicted and in vitro measured fup values. The results of the reduced dataset show a significant reduction in the variation in Cmax predictions and the related predicted:observed Cmax ratios. Within the reduced dataset, the median Cmax values of 34 out of the 44 compounds are predicted within 5-fold, and the Cmax values of 19 compounds are predicted within 2-fold. Ten out of the 12 initially overestimated reference compounds remain overestimated with more than 5-fold. The median Cmax predictions for chlorpromazine and coumarin now fall within 5-fold of the observed Cmax. No underestimations of more than 5-fold occur.  Fig. 3: Differences in predicted:observed Cmax ratios related to the applied input approaches Highlighted are those results for which a more than 3-fold difference in median Cmax predictions occur between the applied input approaches. In case of Figure 3B, S9-derived intrinsic clearance data are included in the comparison of input approaches for compounds marked with an asterisk. For each approach the number of compounds that was predicted within 2-fold of the observed Cmax is provided between brackets.

Fig. 4: Ratios between PBK model-predicted Cmax values and in vivo observed Cmax values observed for 44 reference compounds in human as obtained after removal of the simulations based on input methods that led to significant differences between predicted and observed ratios as described in the main text
3.5 Characteristics of the compounds for which the Cmax is more than fivefold overpredicted Irrespective of the input method that was applied, the majority of the Cmax predictions for bisphenol A, buspirone, curcumin, desipramine, dextromethorphan, fluvastatin, genistein, naloxone, resveratrol, and rosuvastatin were more than 5-fold overpredicted. Underpredictions of more than 5-fold did not occur. To obtain insights in the common causes for the overpredictions, a literature search was performed to explore available chemical-specific PBK models for all model compounds. The overpredictions either indicate that the bioavailability of the compounds is predicted to be too high or that the volume of distribution is predicted to be too low. This means that additional kinetic processes need to be considered that are related to either a reduced bioavailability or an increased volume of distribution. The PBK model predictions of the present study were performed based on a minimal PBK model that includes the intrinsic liver clearance, partition coefficients and the fraction unbound in plasma as main input parameters, lacking description of kinetic processes that could reduce the bioavailability like, for example, extrahepatic metabolism or active efflux transport. In addition, parameters that affect the volume of distribution, like the fup, partition coefficients or tissue uptake transport, might be different in the current model compared with the chemical-specific PBK models. The key differences between the chemical-specific PBK models for the different compounds available in literature and the current generic model were explored. Table 2 reveals the results of this comparison and shows that phase I metabolism or glucuronidation by intestinal epithelial cells is most frequently considered in the chemical-specific PBK models. This is particularly true for the compounds that were substantially overpredicted in the present study, but also for some of the compounds that were predicted within 5-fold (i.e., lorazepam, midazolam, sildenafil, and verapamil). Apart from intestinal metabolism, hepatic uptake and biliary excretion by active transport have been considered for the kinetic models of bosentan, bisphenol A, fluvastatin, rosuvastatin, and S-warfarin. In this context, the inclusion of active hepatic uptake in the PBK model may increase the liver concentration, resulting in enhanced metabolism and/or biliary excretion and reduced (initial) peak plasma concentrations. Biliary excretion might increase plasma concentrations again at later time points due to enterohepatic circulation. Intestinal efflux transport and renal active resorption and excretion have only been considered for a few compounds in the identified PBK models.

Discussion
Adequate predictions of internal dosimetry, such as Cmax, are crucial in the transition towards next generation (animal-free) testing strategies for chemical safety evaluations to convert in vitro toxicity data into in vivo dose−response or at least potency information (e.g., Fabian et al., 2019;Louisse et al., 2017;Punt et al., 2019;Wetmore et al., 2015). To facilitate the application of PBK models in next generation risk assessment, it is important to gain confidence in the predictive performance of these models without case-by-case validation against in vivo data (e.g., plasma concentrations). The goal of the present study was to assess the predictive performance of a minimal generic human PBK model based on in vitro and in silico input data to predict the Cmax. Two cut-off values (2-fold and 5-fold) were used as performance indicators in the current study. Discussions are presently still ongoing on what level of deviation between predicted and observed kinetics is acceptable within a regulatory context (Shebly et al., 2018). The 2-fold cut off value is frequently requested within a regulatory context to demonstrate that the proposed model is fit for purpose (Peters and Dolgos, 2019;Sager et al., 2015). A key challenge with this 2-fold is that the differences between in vivo studies, to which the model predictions are compared, tend to be higher than 2-fold themselves, possibly related to inter-individual differences in biology or technical aspects (Shebly et al., 2018). This was also observed in the present study, revealing that the variation in normalized in vivo Cmax values generally range between 1.3 to 5-fold, but can  Li et al., 2015;Yang et al., 2020;Posada et al., 2020 Bisphenol A -Intestinal glucuronidation -Hepatic uptake transport Kawamoto et al., 2007;Teeguarden et al., 2005 Lorazepam -Intestinal glucuronidation Docci et al., 2020 Buspirone -Intestinal phase I metabolism (CYP3A4) Heikkinen et al., 2012;Gertz et al., 2011;Karlsson et al., 2013 Midazolam -Intestinal phase I metabolism (CYP3A4) Karlsson et al., 2013;Heikkinen et al., 2012;Cao and Jusko, 2012;Gertz et al., 2011;Nguyen et al., 2016;Jamei et al., 2009  also be as high as 16-fold. Based on these results a 5-fold cut-off value was selected to determine whether a compound fits the applicability domain of the applied minimal PBK model Overall, the medians of the predicted Cmax values for 34 out of the 44 compounds, corresponding to 77%, were within 5-fold of the observed Cmax values, whereas 19 compounds could be predicted within 2-fold, corresponding to 43%. The medians of the predicted Cmax values of 10 out of 44 compounds (23%) were more than 5-fold overestimated. Underestimations of the median Cmax (higher than 5-fold) did not occur. A bias towards overestimations has been reported before for PBK model predictions based on minimal input (i.e. liver clearance, partition coefficients, and passive intestinal uptake) and might indicate that such minimal PBK model predictions represent a worst case approach for predicting Cmax values (Wambaugh et al., 2018). The choice of the calculation method for tissue:plasma partition coefficients had a high impact on the Cmax predictions. In addition, occasional differences in Cmax predictions were observed as a result of differences between input approaches (in silico, in vitro) for the fraction unbound in plasma and the different input approaches (in silico, primary hepatocytes, or S9) for intrinsic hepatic clearance. In case of the calculation methods for the partition coefficients, particularly the calculation method of Berezhkovskiy led to underpredictions of the Cmax of acidic compounds (pKa <6). For these compounds, the predicted partition coeffects for the different organs (except for the adipose tissue) were higher than with the method of Rodgers and Rowland or the method of Schmitt, resulting in lower plasma concentrations. The method of Berezhkovskiy takes the charge of the molecules into account for the calculation of the partition coefficient of the adipose tissues, but not for the other organs. This might explain the differences between the method of Berezhkovskiy and the calculations methods of Rodgers and Rowland and Schmitt, with the latter specifically focusing on the impact of drug ionization on partitioning (Utsey et al., 2020). The calculation method of Schmitt appeared to work less well for lipophilic compounds. This was previously also observed by Punt et al. (2022), and might be due to the fact that the method of Schmitt largely depends on the membrane affinity as input, which is difficult to obtain for different chemicals. The applied QSAR from Yun and Edginton (2013) to calculate the membrane affinity, might not be applicable to lipophilic compounds. Overall, the calculation method of Rodgers and Rowland performed best in the current study, resulting in the highest number of compounds for which the Cmax predictions were within 2-fold (and 5-fold) of the observed Cmax.
The in silico input approaches that were used in the present study for the calculation of the apparent Caco-2 permeability, fraction unbound in plasma (particularly the method of Lobell and Sivarajah), and intrinsic hepatic clearance, provided relevant Cmax predictions for the majority of the compounds. However, particularly in case of the in silico calculated intrinsic clearance values and the in silico calculated Caco-2 permeability, care should still be taken with using these in silico approaches for the parameterization of PBK models as direct comparisons of the in silico and in vitro measured values revealed poor correlations between the in vitro and in silico data ( Figure S3). It should also be noted that the in silico intrinsic clearance results that were obtained with the ADMET Predictor 5 only cover the cytochrome P450-dependent human hepatic clearance, and not sulfation or glucuronidation reactions. In general, in vitro input parameters remain therefore the preferred input approach for these input parameters. In contrast, the values of the in silico calculated fractions unbound in plasma were found to correlate better with the in vitro measured ones (R 2 = 0.6) and even provided better Cmax predictions in the present study than the in vitro measured values. This might be due to challenges with measuring this parameter in vitro, particularly for highly lipophilic or chemically unstable compounds (Bowman et al., 2021;Emami Riedmaier et al., 2016;Wambaugh et al., 2019).
Given the importance of the in vitro measured input data for PBK model development the quality of these measurements is crucial. This is true for the measurements of the fraction unbound as described above as well as for in vitro intrinsic hepatic clearance measurements and Caco-2 permeability studies. Recently, Louisse et al. (2020) evaluated the influence of experimental conditions of clearance studies on intrinsic clearance (CLint) values obtained from literature data. The CLint values for the majority of compounds differed by more than one order of magnitude, which is expected to partly depend on the in vitro protocol that was used. Such variation can have a large impact on the Cmax predictions. These results highlight the importance of obtaining harmonized in vitro approaches (Gouliarmou et al., 2018;Louisse et al., 2020;Paini et al., 2019) to parameterize PBK models.
The frequently observed overpredictions, indicate that particularly kinetic processes are missing in the minimal PBK model that are related to reduced predicted bioavailability or increased predicted volume of distribution. Comparison of the current model structure with chemical-specific PBK models reported in literature revealed that particularly intestinal first pass metabolism (CYP3A4-mediated oxidation or glucuronidation) is frequently considered as additional kinetic process, which was not included in the current generic PBK model structure. In addition, active hepatic uptake and biliary excretion have been included in literature available chemical-specific PBK models of some of the compounds of the present study (Emami Riedmaier et al., 2016). Intestinal efflux and active renal excretion and/or resorption are less frequently considered, yet might also be important. Solubility and dissolution issues might also be an issue, resulting in a poorer intestinal fraction absorbed. New approaches to describe passive renal excretion based on Caco-2 permeability data might also help in obtaining better estimates of plasma concentrations (Scotcher et al., 2016). Further work will be needed to determine whether the Cmax predictions can be improved by inclusion of some of these additional kinetic processes.
The current study focused on Cmax predictions with the minimal PBK model. This is a critical kinetic parameter within PBK modelling-based reverse dosimetry approaches to predict in vivo toxicity from in vitro data (Louisse et al., 2017). Further work will, however, be needed to also evaluate the predictive performance with respect to other relevant kinetic descriptors like the area under the plasma concentration-time curve (AUC), steady state plasma concentrations (Css). Together these different parameters describe the full kinetic profile of a chemical, and it remains to be elucidated to what extent the optimized selection of input parameters, as obtained in the present study for Cmax, will also work for other biokinetic parameters. In addition, it should be noted that the 44 reference compounds only represent a limited chemical space. Even though the chemicals within the dataset are diverse (e.g., log P between -0.55 and 5.3, MW between 179 and 552, and different ionization characteristics), the observations made and conclusions drawn within the present study cannot be directly extrapolated to a larger dataset. Overall, the results of the current study provide relevant insights into the predictive performance of a minimal PBK model and the influence of different input approaches on the model predictions. Further work will be needed in particular to find ways to determine when and which additional kinetic processes (like liver metabolism or transporter-mediated processes) need to be added in the absence of prior knowledge on the chemical's in vivo toxicokinetics.