Skip to contents

Fluconazole (Han 2013)

Source: vignettes/articles/Han_2013_fluconazole.Rmd
Han_2013_fluconazole.Rmd

Model and source 

mod_meta <- nlmixr2est::nlmixr(readModelDb("Han_2013_fluconazole"))$meta
#> ℹ parameter labels from comments will be replaced by 'label()'
  • Citation: Han S, Kim J, Yim H, Hur J, Song W, Lee J, Jeon S, Hong T, Woo H, Yim DS. Population pharmacokinetic analysis of fluconazole to predict therapeutic outcome in burn patients with Candida infection. Antimicrob Agents Chemother. 2013;57(2):1006-1011. doi:10.1128/AAC.01372-12
  • Description: One-compartment IV population PK model for fluconazole in adult burn-ICU patients with suspected or confirmed Candida infection, with a piecewise CL covariate model that switches between a fixed CRRT-cohort CL and a Cockcroft-Gault-CrCl / postburn-recency / sepsis-shifted non-CRRT CL plus an additive WT / edema / postburn-recency model on volume (Han 2013)
  • Article (DOI): https://doi.org/10.1128/AAC.01372-12

This vignette validates the packaged Han_2013_fluconazole model – a one-compartment IV population PK model for fluconazole in 60 adult burn-ICU patients with suspected or confirmed Candida infection – against the source publication’s Figure 2 (steady-state VPC by 100, 200, and 400 mg daily regimens) and Figure 3 (probability of fAUC/MIC target attainment across patient subtypes at 200 and 400 mg/day). The packaged model encodes the piecewise CL covariate structure of Han 2013 Table 3, in which a fixed CRRT-cohort CL (theta3 = 1.85 L/h) substitutes for the additive non-CRRT formula (theta1 + (CLCR/120)*theta4 + TBI*theta7 + SEPS*theta9) whenever RRT_CRRT_STATUS = 1, and an additive WT / edema / postburn-recency formula on V (theta5*WT/65 + theta6*EDEM + theta8*TBI) without the not-estimable theta2 intercept.

Population

Han 2013 enrolled 60 burn-ICU patients (50 male / 10 female) at Hangang Sacred Heart Hospital between December 2008 and May 2010, with burns ranging from 11 to 95 percent of total body surface area (mean TBSA 50.3 percent). Patients were intravenously infused with 100-400 mg fluconazole over 10-40 minutes every 24 h for treatment of suspected (47/60) or confirmed (13/60) fungal infection, and PK sampling was performed at steady state after 5-17 days (mean 6.5 days) of fluconazole therapy. The cohort had mean age 50.3 years (range 20-82), mean weight 65.6 kg (range 40-90 kg), mean Cockcroft-Gault creatinine clearance 123.5 mL/min (range 21.6-282.7), and APACHE II score mean 7.6 (range 2-26). Twenty-nine of 60 had clinically diagnosed sepsis, 20/60 had clinical edema (puffy face and pitting peripheral edema), 15/60 were on continuous renal replacement therapy, and mean time from burn injury to fluconazole administration was 23 days (range 7-88) (Han 2013 Table 1).

The same information is available programmatically via the model’s population metadata:

