Skip to contents

Spectinamide 1810 in mouse TB (Wagh 2021)

Source: vignettes/articles/Wagh_2021_spectinamide_1810_mouse.Rmd
Wagh_2021_spectinamide_1810_mouse.Rmd

Model and source

  • Citation: Wagh S, Rathi C, Lukka PB, Parmar K, Temrikar Z, Liu J, Scherman MS, Lee RE, Robertson GT, Lenaerts AJ, Meibohm B. (2021). Model-Based Exposure-Response Assessment for Spectinamide 1810 in a Mouse Model of Tuberculosis. Antimicrob Agents Chemother 65(11):e01744-20. doi:10.1128/AAC.01744-20. PMID 34424046.
  • Description: Population PK + PK/PD model with a novel post-antibiotic-effect (PAE) compartment for subcutaneous spectinamide 1810 in BALB/c mice chronically infected with Mycobacterium tuberculosis (Erdman strain).
  • Article: https://doi.org/10.1128/AAC.01744-20
  • Supplemental material (PRIOR NWPRI methodology, GoF + VPC figures): linked from the article landing page.

Population

A pooled dataset of 196 Mtb-infected female BALB/c mice (147 from study 1 + 49 from study 2; 8 weeks old, 18-22 g body weight) plus 84 healthy female BALB/c mice that provided dense-sampling PK to seed the Bayesian priors in the infected-animal sparse-sampling fit. Mice were infected via low-dose aerosol exposure (Mtb ATCC 35801 / TMC 107 / Erdman; MIC 1.6 mg/L) at day 0, with treatment initiated at day 34 post-infection and continued for 4 weeks with weekend drug holidays (5 dosing days per week). 29 dose-group combinations spanned 20-800 mg/kg per dose, weekly totals 20-4000 mg/kg, at frequencies QW / BIW / TIW / QD / BID (Wagh 2021 Tables 4 and 5). Study 1 and study 2 differed materially in baseline bacterial load at the start of treatment (7.08 vs 5.37 log CFU), which the integrated PK/PD model captures via a multiplicative 1.15 study-2 effect on the maximum kill rate.

The model’s population metadata also exposes the same information programmatically:

pop <- nlmixr2lib::readModelDb("Wagh_2021_spectinamide_1810_mouse")()$meta$population
str(pop)
#> List of 11
#>  $ species       : chr "mouse (BALB/c, female)"
#>  $ n_subjects    : int 376
#>  $ n_studies     : int 2
#>  $ age_range     : chr "8 weeks at study start (acclimatized 72 h prior to dosing)"
#>  $ weight_range  : chr "18-22 g body weight (Wagh 2021 Methods)"
#>  $ sex_female_pct: num 100
#>  $ race_ethnicity: chr NA
#>  $ disease_state : chr "Females BALB/c mice chronically infected with Mycobacterium tuberculosis (ATCC 35801 / TMC 107 / Erdman; MIC of"| __truncated__
#>  $ dose_range    : chr "Spectinamide 1810 subcutaneous, 20-800 mg/kg per dose, weekly totals 20-4000 mg/kg, frequencies QW / BIW / TIW "| __truncated__
#>  $ regions       : chr "USA (Univ. of Tennessee Health Science Center for healthy-animal PK; Colorado State University BSL-3 for infect"| __truncated__
#>  $ notes         : chr "Total infected mice n = 196 (147 with evaluable PK + log CFU in study 1 + 49 from study 2). 84 healthy mice pro"| __truncated__

Source trace

Every ini() parameter in the model file carries an in-file comment pointing to its source table. The table below collects the structural references in one place for review.

