Comparisons of Four Protein-Binding Models Characterizing the Pharmacokinetics of Unbound Phenytoin in Adult Patients Using Non-Linear Mixed-Effects Modeling.

BACKGROUND AND
OBJECTIVE: Phenytoin is extensively protein bound with a narrow therapeutic range. The unbound phenytoin is pharmacologically active, but total concentrations are routinely measured in clinical practice. The relationship between free and total phenytoin has been described by various binding models with inconsistent findings. Systematic comparison of these binding models in a single experimental setting is warranted to determine the optimal binding behaviors.
METHODS: Non-linear mixed-effects modeling was conducted on retrospectively collected data (n = 37 adults receiving oral or intravenous phenytoin) using a stochastic approximation expectation-maximization algorithm in MonolixSuite-2019R2. The optimal base structural model was initially developed and utilized to compare four binding models: Winter-Tozer, linear binding, non-linear single-binding site, and non-linear multiple-binding site. Each binding model was subjected to error and covariate modeling. The final model was evaluated using relative standard errors (RSEs), goodness-of-fit plots, visual predictive check, and bootstrapping.
RESULTS: A one-compartment, first-order absorption, Michaelis-Menten elimination, and linear protein-binding model best described the population pharmacokinetics of free phenytoin at typical clinical concentrations. The non-linear single-binding-site model also adequately described phenytoin binding but generated larger RSEs. The non-linear multiple-binding-site model performed the worst, with no identified covariates. The optimal linear binding model suggested a relatively high binding capacity using a single albumin site. Covariate modeling indicated a positive relationship between albumin concentration and the binding proportionality constant.
CONCLUSIONS: The linear binding model best described the population pharmacokinetics of unbound phenytoin in adult subjects and may be used to improve the prediction of free phenytoin concentrations.

Key Points

Introduction

Phenytoin has been one of the most widely prescribed anticonvulsants since its initial approval in the 1950s, although its usage has been slowly declining with the availability of newer-generation antiepileptic drugs with lower risks of drug–drug interactions [1-3]. Phenytoin is highly bound to plasma albumin (~ 90%) in healthy adults [4, 5], where unbound phenytoin is primarily eliminated by hepatic intrinsic clearance mediated by cytochrome P450 (CYP)-2C9 and CYP2C19 [5-8]. Only minimal amounts of phenytoin (~ 1.2%) are excreted unchanged into the urine [9]. The CYP2C9/19-mediated oxidation of phenytoin is usually saturated at the typical phenytoin plasma concentrations observed clinically, leading to non-linear (i.e., Michaelis–Menten) elimination [6, 10]. Moreover, phenytoin is a broad spectrum inducer of metabolism enzymes [11, 12], which can increase the intrinsic clearance of coadministered drugs.

Therapeutic drug monitoring of phenytoin is recommended because of its narrow therapeutic range (i.e., 10–20 mg/L based on total concentrations), non-linear kinetics, and potential for drug–drug interactions [4, 13]. Loss of antiepileptic effect is typically evident at total phenytoin serum concentrations < 10 mg/L [14], whereas phenytoin-induced neurotoxicity (e.g., nystagmus, fatigue) can usually be observed in patients with total serum concentrations > 20 mg/L [15, 16]. Hence, a delicate balance within this narrow therapeutic window is required. However, the attainment of phenytoin concentrations within the narrow target can be complicated by non-linear kinetics and other factors known to influence phenytoin’s protein binding or intrinsic clearance (e.g., [4, 10, 13, 17]).

Pharmacologically active unbound phenytoin is the ideal marker for therapeutic drug monitoring. However, free phenytoin concentrations are not routinely measured because of the complexity of the analysis and the higher costs [13, 18]. Therefore, phenytoin therapeutic drug monitoring is commonly conducted with the total concentration, and free concentrations are only estimated in patients with renal dysfunction and/or reduced albumin [13, 18]. Several approaches have been proposed to describe the relationship between free and total phenytoin by characterizing the binding of free phenytoin with serum albumin. The most widely used approach is the Winter–Tozer equation, which has been modified by many others to better tailor to individual patient populations [19-23]. Moreover, non-linear mixed-effects modeling has been used to predict unbound phenytoin concentrations and characterize the nature of protein binding. For example, Hennig et al. [24] used a proportional linear binding model in critically ill pediatric subjects, ter Heine et al. [25] validated a non-linear binding model in various populations, and various investigative groups examined clinical covariates influencing the protein binding or free fraction of phenytoin [24-27]. However, to our knowledge, individual phenytoin protein-binding models have never been directly compared in a single experimental setting, so the optimal binding model remains equivocal. The objective of this study was to systematically compare four individual protein-binding models (i.e., Winter–Tozer, linear binding, non-linear single-binding site, and non-linear multiple-binding site) within a fully developed population pharmacokinetic model to identify the optimal albumin-binding characteristics of phenytoin and the influence of clinical covariates.

Methods

Patients and Study Design