str(mod_meta$population)
#> List of 20
#>  $ species       : chr "human"
#>  $ n_subjects    : int 60
#>  $ n_studies     : int 1
#>  $ age_range     : chr "20-82 years"
#>  $ age_median    : chr "mean 50.3 years (range 20-82)"
#>  $ weight_range  : chr "40.0-90.0 kg"
#>  $ weight_median : chr "mean 65.6 kg (range 40.0-90.0)"
#>  $ sex_female_pct: num 17
#>  $ race_ethnicity: chr "Not reported (single-centre Korean burn-ICU cohort)"
#>  $ disease_state : chr "Adult burn-ICU patients (burns 11-95% of total body surface area; mean TBSA 50.3%) with suspected (47/60) or co"| __truncated__
#>  $ dose_range    : chr "Fluconazole 100-400 mg IV infused over 10-40 minutes every 24 h for treatment of suspected or confirmed fungal "| __truncated__
#>  $ regions       : chr "South Korea (single-centre burn ICU at Hangang Sacred Heart Hospital, Hallym University, Seoul)"
#>  $ tbsa_range    : chr "11-95% (mean 50.3% of total body surface area burned)"
#>  $ apache_ii     : chr "mean 7.6 (range 2-26)"
#>  $ renal_function: chr "Cockcroft-Gault creatinine clearance mean 123.5 mL/min (range 21.6-282.7), raw mL/min, not BSA-normalized"
#>  $ rrt_status    : chr "15/60 (25%) on continuous renal replacement therapy (CRRT)"
#>  $ sepsis        : chr "29/60 (48%) septic at PK sampling"
#>  $ edema         : chr "20/60 (33%) with clinical edema (puffy face and pitting peripheral edema)"
#>  $ postburn_time : chr "Mean 23 days from burn injury to fluconazole administration (range 7-88)"
#>  $ notes         : chr "Baseline demographics per Han 2013 Table 1. 60 burn-ICU patients enrolled between December 2008 and May 2010. P"| __truncated__

Source trace

The per-parameter origin is recorded as an in-file comment next to each ini() entry in inst/modeldb/specificDrugs/Han_2013_fluconazole.R. The table below collects them in one place; all values come from Han 2013 Table 3 (Final parameter estimates and bootstrap results).

Parameter / equation Value Source location
lcl (CL intercept, non-CRRT) log(0.693) Table 3 row theta1 = 0.693 L/h
lcl_rrt (CL on CRRT) log(1.85) Table 3 row theta3 = 1.85 L/h
lvc (V proportionality at WT = 65 kg) log(50.3) Table 3 row theta5 = 50.3 L
e_crcl_cl (CL slope per CLCR/120) 0.557 Table 3 row theta4 = 0.557 L/h
e_dis_burn_recent_cl (CL shift if within 30 postburn days) 0.504 Table 3 row theta7 = 0.504 L/h
e_dis_sepsis_cl (CL shift if sepsis) -0.369 Table 3 row theta9 = -0.369 L/h
e_dis_edema_vc (V shift if edema) 13.6 Table 3 row theta6 = 13.6 L
e_dis_burn_recent_vc (V shift if within 30 postburn days) 9.61 Table 3 row theta8 = 9.61 L
etalcl ~ 0.12000 log(0.357^2 + 1) Table 3 row omega1^2 BSV CL non-CRRT = 35.7 %CV
etalvc ~ 0.02983 log(0.174^2 + 1) Table 3 row omega2^2 BSV V = 17.4 %CV
etalcl_rrt ~ 0.02559 log(0.161^2 + 1) Table 3 row omega3^2 BSV CL CRRT = 16.1 %CV
propSd <- 0.3332 sqrt(0.111) Table 3 row sigma1^2 = 0.111 (proportional variance)
cl_norrt <- (exp(lcl) + e_crcl_cl*(CRCL/120) + e_dis_burn_recent_cl*DIS_BURN_RECENT + e_dis_sepsis_cl*DIS_SEPSIS) * exp(etalcl) n/a Table 3 caption: CL = (1 - CRRT)*(theta1 + CLCR*theta4 + TBI*theta7 + SEPS*theta9) + CRRT*theta3; Table 3 theta4 description “Proportionality constant between CL and CLCR/120” fixes the (CLCR/120) normalisation.
cl_rrt <- exp(lcl_rrt) * exp(etalcl_rrt) n/a Table 3 caption: + CRRT*theta3
cl <- (1 - RRT_CRRT_STATUS)*cl_norrt + RRT_CRRT_STATUS*cl_rrt n/a Table 3 caption
vc <- ((WT/65)*exp(lvc) + e_dis_edema_vc*DIS_EDEMA + e_dis_burn_recent_vc*DIS_BURN_RECENT) * exp(etalvc) n/a Table 3 caption: V = theta2 + (WT/65)*theta5 + EDEM*theta6 + TBI*theta8; Table 3 footnote b: theta2 “not successfully estimated” -> dropped
d/dt(central) <- -kel*central n/a Han 2013 Results, “one-compartment first-order elimination model”
Cc ~ prop(propSd) n/a Han 2013 Results, “proportional residual error was chosen as the base PK model”