Quantity Value Source location
Plasma PK: 2-cmt, 1st-order SC structure Wagh 2021 Results ‘Pharmacokinetic analysis’ + Figure 5
K_a 15.4 1/h Table 1, Infected (RSE 8%)
V_c/F 0.217 L/kg Table 1, Infected (RSE 3%)
V_p/F 0.160 L/kg Table 1, Infected (RSE 12%)
CL/F 0.697 L/h/kg Table 1, Infected (RSE 2%)
Q/F 0.0089 L/h/kg Table 1, Infected (RSE 6%)
IIV on CL/F 13.2% CV Table 1, Infected (RSE 24%)
Plasma RUV (LTBS proportional) 60.0% CV Table 1, Infected (RSE 8%)
Bacterial natural growth: logistic structure Methods ‘Bacterial growth kinetic model’, Eq. 2
K_gs 0.0357 1/h (FIXED) Table 2 (RSE 4%)
log_10 N_max 6.11 (FIXED) Table 2 (RSE 1%)
log_10 initial CFU 1.77 (FIXED) Table 2 (RSE 3%)
PAE + sigmoidal-Emax PD: structure high-water-mark + 1st-order decay Methods ‘PK/PD modeling with PAE estimation’, Eq. 4 + Figure 5
K_kill_max (study 1) 0.0374 1/h Table 3 (RSE 5.3%)
EC_50 79.6 mg/L Table 3 (RSE 60.2%)
Hill coefficient g 1.58 Table 3 (RSE 23.1%)
K_PAE 0.0142 1/h Table 3 (RSE 74.6%; PAE half-life ~ 48.8 h)
Study 2 K_kill_max multiplier 1.15 (FIXED) Table 3 (RSE 4.0%)
PD RUV1 (study 1, log10 CFU) 0.284 Table 3 (RSE 12.1%)
PD RUV2 (study 2, log10 CFU) 0.173 Table 3 (RSE 27.2%); see Errata below
MIC (spectinamide 1810 vs Mtb Erdman) 1.6 mg/L Methods ‘Animal experimentation’; also Wagh 2021 ref [6]
Mouse plasma free fraction (1810) 0.31 Methods ‘Identification of exposure drivers for efficacy’; Wagh 2021 ref [6]

Virtual cohort

The simulation reproduces a small subset of the Wagh 2021 dose-fractionation panel that compares equal weekly doses delivered at different frequencies. Per the skill’s cap, no arm exceeds 200 mice; here each arm is 25 mice for the VPC, well under the cap.

# Helper: build one cohort's event table. id_offset shifts subject IDs so
# multiple cohorts bind_rows cleanly without rxSolve silently merging
# duplicate IDs into a single subject.
build_arm <- function(per_dose_mg_per_kg,
                      dose_times_h,
                      pd_time_h,
                      arm_label,
                      study_wagh_2,
                      n_mice,
                      id_offset = 0L) {
  ids <- id_offset + seq_len(n_mice)

  # Doses: every (id, time) gets a dosing row pointing to depot.
  doses <- tidyr::expand_grid(id = ids, time = dose_times_h) |>
    dplyr::mutate(
      evid  = 1L,
      amt   = per_dose_mg_per_kg,
      cmt   = "depot",
      dvid  = NA_integer_
    )

  # PK observation grid: 10 points across the last dosing interval
  # (steady state) to capture Cmax / decline shape. Anchored on the next-
  # to-last dose. Use cmt = "central" (the ODE state); rxode2 returns the
  # algebraic Cc as a column in the output frame.
  last_dose <- max(dose_times_h)
  prev_dose <- sort(dose_times_h, decreasing = TRUE)[2]
  if (is.na(prev_dose)) prev_dose <- last_dose - 24
  pk_times_steady <- seq(prev_dose, last_dose, length.out = 11)
  pk_times_after  <- last_dose + c(0.1, 0.25, 0.5, 1, 2, 4, 8, 16, 24, 48)
  pk_obs <- tidyr::expand_grid(
    id   = ids,
    time = c(pk_times_steady, pk_times_after)
  ) |>
    dplyr::mutate(
      evid  = 0L,
      amt   = NA_real_,
      cmt   = "central",
      dvid  = 1L
    )

  # PD observation: one log_cfu reading per mouse at end of treatment + washout.
  pd_obs <- tibble::tibble(
    id   = ids,
    time = pd_time_h,
    evid = 0L,
    amt  = NA_real_,
    cmt  = "bacteria",
    dvid = 2L
  )

  out <- dplyr::bind_rows(doses, pk_obs, pd_obs) |>
    dplyr::arrange(id, time, dplyr::desc(evid))
  out$STUDY_WAGH_2 <- study_wagh_2
  out$arm          <- arm_label
  out
}

