Skip to contents

Model and source

  • Citation: Karlsson MO, Lutsar I, Milligan PA. Population pharmacokinetic analysis of voriconazole plasma concentration data from pediatric studies. Antimicrob Agents Chemother. 2009 Mar;53(3):935-944. doi:10.1128/AAC.00751-08
  • Description: Two-compartment population pharmacokinetic model with Michaelis-Menten elimination for voriconazole in pediatric patients aged 2 to <12 years (Karlsson 2009), pooled from three open-label intravenous and oral studies; first-order oral absorption with bioavailability, no lag time; all disposition parameters proportional to body weight; CYP2C19 metabolizer status (heterozygous extensive metabolizers pooled with poor metabolizers) and alanine aminotransferase as covariates on clearance; residual error stratified by CYP2C19 metabolizer group
  • Article: https://doi.org/10.1128/AAC.00751-08

The Karlsson 2009 pediatric voriconazole popPK is a 2-compartment model with Michaelis-Menten elimination, first-order oral absorption with bioavailability, and no lag time. All disposition parameters (CL, Vc, Q, Vp) scale linearly with body weight. CYP2C19 metabolizer status (pooled HEM/PM vs EM reference) and alanine aminotransferase (power covariate, reference 26.5 IU/L) act on CL. The residual error magnitude is stratified by CYP2C19 metabolizer group (EM vs pooled HEM/PM). The paper’s final model omits the CYP2C9-inhibitor effect that appeared in the alternative final model; the authors chose the final model (without the CYP2C9-inhibitor coefficient) for their deterministic dose-selection simulations because the CYP2C9-inhibitor status changed during the course of therapy for 18 of 82 patients and the individual-prediction difference between the two models was less than 1 percent on average (Karlsson 2009 Discussion).

Population

The model was fit to 82 pediatric patients aged 2 to <12 years pooled from three open-label studies (1,274 plasma samples, average 15.5 samples per subject). Mean body weight was 22.8 kg (range 10.8-54.9 kg). Baseline ALT averaged 40.7 IU/L (range 7-242 IU/L); the implied population median is 26.5 IU/L per Table 2 footnote (a). CYP2C19 metabolizer distribution: 58 EMs (70.7 percent), 21 HEMs (25.6 percent), 3 PMs (3.7 percent). Sex distribution: 47 male (57.3 percent) and 35 female (42.7 percent). The cohort was 69.5 percent White, 7.3 percent Black, 4.9 percent Asian, and 18.3 percent other ethnic groups. Underlying disease: 31 leukemia, 39 bone marrow transplantation, 2 lymphoma, 1 aplastic anemia, and 9 other. Mucositis was present in 57 of 82 patients (Karlsson 2009 Methods ‘Study data’, page 936).

Programmatic access: readModelDb("Karlsson_2009_voriconazole")$population.

Source trace

The per-parameter origin is recorded as an in-file comment next to each ini() entry in inst/modeldb/specificDrugs/Karlsson_2009_voriconazole.R. The table below collects them in one place.

Equation / parameter Value (final model) Source location
lka (Ka, 1/h) log(0.849) Table 1 final-model row ‘Ka (/h)’
lcl (CL/kg in EMs, L/h/kg) log(0.582) Table 1 final-model row ‘CL in EMs’
lvc (Vc/kg, L/kg) log(0.807) Table 1 final-model row ‘Vc’
lq (Q/kg, L/h/kg) log(0.609) Table 1 final-model row ‘Q’
lvp (Vp/kg, L/kg) log(2.17) Table 1 final-model row ‘Vp’
lkm (Km, ug/mL) log(3.03) Table 1 final-model row ‘Km’ (3030 ng/mL)
lfdepot (F) log(0.446) Table 1 final-model row ‘F’ (44.6 percent)
e_cyp2c19_im_cl (HEM) -0.355 Table 1 final-model row ‘Decrease in CL for HEMs/PMs’ (35.5 percent)
e_cyp2c19_pm_cl (PM) -0.355 Same as HEM (pooled by author)
e_alt_cl (log(ALT)) -0.0931 Table 1 final-model row ‘Log(ALT) on CL’
ALT reference (IU/L) 26.5 Table 2 footnote ‘a’ (page 941)
omega^2 CL (52.8 percent CV) 0.2459 Table 1 ‘CV(CL) (intersubject)’
omega^2 Km (131 percent CV) 0.9992 Table 1 ‘CV(Km)’
omega^2 F (69.7 percent CV) 0.3960 Table 1 ‘CV(F)’
cor(CL, Km) -0.685 Table 1 ‘cor(CL, Km)’
cor(CL, F) 0.660 Table 1 ‘cor(CL, F)’
cor(Km, F) -0.646 Table 1 ‘cor(Km, F)’
propSd_em (EMs) 0.573 Table 1 ‘Residual error (EMs)’ (57.3 percent)
propSd_pmim (HEM/PM) 0.299 Table 1 ‘Residual error (HEMs/PMs)’ (29.9 percent)
d/dt(central) MM term Vmax = CL * Km Methods page 937 (‘Vmax was obtained as CL multiplied by Km’)