The study protocol was approved by the University of British Columbia (H18-02215) and the University of Alberta (Pro00100357) research ethics boards. Non-linear mixed-effects modeling was conducted on retrospectively collected data from Vancouver General Hospital in Vancouver, British Columbia, Canada. Data from patients aged > 18 years receiving oral or intravenous phenytoin were included. Paired measurements of free and total phenytoin concentrations were extracted from the hospital electronic database between 31 August 2014 and 31 August 2018. All included patients received at least 5 days of the same phenytoin regimen. Patients were excluded if they were admitted to the emergency department at the time of receiving phenytoin or lacked clinical laboratory values. A total of 37 patients with 37 paired free and total phenytoin concentrations were included in this study. The following data were collected: sex, critical care status (i.e., admission to the general or neurosurgical intensive care units), age, actual body weight, serum albumin, alanine transaminase (ALT), aspartate transaminase (AST), total bilirubin, international normalized ratio (INR), serum creatinine (SCr), hemodialysis status, and concomitant medications (i.e., low-dose aspirin, carbamazepine, heparin, phenobarbital, sulfonamides, valproic acid, or warfarin). The majority of laboratory values were gathered on the same day, and all values were obtained within 1 week of phenytoin concentration collection. Glomerular filtration rate (GFR) was estimated using the chronic kidney disease epidemiology collaboration equation [28]. Phenytoin plasma samples were collected 3–24 h after dose administration. Free and total phenytoin plasma concentrations were analyzed at the Vancouver General Hospital clinical laboratory as part of routine care using a validated immunoassay [29]. Throughout the entire manuscript, all phenytoin concentrations are designated as the free (unbound) concentration unless otherwise specified (i.e., total, bound).

Model Development

Software

Modeling was conducted using a stochastic approximation expectation–maximization algorithm in MonolixSuite-2019R2 (Lixoft SAS, Antony, France). For covariate and bootstrapping analyses, we used the Rsmlx package (Rsmlx, version 2.0.2, Xpop Inria Team) with R software v3.6.3 (R Core Team, Vienna, Austria).

Base Model Selection

First-order extravascular absorption with and without lag time or transit compartments were tested to describe phenytoin oral absorption. One- or two-compartment structural models were evaluated to obtain the best fit for free phenytoin concentration–time profiles. To describe the elimination of free phenytoin, we examined linear and Michaelis–Menten models. Based on the literature [19–25, 30–32], protein-binding relationships of phenytoin could be characterized by four types of binding-dissociation models, described as follows:

(1) Empiric Winter–Tozer equation [19] with an albumin coefficient [20-23],where (mg/L) is the free phenytoin concentration, (mg/L) is the total concentration, is the albumin coefficient, and is expressed in g/dL.

(2) Linear binding model with the assumption that is directly proportional to with no saturable mechanism [24, 33],where (mg/L) is the bound concentration and is the binding proportionality constant (i.e., Cbound/Cfree).

(3) Non-linear binding model with the assumption that phenytoin binds to a single site on albumin [25, 31],where (mg/L) represents the maximal binding capacity and (mg/L) defines the equilibrium dissociation constant (i.e., free concentration at which 50% of the binding site is occupied).

(4) Non-linear binding model with the assumption that phenytoin binds to more than one saturable binding site [30-32],where and are the binding-dissociation parameters of the th binding site of albumin.

Combined, proportional, or constant error models with normal, lognormal, or logit-normal distributions were tested. The most suitable structural, protein binding, and error models were selected based on biological plausibility, objective function values (OFVs), Akaike information criterion (AIC), Bayesian information criterion (BIC), relative standard errors (RSEs), and η-shrinkage values. A reduction in OFVs by more than 3.84 (i.e., p < 0.05) was considered a significant improvement in modeling.

Covariate Model Selection

The four protein-binding models were further subjected to covariate analyses under identical pre-conditions (i.e., same structural model, error model, and initial parameter estimates) to identify the best final protein-binding model. The effects of an oral formulation (extended-release capsule or suspension) and clinical factors listed in Sect. 2.1 were assessed. For continuous factors, both original and natural log-transformed values were tested. The relationship between population pharmacokinetic parameters and covariates are provided in electronic supplementary material (ESM) 1 [34]. Clinical factors were identified as potential covariates by Pearson’s correlation test (p < 0.05), the Wald test (p < 0.05), or stepwise forward inclusion using the R software (e.g., [35, 36]). All potential covariates were included in the base model and further subjected to stepwise backward elimination. The biological plausibility of all included covariates was assessed [4, 5, 13, 37, 38].

Model Validation

The final model was internally validated using goodness-of-fit plots, prediction-corrected visual predictive checks (pcVPC), and bootstrap analysis, as described previously [35, 36, 39, 40] (summarized in ESM 1).

Results

Sample Population

Data from 37 adult patients (mean age 61.1 ± standard deviation 17.9 years) were available for analysis (Table 1). The majority (73.0%) of the patients were male and classified as non-critical care (56.8%). The average weight, AST, total bilirubin, INR, and SCr values (Table 1) were within the normal reference ranges [41, 42]. However, a substantial portion of patients (32.4%) had increased ALT (i.e., > 53 IU/L), and a small portion of patients (5.4%) had severe renal dysfunction (i.e., GFR < 30 mL/min/1.73 m2). Moreover, the majority of patients (94.6%) had mild hypoalbuminemia (i.e., serum albumin < 3.5 g/dL), but the unbound fraction (12.4 ± 3.1%) was within the expected range [13]. On average, patients received 378.5 ± 148.3 mg/day of phenytoin in one to three divided doses, corresponding to a standard maintenance dose of 200–600 mg/day [43]. The reported total and free plasma phenytoin concentrations were also within the proposed monitoring targets (Table 1) [4, 18]. A small number of patients received co-medications (Table 1).

Table 1

