Theoretical Biology and Medical

Background: Drug-drug interactions resulting from the inhibition of an enzymatic process can have serious implications for clinical drug therapy. Quantification of the drugs internal exposure increase upon administration with an inhibitor requires understanding to avoid the drug reaching toxic thresholds. In this study, we aim to predict the effect of the CYP3A4 inhibitors, itraconazole (ITZ) and its primary metabolite, hydroxyitraconazole (OH-ITZ) on the pharmacokinetics of the anesthetic, midazolam (MDZ) and its metabolites, 1' hydroxymidazolam (1OH-MDZ) and 1' hydroxymidazolam glucuronide (1OH-MDZ-Glu) using mechanistic whole body physiologically-based pharmacokinetic simulation models. The model is build on MDZ, 1OH-MDZ and 1OH-MDZ-Glu plasma concentration time data experimentally determined in 19 CYP3A5 genotyped adult male individuals, who received MDZ intravenously in a basal state. The model is then used to predict MDZ, 1OH-MDZ and 1OH-MDZ-Glu concentrations in an CYP3A-inhibited state following ITZ administration.


Background
Drug-drug interactions resulting from the inhibition of an enzymatic process can have serious implications for clinical therapy. Quantifying the magnitude of the inhibitor effect in vivo is an active area of study although methods of quantifying the exposure increase of a drug concomitantly administered with an inhibitor have focused on, until now, simplistic, static models [1][2][3]. These approaches assume that there is a proportional increase in exposure at high inhibitor concentrations and do not account for the time course of inhibitor concentrations. The approach that is taken is one which dynamically links inhibitor and drug models using whole-body physiologically-based pharmacokinetic models (WB-PBPK) to quantify, under any administration time and dose regimen, the changes that occur in parent compound exposure as well as the dynamic changes in the respective metabolite exposures. This has been done for midazolam (MDZ), and its two major metabolites 1' hydroxymidazolam (1OH-MDZ) and the glucuronide of 1' hydroxymidazolam (1OH-MDZ-Glu), in the presence of the CYP3A4 inhibitors itraconazole (ITZ) and its major metabolite hydroxyitraconazole (OH-ITZ). This example was used because of the importance of CYP3A4 to drug metabolism and the availability of a full clinical data set for MDZ given in the basal and ITZ/OH-ITZ inhibited state [4].
Cytochrome P450 (P450) enzymes play an important role in the metabolism of exogenous and endogenous molecules. In humans, CYP3A represents one of the most important subfamilies of the P450 superfamily. CYP3A4 is the major P450 in the liver and intestine and has been reported to be involved in the metabolism of more than 60% of all medically relevant drugs [5]. The expression of CYP3A5 is highly polymorphic, due to a single nucleotide polymorphism, which is designated CYP3A5*3 [6]. Population frequencies for CYP3A5 variants in mixed American and Korean individuals are 61-77% for CYP3A5*3/ *3, 22-33% for CYP3A5*1/*3 and 1-5% for CYP3A5*1/ *1 [4,7] with CYP3A5*3/*3, CYP3A5*1/*3 and CYP3A5*1/*1 constituting 5%, 50% and 76% of the total CYP3A concentration, respectively. Total CYP3A content was more than 2-fold higher for livers with at least one CYP3A5*1 allele compared with CYP3A5*3/*3 livers [7]. Because CYP3A5 exhibits an overlapping substrate specificity with that of CYP3A4, it may contribute significantly to the metabolic elimination of CYP3A substrates in people carrying the wild-type CYP3A5*1 allele, although in vivo data as well as in vitro evidence are conflicting [4,8]. Because CYP3A is significantly involved in drug biotransformation, drug-drug interactions resulting from the inhibition of CYP3A-mediated metabolism by a coadministered therapeutic agent are of clinical importance.
Itraconazole (ITZ) is an orally active triazole antimycotic agent, which is active against a broad spectrum of fungal species. ITZ is extensively metabolized in humans, yielding over 30 metabolites, including its primary active metabolite hydroxy-itraconazole (OH-ITZ). ITZ and its subsequent sequential metabolites [OH-ITZ, keto-itraconazole (keto-ITZ) and N-desalkyl-itraconazole (ND-ITZ)] are all high affinity ligands and substrates of CYP3A4 [13]. ITZ and OH-ITZ are also competitive inhibitors of CYP3A4. Keto-ITZ and ND-ITZ may also contribute to CYP3A4 inhibition in vivo following ITZ therapy, but their concentration following ITZ administration is significantly lower resulting in a low inhibitory influence [13].
In this study, we aim to predict the effect of ITZ and OH-ITZ CYP3A4 inhibition following oral ITZ administration on the pharmacokinetics of intravenously administered MDZ using a mechanistic WB-PBPK simulation model. WB-PBPK modeling allows for the simulation of the fate of xenobiotics in the human body on the basis of individual physiological characteristics [14]. The model is build on MDZ, 1OH-MDZ and 1OH-MDZ-Glu plasma concentration time data experimentally determined in 19 CYP3A5 genotyped adult male individuals, who received MDZ intravenously in basal and ITZ-inhibited CYP3A metabolic states [4]. At first, three WB-PBPK models, one for MDZ, one for 1OH-MDZ and one for 1OH-MDZ-Glu, will be coupled per individual to dynamically simulate the kinetics of MDZ hydroxylation and glucuronidation for each study volunteer. In a second step, the MDZ models will be extended by linking WB-PBPK models of ITZ and OH-ITZ in order to predict the interaction between MDZ and ITZ via ITZ and OH-ITZ-mediated CYP3A4inhibition in a time-dependent manner.

Results
The simulated plasma concentration time curves following intravenous administration of MDZ and 1OH-MDZ adequately represented the corresponding in vivo time course data reported by Mandema et al. [15] (Figure 2).  Figure 3 presents the mean (± standard deviation) of the six individual optimized plasma concentration time profiles for MDZ and 1OH-MDZ (lines) along with the mean (± standard deviation) of the experimental data (symbols) [11]. The elimination optimized curves matched the experimental data well. Table 1 presents the minimum, mean and maximum velocity-rate constants (k 1,-1,2 ), K m , V max, V max /K m and CL(1OH) values which were the results of the individual optimizations for the six individuals.
The optimized plasma concentration time profiles for MDZ, 1OH-MDZ and 1OH-MDZ-Glu (lines) are presented in Figure 4 together with the experimental data (symbols) [16]. The elimination optimized curves well represented the experimental data. Numerical results are shown in Table 1. Figure 5 presents the mean (± standard deviation) optimized plasma concentration time profiles for MDZ and the sum of 1OH-MDZ and 1OH-MDZ-Glu (lines) along with the mean experimental data (symbols) separated by the CYP3A5 genotype group [4]. The elimination optimized curves well represented the experimental data.  [5.1-34.9].
In Figure 6, the individuals from the Yu et al [4] are ranked according to their K m and V max values for the hydroxylation, normalized to the mean. There is no correlation between the CYP3A5 genotype and the affinity of MDZ to CYP3A (K m ) or the rate of hydroxylation (V max ) in the basal state. Figure 7 presents the simulated ITZ plasma concentration time curve following 200 mg dosing of ITZ in the fed state as compared to experimental data from Barone et al. [17] for the first day and Hardin et al. [18] for the trough concentrations on subsequent days. The fraction of remaining CYP3A activity resulting from the simulated unbound intracellular concentrations of ITZ and OH-ITZ in the liver over time is also presented in Figure 7. rated by the CYP3A5 genotype group following ITZ administration [4]. MDZ exposure increases with a consequent decrease in initial 1OH-MDZ+1OH-MDZ-Glu concentrations and increase at later time points. While the shape of the curve is generally well described, the simulated concentrations tend to exceed the experimental data between approximately 6 to 12 hours. This deviations appears to be greater in the CYP3A5 *1/*1 and CYP3A5*1/*3 group than in the CYP3A5 *3/*3 group.

Discussion
A dynamic, coupled WB-PBPK model of MDZ and its metabolites was developed to predict the time-dependent interaction of MDZ when co-administered with the CYP3A4 inhibitor ITZ. This novel method differs greatly from previously used models to identify the potential for drug-drug interactions. Static models generally prevail where competitive and irreversible inhibition is defined by the kinetic parameters of the inhibitor and are only dose dependent in relation to the inhibitor [2,3]. Using the static model, a ratio of basal state area under the curve (AUC) to inhibited state AUC is calculated, which provides a means of estimating the mean extent of inhibition in a given time period. Identifying parent compound and metabolite curve shape changes along with the generation of any metabolite profiles is not possible using this method. Since the static nature of these models limits their flexibility, dynamically interacting WB-PBPK models were seen as the method of choice for understanding how drug-drug interactions affect pharmacokinetic profiles of both the parent compound and its metabolites. Recently, Özdemir et al. [19] presented an attempt to predict the drug-drug interaction of MDZ with another antifungal agent, ketoconazole by means of the AUC increase of MDZ when co-administered with KTZ, and concluded that comprehensive PBPK models are required to correctly predict the inhibitory effect of ketoconazole.
The study of Yu et al. [4] provided an excellent data set of anthropometric (weight, height, age) data, CYP3A5 genotype and plasma concentration profiles for 19 individuals given MDZ alone and in combination with ITZ. Using this experimental data, each basal state individual MDZ and metabolite profile was optimized and, followed by the addition of the ITZ data and models, was used to test the Mean of the predicted plasma concentration time curves (line) of MDZ and 1OH-MDZ following an intravenous administration of 0.15 mg/kg as compared to the mean (n = 6) of the experimental data Figure 3 Mean of the predicted plasma concentration time curves (line) of MDZ and 1OH-MDZ following an intravenous administration of 0.15 mg/kg as compared to the mean (n = 6) of the experimental data. Experimental data was taken from Heizmann et al. [11]. model's ability to predict the MDZ-ITZ drug-drug interaction potential.
In a first step, the MDZ plus metabolite model was developed. Distribution properties of MDZ, 1OH-MDZ and 1OH-MDZ-Glu in this study were based on previously established mechanistic models to estimate organ partitioning parameters [14,20,21]. Four out of 17 organ partition coefficients were ultimately modified using experimental partition coefficients for MDZ scaled from rats [22]. For MDZ and 1OH-MDZ, this modification accurately led to simulated plasma concentration time profiles that mirrored observed profiles following IV administration.
The subsequent coupling of MDZ, 1OH-MDZ and 1OH-MDZ-Glu was achieved using velocity rate constants for metabolism parameterization. Using a step wise process inclusive of different experimental data sets [11,16], the velocity rate constants that defined the basal state kinetics of MDZ, 1OH-MDZ and 1OH-MDZ-Glu, were optimized for our study data [4]. In order to evaluate the optimization results, 1OH_K m values were compared to in vitro findings. The mean ± SD 1OH_K m values for the Heizmann et al. [11], Kharasch et al. [16] and Yu et al. [4] data sets were 2.2 ± 0.99 umol/L, 5.3 umol/L and 2.1 ± 1.17 umol/L, respectively. These were all similar to each other and also to the in vitro mean ± SD 1OH_K m = 3.9 ± 3.1 umol/L reported by Patki et al. [23]. The coefficient of variability (CV%) for 1OH_K m and 1OH_V max were also compared to in vitro findings to evaluate the reasonability of the generated inter-individual variation. The in vitro CV% for 1OH_K m and 1OH_V max generated from human liver microsomes (HLMs) of twelve individuals being 80% and 58%, respectively [23], were similar to the results of the optimization of the Yu et al. [4] data, namely 1OH_K m CV% = 57% and 1OH_V max CV% = 47%. Therefore, there is a correlation between in vitro and in vivo 1OH_K m and 1OH_V max for the hydroxylation of MDZ via CYP3A4/5. Furthermore, the individual 1OH_V max variability was linked to the underlying process by comparing the CV% of 1OH_V max with that of in vitro hepatic CYP3A4/5 concentration. The CV% for the CYP3A4/5 concentration in a group of Japanese, Caucasian and mixed Americans was 160% [23,24] whereas that for the Japanese population alone was 35% [24]. The CV% of V max for the Koreans of the Yu et al. [4] data, 47%, is between both values and closer to that for the Japanese population. This suggests that our model was able to capture the inter-individual variability in hepatic CYP3A concentration using 1OH_V max as a surrogate.
The enzyme responsible for the glucuronidation of 1OH-MDZ is unknown. In order to estimate the relevant UDPglucuronosyltransferase (UGT) concentration, the findings by Reinach et al. [25] were used. They observed that the metabolism of MDZ was 0.080 nmol/min/10 6 cells Predicted plasma concentration time curves (lines) of MDZ, 1OH-MDZ and 1OH-MDZ-Glu following an intravenous bolus administration of 1 mg MDZ as compared to experi-mental data from a typical individual (symbols) Figure 4 Predicted plasma concentration time curves (lines) of MDZ, 1OH-MDZ and 1OH-MDZ-Glu following an intravenous bolus administration of 1 mg MDZ as compared to experimental data from a typical individual (symbols). Experimental data was taken from Kharasch et al. [16]. Mean ± standard deviation plots of the experimental observed plasma concentration time profiles for MDZ and the sum of its metabolites, 1OH-MDZ and 1OH-MDZ-Glu (symbols: squares) and the corresponding mean elimination optimized curves (solid lines) for the CYP3A5 genotypes, CYP3A5 *1/*1 (n = 6), CYP3A5 *1/*3 (n = 6) and CYP3A5*3/*3 (n = 7) in the basal state Figure 5 Mean ± standard deviation plots of the experimental observed plasma concentration time profiles for MDZ and the sum of its metabolites, 1OH-MDZ and 1OH-MDZ-Glu (symbols: squares) and the corresponding mean elimination optimized curves (solid lines) for the CYP3A5 genotypes, CYP3A5 *1/*1 (n = 6), CYP3A5 *1/*3 (n = 6) and CYP3A5*3/*3 (n = 7) in the basal state. Also presented is the mean ± standard deviation plots of the experimental observed plasma concentration time profiles for MDZ and the sum of its metabolites, 1OH-MDZ and 1OH-MDZ-Glu (symbols: circles) and the corresponding mean predicted curves (dotted lines) for the CYP3A5 genotypes in the CYP3A inhibited state resulting from ITZ administration. Plasma concentration time curves in the inhibited state were graphed starting at time = 0 to allow for direct comparison with the basal state curves. Experimental data was taken from Yu et al. [4].
while the metabolism of 1OH-MDZ was 0.020 nmol/ min/10 6 cells in human hepatocytes. As there was no additional information available, for our study, it was assumed that the difference between the metabolism of MDZ and 1OH-MDZ was only dependent on differing enzyme concentration. Under this assumption, the UGT concentration was estimated, namely [UGT] = 0.25 [CYP3A]. The absolute Glu_K m and Glu_V max values are therefore provisional. Nevertheless, the optimized curves will not change when the actual UGT concentration is known, as the product of the velocity-rate constants and the enzyme concentration were optimized. Thus, if the enzyme concentration changes, then k 1 , k -1 and k 2 can change accordingly (see equation 1) therefore changing the calculated K m (see equation 2, no steady state assumption is made) and V max values (see equation 3). Despite that the absolute K m and V max values are dependent on the enzyme concentration, inter-individual variability is independent of the absolute values. It is remarkable that the CV% of Glu_K m was relative high with 191%. For Glu_V max , the CV% equaled 299%. This was substantially higher than the CV% for 1OH_V max . As the CV% for V max can be linked to the underlying process, or more precisely the enzyme concentration, it can be inferred that the inter-individual variability of the corresponding UGT concentration might be higher than that of the CYP3A concentration. This hypothesis cannot be confirmed nor rejected because the UGT(s) catalyzing the glucuronidation of 1OH-MDZ or similar compounds are unknown. Once the UGT(s) catalyzing the glucuronidation of 1OH-MDZ and its (their) concentration is known, Glu_K m and Glu_V max can be recalculated and further interpreted.
The information about the contribution of CYP3A5 to overall MDZ metabolism is controversial. According to Huang et al. [8], CYP3A5 accounts for 27% of the total product formation catalyzed by microsomes with at least one CYP3A5*1 allele, whereas Williams et al. [26] reports that CYP3A4 and CYP3A5 contribute equally to the MDZ elimination. Other in vitro studies demonstrate that CYP3A5 has a higher catalytic activity but a lower affinity than CYP3A4 for MDZ 1-hydroxylation [8,23,26]. Despite the fact that the in vitro studies consistently suggested an influence of the CYP3A5 concentration on the MDZ elimination, in vivo findings have been inconsistent. A higher MDZ elimination for white cancer patients [27] and healthy Asian subjects [28] having at least one wild-type CYP3A5*1 allele compared to CYP3A5*3/*3 patients have been reported. In contrast, Floyd et al. [29] and Shih & Huang [30] found no effect of the CYP3A5 genotype on the elimination of MDZ in a mixed population of healthy white and African American adults and Chinese volun- Mean experimental data (symbols: squares) for the first day was taken from Barone et al. [17] and experimental data (mean ± SD, symbols: circles) was taken from Hardin et al. [18]. The inset graph presents the estimated fraction of remaining CYP3A activity over time as a result of ITZ and OH-ITZ inhibition.
Rank-order plots for the optimized 1OH_K m and 1OH_V max values normalized to the mean Figure 6 Rank-order plots for the optimized 1OH_K m and 1OH_V max values normalized to the mean.
teers, respectively. Figure 6 represents the rank-order plot for the optimized 1OH_K m and 1OH_V max values for the Yu et al. [4] data. The result of this figure supports the hypothesis that CYP3A5 plays little to no role in MDZ hydroxylation and is the same result presented in the Yu et al. [4] study in the basal state.
In a second step, a WB-PBPK model for ITZ and OH-ITZ was developed and linked to in vitro data of their CYP3A4 inhibition potential. Through dynamically linking the inhibitor models with MDZ, and thus the MDZ-metabolite models, a prediction of the effect of ITZ given at 200 mg per day over 5 days could be made for each individual from the Yu et al. [4] study. One assumption of our model was that ITZ inhibits total CYP3A activity whereas it is actually specific to the inhibition of CYP3A4. The potential consequences of this were observed when MDZ was modeled in the presence of ITZ and OH-ITZ. There was a slight overprediction of the mean predicted plasma concentrations in relation to those observed between 6 and 12 hours, which was greater in individuals with at least one CYP3A5*1 allele. Since these individuals have a greater proportion of CYP3A5 relative to their total CYP3A content, this suggests that CYP3A5 remains uninhibited by ITZ and OH-ITZ and thus can contribute to MDZ elimination. By comparing the observed and predicted MDZ pharmacokinetic profiles in the basal and inhibited states, the relative contribution of CYP3A5 to overall MDZ elimination in vivo is low regardless of the CYP3A5 genotype. Even in light of this model simplification, the resulting simulations adequately predicted the changes in the MDZ profiles as well as predicting the completely different MDZ metabolite profiles in the inhibited state.
While dynamic simulations of drug-drug interactions are desirable, their conceptualization and implementation are complex compared to the static AUC ratio method, as described above. Major limitations are the availability of pharmacokinetic data for the compound of interest and the inhibitor (especially if multiple metabolites are of concern) as well as a lack of thoroughly studied in vitro inhibition data. Simplification of the model is necessary when this data is not available, which may limit the flexibility and application of the model. In our case, the simplification of the CYP3A complex is one such example. The data that was required to incorporate both a CYP3A4 and CYP3A5 enzyme into the liver was limited by a lack of data on how much each enzyme contributed to 1OH-MDZ plasma concentrations. Further, no individual ITZ plasma concentration time data was available and thus a mean curve had to be generated and used for each individual. Knowing these limitations ahead of gathering the experimental data, may help to inform the clinical drugdrug interaction trial by guiding the planning stage.

Conclusion
This study represents the first attempt to dynamically simulate metabolic enzymatic drug-drug interactions via coupled WB-PBPK models. The workflow described herein, basal state optimization followed by inhibition prediction, is novel and will provide a basis for the development of other inhibitor models that can be used to guide, interpret, and potentially replace clinical drug-drug interaction trials.

WB-PBPK model development and parameterization for defining MDZ CYP3A-mediated metabolism
PK-Sim ® (ver 3.0, Bayer Technology Services GmbH, Leverkusen, Germany) was used to generate the individual WB-PBPK models for MDZ, 1OH-MDZ and 1OH-MDZ-Glu. PK-Sim ® is a commercially available software tool for WB-PBPK modeling of drugs in laboratory animals and humans. It uses validated physiological models to estimate substance-specific absorption [21] and distribution parameters, such as organ/plasma partition coefficients and permeability coefficients, from physicochemical properties of a compound such as lipophilicity, plasma protein binding, molecular weight and solubility [14,20,21]. Physiological databases are included in the software that incorporate the dependencies of organ weights, organ blood flows and intestinal parameters [gastrointestinal length, radius of each section, intestinal surface area [21]] with the weight and height of the individual [31]. For a detailed description of the WB-PBPK model structure implemented in PK-Sim ® , see Willmann et al. [14,20,21].
Individual coupled WB-PBPK models, inclusive of MDZ and its metabolites, were set up. First, the MDZ WB-PBPK model was parameterized by inputting the physico-chemical properties of MDZ (Table 2), the weight of the individual, a mean height and the application regime into PK-Sim ® . The MDZ WB-PBPK model contains a hepatic elimination process, defined in the liver intracellular space, describing the CYP3A mediated elimination of MDZ to its metabolite, 1OH-MDZ. This metabolism process was defined in the model assuming Michaelis-Menten kinetics (Figure 1b) The maximum velocity of the reaction (V max ) depends on the total enzyme concentration [E] total and the rate of the second reaction step: The intrinsic clearance (CL) is the rate constant of metabolization and is defined by the Michaelis-Menten equation: Under the condition that the substrate concentration [S] is much smaller than K m , the intrinsic clearance approaches the ratio of V max and K m : Furthermore, the simplified model includes a CYP3A complex, which represents both iso-enzymes (CYP3A4 and CYP3A5) in the human liver. It is not possible to model the reactions catalyzed by CYP3A4 and CYP3A5 separately, as it is not clearly understood how much each contributes to the production of 1OH-MDZ. Therefore, an assumption has to be made, namely that the CYP3A complex is sufficient to model MDZ-ITZ interactions.
Following generation of a MDZ WB-PBPK model, the WB-PBPK models for 1OH-MDZ and 1OH-MDZ-Glu were parameterized by inputting 1OH-MDZ and 1OH-MDZ-Glu physico-chemical data ( Table 2). All physiological model parameters such as age, weight, height and all depending parameters such as organ weights and bloodflow rates remained the same as in the MDZ model of the same individual. Since 1OH-MDZ is cleared in the liver via glucuronidation, the 1OH-MDZ WB-PBPK model contains a hepatic elimination of 1OH-MDZ to its glucuronide, 1OH-MDZ-Glu via UGT, defined by velocity rate constants (k 1,-1,2 ) as described above in Eq.(1). Elimination of 1OH-MDZ-Glu [CL(Glu)] from the circulation was defined by a first order intrinsic elimination process in the kidney. Coupling of the three WB-PBPK models was done such that the source function generating 1OH-MDZ in the liver intracellular volume was equal to the output of the CYP3A mediated elimination of MDZ (on a molar basis). In the same way, the source function of 1OH-MDZ-Glu was equal to the output of the UGT mediated elimination of 1OH-MDZ within the liver intracellular space.

Enzyme concentration
Hepatic enzyme concentrations in the individuals from where the study data was derived were unknown. Thus, typical enzyme concentrations had to be taken from the literature. The median CYP3A complex concentration of a mixed and a Japanese population was used as the enzyme concentration for the hydroxylation step in this study, such that [E 0 ] equaled 70 pmol/mg microsomal protein [7,23,24]. The value 40 mg protein/g liver [32][33][34]

Volume of distribution
The volume of distribution is an important parameter in pharmacokinetic studies. Distribution volumes in PK-Sim ® are estimated based on physico-chemical data ( Table  2), as described in Willmann et al. [20]. Membrane affinity, as a lipophilicity value, and plasma unbound fraction were assessed for MDZ and 1OH-MDZ using the Nimbus Technology (Nimbus Biotechnology, Leipzig, Germany) [35,36]. 1OH-MDZ-Glu is not commercially available and there exists no published lipophilicity or unbound fraction in plasma value. Therefore, the lipophilicity of 1OH-MDZ-Glu was estimated by using the mean logarithmic difference of the lipophilicities of other hydroxylated compounds and their corresponding glucuronides that were taken from the literature [37] ( Table 2). The glucuronide of valproic acid has an ex vivo unbound fraction in serum that is half that of the parent compound [38]. In another study that examined the protein binding of three compounds and their glucuronides, the unbound fraction of the glucuronide increased by factors of 2.2, 1.3 and 2.6 [39]. For the purposes of the present study, the fraction unbound in plasma for 1OH-MDZ-Glu was kept the same as that of 1OH-MDZ (Table 2).
Plasma concentration time data was located for individuals who received either MDZ or 1OH-MDZ intravenously [15] and experimental data were compared to the curves predicted by the WB-PBPK model using estimated partition coefficients. Because there was a slight deviation of predicted to experimental plasma concentration profiles, the partition coefficients estimated from physico-chemistry for MDZ were compared to experimental partition coefficients, observed in rats and converted for use in humans, for MDZ [22]. With the carcass value being taken as that for bone, the correlation revealed that all predicted partition coefficients agreed with the experimental data within a factor of three with the exception of blood cells, bone, fat and lung, (Figure 8). For these tissues, the experimental partition coefficients were used in the MDZ simulations and the same predicted to experimental conversion factors were used for these tissues in the WB-PBPK models for 1OH-MDZ and 1OH-MDZ-Glu. An independent validation of the 1OH-MDZ-Glu distribution volume was not possible due to a lack of plasma concentration time data following 1OH-MDZ-Glu administration.

Optimization of velocity-rate Constants for the hydroxylation of Midazolam
All relevant model parameters are now defined as described above with the exception of the biotransformation parameters. To generate velocity rate constants for the elimination of MDZ to 1OH-MDZ, individual plasma concentration time data for MDZ and its metabolite 1OH-MDZ following intravenous administration were gathered from Heizmann et al. [11]. Six normal healthy volunteers with reported body weight and age were enrolled in this study.
Comparison of observed and predicted tissue/plasma partition coefficients for the various organs (left: log-log plot, right: enlarged linear plot of the interval [0,3.5])

Figure 8
Comparison of observed and predicted tissue/plasma partition coefficients for the various organs (left: log-log plot, right: enlarged linear plot of the interval [0,3.5]). Small symbols indicate individual data reported by Björkmann et al. [22], large symbols denote the mean value for each organ. The solid line represents the identity, the dotted lines mark the region 3-fold-off the identity.
For each individual from the Heizmann et al. [11] study, WB-PBPK models were generated for MDZ and 1OH-MDZ. These two models were coupled and consequently contained four unknown parameters, 1OH_k 1 , 1OH_k -1 and 1OH_k 2 , which describe MDZ hydroxylation by CYP3A4/5 as in Eq. (1), and a first order elimination rate of OH-MDZ via an intrinsic clearance CL(1OH), which describes the glucuronidation of 1OH-MDZ. These four parameters were simultaneously optimized per individual, without constraints, in order to obtain the best fit of the simulated concentration time curves for MDZ and 1OH-MDZ in plasma to the respective experimental data. The objective function was the root mean squared error (RMSE) measured as the root of the squared difference between the predicted and observed plasma concentrations for both compounds. Objective function minimization was done using the fminsearch optimization routine of MatLab [Version 7 (R14), The MathWorks Inc., Natick, MA, USA].

Optimization of velocity-rate constants for the Glucuronidation of 1-hydroxymidazolam
Plasma concentration time data for MDZ, its metabolite 1OH-MDZ and its conjugated metabolite 1OH-MDZ-Glu following intravenous bolus MDZ administration were gathered from Kharasch et al. [16]. Data from one typical representative was available and used here. The coupled WB-PBPK model for the Kharasch et al. [16] data contained seven unknown parameters, 1OH_k 1 , 1OH_k -1 and 1OH_k 2 , Glu_k 1 , Glu_k -1 and Glu_k 2 , and CL(Glu). These seven parameters were simultaneously optimized in order to obtain the best fit of simulated concentration time curves for MDZ, 1OH-MDZ and 1OH-MDZ-Glu to the corresponding experimental data. For the optimization of the elimination of MDZ to its metabolite 1OH-MDZ, the results of the previous optimization were used as constraints by setting the min and max constraints for 1OH_ k 1,-1,2 as half of and two times the min and max values, respectively, as determined during the optimization of the six Heizmann et al. [11] individuals. The mean values of the optimization of the Heizmann et al. [11] data were used as initial parameters for 1OH_k 1,-1,2 . The optimization of the velocity-rate constants for the metabolization of 1OH-MDZ to its metabolite 1OH-MDZ-Glu were unconstrained, as was the first order elimination of 1OH-MDZ-Glu [CL(Glu)].

Optimization of the velocity rate constants for the Yu et al (2004) study data in the Basal state
Individual plasma concentration time data for MDZ and the sum of its metabolites, 1OH-MDZ and 1OH-MDZ-Glu, in the basal state were taken from Yu et al. [4]. In this study, plasma concentrations of MDZ and 1OH-MDZ were determined (after addition of glucuronidase to the plasma samples, i. e. 1OH-MDZ concentrations represent the sum of 1OH-MDZ and 1OH-MDZ-Glu) for nineteen healthy adults having either the CYP3A5 *1/*1, *1/*3 or *3/*3 genotype. Individual values for height, body weight, sex, age and genotype were also provided [4].
A coupled WB-PBPK model was set up for each individual in the Yu et al. [4] study. Each model contained seven unknown parameters, 1OH_k 1 , 1OH_k -1 and 1OH_k 2 , Glu_k 1 , Glu_k -1 and Glu_k 2 and CL(Glu). These seven parameters were again simultaneously optimized as described above in order to obtain the best fit of simulated concentration time curves for MDZ and the sum of 1OH-MDZ and 1OH-MDZ-Glu to the corresponding experimental data. The same constraints as above were used for the optimization of the elimination of MDZ to 1OH-MDZ. The initial values for optimization for Glu_ k 1,-1,2 are the results of the optimization of the Kharasch et al. data [16]. As Kharasch et al. [16] presented data from only one individual, there were no minima and maxima which could be used as constraints, therefore no constraints were used for the 1OH-MDZ-Glu parameters.

Itraconazole CYP3A4 inhibition model
Itraconazole (ITZ) and its metabolites competitively inhibit CYP3A4 [13]. The major, and only significant, metabolite in vivo is hydroxyitraconazole (OH-ITZ) which is found at concentrations higher than ITZ in plasma [17].
First, a WB-PBPK model for ITZ was developed using the physicochemical properties as presented in Table 2. Elimination, representing CYP3A metabolism, was defined as a Michaelis-Menten process in the liver with a K m of 3.9 ηM [13] and a V max to fit the terminal phase of the plasma concentration profile from Barone et al. [17].
According to the study protocol in Yu et al. [4], the volunteers received a single oral dose of a commercial ITZ tablet, Sporanox ® (Janssen Pharmaceutica), in the fasted state over five days. Since the pharmacokinetics of ITZ and OH-ITZ were not assessed by Yu et al., the literature was searched for similar pharmacokinetic data. Barone et al. [17] provided plasma concentration time data for Sporanox ® tablets, but administered in the fed state. ITZ is a basic compound with a pKa of 3.7, and it is well known that administration of ITZ together with food increases the bioavailability due to a longer retention of the tablet in the stomach. Consequently, a modified dissolution profile had to be generated for the Sporanox ® tablet in the fasted state. Sporanox ® tablets are well dissolved in acidic solutions and, within 45 minutes, 50% are dissolved in an artificial stomach solution [40]. In contrast, in Fed State Simulating Intestinal Fluid (FeSSIF), Sporanox ® tablets reach a maximum fraction dissolved of 2% in 30 minutes [40]. In the WB-PBPK model, dissolution was set as 9, 35, 50% dissolved in 15, 30 and 45 minutes (the average length of time the tablets reside in the stomach in the fed state) to 52% in 75 minutes (2% being the maximum dissolution in the intestines). Further, ITZ is metabolized by CYP3A in the intestine. To account for gut wall metabolism, in each intracellular section of the gastrointestinal tract (duodenum, jejunum, ileum), a Michaelis-Menten process in the gut wall was assumed with a K m the same as that for hepatic metabolism (3.9 nM). A V max value for gut wall metabolism was fit to best represent the absorption phase of the Barone et al. [17] data. Only one V max required fitting since the proportions of CYP3A in each intestinal section are known (duodenum: jejunum: ileum = 9.5: 6.3: 1) [41]. The resulting single dose administration simulation was fit to experimental data from Barone et al. [17].
ITZ and OH-ITZ are primarily metabolized by CYP3A4.
Since both have a high affinity for CYP3A4, saturation of the enzymes can lead to non-linear kinetics. Following multiple administrations, ITZ and OH-ITZ clearances are reduced and the concentrations in plasma of the two substances increase as well as the half-lives. In the Yu et al. [4] study, ITZ was administered once-daily for five days and on the fifth day given in conjunction with an IV midazolam bolus. To simulate the likely ITZ concentrations over 6 days, the time dependence of ITZ elimination was assessed using the data of Hardin et al. [18]. In this study, trough plasma concentrations were taken at days 2, 7, 11 and 13 during 200 mg once-daily administrations of ITZ in the fed-state. V max at time 0 was the same as that previously defined from the single dose study of Barone et al. [17]. The V max was reduced daily to fit the trough concentrations from Hardin et al. [18] (Figure 7). The next step was to generate the fasted state ITZ and OH-ITZ liver unbound concentration time profiles that will be used to calculate the time specific CYP3A activity reduction (see next section). Since most literature data providing pharmacokinetic profiles of ITZ and OH-ITZ were those following ITZ administration in the fed state, from the multiple dose ITZ fed-state profiles (as described above), both the ITZ and OH-ITZ fasted-state profiles needed to be generated. Since the difference between the fed and fasted state plasma profiles are mainly due to differences in tablet dissolution, concentration differences over 24 hours are primarily the ratio of AUCs (AUC fed = 1.7*AUC fasted ) [17]. This ratio remained the same throughout a 15 day BID administration regimen [17]. Therefore, the unbound concentrations in the liver intracellular space were taken from the WB-PBPK model for ITZ fed and divided by 1.7 to generate the relevant concentrations for ITZ fasted . For OH-ITZ fasted , the AUC ratio of ITZ fed to OH-ITZ fasted was used to generate concentrations where AUC OH-ITZ,fasted = 1.6* AUC ITZ,fed [17]. The relevant liver concentrations were therefore the ITZ fed unbound liver concentrations multiplied by 1.6 and divided by the ratio of unbound fractions (f u,ITZ /f u,OH-ITZ = 0.8). The unbound fraction of OH-ITZ is not available in the literature therefore the same ratio of parent to hydroxylated f u , as was determined for MDZ (Table 2), was used.
The effect of the unbound ITZ and OH-ITZ concentrations on CYP3A4 activity in vitro using MDZ hydroxylation as an endpoint has been assessed by Isoherranen et al. [13] and is presented in Figure 9. The unbound ITZ and OH-ITZ concentrations in the intracellular space of the liver, as obtained by the WB-PBPK model, were used in conjunction with Figure 9 to generate a remaining percentage of CYP3A activity. At each time step, a new CYP3A concentration ([CYP3A i ]) was calculated using the following equation: [CYP3A i ] = [CYP3A]*ITZ frac·act *OH -ITZ frac·act (6) Here, ITZ frac·act and OH-ITZ frac·act are the fraction of the initial CYP3A activity based on ITZ and OH-ITZ concentrations from Figure 9 (constrained to values between 0 and 1).

Prediction of MDZ and MDZ metabolite pharmacokinetic profiles following ITZ administration
Coupling of the previously established MDZ and ITZ models (a total of 5 WB-PBPK models) finally facilitates Fraction of remaining CYP3A activity in pooled human liver microsomes in relation to the unbound inhibitor concentra-tions of ITZ and OH-ITZ Figure 9 Fraction of remaining CYP3A activity in pooled human liver microsomes in relation to the unbound inhibitor concentrations of ITZ and OH-ITZ. Data points, as taken from Isoherranen et al. [13], were fit to a logistic, 3 parameter model and the equations were used to generate a percentage of CYP3A activity remaining, constrained between 0 and 1.
the prediction of the effect of ITZ and OH-ITZ CYP3A inhibition in the pharmacokinetics of MDZ and its metabolites. The coupling was done such that the clearance of MDZ was affected by the changing concentrations of free CYP3A [E o ] as defined from the ITZ inhibition model ( Figure 10). WB-PBPK Simulations were performed for each of the individuals in the Yu et al. [4] study. All velocity-rate constants as well as CL(Glu) were fixed to the optimized values determined previously in the simulations for the basal state for each individual. Results were compared to the MDZ and 1OH-MDZ + 1OH-MDZ-Glu plasma concentration profiles. In summary, the assumptions of this combined model were: 1/ MDZ and metabolite elimination kinetics can be adequately described using Michaelis-Menten kinetics.
2/ The 1'-hydroxylation step is the only significant route of elimination for midazolam from the body. 4OH-MDZ and 1,4-di-OH-MDZ plasma pharmacokinetic profiles are not published since they are not found in significant amounts in plasma [15]. Furthermore, as the percentage of administered drug excreted in the urine as 4OH-MDZ and 1,4di-OH-MDZ (or their glucuronides) is less than 4%, suggesting that these MDZ elimination pathways are not significant.
3/ The physico-chemical properties of the compounds were adequate to describe the compound's distribution volume. This was independently evaluated for MDZ, 1OH-MDZ and ITZ using experimental plasma concentration time data following IV administration and using some MDZ experimental partition coefficients. Vmax Maximum velocity of an enzymatic reaction WB-PBPK Whole body physiologically-based pharmacokinetic system to perform the simulations, CN and SW provided valuable inputs into the conceptual design of the study and helped to draft the manuscript, IJ provided experimental data and data interpretation for the study, and ANE conceived of the study and designed and performed the inhibited state MDZ-ITZ simulations. All authors read and approved the final manuscript.