The discussion paragraph below Table 3 provides two structural-model sanity values that the encoded equations must reproduce:

  • “When a patient was in hypermetabolic status with a CLCR of 120, the typical value of CL was expected to be 1.75 liters/h (sum of theta1, theta4, and theta7).” At CRRT = 0, CLCR = 120, TBI = 1, SEPS = 0: CL = 0.693 + 0.557*1 + 0.504*1 + 0 = 1.754 L/h. Matches.
  • “For a patient within 30 days after injury, weighing 65 kg, the typical value of V was estimated to be 59.9 liters (sum of theta5 and theta8) when the patient had no edema.” At WT = 65, EDEM = 0, TBI = 1: V = (65/65)*50.3 + 0 + 9.61 = 59.91 L. Matches.

Virtual cohort

Original observed concentrations are not publicly available. The cohort below reproduces the steady-state VPC design described in Han 2013 Methods (1,000 virtual patients with the covariate composition of the original data set; concentrations simulated at the seven sampling time points 3, 5, 9, 24, 27, 48, and 51 h after the start of infusion). The vignette uses n = 200 per dose regimen (reduced from 1,000 to stay inside the 5-minute pkgdown render budget; the median, 5th, and 95th percentiles remain stable). The covariate composition is sampled to match the Han 2013 Table 1 marginal distributions: CRRT 25%, sepsis 48%, edema 33%, recent-postburn (< 30 days) approximately 65% (consistent with mean postburn time 23 days and majority of subjects in the acute hypermetabolic window), weight log-normally distributed around mean 65.6 kg, and CLCR log-normally distributed around mean 123.5 mL/min (truncated at the Table 1 range 21.6-282.7).

set.seed(20260609)

n_per_dose <- 200L
dose_levels <- c(100, 200, 400)    # mg q24h
infusion_min <- 25                 # minutes (Han 2013 Methods: 10-40 min, midpoint 25 min)
infusion_h <- infusion_min / 60
n_doses    <- 7L                   # 7 days q24h to reach steady state
tau_h      <- 24

# Sampling schedule matches Han 2013 Methods: 0, 3, 5, 9, 24, 27, 48, 51 h
# relative to the start of an arbitrary dosing interval at steady state.
# We pick the final 48 h window (interval (n_doses-1) and (n_doses)) so the
# trough at 48 h lands cleanly between two doses for SS interval comparison.
last_dose_time <- (n_doses - 1) * tau_h                 # 144 h
ss_window_start <- last_dose_time - tau_h               # 120 h (start of penultimate interval)
han_sample_offsets <- c(0, 3, 5, 9, 24, 27, 48, 51)     # h within penultimate interval
ss_sample_times <- ss_window_start + han_sample_offsets

# A denser grid for visual VPC overlay (covers two full intervals after
# steady state).
plot_sample_times <- seq(ss_window_start, ss_window_start + 51, by = 0.5)

sample_times <- sort(unique(c(ss_sample_times, plot_sample_times)))

# Helper draws one virtual cohort matched to Han 2013 Table 1 marginals.
draw_covariates <- function(n, id_offset = 0L) {
  # WT: truncated log-normal around the Table 1 mean 65.6 kg
  wt <- pmin(pmax(rlnorm(n, log(65.6) - 0.5 * 0.15^2, 0.15), 40), 90)
  # CLCR: truncated log-normal around the Table 1 mean 123.5 mL/min
  crcl <- pmin(pmax(rlnorm(n, log(123.5) - 0.5 * 0.45^2, 0.45),
                    21.6), 282.7)
  # Binary covariates per Table 1 marginal prevalences
  sepsis <- rbinom(n, 1, 29 / 60)
  edema  <- rbinom(n, 1, 20 / 60)
  crrt   <- rbinom(n, 1, 15 / 60)
  burn_recent <- rbinom(n, 1, 0.65)  # see Assumptions section
  data.frame(
    id              = id_offset + seq_len(n),
    WT              = wt,
    CRCL            = crcl,
    RRT_CRRT_STATUS     = crrt,
    DIS_SEPSIS      = sepsis,
    DIS_EDEMA       = edema,
    DIS_BURN_RECENT = burn_recent
  )
}