# Treatment week schedule: M-F dosing weeks (drug holidays Sat / Sun).
# Each treatment week occupies 168 h; doses occur on relative days 0,1,2,3,4
# (M-F) at the per-arm intra-day pattern. 4 weeks total (Wagh 2021 design).
# Times are RELATIVE TO START OF TREATMENT (we override bacteria(0) per arm
# to the observed initial-treatment value so we don't have to simulate the
# 34-day natural-growth lead-in).
weeks <- 4L
mf_days <- c(0, 1, 2, 3, 4)
# QD M-F (one dose per dosing day at 0h-of-day)
times_QD  <- as.vector(outer(mf_days * 24, 168 * (0:(weeks - 1)), "+"))
# TIW (M, W, F) at 0h-of-day; days 0, 2, 4 each week
times_TIW <- as.vector(outer(c(0, 2, 4) * 24, 168 * (0:(weeks - 1)), "+"))
# BID M-F at 0h-of-day and 12h-of-day
times_BID <- as.vector(outer(c(mf_days * 24, mf_days * 24 + 12),
                             168 * (0:(weeks - 1)), "+"))

# Two-day washout after the last treatment day (Friday of week 4); PD reading
# 48 h after the last dose, mirroring the paper's design.
pd_time_h <- 168 * (weeks - 1) + 4 * 24 + 48

n_per <- 25L

events <- dplyr::bind_rows(
  build_arm(per_dose_mg_per_kg = 100, dose_times_h = times_QD,
            pd_time_h = pd_time_h, arm_label = "100 mg/kg QD (S1)",
            study_wagh_2 = 0, n_mice = n_per, id_offset =   0L),
  build_arm(per_dose_mg_per_kg = 200, dose_times_h = times_TIW,
            pd_time_h = pd_time_h, arm_label = "200 mg/kg TIW (S1)",
            study_wagh_2 = 0, n_mice = n_per, id_offset =  25L),
  build_arm(per_dose_mg_per_kg = 50,  dose_times_h = times_BID,
            pd_time_h = pd_time_h, arm_label = "50 mg/kg BID (S1)",
            study_wagh_2 = 0, n_mice = n_per, id_offset =  50L),
  build_arm(per_dose_mg_per_kg = 100, dose_times_h = times_QD,
            pd_time_h = pd_time_h, arm_label = "100 mg/kg QD (S2)",
            study_wagh_2 = 1, n_mice = n_per, id_offset =  75L)
)

# Disjoint-ID guard (mandatory for multi-cohort sims, per skill template):
stopifnot(!anyDuplicated(unique(events[, c("id", "time", "evid")])))

Simulation

The model’s bacterial-state initial condition (bacteria(0) <- rbase = 10^1.77 = 58.9 CFU) corresponds to the aerosol-infection time (day 1 post-infection). To avoid simulating the 34-day natural-growth lead-in, we override the bacterial initial condition per arm to the observed bacterial load at the start of treatment (7.08 log10 CFU in study 1; 5.37 log10 CFU in study 2; Wagh 2021 Discussion paragraph 4).

mod <- nlmixr2lib::readModelDb("Wagh_2021_spectinamide_1810_mouse")()

# Override bacteria(0) per arm; rxSolve does not directly accept per-id
# init overrides via `inits`, so we split-simulate by arm and then combine.
sim_one <- function(arm_label, init_bacteria, ev) {
  ev_arm <- dplyr::filter(ev, arm == arm_label)
  rxode2::rxSolve(
    mod,
    events    = ev_arm,
    inits     = c(bacteria = init_bacteria),
    keep      = c("arm", "STUDY_WAGH_2"),
    returnType = "data.frame"
  ) |>
    dplyr::mutate(arm_factor = factor(arm_label,
                                      levels = c("100 mg/kg QD (S1)",
                                                 "200 mg/kg TIW (S1)",
                                                 "50 mg/kg BID (S1)",
                                                 "100 mg/kg QD (S2)")))
}

# Study 1 arms: initial bacterial load 10^7.08 = 1.20e7 CFU/lung
# Study 2 arm:  initial bacterial load 10^5.37 = 2.34e5 CFU/lung
sim <- dplyr::bind_rows(
  sim_one("100 mg/kg QD (S1)",  init_bacteria = 10^7.08, ev = events),
  sim_one("200 mg/kg TIW (S1)", init_bacteria = 10^7.08, ev = events),
  sim_one("50 mg/kg BID (S1)",  init_bacteria = 10^7.08, ev = events),
  sim_one("100 mg/kg QD (S2)",  init_bacteria = 10^5.37, ev = events)
)