Patient demographics (n = 37)

Parameter	Mean ± SD or count	Median (range)
Sex (female/male)	10/27
Critical care (Y/N)a	16/21
Age (years)	61.1 ± 17.9	62 (20–93)
Weight (kg)	68.5 ± 15.6	70 (30–102)
Albumin (g/dL)	2.64 ± 0.52	2.7 (1.7–3.6)
ALT (IU/L)	71.6 ± 116.2	36 (13–693)
AST (IU/L)	42.8 ± 37.1	28 (13–161)
Total bilirubin (μmol/L)	14.6 ± 39.9	6 (2–249)
INR	1.1 ± 0.1	1.0 (0.8–1.5)
SCr (μmol/L)	96.8 ± 91.6	78 (29–570)
Hemodialysis (Y/N)b	1/34
Current medicationsb
Aspirin (Y/N)c	10/26
Carbamazepine (Y/N)	1/23
Heparin (Y/N)	1/34
Phenobarbital (Y/N)	1/34
Sulfonamides (Y/N)	1/34
Valproic acid (Y/N)	2/33
Warfarin (Y/N)	1/34
Phenytoin dosage and measurements
Phenytoin dose (mg/day)	378.5 ± 148.3	300 (120–850)
Administration route (IV/PO)	16/21
Time between dose/test (h)	12.9 ± 7.3	10.5 (3.0–24.0)
Total phenytoin concentration (mg/L)	11.4 ± 5.3	9.8 (5.3–23.5)
Free phenytoin concentration (mg/L)	1.4 ± 0.7	1.1 (0.5–3.5)
Unbound fraction (%)	12.4 ± 3.1	11.8 (7.1–19.4)

ALT alanine transaminase, AST aspartate transaminase, INR international normalized ratio, IV intravenous, N no, PO oral, SCr serum creatinine, SD standard deviation, Y yes

aPatients were considered under the category of “critical care” when admitted to either the general or the neurosurgical intensive care unit

bRecords were missing for some patients

cDose confirmed as 80 mg orally daily in nine of ten subjects

Base Model Selection

Several combinations of structural, protein binding, and error models were tested (Table 2). A systematic approach was implemented to identify the most suitable base model: (1) initial search of the best model without considering protein binding (i.e., models 1–12 in Table 2), (2) screening of four different binding models that best complemented the identified structural model (i.e., models 13–16), and (3) testing for the optimal residual error models (i.e., models 13–28). For the structural model, the Michaelis–Menten elimination model was used because the free phenytoin concentrations in our dataset (1.4 ± 0.7 mg/L) exceeded the literature-estimated (0.40–0.73 mg/L [26, 44]) and model-estimated Km values (Table 3). Furthermore, a one-compartment structural model with first-order absorption without lag time best complemented the Michaelis–Menten elimination model (model 7) in our dataset based on likelihood estimations. Proportional errors were superior to both constant and combined models based on likelihood estimations (models 13–28) and goodness-of-fit plots (data not shown).

Table 2

Base model selection history

Model	Structural model			Protein-binding model	Error model	OFV	AIC	BIC
	Elimination	Distribution	Absorption
1	Linear	One-compartment	First-order with no delay	NA	Proportional	72.08	82.08	90.14
2	Linear	One-compartment	First-order with lag time	NA	Proportional	72.12	86.12	97.39
3	Linear	One-compartment	Transit-compartment	NA	Proportional	70.84a	88.84a	103.33a
4	Linear	Two-compartment	First-order with no delay	NA	Proportional	71.58	89.58	104.08
5	Linear	Two-compartment	First-order with lag time	NA	Proportional	71.51a	93.51a	111.23a
6	Linear	Two-compartment	Transit-compartment	NA	Proportional	70.95a	96.95a	117.89a
7	Michaelis–Menten	One-compartment	First-order with no delay	NA	Proportional	89.21	103.21	114.48
8	Michaelis–Menten	One-compartment	First-order with lag time	NA	Proportional	87.02a	105.02a	119.52a
9	Michaelis–Menten	One-compartment	Transit-compartment	NA	Proportional	87.40a	109.40a	127.12a
10	Michaelis–Menten	Two-compartment	First-order with no delay	NA	Proportional	89.24a	111.24a	128.96a
11	Michaelis–Menten	Two-compartment	First-order with lag time	NA	Proportional	87.39a	113.39a	134.33a
12	Michaelis–Menten	Two-compartment	Transit-compartment	NA	Proportional	86.67a	116.67a	140.84a
13	Michaelis–Menten	One-compartment	First-order with no delay	Winter–Tozer	Proportional	244.00	264.00	280.11
14	Michaelis–Menten	One-compartment	First-order with no delay	Linear binding	Proportional	259.54	279.54	295.65
15	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear single-binding site	Proportional	258.31	282.31	301.64
16	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear multiple-binding site	Proportional	256.21	280.21	299.54
17	Michaelis–Menten	One-compartment	First-order with no delay	Winter–Tozer	Combined 1	238.02a	264.02a	284.96a
18	Michaelis–Menten	One-compartment	First-order with no delay	Linear binding	Combined 1	262.07a	286.07a	305.41a
19	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear single-binding site	Combined 1	253.33a	281.33a	303.88a
20	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear multiple-binding site	Combined 1	254.20a	282.20a	304.75a
21	Michaelis–Menten	One-compartment	First-order with no delay	Winter–Tozer	Combined 2	242.40a	268.40a	289.34a
22	Michaelis–Menten	One-compartment	First-order with no delay	Linear binding	Combined 2	263.18a	287.18a	306.51a
23	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear single-binding site	Combined 2	254.96a	282.96a	305.51a
24	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear multiple-binding site	Combined 2	260.93a	288.93a	311.48a
25	Michaelis–Menten	One-compartment	First-order with no delay	Winter–Tozer	Constant	245.84	267.84	285.56
26	Michaelis–Menten	One-compartment	First-order with no delay	Linear binding	Constant	260.23	280.23	296.34
27	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear single-binding site	Constant	257.24a	281.24a	300.58a
28	Michaelis–Menten	One-compartment	First-order with no delay	Non-linear multiple-binding site	Constant	258.35a	282.35a	301.68a