make_cohort <- function(dose_mg, id_offset) {
  covs <- draw_covariates(n = n_per_dose, id_offset = id_offset)
  dose_rows <- data.frame(
    id   = rep(covs$id, each = n_doses),
    time = rep(seq(0, by = tau_h, length.out = n_doses), times = n_per_dose),
    evid = 1L,
    amt  = dose_mg,
    rate = dose_mg / infusion_h,
    dv   = NA_real_
  )
  obs_rows <- data.frame(
    id   = rep(covs$id, each = length(sample_times)),
    time = rep(sample_times, times = n_per_dose),
    evid = 0L,
    amt  = NA_real_,
    rate = NA_real_,
    dv   = NA_real_
  )
  out <- rbind(dose_rows, obs_rows)
  out <- merge(out, covs, by = "id", all.x = TRUE)
  out$cohort  <- sprintf("%s mg q24h", dose_mg)
  out$dose_mg <- dose_mg
  out[order(out$id, out$time, -out$evid), ]
}

events_list <- vector("list", length(dose_levels))
for (i in seq_along(dose_levels)) {
  events_list[[i]] <- make_cohort(
    dose_mg   = dose_levels[i],
    id_offset = (i - 1L) * n_per_dose
  )
}
events <- dplyr::bind_rows(events_list)

stopifnot(!anyDuplicated(unique(events[, c("id", "time", "evid")])))

Simulation 

mod <- readModelDb("Han_2013_fluconazole")

sim <- rxode2::rxSolve(
  object = mod, events = events,
  keep   = c("cohort", "dose_mg", "WT", "CRCL", "RRT_CRRT_STATUS",
             "DIS_SEPSIS", "DIS_EDEMA", "DIS_BURN_RECENT")
) |>
  as.data.frame()
#> ℹ parameter labels from comments will be replaced by 'label()'

Replicate published figures

Figure 2 – steady-state VPC by daily dose regimen

Han 2013 Figure 2 plots simulated 90% confidence intervals (5th-95th percentiles) of fluconazole plasma concentrations grouped by the 100, 200, and 400 mg/day dosing regimens, with concentrations on a log scale. The figure shows a single 48-h steady-state window. The vignette reproduces the per-dose VPC ribbons and median curve at the matching sampling time points.

ss_vpc <- sim |>
  dplyr::filter(time %in% plot_sample_times) |>
  dplyr::mutate(time_h_in_ss_window = time - ss_window_start)

ss_summary <- ss_vpc |>
  dplyr::group_by(dose_mg, cohort, time_h_in_ss_window) |>
  dplyr::summarise(
    p05 = quantile(Cc, 0.05, na.rm = TRUE),
    p50 = quantile(Cc, 0.50, na.rm = TRUE),
    p95 = quantile(Cc, 0.95, na.rm = TRUE),
    .groups = "drop"
  ) |>
  dplyr::mutate(
    cohort = factor(cohort,
                    levels = sprintf("%s mg q24h", dose_levels))
  )

ggplot(ss_summary,
       aes(time_h_in_ss_window, p50, group = cohort)) +
  geom_ribbon(aes(ymin = p05, ymax = p95), fill = "steelblue", alpha = 0.25) +
  geom_line(linewidth = 0.7, colour = "steelblue") +
  facet_wrap(~ cohort) +
  scale_y_log10() +
  labs(
    x        = "Time within steady-state window (h)",
    y        = "Fluconazole plasma concentration (ug/mL, log scale)",
    title    = "Figure 2 -- steady-state VPC by daily dose",
    subtitle = paste0("Replicates Han 2013 Figure 2 (100, 200, 400 mg/day); ",
                      "n = ", n_per_dose, " per dose; ",
                      "median (line) and 5th-95th percentile (ribbon).")
  ) +
  theme_minimal()