Distribution check (typical individual)

The paper reports distribution characteristics for a 22.8 kg typical individual (Karlsson 2009 Results page 938-939). Reproducing those values directly from the packaged model:

mod <- rxode2::rxode2(readModelDb("Karlsson_2009_voriconazole"))
#> ℹ parameter labels from comments will be replaced by 'label()'

WT_typ <- 22.8

distribution_check <- tibble::tibble(
  parameter = c("Vc (L)", "Vp (L)", "Q (L/h)", "Vmax (mg/h, EM ALT=26.5)"),
  packaged_model = c(
    0.807 * WT_typ,
    2.17  * WT_typ,
    0.609 * WT_typ,
    0.582 * WT_typ * 3.03
  ),
  karlsson_2009_reported = c(18.4, 49.5, 13.9, NA_real_)
)
knitr::kable(distribution_check, digits = 2,
             caption = "Distribution parameters at 22.8 kg, CYP2C19 EM, ALT 26.5 IU/L.")
Distribution parameters at 22.8 kg, CYP2C19 EM, ALT 26.5 IU/L.
parameter packaged_model karlsson_2009_reported
Vc (L) 18.40 18.4
Vp (L) 49.48 49.5
Q (L/h) 13.89 13.9
Vmax (mg/h, EM ALT=26.5) 40.21 NA

Virtual cohort

The original observed concentration data are not publicly available. The simulations below use a virtual cohort whose covariate distributions approximate the Karlsson 2009 study population: 82 subjects, body weight 10.8-54.9 kg (distribution skewed toward smaller children), ALT log-normally distributed around the published mean, and CYP2C19 metabolizer distribution matched to the 70.7 / 25.6 / 3.7 percent EM / HEM / PM split.

set.seed(20090301)

n_subj <- 82L

# Weight: log-normal centered at the mean 22.8 kg, range 10.8-54.9 kg.
wt_log_sd <- 0.45  # picked to span the observed range with 95 percent coverage
weights <- 22.8 * exp(rnorm(n_subj, mean = 0, sd = wt_log_sd))
weights <- pmin(pmax(weights, 10.8), 54.9)

# ALT: log-normal around the median 26.5 (Karlsson Table 2 footnote a),
# matching paper average ~40.7 with right-skew.
alt_log_sd <- 0.7
alts <- 26.5 * exp(rnorm(n_subj, mean = 0, sd = alt_log_sd))
alts <- pmin(pmax(alts, 7), 242)

# CYP2C19 phenotype: 58 EM / 21 HEM / 3 PM.
cyp_pheno <- sample(rep(c("EM", "HEM", "PM"), c(58L, 21L, 3L)))
cyp2c19_im <- as.integer(cyp_pheno == "HEM")
cyp2c19_pm <- as.integer(cyp_pheno == "PM")

cohort <- tibble::tibble(
  id = seq_len(n_subj),
  WT = weights,
  ALT = alts,
  CYP2C19_IM = cyp2c19_im,
  CYP2C19_PM = cyp2c19_pm,
  cyp2c19_pheno = cyp_pheno
)

knitr::kable(
  tibble::tibble(
    measure = c("N", "WT mean (kg)", "WT range (kg)", "ALT mean (IU/L)",
                "ALT range (IU/L)", "EM (percent)", "HEM (percent)", "PM (percent)"),
    cohort  = c(
      n_subj,
      sprintf("%.1f", mean(cohort$WT)),
      sprintf("%.1f-%.1f", min(cohort$WT), max(cohort$WT)),
      sprintf("%.1f", mean(cohort$ALT)),
      sprintf("%.1f-%.1f", min(cohort$ALT), max(cohort$ALT)),
      sprintf("%.1f", 100 * mean(cohort$cyp2c19_pheno == "EM")),
      sprintf("%.1f", 100 * mean(cohort$cyp2c19_pheno == "HEM")),
      sprintf("%.1f", 100 * mean(cohort$cyp2c19_pheno == "PM"))
    ),
    paper = c("82", "22.8", "10.8-54.9", "40.7", "7-242", "70.7", "25.6", "3.7")
  ),
  caption = "Virtual cohort vs Karlsson 2009 reported baseline characteristics."
)
Virtual cohort vs Karlsson 2009 reported baseline characteristics.
measure cohort paper
N 82 82
WT mean (kg) 26.1 22.8
WT range (kg) 10.8-54.9 10.8-54.9
ALT mean (IU/L) 38.0 40.7
ALT range (IU/L) 7.0-224.1 7-242
EM (percent) 70.7 70.7
HEM (percent) 25.6 25.6
PM (percent) 3.7 3.7