AIC Akaike information criterion, BIC Bayesian information criterion, NA not applicable, OFV objective function value

aThe fisher information matrix was not correctly estimated (i.e., the models failed to converge using stochastic approximation expectation–maximization algorithm)

Table 3

Population pharmacokinetic parameters of the selected (Winter–Tozer and linear binding) base and final models

Models	Winter–Tozer model						Linear binding model
Parameters	Base modela			Final model			Base modelb			Final model
	Estimated mean value (RSE%)	η (%)	Bootstrap mean(95% CI)	Estimated mean value (RSE%)	η (%)	Bootstrap mean (95% CI)	Estimated mean value (RSE%)	η (%)	Bootstrap mean (95% CI)	Estimated mean value (RSE%)	η (%)	Bootstrap mean(95% CI)
Fixed effects
F	0.859 fix			0.859 fix			0.859 fix			0.859 fix
ka (h−1)	0.225 fix			0.225 fix			0.225 fix			0.225 fix
Vd (L)	1620 (35.9)		1620 (1240–2190)	1650 (16.6)		1650 (1200–2140)	1700 (34.2)		1700 (1270–2150)	1650 (15.0)		1650 (1340–2750)
Vmax (mg/h)	10.8 (85.5)		10.8 (9.55–15.1)	11.7 (22.0)		11.7 (9.67–15.2)	12.6 (66.9)		12.6 (9.86–15.2)	15.3 (33.3)		15.3 (12.1–27.5)
βweight_Vmax	NA			2.11 (37.4)		2.11 (0.445–5.72)	NA			2.04 (29.0)		2.04 (0.608–8.47)
βbilirubin_Vmax	NA			NA			NA			− 0.0355 (11.3)		− 0.0355 (− 0.456 to − 0.0136)
Km (mg/L)	0.541 (144)		0.541 (0.516–0.595)	0.589 (91.3)		0.589 (0.508–0.597)	0.581 (301)		0.581 (0.512–0.600)	0.532 (78.6)		0.532 (0.496–0.569)
CF	0.281 (4.15)		0.281 (0.259–0.303)	0.276 (4.69)		0.276 (0.257–0.302)	NA			NA
Kb	NA			NA			7.26 (6.81)		7.26 (6.49–7.96)	7.66 (4.82)		7.66 (6.87–8.12)
βalbumin_Kb	NA			NA			NA			0.799 (19.0)		0.799 (0.358–1.07)
Bmax (mg/L)	NA			NA			NA			NA
βalbumin_Bmax	NA			NA			NA			NA
KD (mg/L)	NA			NA			NA			NA
Interindividual variability
ω_Vd	0.559 (28.3)	− 7.02	0.559 (0.122–0.652)	0.513 (33.0)	− 7.49	0.513 (0.127–0.631)	0.517 (65.6)	1.19	0.517 (0.0806–0.652)	0.417 (35.5)	− 10.2	0.417 (0.143–0.535)
ω_Vmax	0.450 (75.2)	14.3	0.450 (0.124–0.903)	0.263 (61.1)	4.42	0.263 (0.0948–0.599)	0.385 (87.8)	18.8	0.385 (0.122–0.883)	0.377 (64.3)	14.9	0.377 (0.0450–0.805)
ω_Km	0.890 (68.7)	9.87	0.890 (0.461–5.21)	0.847 (152)	− 7.03	0.847 (0.466–5.99)	1.57 (293)	4.21	1.57 (0.537–6.52)	0.985 (57.4)	− 2.35	0.985 (0.488–21.1)
ω_CF	0.156 (50.5)	12.6	0.156 (0.0545–0.195)	0.134 (85.7)	− 2.35	0.134 (0.0561–0.194)	NA			NA
ω_Kb	NA			NA			0.192 (78.5)	− 4.32	0.192 (0.0608–0.247)	0.0671 (67.8)	− 0.0426	0.0671 (0.0365–0.186)
ω_Bmax	NA			NA			NA
ω_KD	NA			NA			NA
Residual variability
b1	0.108 (72.4)		0.108 (0.0355–0.157)	0.0883 (70.1)		0.0883 (0.0340–0.157)	0.0575 (80.7)		0.0575 (0.0332–0.156)	0.147 (30.5)		0.147 (0.0317–0.161)
b2	0.109 (118)		0.109 (0.0499–0.217)	0.146 (65.5)		0.146 (0.0465–0.223)	0.160 (74.3)		0.160 (0.0600–0.256)	0.122 (42.4)		0.122 (0.0509–0.213)
Likelihood estimations
OFV	244.00			237.16			259.54			232.48
AIC	264.00			259.16			279.54			258.48
BIC	280.11			276.88			295.65			279.42