For deterministic typical-value plots (no between-subject IIV), set zeroRe() first; we keep IIV here for the VPC-style 5-95% ribbons.

Replicate published figures

Figure 1 – peak plasma concentration at steady state by dosing regimen

The paper’s Figure 1 shows peak (0.25 h post-dose) plasma concentrations by dose / regimen with the 95% prediction interval. Below we render the VPC ribbon and the typical-value line over the final steady-state dosing interval (week 4) of each arm.

last_dose <- max(c(times_QD, times_TIW, times_BID))
window_t  <- last_dose + seq(0, 12, by = 0.25)

# Pull simulated Cc over the window from the stochastic VPC simulation
pk_window <- sim |>
  dplyr::filter(time >= last_dose - 0.001, time <= last_dose + 12.001) |>
  dplyr::mutate(t_after_last_dose = time - last_dose)

pk_summary <- pk_window |>
  dplyr::group_by(arm_factor, t_after_last_dose) |>
  dplyr::summarise(
    Q05 = quantile(Cc, 0.05, na.rm = TRUE),
    Q50 = quantile(Cc, 0.50, na.rm = TRUE),
    Q95 = quantile(Cc, 0.95, na.rm = TRUE),
    .groups = "drop"
  )

ggplot(pk_summary, aes(t_after_last_dose, Q50)) +
  geom_ribbon(aes(ymin = Q05, ymax = Q95), alpha = 0.25, fill = "steelblue") +
  geom_line(colour = "steelblue") +
  facet_wrap(~arm_factor) +
  scale_y_log10() +
  labs(x = "Time after last dose (h)",
       y = "Plasma spectinamide 1810 (mg/L)",
       title = "Simulated plasma PK over the final steady-state dosing interval",
       caption = "VPC ribbon = simulated 5%-95% range; line = simulated median. Replicates the design of Figure 1 of Wagh 2021.")

Figure 2 – lung bacterial burden after 4 weeks of therapy

The paper’s Figure 2 plots lung log CFU after 4 weeks of treatment by dose group, showing the dose-fractionation effect (less frequent / larger doses beat more frequent / smaller doses at equal weekly totals).

pd <- sim |>
  dplyr::filter(time == pd_time_h) |>
  dplyr::distinct(arm_factor, STUDY_WAGH_2, id, log_cfu)

# Weekly dose per arm (mg/kg/week)
weekly_dose <- c(
  "100 mg/kg QD (S1)"  = 100 * 5,
  "200 mg/kg TIW (S1)" = 200 * 3,
  "50 mg/kg BID (S1)"  = 50 * 2 * 5,
  "100 mg/kg QD (S2)"  = 100 * 5
)
pd$weekly_dose_mg_kg <- weekly_dose[as.character(pd$arm_factor)]

ggplot(pd, aes(arm_factor, log_cfu)) +
  geom_boxplot(width = 0.45, fill = "grey92") +
  geom_jitter(width = 0.06, alpha = 0.45, size = 1.2) +
  labs(x = NULL,
       y = "Lung log10 CFU at end of treatment + 2-day washout",
       title = "Bacterial burden by dosing regimen after 4 weeks of therapy",
       caption = paste0(
         "Initial bacterial load before treatment: study 1 = 7.08 log CFU, ",
         "study 2 = 5.37 log CFU. Replicates the design of Figure 2 of Wagh 2021."
       )) +
  theme(axis.text.x = element_text(angle = 20, hjust = 1))

PKNCA validation

Compute Cmax, Tmax, AUC over the last dosing interval (steady-state approximation) for each arm using PKNCA. Per the skill template, the filter is !is.na(Cc) only – adding time > 0 or Cc > 0 would drop the pre-dose anchor row that PKNCA needs.

last_dose <- max(c(times_QD, times_TIW, times_BID))
# Use the pre-last-dose anchor through 48 h after the last dose so PKNCA
# can compute Cmax / Tmax / AUC for the final interval.
sim_for_nca <- sim |>
  dplyr::filter(time >= (last_dose - 0.001), time <= last_dose + 48.001) |>
  dplyr::mutate(t_rel = time - last_dose) |>
  dplyr::filter(!is.na(Cc)) |>
  dplyr::select(id, arm = arm_factor, t_rel, Cc)