Figure 3 – probability of fAUC/MIC target attainment

Han 2013 Figure 3 plots the probability of target attainment (PTA) for fAUC/MIC targets of 25 (upper panels) and 50 (lower panels) versus MIC across patient subtypes that affect fluconazole CL (CRRT, sepsis, and burn-recency status). The fAUC is computed as 0.88 * AUC_ss because fluconazole is approximately 12% protein bound (Han 2013 Methods, citing EUCAST 2007 reference 14); AUC_ss = daily dose / CL because fluconazole PK is linear in the dose range studied. The replication below tabulates PTA at MIC = 0.5, 1, and 2 ug/mL for 200 and 400 mg/day at three patient profiles: CRRT, no-CRRT/sepsis/recent burn, and no-CRRT/no-sepsis/resolved (no recent burn). All profiles use the median CLCR (123.5 mL/min) and median WT (65 kg).

# Per Han 2013 simulation: AUC_ss = daily_dose / CL at the typical-value level
# is shifted by IIV via the standard exp(eta) on CL. For PTA simulation,
# match the paper's approach: draw CL from the typical structural model under
# each patient profile and compute AUC = daily_dose / CL_i for each draw.

draw_cl <- function(n, profile, dose_mg) {
  # profile is a list with WT, CRCL, RRT_CRRT_STATUS, DIS_SEPSIS, DIS_EDEMA, DIS_BURN_RECENT
  if (profile$RRT_CRRT_STATUS == 1L) {
    eta <- rnorm(n, 0, sqrt(0.02559))
    cl_typ <- 1.85
  } else {
    eta <- rnorm(n, 0, sqrt(0.12000))
    cl_typ <- 0.693 +
              0.557 * (profile$CRCL / 120) +
              0.504 * profile$DIS_BURN_RECENT +
              (-0.369) * profile$DIS_SEPSIS
  }
  cl_typ * exp(eta)
}

pta_one <- function(dose_mg, profile, target, mics, n_draws = 10000L) {
  cl <- draw_cl(n_draws, profile, dose_mg)
  auc <- dose_mg / cl                          # ug*h/mL = mg*h/L (since 1 mg/L = 1 ug/mL)
  fauc <- 0.88 * auc                           # 12% protein binding -> 88% free
  sapply(mics, function(m) mean(fauc / m >= target))
}

profiles <- list(
  "CRRT (1g equivalent)" = list(
    WT = 65, CRCL = 0, RRT_CRRT_STATUS = 1L,
    DIS_SEPSIS = 0L, DIS_EDEMA = 0L, DIS_BURN_RECENT = 0L
  ),
  "non-CRRT, septic, recent burn" = list(
    WT = 65, CRCL = 123.5, RRT_CRRT_STATUS = 0L,
    DIS_SEPSIS = 1L, DIS_EDEMA = 0L, DIS_BURN_RECENT = 1L
  ),
  "non-CRRT, no sepsis, resolved" = list(
    WT = 65, CRCL = 123.5, RRT_CRRT_STATUS = 0L,
    DIS_SEPSIS = 0L, DIS_EDEMA = 0L, DIS_BURN_RECENT = 0L
  )
)

mic_grid <- c(0.5, 1, 2)
pta_grid <- expand.grid(
  profile = names(profiles),
  dose_mg = c(200, 400),
  target  = c(25, 50),
  stringsAsFactors = FALSE
)

pta_rows <- vector("list", nrow(pta_grid))
for (i in seq_len(nrow(pta_grid))) {
  row <- pta_grid[i, ]
  vals <- pta_one(
    dose_mg = row$dose_mg,
    profile = profiles[[row$profile]],
    target  = row$target,
    mics    = mic_grid
  )
  pta_rows[[i]] <- data.frame(
    profile = row$profile,
    dose_mg = row$dose_mg,
    target  = row$target,
    mic     = mic_grid,
    pta     = round(vals, 3)
  )
}
pta_long <- dplyr::bind_rows(pta_rows)

pta_wide <- pta_long |>
  tidyr::pivot_wider(
    names_from  = mic,
    values_from = pta,
    names_glue  = "PTA at MIC = {mic} ug/mL"
  ) |>
  dplyr::arrange(profile, target, dose_mg)