AIC Akaike information criterion, BIC Bayesian information criterion, B maximal binding capacity, b1 proportional error for total phenytoin concentrations, b2 proportional error for free phenytoin concentrations, β covariate parameter estimate, CF albumin coefficient, CI confidence interval, F bioavailability, η η-shrinkage, k absorption rate constant, K binding proportionality constant, K equilibrium dissociation constant, K Michaelis–Menten constant, NA not applicable, OFV objective function value, RSE relative standard error, V volume of distribution, V maximal elimination rate, ω interindividual variability

aModel 13 in Table 2

bModel 14 in Table 2

The structural model (i.e., model 7 in Table 2) was parameterized as bioavailability (F), absorption rate constant (ka, h−1), volume of distribution (Vd, L), maximal elimination rate (Vmax, mg/h), and Michaelis–Menten constant (i.e., free phenytoin concentration at which the elimination rate reaches 50% of Vmax) (Km, mg/L). Moreover, additional parameters tailored to each protein-binding model were considered: CF for the Winter–Tozer model (model 13); Kb for the linear binding model (model 14); and Bmax and KD for the non-linear binding models (models 15 and 16). Neither F nor ka could be precisely estimated in our dataset because of the limited sample size, and these values were fixed to 0.859 and 0.225 h−1, respectively (e.g., [24, 45]). To obtain biologically reasonable parameter estimates, Vmax and Km were initially fixed to 25 ± 0.2 mg/h and 0.5 ± 0.1 mg/L, respectively, using reference data collected in prior studies [26, 44]. Because of the computational complexity of the non-linear multiple-binding-site model, two KD values were also fixed to 38.8 and 24.3 mg/L [30].

Covariate Model Selection

Covariate analyses were conducted on all four protein-binding models. For the linear binding model, log-transformed weight (positive relationship) and total bilirubin (negative relationship) independently affected Vmax, whereas log-transformed albumin positively affected Kb. The inclusion of these covariates reduced the overall OFV by 27.06 compared with the base model. OFV was increased by 14.08 (p < 0.001), 6.21 (p < 0.05), and 6.77 (p < 0.01) when albumin, weight, and total bilirubin were sequentially eliminated. Vmax and Kb can be estimated as follows:where and are the estimated Vmax and Kb of the th individual with (kg), (μmol/L), and (g/dL). represents the random effects of the th individual. The results of covariate analyses conducted on other protein-binding models are provided in ESM 1 [21-23].

Final Model

The population pharmacokinetic parameters from each model are summarized in Tables 3 and 4. Significantly higher OFV, AIC, and BIC values observed in the final non-linear multiple-binding site model indicated poor fitting. AIC and BIC were comparable in other protein-binding models, but the linear binding model showed significantly lower OFV (ΔOFV =  − 4.68, p < 0.05) than the Winter–Tozer model and relatively smaller RSEs in fixed-effects estimates compared with the non-linear single-binding-site model. Overall assessments of likelihood estimations, physiological relevance, and RSEs indicated that the linear binding model best described the protein binding of phenytoin in this adult population (please see Sect. 4). This result is supported by the linear relationship between actual free and total phenytoin concentrations (R2 = 0.95, p < 0.0001; Fig. 1). Figure 2 represents the final structural model incorporating the linear binding model.

Table 4

Population pharmacokinetic parameters of the selected (non-linear) base and final models

Models	Non-linear single-binding site model						Non-linear multiple-binding site model
	Base modela			Final model			Base and final modelb,c
Parameters	Estimated mean value (RSE%)	η (%)	Bootstrap mean (95% CI)	Estimated mean value (RSE%)	η (%)	Bootstrap mean (95% CI)	Estimated mean value (RSE%)	η (%)	Bootstrap mean (95% CI)
Fixed effects
F	0.859 fix			0.859 fix			0.859 fix
ka (h−1)	0.225 fix			0.225 fix			0.225 fix
Vd (L)	1440 (30.7)		1440 (1160–2290)	1440 (16.9)		1440 (1100–2170)	1490 (23.4)		1490 (1240–2230)
Vmax (mg/h)	13.0 (38.0)		13.0 (9.35–15.1)	13.2 (37.0)		13.2 (9.76–15.4)	13.0 (42.3)		13.0 (9.52–15.2)
βweight_Vmax	NA			1.21 (33.8)		1.21 (0.390–5.25)	NA
βbilirubin_Vmax	NA			NA			NA
Km (mg/L)	0.545 (143)		0.545 (0.505–0.597)	0.560 (105)		0.560 (0.506–0.601)	0.603 (113)		0.603 (0.510–0.595)
CF	NA			NA			NA
Kb	NA			NA			NA
βalbumin_Kb	NA			NA			NA
Bmax (mg/L)	127 (9.75)		127 (53.1–161)	94.1 (8.23)		94.1 (54.7–162)	191 (20.5) 61.5 (38.6)		191 (40.3–246) 61.5 (26.6–160)
βalbumin_Bmax	NA			0.823 (20.0)		0.823 (0.416–1.05)	NA
KD (mg/L)	16.0 (11.2)		16.0 (5.61–22.0)	11.1 (10.5)		11.1 (5.84–19.7)	38.8 fix 24.3 fix
Interindividual variability
ω_Vd	0.518 (48.2)	− 5.21	0.518 (0.102–0.693)	0.605 (17.7)	− 10.7	0.605 (0.115–0.687)	0.619 (19.2)	4.14	0.619 (0.0878–0.674)
ω_Vmax	0.447 (105)	10.6	0.447 (0.117–0.956)	0.177 (43.8)	6.63	0.177 (0.0759–0.629)	0.255 (45.9)	1.27	0.255 (0.127–0.842)
ω_Km	1.10 (58.3)	8.65	1.10 (0.425–13.8)	0.752 (46.7)	0.346	0.752 (0.284–5.72)	0.907 (55.8)	10.2	0.907 (0.425–9.43)
ω_CF	NA			NA			NA
ω_Kb	NA			NA			NA
ω_Bmax	0.101 (98.9)	7.66	0.101 (0.0494–0.218)	0.107 (155)	1.62	0.107 (0.0325–0.161)	0.179 (39.3) 0.336 (210)	− 4.38  − 10.1	0.179 (0.0906–0.577) 0.336 (0.0948–0.520)
ω_KD	0.193 (39.0)	− 5.79	0.193 (0.0660–0.245)	0.122 (79.3)	5.57	0.122 (0.0534–0.199)	NA
Residual variability
b1	0.124 (52.7)		0.124 (0.0287–0.142)	0.0551 (174)		0.0551 (0.0257–0.138)	0.0909 (44.6)		0.0909 (0.0394–0.151)
b2	0.0765 (104)		0.0765 (0.0425–0.229)	0.124 (31.4)		0.124 (0.0356–0.190)	0.163 (61.6)		0.163 (0.0456–0.251)
Likelihood estimations
OFV	258.31			232.02			256.21
AIC	282.31			260.02			280.21
BIC	301.64			282.58			299.54