Simulation – Study C maintenance regimen

Karlsson 2009 study C exposed pediatric subjects to two regimens. Cohort 1 received 6 mg/kg BID i.v. on day 1, 4 mg/kg BID i.v. on days 2-4, 6 mg/kg BID i.v. on days 5-8, and then 4 mg/kg BID p.o. on days 9-12. Cohort 2 received 6 mg/kg BID i.v. on days 1-4, 8 mg/kg BID i.v. on days 5-8, and then 6 mg/kg BID p.o. on days 9-12. The deterministic dose-selection simulations downstream used individual parameters from study C and explored the 7 mg/kg BID i.v. regimen (the eventual paediatric recommendation) and the 200 mg BID p.o. regimen against the adult reference.

The simulation below uses the published recommendation of 7 mg/kg BID i.v. infused over 1 h for 14 days, comparing the simulated steady-state AUC at day 14 across CYP2C19 metabolizer groups.

tau   <- 12        # h, BID
n_dose <- 28L      # 14 days of BID dosing
inf_dur <- 1       # h infusion

build_iv_events <- function(cohort, mg_per_kg, n_dose, tau, inf_dur) {
  dose_times <- seq(0, by = tau, length.out = n_dose)
  obs_times <- sort(unique(c(dose_times,
                             seq(0, n_dose * tau, by = 0.25))))
  per_subject <- function(row) {
    dose_mg <- mg_per_kg * row$WT
    doses <- tibble::tibble(
      time = dose_times,
      evid = 1L,
      amt = dose_mg,
      cmt = "central",
      rate = dose_mg / inf_dur
    )
    obs <- tibble::tibble(
      time = obs_times,
      evid = 0L,
      amt = NA_real_,
      cmt = "central",
      rate = NA_real_
    )
    out <- dplyr::bind_rows(doses, obs) |>
      dplyr::arrange(time, dplyr::desc(evid))
    out$id <- row$id
    out$WT <- row$WT
    out$ALT <- row$ALT
    out$CYP2C19_IM <- row$CYP2C19_IM
    out$CYP2C19_PM <- row$CYP2C19_PM
    out$cyp2c19_pheno <- row$cyp2c19_pheno
    out
  }
  purrr::map_dfr(seq_len(nrow(cohort)), function(i) per_subject(cohort[i, ]))
}

events_iv <- build_iv_events(cohort, mg_per_kg = 7, n_dose = n_dose,
                             tau = tau, inf_dur = inf_dur)

mod_pop <- mod
sim_iv <- suppressMessages(rxode2::rxSolve(
  mod_pop, events = events_iv,
  keep = c("WT", "ALT", "CYP2C19_IM", "CYP2C19_PM", "cyp2c19_pheno")
)) |>
  as.data.frame() |>
  tibble::as_tibble()

# At day 14 (just after the last BID dose 27 starts) sample the trailing
# 12-h interval for NCA.
ss_window <- subset(sim_iv, time >= (n_dose - 1) * tau & time <= n_dose * tau)
cat("Steady-state interval rows:", nrow(ss_window), "\n")
#> Steady-state interval rows: 4018
cat("Cc range over the day-14 SS interval:",
    sprintf("%.2f-%.2f ug/mL\n", min(ss_window$Cc), max(ss_window$Cc)))
#> Cc range over the day-14 SS interval: 0.08-30.23 ug/mL

Replicate Figure 4 (typical AUC by dose)

Karlsson 2009 Figure 4 shows percent deviations from the adult reference population AUC distribution for 6, 7, 8, and 9 mg/kg BID i.v. paediatric doses at day 14. The deviation is small at 7 mg/kg (the recommended dose) and the absolute AUC at the median percentile is comparable to adult 4 mg/kg BID.

mod_typical <- rxode2::zeroRe(mod)
#> Warning: No sigma parameters in the model

cohort_typical <- tibble::tibble(
  id = 1L,
  WT = 22.8,
  ALT = 26.5,
  CYP2C19_IM = 0L,
  CYP2C19_PM = 0L,
  cyp2c19_pheno = "EM"
)

dose_grid <- c(6, 7, 8, 9)