knitr::kable(
  pta_wide,
  caption = paste0("Probability of fAUC/MIC target attainment (PTA) ",
                   "for fluconazole 200 and 400 mg/day across three ",
                   "patient profiles, computed from 10,000 Monte-Carlo ",
                   "draws of CL under the Han 2013 covariate model. ",
                   "Compare against Han 2013 Figure 3.")
)
Probability of fAUC/MIC target attainment (PTA) for fluconazole 200 and 400 mg/day across three patient profiles, computed from 10,000 Monte-Carlo draws of CL under the Han 2013 covariate model. Compare against Han 2013 Figure 3.
profile dose_mg target PTA at MIC = 0.5 ug/mL PTA at MIC = 1 ug/mL PTA at MIC = 2 ug/mL
CRRT (1g equivalent) 200 25 1 1.000 1.000
CRRT (1g equivalent) 400 25 1 1.000 1.000
CRRT (1g equivalent) 200 50 1 1.000 0.370
CRRT (1g equivalent) 400 50 1 1.000 1.000
non-CRRT, no sepsis, resolved 200 25 1 1.000 0.999
non-CRRT, no sepsis, resolved 400 25 1 1.000 1.000
non-CRRT, no sepsis, resolved 200 50 1 0.999 0.839
non-CRRT, no sepsis, resolved 400 50 1 1.000 0.999
non-CRRT, septic, recent burn 200 25 1 1.000 0.996
non-CRRT, septic, recent burn 400 25 1 1.000 1.000
non-CRRT, septic, recent burn 200 50 1 0.997 0.747
non-CRRT, septic, recent burn 400 50 1 1.000 0.995

Han 2013 Figure 3 reports near-100 % PTA for 400 mg/day at a target of 25 across all profiles when MIC <= 2 ug/mL, and shows the CRRT / recent- burn profiles producing markedly lower PTA at the higher target (50) because their CL is higher. The simulated table reproduces this qualitative pattern: the recent-burn / non-CRRT / septic profile and the CRRT profile have lower PTA than the resolved-phase reference at the stricter target of 50.

PKNCA validation

PKNCA computes the steady-state AUC0-tau on the penultimate dosing interval (ss_window_start to ss_window_start + tau_h) across the three dose regimens. The treatment grouping in the formula is the cohort label so per-dose Cmax, Cmin, Cavg, and AUC0-tau are summarised together.

pk_conc <- sim |>
  dplyr::filter(time >= ss_window_start, time <= ss_window_start + tau_h,
                !is.na(Cc), Cc > 0) |>
  dplyr::select(id, time, Cc, cohort)

pk_dose <- events |>
  dplyr::filter(evid == 1L,
                time >= ss_window_start, time < ss_window_start + tau_h) |>
  dplyr::select(id, time, amt, cohort)

conc_obj <- PKNCA::PKNCAconc(
  data    = pk_conc,
  formula = Cc ~ time | cohort + id,
  concu   = "ug/mL",
  timeu   = "hr"
)
dose_obj <- PKNCA::PKNCAdose(
  data    = pk_dose,
  formula = amt ~ time | cohort + id,
  doseu   = "mg"
)

intervals <- data.frame(
  start   = ss_window_start,
  end     = ss_window_start + tau_h,
  cmax    = TRUE,
  tmax    = TRUE,
  cmin    = TRUE,
  cav     = TRUE,
  auclast = TRUE
)

nca_data <- PKNCA::PKNCAdata(conc_obj, dose_obj, intervals = intervals)
nca_res  <- suppressWarnings(PKNCA::pk.nca(nca_data))