AIC Akaike information criterion, BIC Bayesian information criterion, B maximal binding capacity, b1 proportional error for total phenytoin concentrations, b2 proportional error for free phenytoin concentrations, β covariate parameter estimate, CF albumin coefficient, CI confidence interval, F bioavailability, η η-shrinkage, k absorption rate constant, K binding proportionality constant, K equilibrium dissociation constant, K Michaelis–Menten constant, NA not applicable, OFV objective function value, RSE relative standard error, V volume of distribution, V maximal elimination rate, ω interindividual variability

aModel 15 in Table 2

bModel 16 in Table 2

cFinal model is the same as the base model

Fig. 1

Linear correlation between actual free and total phenytoin concentrations observed in 37 subjects

Fig. 2

Final population pharmacokinetic structural model of phenytoin incorporating the linear binding model. Cfree free phenytoin concentration, Ctotal total phenytoin concentration, F bioavailability, IV intravenous, ka absorption rate constant, Kb binding proportionality constant, Km Michaelis–Menten constant, PO oral, Vd volume of distribution, Vmax maximal elimination rate

Model Evaluation

Diagnostic plots for the final model with linear protein binding are shown in Figs. 3, 4, 5. The plots of observed free plasma concentrations versus population- or individual-predicted free concentrations showed uniformly scattered data along the line of identity, indicating the lack of significant bias and relatively good precision (Fig. 3) [40, 46]. Except for slight positive bias, the data generally centered around the line of zero identity without major systematic bias in the population- and individual-weighted residual plots of free phenytoin (Fig. 4), confirming the lack of obvious defects in the structural (Fig. 4a,c) and error (Fig. 4b,d) models [46]. The pcVPC plot of free phenytoin (Fig. 5) indicated the lack of apparent model misspecifications [40, 46], and all population parameter estimates were within the 95% confidence interval generated by bootstrap analysis (Table 3), verifying the robustness of the final model [40, 46]. The η-shrinkage values were < 15% in the final linear binding model (Table 3), confirming adequate model parameterization [46, 47].

Fig. 3

a Observed plasma concentration of free phenytoin vs. population predicted free concentration; b observed plasma concentration of free phenytoin vs. individual-predicted free concentration in the final population pharmacokinetic structural model of phenytoin incorporating the linear binding model

Fig. 4

a Population-weighted residuals of free phenytoin (PWRES) vs. time; b PWRES vs. predicted plasma concentration of free phenytoin; c Individual-weighted residuals of free phenytoin (IWRES) vs. time; and d IWRES vs. predicted plasma concentration of free phenytoin in the final population pharmacokinetic structural model of phenytoin incorporating the linear binding model

Fig. 5

Prediction-corrected visual predictive check for free phenytoin based on 1000 simulations in the final population pharmacokinetic structural model incorporating the linear binding model. Individual plasma concentrations of free phenytoin (·); 5th, median, and 95th empirical percentiles of plasma concentrations of free phenytoin (─); 5th and 95th percentiles (blue) or median (pink) prediction interval areas

Discussion