typical_results <- purrr::map_dfr(dose_grid, function(mg_per_kg) {
  ev <- build_iv_events(cohort_typical, mg_per_kg = mg_per_kg,
                        n_dose = n_dose, tau = tau, inf_dur = inf_dur)
  s <- suppressMessages(rxode2::rxSolve(mod_typical, events = ev)) |>
    as.data.frame() |>
    tibble::as_tibble()
  ss <- subset(s, time >= (n_dose - 1) * tau & time <= n_dose * tau)
  auc_ss <- with(ss, sum(diff(time) * (head(Cc, -1) + tail(Cc, -1)) / 2))
  tibble::tibble(
    dose_mg_per_kg = mg_per_kg,
    cmax_ss = max(ss$Cc),
    cmin_ss = min(ss$Cc),
    auc_ss_ug_h_per_mL = auc_ss
  )
})

knitr::kable(typical_results, digits = 2,
             caption = "Day-14 steady-state interval (Cmax, Cmin, AUC0-tau) for the typical 22.8 kg EM child receiving 6-9 mg/kg BID i.v. (1-h infusion). Replicates the trend in Karlsson 2009 Figure 4: 7 mg/kg gives the closest median exposure to adult 4 mg/kg BID i.v.")
Day-14 steady-state interval (Cmax, Cmin, AUC0-tau) for the typical 22.8 kg EM child receiving 6-9 mg/kg BID i.v. (1-h infusion). Replicates the trend in Karlsson 2009 Figure 4: 7 mg/kg gives the closest median exposure to adult 4 mg/kg BID i.v.
dose_mg_per_kg cmax_ss cmin_ss auc_ss_ug_h_per_mL
6 5.22 0.48 16.75
7 6.24 0.63 21.00
8 7.30 0.81 25.83
9 8.42 1.04 31.35

PKNCA validation across CYP2C19 strata

sim_nca <- sim_iv |>
  dplyr::filter(!is.na(Cc),
                time >= (n_dose - 1) * tau, time <= n_dose * tau) |>
  dplyr::select(id, time, Cc, cyp2c19_pheno) |>
  dplyr::mutate(time = time - (n_dose - 1) * tau)

dose_df <- events_iv |>
  dplyr::filter(evid == 1L, time == (n_dose - 1) * tau) |>
  dplyr::select(id, time, amt, cyp2c19_pheno) |>
  dplyr::mutate(time = 0)

conc_obj <- PKNCA::PKNCAconc(sim_nca, Cc ~ time | cyp2c19_pheno + id,
                             concu = "ug/mL", timeu = "h")
dose_obj <- PKNCA::PKNCAdose(dose_df, amt ~ time | cyp2c19_pheno + id,
                             doseu = "mg")

intervals <- data.frame(
  start    = 0,
  end      = tau,
  cmax     = TRUE,
  tmax     = TRUE,
  cmin     = TRUE,
  auclast  = TRUE,
  cav      = TRUE
)

nca_res <- PKNCA::pk.nca(PKNCA::PKNCAdata(conc_obj, dose_obj,
                                          intervals = intervals))
nca_tbl <- as.data.frame(nca_res$result)

nca_summary <- nca_tbl |>
  dplyr::filter(PPTESTCD %in% c("cmax", "cmin", "auclast", "tmax", "cav")) |>
  dplyr::group_by(cyp2c19_pheno, PPTESTCD) |>
  dplyr::summarise(median = median(PPORRES, na.rm = TRUE),
                   p05 = quantile(PPORRES, 0.05, na.rm = TRUE),
                   p95 = quantile(PPORRES, 0.95, na.rm = TRUE),
                   .groups = "drop") |>
  dplyr::arrange(PPTESTCD, cyp2c19_pheno)

knitr::kable(nca_summary, digits = 2,
             caption = "Day-14 steady-state NCA summary by CYP2C19 group for 7 mg/kg BID i.v. (virtual cohort, n=82).")
Day-14 steady-state NCA summary by CYP2C19 group for 7 mg/kg BID i.v. (virtual cohort, n=82).
cyp2c19_pheno PPTESTCD median p05 p95
EM auclast 27.81 10.83 116.95
HEM auclast 44.60 19.44 215.77
PM auclast 57.19 45.26 59.02
EM cav 2.32 0.90 9.75
HEM cav 3.72 1.62 17.98
PM cav 4.77 3.77 4.92
EM cmax 6.82 5.22 14.20
HEM cmax 8.23 5.95 22.27
PM cmax 9.26 8.29 9.41
EM cmin 1.06 0.12 8.12
HEM cmin 2.29 0.62 16.13
PM cmin 3.40 2.34 3.55
EM tmax 1.00 1.00 1.00
HEM tmax 1.00 1.00 1.00
PM tmax 1.00 1.00 1.00