knitr::kable(
  summary(nca_res),
  caption = paste0("Steady-state NCA parameters (Cmax,ss, Cmin,ss, ",
                   "Cavg,ss, AUC0-tau,ss) by dose regimen at ",
                   "the penultimate dosing interval; n = ",
                   n_per_dose, " per regimen.")
)
Steady-state NCA parameters (Cmax,ss, Cmin,ss, Cavg,ss, AUC0-tau,ss) by dose regimen at the penultimate dosing interval; n = 200 per regimen.
Interval Start Interval End cohort N AUClast (hr*ug/mL) Cmax (ug/mL) Cmin (ug/mL) Tmax (hr) Cav (ug/mL)
120 144 100 mg q24h 200 63.7 [41.3] 3.59 [31.2] 1.82 [56.5] 0.500 [0.500, 0.500] 2.65 [41.3]
120 144 200 mg q24h 200 128 [36.7] 7.23 [27.6] 3.68 [51.3] 0.500 [0.500, 0.500] 5.35 [36.7]
120 144 400 mg q24h 200 251 [39.4] 14.1 [30.4] 7.26 [51.3] 0.500 [0.500, 0.500] 10.5 [39.4]

Comparison against derived AUC_ss

Han 2013 did not publish an NCA table directly. The Methods section, however, states that the simulated 24-h steady-state AUC was computed as daily_dose / CL because fluconazole shows linear PK in the dose range studied. For the resolved-phase typical patient (CRRT = 0, CLCR = 120, TBI = 0, SEPS = 0), CL = 0.693 + 0.557 = 1.25 L/h (Han 2013 Discussion explicit value), so:

Daily dose Expected AUC_ss (mgh/L = ugh/mL)
100 mg 100 / 1.25 = 80
200 mg 200 / 1.25 = 160
400 mg 400 / 1.25 = 320

These expected AUC_ss values are the typical-value reference; the simulated PKNCA AUC0-tau median should be close to the expected value at the resolved-phase typical patient. The simulated medians above include the full cohort with CRRT (CL = 1.85, AUC reduced by 1.25/1.85), recent burn (CL increased), and sepsis (CL decreased) variation, so the cohort median AUC0-tau will be slightly below the resolved-phase typical-patient value of 80 / 160 / 320 ug*h/mL.