# Time-zero anchor per (id, arm) -- pre-dose concentration is the
# steady-state trough; here we use the actual simulated value at t_rel = 0.
# (For the multi-cohort case, time = last_dose is in the simulation grid by
# construction.)
sim_for_nca <- dplyr::bind_rows(
  sim_for_nca,
  sim_for_nca |>
    dplyr::group_by(id, arm) |>
    dplyr::slice(1) |>
    dplyr::ungroup() |>
    dplyr::mutate(t_rel = 0)
) |>
  dplyr::distinct(id, arm, t_rel, .keep_all = TRUE) |>
  dplyr::arrange(arm, id, t_rel)

conc_obj <- PKNCA::PKNCAconc(sim_for_nca, Cc ~ t_rel | arm + id)

# Dose record: last per-mouse dose at t_rel = 0
dose_df <- events |>
  dplyr::filter(evid == 1, time == last_dose) |>
  dplyr::transmute(id, t_rel = 0, amt, arm = arm)
dose_obj <- PKNCA::PKNCAdose(dose_df, amt ~ t_rel | arm + id)

intervals <- data.frame(
  start      = 0,
  end        = 24,
  cmax       = TRUE,
  tmax       = TRUE,
  aucinf.obs = TRUE,
  half.life  = TRUE
)

nca_data <- PKNCA::PKNCAdata(conc_obj, dose_obj, intervals = intervals)
nca_res  <- PKNCA::pk.nca(nca_data)
#> Warning: Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)
#> Too few points for half-life calculation (min.hl.points=3 with only 2 points)

Comparison against published PK behaviour

Wagh 2021 does not tabulate per-dose-group NCA, but reports descriptive PK behaviour in the Results paragraph ‘Pharmacokinetic analysis’: an absorption half-life of about 3 min, an alpha-phase half-life of about 0.20 h at therapeutic concentrations, a total V_d around 0.35 L/kg, and a terminal half-life of 9.7 h at very low concentrations. The simulated typical-value profile reflects these qualitative features (rapid SC absorption, biphasic decline).

# Per-arm NCA summary (median across the cohort)
nca_summary <- as.data.frame(nca_res$result) |>
  dplyr::filter(PPTESTCD %in% c("cmax", "tmax", "aucinf.obs", "half.life")) |>
  dplyr::group_by(arm, `NCA parameter` = PPTESTCD) |>
  dplyr::summarise(
    Median = round(median(PPORRES, na.rm = TRUE), 3),
    `5%`   = round(quantile(PPORRES, 0.05, na.rm = TRUE), 3),
    `95%`  = round(quantile(PPORRES, 0.95, na.rm = TRUE), 3),
    .groups = "drop"
  ) |>
  dplyr::mutate(
    `NCA parameter` = dplyr::recode(`NCA parameter`,
                                    cmax       = "Cmax (mg/L)",
                                    tmax       = "Tmax (h)",
                                    aucinf.obs = "AUC0-inf (mg*h/L)",
                                    half.life  = "t1/2 terminal (h)")
  )

knitr::kable(
  nca_summary,
  caption = "Simulated NCA summary by arm (steady-state interval after last dose).",
  align   = c("l", "l", "r", "r", "r")
)
Simulated NCA summary by arm (steady-state interval after last dose).
arm NCA parameter Median 5% 95%
100 mg/kg QD (S1) AUC0-inf (mg*h/L) NA NA NA
100 mg/kg QD (S1) Cmax (mg/L) 0.055 0.045 0.084
100 mg/kg QD (S1) t1/2 terminal (h) NA NA NA
100 mg/kg QD (S1) Tmax (h) 0.000 0.000 0.000
100 mg/kg QD (S2) AUC0-inf (mg*h/L) NA NA NA
100 mg/kg QD (S2) Cmax (mg/L) 0.057 0.033 0.079
100 mg/kg QD (S2) t1/2 terminal (h) NA NA NA
100 mg/kg QD (S2) Tmax (h) 0.000 0.000 0.000
200 mg/kg TIW (S1) AUC0-inf (mg*h/L) NA NA NA
200 mg/kg TIW (S1) Cmax (mg/L) 0.092 0.064 0.131
200 mg/kg TIW (S1) t1/2 terminal (h) NA NA NA
200 mg/kg TIW (S1) Tmax (h) 0.000 0.000 0.000
50 mg/kg BID (S1) AUC0-inf (mg*h/L) 69.499 52.779 91.463
50 mg/kg BID (S1) Cmax (mg/L) 148.581 140.595 155.207
50 mg/kg BID (S1) t1/2 terminal (h) 12.586 12.555 12.672
50 mg/kg BID (S1) Tmax (h) 0.100 0.100 0.100