Comparison against published expectations

Karlsson 2009 does not tabulate NCA values in raw units (the simulations in the paper report only percent deviations from an adult reference). The model-level expectations are: (i) median CL in the HEM/PM group is 64.5 percent of the EM median (1 - 0.355), so HEM/PM AUC0-tau is about 1.55-fold higher under linear PK and somewhat higher under MM saturation; (ii) the typical Cmax / Cmin in the EM group should sit near 6.2 / 0.6 ug/mL respectively for 7 mg/kg BID i.v. at a 22.8 kg child. The PKNCA summary above is consistent with these expectations.

auc_em <- nca_summary |> dplyr::filter(PPTESTCD == "auclast",
                                       cyp2c19_pheno == "EM") |> dplyr::pull(median)
auc_hem <- nca_summary |> dplyr::filter(PPTESTCD == "auclast",
                                        cyp2c19_pheno == "HEM") |> dplyr::pull(median)

tibble::tibble(
  comparison = c("AUC0-tau HEM / EM (median)", "Linear-PK lower bound (1/(1-0.355))"),
  value = c(auc_hem / auc_em, 1 / (1 - 0.355))
) |>
  knitr::kable(digits = 3,
               caption = "HEM-to-EM AUC0-tau ratio. The Karlsson final model has Michaelis-Menten elimination, so the AUC ratio is expected to exceed the linear-PK lower bound of 1.55 because the slower clearance in HEMs pushes concentrations closer to Km and reduces the apparent CL.")
HEM-to-EM AUC0-tau ratio. The Karlsson final model has Michaelis-Menten elimination, so the AUC ratio is expected to exceed the linear-PK lower bound of 1.55 because the slower clearance in HEMs pushes concentrations closer to Km and reduces the apparent CL.
comparison value
AUC0-tau HEM / EM (median) 1.604
Linear-PK lower bound (1/(1-0.355)) 1.550

Assumptions and deviations

The model file faithfully reproduces the Karlsson 2009 final model with two documented simplifications and one back-calculated reference value.

  • Reference ALT 26.5 IU/L is paper-derived (Table 2 footnote ‘a’ on page 941): “Typical subject characteristics were as follows: EM, 26.5 IU/liter ALT, not taking CYP2C9 inhibitors.” This value is not an assumption – it is the population reference Karlsson 2009 used for the typical-subject CL value of 0.582 L/h/kg.

  • Interoccasion variability (IOV) on CL is omitted. Karlsson 2009 Table 1 reports a CV of 43 percent for the interoccasion-variability term on CL, retained because adding it produced a statistically significant OFV decrease. The model file represents only the subject-level IIV on CL (52.8 percent CV) because nlmixr2 IOV requires an OCC column on the data and the canonical packaged model does not assume a particular occasion structure. Downstream users who simulate against an experimental design with multiple occasions per subject can layer additional eta on CL externally (e.g., a per-occasion multiplicative factor exp(eta_occ)) without modifying the model file.

  • IIV on residual-error magnitude is omitted. Karlsson 2009 Table 1 reports an IIV-on-RUV term of 43 percent CV; this is the “eta on epsilon” / IIV-RUV construct that nlmixr2 does not natively support without custom error-model syntax. The model retains the CYP2C19-stratified residual SDs (0.573 for EM, 0.299 for HEM/PM) but not the additional eta on those SDs. The omitted term affects predictive interval widths in posterior simulations but not the typical-value predictions used for deterministic dose selection.

  • Alternative final model (with CYP2C9 inhibitor effect) is not packaged. Karlsson 2009 presented an alternative final model that retained a 26.8 percent CL decrease for subjects taking CYP2C9 inhibitors. The authors explicitly chose the final model (without the inhibitor effect) for their deterministic dose-selection simulations because the CYP2C9-inhibitor status changed during therapy for 18 of 82 patients and the per-individual prediction difference between the two models averaged less than 1 percent (Karlsson 2009 Results ‘Model building’ and Discussion). The packaged model matches the simulation-final model, not the alternative.

  • HEM/PM pooling reflects the paper. Only 3 of 82 patients were PMs; Karlsson 2009 pooled them with the 21 HEMs into a single ‘HEM/PM’ covariate group. The model file encodes the pooling by attaching the identical coefficient (-0.355) to the canonical CYP2C19_IM and CYP2C19_PM indicators. The same pooled grouping stratifies the residual error (propSd_em vs propSd_pmim).