Assumptions and deviations

  • theta2 (V intercept) dropped per Han 2013 Table 3 footnote b. The Table 3 row for theta2 records “NE” (not estimated) with the footnote “In the final model, when the covariates were fully included, the slope of the proportional relationship (theta2) was not successfully estimated.” The packaged model encodes the V equation as V = (WT/65)*exp(lvc) + e_dis_edema_vc*EDEM + e_dis_burn_recent_vc*TBI with theta2 = 0, consistent with the Han 2013 Discussion line “For a patient within 30 days after injury, weighing 65 kg, the typical value of V was estimated to be 59.9 liters (sum of theta5 and theta8)” – no theta2 contribution is included in the published reference value.

  • CLCR normalisation factor 120 mL/min from Han 2013 Table 3 parameter description. Han 2013 Table 2 writes the model-development equations as CL_RENAL = (theta1 + CLCR * theta4) * exp(eta1) (CLCR as bare mL/min), but Han 2013 Table 3 describes theta4 as the “Proportionality constant between CL and CLCR/120” and Han 2013 Discussion reports that at CLCR = 120 the CL contribution from this term equals theta4 = 0.557 L/h (which is consistent only if (CLCR/120) * theta4 is used, not CLCR * theta4). The packaged model uses the Table 3 / Discussion form (CLCR / 120) * theta4; the Table 2 form appears to be a typographical contraction of the normalised version.

  • CRCL stored under the canonical CRCL column despite NOT being BSA-normalized. Han 2013 uses raw Cockcroft-Gault CrCL (mL/min, not BSA-normalized). Following the precedent of Delattre_2010_amikacin.R and Shekar_2014_meropenem.R, the model stores the source CLCR column under CRCL, with the raw / non-BSA-normalized status documented in the per-model covariateData[[CRCL]]$units and notes fields. Do not compare the magnitude of e_crcl_cl = 0.557 against BSA-normalized reference values listed in the canonical entry – the units differ.

  • Three new canonical covariates ratified alongside this extraction. DIS_SEPSIS, DIS_EDEMA, and DIS_BURN_RECENT are added to inst/references/covariate-columns.md. The Han 2013 source column names SEPS, EDEM, and TBI are recorded as aliases in each respective entry. The TBI alias carries an explicit note that this Han 2013 column does NOT denote traumatic brain injury – it is a recoded binary indicator at the 30-day postburn threshold.

  • Conditional IIV split by CRRT status; etalcl_rrt declared as paper-specific. Han 2013 reports two distinct CL BSV estimates (omega1^2 = 35.7 %CV for the non-RRT arm, omega3^2 = 16.1 %CV for the RRT arm). The packaged model encodes this by pairing each arm with its own multiplicative eta (etalcl for non-RRT, etalcl_rrt for RRT). etalcl_rrt is declared in paper_specific_etas at the top of the model function because the conditional-IIV split-by-CRRT structure does not match the 1-to-1 eta-to-lparam pairing rule enforced by checkModelConventions(). The non-applicable arm’s eta evaluates with the indicator-zeroed multiplier in model(), so each subject’s per-time-point CL only carries the IIV from the arm matching their CRRT status.

  • Additive covariate model implemented in linear space. Han 2013 reports CL = theta1 + ... and V = theta5*WT/65 + ... as additive formulas, not multiplicative. The packaged model preserves the additive structure by computing the linear sum first and then multiplying by exp(eta); the structural log-transformed lcl, lcl_rrt, and lvc are back-transformed via exp(...) inside model() so only positive intercepts enter the additive sum. The additive covariate coefficients (e_crcl_cl, e_dis_burn_recent_cl, e_dis_sepsis_cl, e_dis_edema_vc, e_dis_burn_recent_vc) are kept in linear space because they can take negative values (the sepsis coefficient is -0.369 L/h).

  • propSd <- sqrt(0.111) = 0.3332 not 0.111. Han 2013 Table 3 reports the proportional residual error as sigma1^2 = 0.111, the NONMEM $SIGMA value (variance on the proportional residual scale, not SD). For nlmixr2’s Cc ~ prop(propSd) the propSd argument is the SD of the proportional residual, so propSd = sqrt(0.111) = 0.3332.

  • Recent-postburn (DIS_BURN_RECENT = 1) prevalence in the virtual cohort set to ~65 %. Han 2013 Table 1 reports the mean time from burn injury was 23 days (range 7-88), which places the majority of the cohort in the acute hypermetabolic window (within 30 days), but the paper does not tabulate the exact TBI prevalence. The vignette uses 65 % as a reasonable approximation (the mean of 23 days suggests well over half the cohort had time from burn injury < 30 days); changing this rate affects only the marginal cohort statistics in Figure 2, not the per-subject typical-value predictions, and not the PTA table (which fixes the patient profile rather than sampling it).

  • Infusion duration set to 25 min (midpoint of Han 2013 Methods). Han 2013 Methods reports fluconazole was “intravenously infused with 100 to 400 mg of fluconazole for 10 to 40 min every 24 h”; the vignette uses 25 min as the midpoint. The exact infusion duration affects the very early time-course (Cmax shape during infusion) but has negligible effect on AUC and on Cmin / Cavg / Ctau at steady state, which are the quantities the validation focuses on.

  • Steady state defined at 7 days of q24h dosing. Han 2013 sampled patients who had received 5-17 days (mean 6.5 days) of fluconazole therapy at the start of PK sampling. The vignette uses 7 doses (7 days) before sampling the penultimate interval at days 5-7; with the typical-patient terminal half-life of ln(2) * 50.3 / 1.25 = 27.9 h, steady state is reached within 5-6 days, so 7 days is sufficient.

  • Race / ethnicity not modeled. Han 2013 enrolled a single-centre Korean cohort. The paper does not test race as a covariate; the model does not include a race effect.

  • No additive residual error in the final model. Han 2013 Methods describes the candidate residual-error structure as a combined additive plus proportional model. Han 2013 Results section (“PK model”) states: “A one-compartment first-order elimination model with proportional residual error was chosen as the base PK model.” Only the proportional component is reported in Table 3 (sigma1^2 = 0.111); no additive sigma is reported. The packaged model uses Cc ~ prop(propSd) (no additive component).