Therapeutic drug monitoring of phenytoin is typically conducted by measuring total concentrations, whereas free concentrations are usually estimated with predictive equations [13, 18]. The original Winter–Tozer equation provided reliable estimations of free phenytoin for patients with normal serum albumin [19] but exhibited bias and imprecision in patients with compromised renal function or hypoalbuminemia [48]. Modified Winter–Tozer equations with revised albumin coefficients have been proposed to improve predictability in various patient populations [20-23]. Population pharmacokinetic modeling has also been used to characterize the binding of phenytoin, with apparently improved prediction of free concentrations compared with the original Winter–Tozer equation. Hennig et al. [24] utilized a linear partition coefficient for protein binding and reported enhanced predictive performance in critically ill children. ter Heine et al. [25] used a non-linear single-site-specific binding model and also found increased accuracy and precision in a variety of patient populations. In addition, clinical covariates influencing the free fraction or binding of phenytoin have been reported in a patient cohort with a large age distribution (e.g., albumin, severe renal dysfunction, valproic acid) [25], in adult patients with cancer (e.g., albumin, carbamazepine, valproic acid) [26], in pediatric patients (e.g., albumin, serum urea, valproic acid) [27], and in critically ill children (e.g., albumin) [24] using non-linear mixed-effects modeling. Based on these findings, it remains equivocal which binding model best described the interaction between phenytoin and albumin and whether additional covariates can influence the binding process.

A one-compartment, first-order absorption, and Michaelis–Menten elimination structural model provided the best fit for free phenytoin pharmacokinetic profiles in our dataset. Although the majority of the population pharmacokinetic studies for phenytoin have only focused on estimating Vmax and Km without providing additional structural components (e.g., [44, 49, 50]), our findings are consistent with those that have employed the one-compartment (e.g., [24, 51–53]) and first-order absorption models (e.g., [24, 51, 52]). However, because insufficient data were collected during early post-dose periods in our study, it was not possible to precisely estimate the absorption process in our model; therefore, the bioavailability and absorption rate constant values had to be fixed to literature data (e.g., [24, 45]). This is a general limitation for other published phenytoin models (e.g., [24, 51, 53]). Moreover, phenytoin clearance was described by Michaelis–Menten kinetics in our model, which is consistent with capacity-limiting metabolism at therapeutic concentrations of phenytoin [4, 6] demonstrated in the literature (e.g., [44, 49, 50, 54, 55]). The only exception to our knowledge is Hennig et al. [24] who used a linear elimination model in critically ill children, possibly because of increased Km values in their patients with head trauma [56]. The population pharmacokinetic parameter estimates for Km (0.532 mg/L; Table 3) and Vmax (5.36 mg/day/kg; converted from data presented in Table 3) obtained in the final linear binding model were generally in agreement with reported ranges (i.e., 0.40–0.73 mg/L and 7.9–9.1 mg/day/kg, respectively [26, 44]) for free phenytoin. However, the estimated Vd for free phenytoin in our model (24.1 L/kg; converted from data presented in Table 3) was larger than previously reported in adult patients with cancer (0.15 L/kg) [26] and children (6.39 L/kg) [24]. Given the lipophilic nature of phenytoin [57], this observation may be attributed to decreased fat tissue distribution in these patient groups [58, 59].

Four protein-binding models were directly compared in this study. The CF of the Winter–Tozer model was estimated to be 0.276 (Table 3), which is consistent with the coefficients reported by Anderson et al. [21] (0.25, hypoalbuminemia), Cheng et al. [22] (0.275, hypoalbuminemia), and Kane et al. [23] (0.29, neurosurgical intensive care) but inconsistent with Soriano et al. [20] (0.20, severe renal impairment with a different assay temperature) and the original Winter–Tozer equation (0.20, normal albumin) [19]. This observation might be explained by similarities of our patient population (i.e. with [21-23]), where the majority of subjects (94.6%) had mild hypoalbuminemia (i.e., albumin < 3.5 g/dL), and a substantial portion (43.2%) was admitted to critical care wards (Table 1). The Winter–Tozer model was not the optimal protein-binding model in this study possibly because it was originally developed with the assumptions of a normal phenytoin free fraction of 10% and serum albumin of 4.4 g/dL [19]. However, as evident by Cheng et al. [22], optimizing the albumin coefficient of the Winter–Tozer equation alone did not fully improve the accuracy and precision of the prediction, indicating additional limitations of the model.

The non-linear single-binding-site model also adequately described phenytoin binding, but it generated relatively larger RSEs in fixed-effects estimates, possibly because only the non-saturable portion of the non-linear binding model was used in our dataset. This is evident mathematically, where the estimated KD value (11.1 mg/L) was larger than the free phenytoin concentrations observed in our subjects (1.4 ± 0.7 mg/L). In this scenario, Eq. (3) (i.e., the non-linear binding model) can be simplified as follows:

Equation (7) is essentially a linear binding model (Eq. 2) where Kb is characterized as Bmax/KD, suggesting that the linear binding model provides a more precise and “parsimonious” characterization of phenytoin binding at the typical free concentrations observed clinically. On the other hand, in patients exhibiting relatively high concentrations of free phenytoin, the non-linear single-binding-site model may be more suitable. This observation may also explain the discrepancy in the protein-binding models utilized by Hennig et al. [24] (i.e., linear binding model) and ter Heine et al. [25] (i.e., non-linear binding model), where the documented free phenytoin concentrations were relatively low (0.88 ± 0.79 mg/L) and high (up to ~ 9 mg/L reported in the study), respectively. Moreover, the Bmax value (94.1 mg/L [RSE 8.23%]) in our non-linear single-binding-site model was lower than that reported in ter Heine et al. [25] (117 mg/L [14%]). The discrepancy may be attributed to the mild hypoalbuminemia observed in the majority of our patients. Supporting this observation, covariate analysis identified albumin as a positive modulator of Bmax in our non-linear single-binding-site model, although the interindividual variability of Bmax and its RSE was not improved with the inclusion of albumin ( Table 4). On the other hand, the non-linear multiple-binding-site model was the worst performing of the four binding models, possibly because of over parametrization in the setting of a relatively small number of subjects and that phenytoin may only bind to a single site on albumin (discussed below in Sect. 4).