Assumptions and deviations

  • Smoothed PAE compartment. The paper’s PAE-compartment definition is piecewise (d C_PAE / dt = d Cc / dt while Cc > C_PAE; otherwise -K_PAE * C_PAE). A literal ifelse implementation creates a discontinuity at each peak crossing that breaks the rxode2 LSODA solver (“excess accuracy requested for precision of machine”). The model file approximates the same dynamic with a fast-load + slow-decay formulation:

    d C_PAE / dt = k_load_pae * weight_up * (Cc - C_PAE) - K_PAE * C_PAE

    where weight_up = 1 / (1 + exp(-50 * (Cc - C_PAE))) is a sigmoidal switch (sharpness 50/mg/L) and k_load_pae = 100/h is much faster than K_a = 15.4/h. Across the Wagh 2021 dose range the simulated peak C_PAE / max(Cc) ratio is > 0.99, the post-peak decay matches the paper’s reported PAE half-life (48.8 h) to within the solver tolerance, and the parameter values reported by Wagh 2021 enter unchanged.

  • Healthy-animal PK parameters not encoded. Wagh 2021 Table 1 reports the healthy-animal popPK estimates alongside the infected-animal estimates and concluded (“no relevant difference”) that the two sets describe the same drug behaviour. The packaged model uses the infected-animal values (the final integrated PK/PD model parameters); the healthy estimates are documented in the population metadata for reference but are not encoded as a separate STUDY covariate.

  • Per-study log10-CFU residual SD. The packaged residual SD addSd_log_cfu = 0.284 is the study 1 RUV1 (the larger / more conservative of the two). Study 2’s RUV2 = 0.173 (Wagh 2021 Table 3) is documented here but not exposed as a separate parameter – a user who needs study-specific stochastic VPCs for study 2 should override addSd_log_cfu directly in params = c(addSd_log_cfu = 0.173) at simulation time.

  • K_PAE precision. K_PAE = 0.0142 1/h is reported with RSE 74.6% and a bootstrap 90% CI of 0.0048 - 0.0254 (Wagh 2021 Table 3). The point estimate is encoded directly; stochastic exploration of the K_PAE uncertainty should use the bootstrap distribution rather than the asymptotic SE.

  • Simulation start time. The packaged model’s bacteria(0) is set to 10^1.77 = 58.9 CFU (the aerosol-infection-time bacterial load fitted by the natural-growth model; Wagh 2021 Table 2). Simulating the full 34-day natural-growth lead-in adds wall-clock time without changing the post-treatment endpoint, so this vignette overrides bacteria(0) per arm to the observed start-of-treatment bacterial load (7.08 log CFU for study 1; 5.37 log CFU for study 2) via the rxSolve(..., inits = ...) argument. Downstream users who do want to simulate the natural-growth lead-in can use the model’s default initial condition and extend the event table back to time = 0 (aerosol-infection time).

  • Bioavailability anchor. Subcutaneous F is not estimated in the paper; all PK parameters are reported as apparent (/F). The model anchors lfdepot <- fixed(log(1)) so the apparent V and CL values match the paper without imposing an unverified F.

  • Plasma free fraction (1810 in mouse plasma) and MIC are reported in the model file but enter only the (post hoc) PK/PD-index analysis in the paper (free f C_max / MIC, f AUC / MIC, f % T_MIC); they are NOT part of the integrated PK/PD with PAE that this model encodes. The numeric values (f_u = 0.31, MIC = 1.6 mg/L) are reproduced in the Source-trace section above for reference.