Overall, our analysis indicated that the linear binding model provided the best description of phenytoin–albumin binding in this patient population (i.e., better diagnostics, improved RSE values for fixed-effects estimates, physiological relevance, and reduced parametrization). In the linear binding model, the relationship between Kb (Cbound/Cfree) and free fraction (fu, Cfree/Ctotal) can be expressed as follows:

The estimated Kb was 7.66 (Table 3) in this study population, with a median fu of 0.118 (Table 1). As phenytoin is extensively bound to albumin (~ 90%) [5, 13], reduced albumin concentrations can decrease the percentage of bound phenytoin, followed by decreased Kb and increased fu. This observation is consistent with that demonstrated by several investigators [24, 26, 27]. In the clinic, reduced albumin concentration leading to increased free fraction may result in reduced total phenytoin concentration, which should be interpreted with caution because the free, pharmacologically active concentration would, in theory, remain the same. In concordance with this, the inclusion of albumin further reduced both the RSE (from 6.81 to 4.82%) and the interindividual variability (from 19.2 to 6.71%) of Kb in the final model compared with the base model (Table 3). In general, linear binding models are typically supported by high binding capacities [60], as evident in our model, with a relatively high Bmax of 373 μM (i.e., based on estimates obtained in the non-linear single-binding-site model; converted to molar units; Table 4). Maximal binding capacity is in turn reflective of the concentration of the binding protein (i.e., albumin) and the number of available binding sites [61]. As the average albumin concentration in our sample population was 397 μM (converted to molar units; Table 1), it may be reasoned that the number of phenytoin binding sites on albumin would be ~ 0.940 (i.e., 373 µM/391 µM). This provides indirect evidence that phenytoin might be primarily bound to only one site on albumin in the clinical setting, which is inconsistent with Chen et al. [30], who suggested the presence of two binding sites using an in vitro model. Further investigations are required to determine the clinical relevance of this information.

In the linear binding model, weight and total bilirubin independently influenced Vmax. Weight was positively associated with Vmax (Table 3, Eq. 5), an effect also observed in studies that estimated Vmax based on the total concentration [52, 53]. Furthermore, increased total bilirubin concentration (representative of compromised hepatic excretory function) was correlated with decreased Vmax, which may be attributed to intrahepatic accumulation of toxic biliary substances, potentially leading to the impairment of hepatic metabolism [62]. On the other hand, clinical factors such as renal function and co-medications were not found to be significant covariates in this dataset. Uremic toxins accumulated in end-stage renal disease may alter protein binding [63, 64] and reduce the intrinsic clearance of phenytoin [65]; however, renal function (measured by SCr) was not identified as a significant covariate in this study, possibly because only a small number of patients had severe renal disease (i.e., GFR < 30 mL/min/1.73 m2). Furthermore, concomitant low-dose aspirin did not impact phenytoin binding in this study, presumably because aspirin’s effect is only evident at relatively higher doses [66]. Carbamazepine, heparin, phenobarbital, sulfonamides, valproic acid, and warfarin also did not affect the pharmacokinetics of phenytoin despite known interacting mechanisms [38, 67–69], which may likely be because of the insufficient study sample size (Table 1).

Despite successfully identifying the optimal phenytoin protein-binding model within a fully developed population pharmacokinetic model, this study has the following limitations. The developed linear binding model can only be applied to patients who have plasma free phenytoin concentrations much less than ~ 11.1 mg/L (i.e., Cfree <  < KD), where the assumption of linearity is still valid. Further validations are also required to ensure the applicability of our model in patient populations with different degrees of hypoalbuminemia, renal dysfunction, or concurrent protein-binding-displacing drugs, as these variables were not characterized systematically in this study. As the study was conducted retrospectively, data on patient ethnicity and genetic polymorphisms (e.g., CYP2C9*3) were not available. Ethnic variations in minor allele frequencies have been reported for CYP2C9 [70], and these polymorphic alleles have been associated with reduced metabolism [71-73] and increased adverse effects of phenytoin [74-76]. Although polymorphic forms of albumin may exhibit altered binding affinity [77], the pharmacogenomics of phenytoin binding was not characterized in this study. Finally, there may be additional binding models that were not investigated in this study.

Conclusion

Systematic comparisons of four common protein-binding models describing the relationships between free and total phenytoin have been conducted using non-linear mixed-effects modeling. We report the novel finding that the linear binding model best characterized phenytoin binding under typical clinical conditions in a process affected by albumin concentration. In patients with higher phenytoin concentrations, the non-linear single-binding-site model might be more suitable. Further validations are required to determine the predictive performances of these models for the estimation of free phenytoin concentrations in the clinical setting.

Below is the link to the electronic supplementary material.

Supplementary file1 (PDF 213 kb)

A systematic comparison of four binding-dissociation models (i.e. Winter-Tozer, linear binding, non-linear single-binding site, and non-linear multiple-binding site) within a fully developed population pharmacokinetic model indicated that phenytoin binding to serum albumin is non-saturable and affected by albumin concentration in adult patients at typical clinical free concentrations.

The final model suggested that phenytoin might be bound to a single site on albumin with a relatively high binding capacity.