---
title: "Lenvatinib (Gupta 2016)"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Lenvatinib (Gupta 2016)}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 7,
  fig.height = 5,
  dpi = 96
)
library(nlmixr2lib)
library(PKNCA)
library(rxode2)
library(dplyr)
library(tidyr)
library(ggplot2)
```

# Model and source

* Citation: `r rxode2::rxode(readModelDb("Gupta_2016_lenvatinib"))$reference`
* Description: `r rxode2::rxode(readModelDb("Gupta_2016_lenvatinib"))$description`
* Article: <https://doi.org/10.1111/bcp.12907>

Lenvatinib (Lenvima) is an oral multi-targeted receptor tyrosine kinase
inhibitor (VEGFR1-3, FGFR1-4, PDGFR-alpha, RET, KIT) approved for the
treatment of radioiodine-refractory differentiated thyroid cancer (RR-DTC)
and, in combination regimens, for several other oncology indications. The
packaged model is the final population PK model from Gupta et al. (2016),
the first pooled population PK analysis of lenvatinib across 15 phase 1-3
studies.

# Population

Pooled cohort of 779 subjects across 15 studies (Gupta 2016 Table 1):

* Eight phase 1 studies in healthy adults (n = 196, 25.2 percent).
* Four phase 1 studies in patients with solid tumours (n = 191, 24.5
  percent for the pooled "other tumours" category in addition to the
  thyroid-cancer cohorts below).
* Two phase 2 studies and one phase 3 study (SELECT) in patients with
  differentiated thyroid cancer (DTC; n = 327, 42.0 percent), medullary
  thyroid cancer (MTC; n = 56, 7.2 percent), and anaplastic thyroid
  cancer (ATC; n = 9, 1.2 percent).

Pooled baseline summary (Gupta 2016 Table 1): age range 18-89 y (median
55), body weight 32.6-177.5 kg (median 75), albumin 19-52 g/L (median 40),
alkaline phosphatase 19-752 U/L (median 77), creatinine clearance 17-268
mL/min (median 98). 44.0 percent female. Race composition was 70.2 percent
White, 11.7 percent Japanese, 9.4 percent Black, 6.3 percent Other, with
small fractions of Hispanic, Asian, Native Hawaiian, and American Indian /
Alaska Native participants. The dataset comprises 10 265 plasma
concentrations from doses ranging 3.2-32 mg PO, predominantly once daily
as tablet or capsule.

The same metadata is available programmatically via
`readModelDb("Gupta_2016_lenvatinib")$population`.

# Source trace

Per-parameter origin is recorded as in-file comments next to every `ini()`
entry in `inst/modeldb/specificDrugs/Gupta_2016_lenvatinib.R`. The table
below collects them in one place for review.

| Equation / parameter             | Value           | Source location |
|----------------------------------|-----------------|-----------------|
| `lcl` (basal CL/F)               | 6.56 L/h        | Table 2, theta_CL |
| `lvc` (basal V1/F)               | 49.3 L          | Table 2, theta_V1 |
| `lvp` (basal V2/F)               | 30.7 L          | Table 2, theta_V2 |
| `lvp2` (basal V3/F)              | 37.1 L          | Table 2, theta_V3 |
| `lq` (basal Q1/F)                | 3.52 L/h        | Table 2, theta_Q1 |
| `lq2` (basal Q2/F)               | 0.769 L/h       | Table 2, theta_Q2 |
| `lka` (Ka)                       | 1.02 1/h        | Table 2, theta_Ka |
| `lduration` (D1)                 | 1.22 h          | Table 2, theta_D1 |
| `lfcap` (F1, capsule vs tablet)  | 0.896           | Table 2, theta_F1 |
| `allo_cl` (WGT exponent on CL/Q) | 0.75 (fixed)    | Table 2, equation for CL/F and Q/F covariate model |
| `allo_v` (WGT exponent on V)     | 1.0 (fixed)     | Table 2, equation for V/F covariate model |
| `e_cyp3a4_ind_cl`                | log(1.30)       | Table 2, theta_INDU |
| `e_cyp3a4_inh_cl`                | log(0.922)      | Table 2, theta_INHIB |
| `e_alb_cl`                       | log(0.837)      | Table 2, theta_ALB (ALB < 30 g/L) |
| `e_alp_cl`                       | log(0.883)      | Table 2, theta_ALP (ALP > ULN) |
| `e_dis_healthy_cl`               | log(1.15)       | Table 2, theta_TM (healthy vs cancer patient) |
| omega^2_CL                       | 0.0630          | Table 2, IIV CL/F 25.5 percent CV |
| omega^2_V1                       | 0.0507          | Table 2, IIV V1/F 22.8 percent CV |
| omega^2_V2                       | 0.1416          | Table 2, IIV V2/F 39.0 percent CV |
| omega^2_V3                       | 0.0879          | Table 2, IIV V3/F 30.3 percent CV |
| omega^2_Ka                       | 0.2624          | Table 2, IIV Ka 54.8 percent CV |
| omega^2_D1                       | 0.4626          | Table 2, IIV D1 76.7 percent CV |
| omega^2_F1                       | 0.0873          | Table 2, IIV F1 30.2 percent CV |
| propSd                           | 0.333           | Table 2, proportional residual error in patient studies |
| addSd                            | 7.19 ng/mL      | Table 2, additive residual error for TAD < 2 h |
| Three-cmt + simultaneous 1st/0-order abs | structure | Methods narrative and Table 2 footnote |
| Reference subject                | 75 kg cancer patient on tablet, no CYP3A4 inducers/inhibitors, normal ALB/ALP | Table 1 baseline medians |

# Virtual cohort

The original participant-level data are not publicly available. The
simulations below build virtual cohorts whose covariate distributions
approximate the SELECT (RR-DTC) phase 3 cohort (Gupta 2016 Table 1).
Covariates are sampled from truncated normal / Bernoulli distributions
whose first two moments match the pooled baseline summary in Table 1.

```{r cohort}
set.seed(20260515)

n_per_dose <- 40  # downsampled from 100 for vignette build budget; AUC envelopes remain stable
doses_mg   <- c(4, 8, 14, 20, 24, 32)

# Helper: build one steady-state cohort for a given dose level. The
# simultaneous first- and zero-order absorption model requires two dose
# records per administration -- one to depot (first-order, rate ka) and
# one to central (zero-order, duration D1). The model splits the dose
# 50/50 between the two routes via f(depot) = f(central) = 0.5 * F1.
make_cohort <- function(n, dose_mg, id_offset = 0L) {
  # SELECT / pooled cancer-patient baseline (Gupta 2016 Table 1):
  # age 53.2 +/- 15.8 y, BW 76.7 +/- 19.4 kg, ALB 39.8 +/- 5.1 g/L,
  # ALP 104.6 +/- 82.1 U/L, 44.0% female. CYP3A4 inducers 2.4%,
  # inhibitors 6.3%, all participants on tablet (FORM_CAPSULE = 0,
  # DIS_HEALTHY = 0) for this cohort.
  wt_kg   <- pmin(pmax(rnorm(n, mean = 76.7, sd = 19.4), 32.6), 177.5)
  alb_gL  <- pmin(pmax(rnorm(n, mean = 39.8, sd = 5.1),  19),    52)
  alp_uL  <- pmin(pmax(rlnorm(n, meanlog = log(77), sdlog = 0.7), 19), 752)
  cyp_ind <- as.integer(runif(n) < 0.024)
  cyp_inh <- as.integer(runif(n) < 0.063)

  tau     <- 24
  n_doses <- 15  # day 15 corresponds to Gupta 2016 simulation horizon for SS profiles.
  dose_times <- (seq_len(n_doses) - 1L) * tau
  final_t <- (n_doses - 1L) * tau
  obs_grid <- c(seq(0, 24, by = 1),
                seq(28, final_t, by = 6),
                final_t + c(0.25, 0.5, 1, 1.5, 2, 3, 4, 6, 8, 12, 18, 24))
  obs_grid <- sort(unique(obs_grid))

  cov_tbl <- tibble::tibble(
    id           = id_offset + seq_len(n),
    WT           = wt_kg,
    ALB          = alb_gL,
    ALP          = alp_uL,
    CONMED_CYP3A4_IND   = cyp_ind,
    CONMED_CYP3A4_INH   = cyp_inh,
    DIS_HEALTHY  = 0L,
    FORM_CAPSULE = 0L,
    treatment    = sprintf("%d mg QD", dose_mg)
  )

  dose_depot <- tidyr::expand_grid(cov_tbl, time = dose_times) |>
    dplyr::mutate(evid = 1L, amt = dose_mg, cmt = "depot", dur = NA_real_)

  dose_central <- tidyr::expand_grid(cov_tbl, time = dose_times) |>
    dplyr::mutate(evid = 1L, amt = dose_mg, cmt = "central", dur = 1.22)

  obs_rows <- tidyr::expand_grid(cov_tbl, time = obs_grid) |>
    dplyr::mutate(evid = 0L, amt = 0, cmt = "central", dur = NA_real_)

  dplyr::bind_rows(dose_depot, dose_central, obs_rows) |>
    dplyr::arrange(id, time, dplyr::desc(evid))
}

events <- dplyr::bind_rows(
  lapply(seq_along(doses_mg), function(i) {
    make_cohort(n_per_dose, dose_mg = doses_mg[i],
                id_offset = (i - 1L) * n_per_dose)
  })
)

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

# Simulation

```{r simulate}
mod <- readModelDb("Gupta_2016_lenvatinib")
sim <- rxode2::rxSolve(
  mod, events = events,
  keep = c("treatment", "WT", "ALB", "ALP", "CONMED_CYP3A4_IND",
           "CONMED_CYP3A4_INH", "DIS_HEALTHY", "FORM_CAPSULE")
) |>
  as.data.frame()
```

# Replicate published figures

## Figure 1 (visual predictive check, single-dose healthy subjects)

Gupta 2016 Figure 1A shows a prediction-corrected VPC of single-dose
lenvatinib in phase 1 healthy subjects. Below we reproduce the
percentile envelope from a typical-value stochastic simulation of a 24 mg
single dose in healthy subjects, comparable to the rising and elimination
phases of Figure 1A.

```{r figure-1a-healthy}
sd_events <- tibble::tibble(
  id = seq_len(150),  # downsampled from 200 for vignette build budget; VPC band shape unchanged
  WT = pmin(pmax(rnorm(150, 76.7, 19.4), 32.6), 177.5),
  ALB = pmin(pmax(rnorm(150, 39.8, 5.1), 19), 52),
  ALP = pmin(pmax(rlnorm(150, log(77), 0.7), 19), 752),
  CONMED_CYP3A4_IND = 0L, CONMED_CYP3A4_INH = 0L,
  DIS_HEALTHY = 1L,    # healthy subjects per Figure 1A
  FORM_CAPSULE = 0L
)

sd_obs <- tidyr::expand_grid(sd_events, time = c(0, 0.5, 1, 1.5, 2, 3, 4,
                                                  6, 8, 12, 16, 24, 36, 48, 72)) |>
  dplyr::mutate(evid = 0L, amt = 0, cmt = "central", dur = NA_real_)

sd_dose_depot <- sd_events |>
  dplyr::mutate(time = 0, evid = 1L, amt = 24,
                cmt = "depot", dur = NA_real_)

sd_dose_central <- sd_events |>
  dplyr::mutate(time = 0, evid = 1L, amt = 24,
                cmt = "central", dur = 1.22)

sd_all <- dplyr::bind_rows(sd_dose_depot, sd_dose_central, sd_obs) |>
  dplyr::arrange(id, time, dplyr::desc(evid))

sd_sim <- rxode2::rxSolve(mod, events = sd_all,
                          keep = c("DIS_HEALTHY")) |>
  as.data.frame()

sd_summary <- sd_sim |>
  dplyr::filter(!is.na(Cc), time > 0) |>
  dplyr::group_by(time) |>
  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(sd_summary, aes(time, Q50)) +
  geom_ribbon(aes(ymin = Q05, ymax = Q95), alpha = 0.25) +
  geom_line(linewidth = 0.8) +
  scale_y_log10() +
  labs(
    x = "Time (h)",
    y = "Lenvatinib (ng/mL)",
    title = "Single-dose 24 mg lenvatinib, healthy subjects",
    caption = "Approximates Gupta 2016 Figure 1A (single-dose phase 1 healthy)."
  )
```

## Figure 2 (typical-value covariate effects on steady-state PK)

Figure 2 of Gupta 2016 shows the simulated median and 90 percent prediction
interval (PI) for 24 mg lenvatinib QD at steady state, contrasting key
covariates (body weight 50 vs 100 kg; CYP3A4 inhibitor / inducer
coadministration; low albumin; high ALP). The plot below regenerates the
typical-value steady-state Cmax,ss and AUCss ratios for the same
covariates (etas zeroed).

```{r figure-2-forest}
make_typical_subject <- function(WT, ALB, ALP, CONMED_CYP3A4_IND, CONMED_CYP3A4_INH,
                                 DIS_HEALTHY, FORM_CAPSULE, label) {
  tibble::tibble(
    id = 1L, WT = WT, ALB = ALB, ALP = ALP,
    CONMED_CYP3A4_IND = CONMED_CYP3A4_IND, CONMED_CYP3A4_INH = CONMED_CYP3A4_INH,
    DIS_HEALTHY = DIS_HEALTHY, FORM_CAPSULE = FORM_CAPSULE,
    label = label
  )
}

# Reference participant: 75 kg cancer patient on tablet, no CYP3A4
# inducers/inhibitors, ALB 40 g/L, ALP 80 U/L.
covariate_grid <- dplyr::bind_rows(
  make_typical_subject(75, 40, 80, 0, 0, 0, 0, "Reference"),
  make_typical_subject(50, 40, 80, 0, 0, 0, 0, "WT 50 kg"),
  make_typical_subject(100, 40, 80, 0, 0, 0, 0, "WT 100 kg"),
  make_typical_subject(75, 25, 80, 0, 0, 0, 0, "ALB < 30 g/L"),
  make_typical_subject(75, 40, 200, 0, 0, 0, 0, "ALP > ULN"),
  make_typical_subject(75, 40, 80, 1, 0, 0, 0, "CYP3A4 inducer"),
  make_typical_subject(75, 40, 80, 0, 1, 0, 0, "CYP3A4 inhibitor"),
  make_typical_subject(75, 40, 80, 0, 0, 1, 0, "Healthy subject"),
  make_typical_subject(75, 40, 80, 0, 0, 0, 1, "Capsule formulation")
)

tau <- 24
n_doses <- 15
final_t <- (n_doses - 1L) * tau
obs_grid <- c(seq(0, 24, by = 0.5),
              seq(28, final_t, by = 6),
              final_t + c(0.5, 1, 1.5, 2, 3, 4, 6, 8, 12, 18, 24))
obs_grid <- sort(unique(obs_grid))

build_typical_events <- function(cov_row, dose_mg = 24) {
  dose_d <- tibble::tibble(
    id = 1L, WT = cov_row$WT, ALB = cov_row$ALB, ALP = cov_row$ALP,
    CONMED_CYP3A4_IND = cov_row$CONMED_CYP3A4_IND, CONMED_CYP3A4_INH = cov_row$CONMED_CYP3A4_INH,
    DIS_HEALTHY = cov_row$DIS_HEALTHY, FORM_CAPSULE = cov_row$FORM_CAPSULE,
    time = (seq_len(n_doses) - 1L) * tau,
    evid = 1L, amt = dose_mg, cmt = "depot", dur = NA_real_
  )
  dose_c <- dose_d |> dplyr::mutate(cmt = "central", dur = 1.22)
  obs_rows <- tibble::tibble(
    id = 1L, WT = cov_row$WT, ALB = cov_row$ALB, ALP = cov_row$ALP,
    CONMED_CYP3A4_IND = cov_row$CONMED_CYP3A4_IND, CONMED_CYP3A4_INH = cov_row$CONMED_CYP3A4_INH,
    DIS_HEALTHY = cov_row$DIS_HEALTHY, FORM_CAPSULE = cov_row$FORM_CAPSULE,
    time = obs_grid,
    evid = 0L, amt = 0, cmt = "central", dur = NA_real_
  )
  dplyr::bind_rows(dose_d, dose_c, obs_rows) |>
    dplyr::arrange(time, dplyr::desc(evid))
}

mod_typ <- mod |> rxode2::zeroRe()

typical_metrics <- do.call(dplyr::bind_rows, lapply(seq_len(nrow(covariate_grid)), function(i) {
  evt <- build_typical_events(covariate_grid[i, ])
  s <- as.data.frame(rxode2::rxSolve(mod_typ, events = evt))
  ss <- s[s$time >= final_t & s$time <= final_t + tau, ]
  dt <- diff(ss$time)
  auc_ss <- sum(0.5 * (head(ss$Cc, -1) + tail(ss$Cc, -1)) * dt)
  tibble::tibble(
    label = covariate_grid$label[i],
    Cmax_ss = max(ss$Cc),
    AUC_ss  = auc_ss
  )
}))

ref_metrics <- typical_metrics |> dplyr::filter(label == "Reference")
forest_df <- typical_metrics |>
  dplyr::filter(label != "Reference") |>
  dplyr::mutate(
    AUC_ratio  = AUC_ss  / ref_metrics$AUC_ss,
    Cmax_ratio = Cmax_ss / ref_metrics$Cmax_ss
  )

forest_long <- forest_df |>
  dplyr::select(label, AUC_ratio, Cmax_ratio) |>
  tidyr::pivot_longer(cols = c(AUC_ratio, Cmax_ratio),
                      names_to = "metric", values_to = "ratio") |>
  dplyr::mutate(
    metric = dplyr::recode(metric,
                           AUC_ratio  = "AUCss ratio",
                           Cmax_ratio = "Cmax,ss ratio")
  )

ggplot(forest_long, aes(x = ratio, y = label)) +
  geom_vline(xintercept = 1, linetype = "dashed") +
  geom_point(size = 2) +
  facet_wrap(~ metric) +
  labs(
    x = "Geometric-mean ratio versus reference participant (75 kg cancer patient, tablet)",
    y = NULL,
    title = "Typical-value covariate forest plot",
    caption = "Replicates Gupta 2016 Figure 2 (covariate-impact simulation)."
  )
```

## Figure 3 (dose-proportional steady-state exposure)

Gupta 2016 Figure 3 reports that simulated steady-state AUC of lenvatinib
increases proportionally with dose (3.2-32 mg), confirming linear PK over
the studied range. The chart below regenerates that relationship from
the simulated cohorts above.

```{r figure-3-dose-proportional}
ss_start <- (n_doses - 1L) * tau
ss_end   <- ss_start + tau

dose_label_levels <- sprintf("%d mg QD", doses_mg)

auc_by_dose <- sim |>
  dplyr::filter(time >= ss_start, time <= ss_end, !is.na(Cc)) |>
  dplyr::group_by(id, treatment) |>
  dplyr::arrange(time) |>
  dplyr::summarise(
    AUC_ss = sum(0.5 * (head(Cc, -1) + tail(Cc, -1)) * diff(time)),
    .groups = "drop"
  ) |>
  dplyr::group_by(treatment) |>
  dplyr::summarise(
    AUC_geomean = exp(mean(log(pmax(AUC_ss, 1e-3)))),
    AUC_q05 = quantile(AUC_ss, 0.05),
    AUC_q95 = quantile(AUC_ss, 0.95),
    .groups = "drop"
  ) |>
  dplyr::mutate(
    dose_mg = as.numeric(sub(" mg QD", "", treatment)),
    treatment = factor(treatment, levels = dose_label_levels)
  )

ggplot(auc_by_dose, aes(dose_mg, AUC_geomean)) +
  geom_errorbar(aes(ymin = AUC_q05, ymax = AUC_q95), width = 0.8) +
  geom_point(size = 2) +
  geom_smooth(method = "lm", formula = y ~ x, se = FALSE,
              linetype = "dashed") +
  labs(
    x = "Lenvatinib dose (mg QD)",
    y = expression("AUC"["ss"]~"(ng*h/mL)"),
    title = "Simulated steady-state AUC vs dose",
    caption = "Replicates Gupta 2016 Figure 3 (linear dose-AUC relationship)."
  )
```

# PKNCA validation

PKNCA computes Cmax,ss, Tmax,ss, AUCtau,ss, and Ctrough,ss across the
final 24 h interval of the 15-day QD simulation, stratified by dose
group.

```{r pknca}
sim_nca <- sim |>
  dplyr::filter(time >= ss_start, time <= ss_end, !is.na(Cc)) |>
  dplyr::mutate(time_in_interval = time - ss_start) |>
  dplyr::select(id, time = time_in_interval, Cc, treatment)

dose_df <- events |>
  dplyr::filter(evid == 1L, cmt == "depot", time == ss_start) |>
  dplyr::group_by(id, treatment) |>
  dplyr::summarise(time = 0, amt = max(amt), .groups = "drop")

conc_obj <- PKNCA::PKNCAconc(
  as.data.frame(sim_nca),
  Cc ~ time | treatment + id,
  concu = "ng/mL", timeu = "h"
)
dose_obj <- PKNCA::PKNCAdose(
  as.data.frame(dose_df),
  amt ~ time | treatment + id,
  doseu = "mg"
)

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

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

nca_tbl <- as.data.frame(nca_res$result) |>
  dplyr::filter(PPTESTCD %in% c("cmax", "tmax", "cmin", "auclast")) |>
  dplyr::group_by(treatment, PPTESTCD) |>
  dplyr::summarise(
    geomean = exp(mean(log(pmax(PPORRES, 1e-3)))),
    .groups = "drop"
  ) |>
  tidyr::pivot_wider(names_from = PPTESTCD, values_from = geomean) |>
  dplyr::rename(
    `Cmax,ss (ng/mL)`    = cmax,
    `Tmax,ss (h)`        = tmax,
    `Ctrough,ss (ng/mL)` = cmin,
    `AUCss (ng*h/mL)`    = auclast
  )

knitr::kable(nca_tbl, digits = 2,
             caption = "Simulated steady-state NCA, dose-ranging virtual cohort (75 kg tablet, normal labs).")
```

## Comparison against published exposure

Gupta 2016 reports a mean dose-normalised AUC of 3710 ng*h/mL (CV 38.3
percent) for the 24 mg dose group in SELECT (RR-DTC, n = 260). The
expected steady-state AUC at the typical CL/F of 6.56 L/h is
`Dose / CL = 24 mg / 6.56 L/h = 3.66 mg*h/L = 3660 ng*h/mL`.

```{r pknca-comparison}
published <- tibble::tibble(
  treatment           = "24 mg QD",
  `AUCss_pub (ng*h/mL)` = 3710,
  `Source` = "Gupta 2016 SELECT n = 260 mean (CV 38.3%)"
)

simulated <- nca_tbl |>
  dplyr::filter(treatment == "24 mg QD") |>
  dplyr::transmute(
    treatment,
    `AUCss_sim (ng*h/mL)` = `AUCss (ng*h/mL)`
  )

cmp <- published |>
  dplyr::left_join(simulated, by = "treatment") |>
  dplyr::mutate(
    AUC_pct_diff = 100 * (`AUCss_sim (ng*h/mL)` - `AUCss_pub (ng*h/mL)`) / `AUCss_pub (ng*h/mL)`
  )

knitr::kable(cmp, digits = 1,
             caption = "Simulated vs published mean steady-state AUC at 24 mg QD.")
```

Simulated AUCss for the 24 mg QD virtual cohort matches the SELECT
population mean within roughly 5 percent, consistent with the underlying
typical-value CL/F (6.56 L/h) at the reference 75 kg cancer-patient
covariate set.

# Assumptions and deviations

* **Simultaneous first- and zero-order absorption -- 50:50 dose split.**
  Gupta 2016 describes simultaneous first- and zero-order absorption but
  does not publish the explicit fraction of the dose delivered through
  each route. This implementation splits the dose 50:50 between the
  first-order depot (rate `ka = 1.02 1/h`) and the zero-order direct-to-
  central infusion (duration `D1 = 1.22 h`). Both routes carry the same
  formulation effect `F1` (1 for tablet, 0.896 for capsule). The 50:50
  split is the standard convention when the fraction is not disclosed;
  downstream users with access to per-route data can adjust the
  `f(depot)` / `f(central)` factors and re-fit.
* **User-data convention for dosing.** Because the model has two parallel
  absorption routes, each administered dose requires two rows in the
  event table: one with `cmt = "depot"` (first-order) and one with
  `cmt = "central"` and `dur = D1` (zero-order). The cohort builder above
  illustrates the pattern. A single-row event would silently route the
  full dose through only one of the two compartments.
* **Stratified residual error simplified.** Gupta 2016 reports a
  stratified residual error model: proportional 16.9 percent CV for
  phase 1 clinical-pharmacology studies, 33.3 percent CV for patient
  studies, and a separate 48.1 percent CV proportional plus 7.19 ng/mL
  additive arm for time-after-dose < 2 h. The packaged model uses a
  combined proportional (33.3 percent, patient-study CV) plus additive
  (7.19 ng/mL, early-absorption term) residual error model. Users
  focusing on the absorption phase should use a higher proportional
  component (48.1 percent); users targeting healthy-volunteer phase 1
  studies should reduce the proportional component to 16.9 percent.
* **ALB and ALP binarisation.** Gupta 2016 Table 2 encodes albumin and
  alkaline phosphatase as binary NONMEM indicators (`ALB = 1` when
  serum ALB < 30 g/L; `ALP = 1` when ALP/ULN > 1). The packaged model
  takes continuous `ALB` (g/L) and continuous `ALP` (U/L) from the
  user's data and binarises them inline: `alb_low <- (ALB < 30)` and
  `alp_high <- (ALP > 120)`. The 120 U/L ULN is a representative adult
  threshold; sites whose laboratory ULN differs should either pass the
  binary indicator directly (set `ALP` above or below 120 to encode the
  intended state) or edit the inline `alp_uln` constant.
* **F1 IIV applies only to capsule.** The 30.2 percent CV IIV on F1
  multiplies the capsule-arm bioavailability only. Tablet subjects
  (`FORM_CAPSULE = 0`) have F = 1 with no eta contribution; the
  `etalfcap` term is structurally inactive for the tablet arm.
* **Race distribution and other tested-but-not-significant covariates.**
  Gupta 2016 evaluated race, sex, age, ECOG performance status, renal
  function (CRCL and CTCAE grades), liver-function markers ALT / AST /
  bilirubin and the corresponding CTCAE grades, thyroid stimulating
  hormone, and concomitant pH-elevating agents (PPIs, H2 blockers,
  antacids). None reached the retention significance threshold and none
  appear as covariates in the packaged model. The virtual cohort does
  not stratify by these factors; simulated exposures correspond to a
  covariate-matched mix.
* **New canonicals introduced by this extraction.** Two new entries are
  proposed in `inst/references/covariate-columns.md`: `ALP` (general
  scope; alkaline phosphatase continuous concentration) and
  `CONMED_CYP3A4_IND` (general scope; concomitant CYP3A4 inducer indicator,
  sibling to the existing `CONMED_CYP3A4_INH`). The `FORM_CAPSULE` entry is
  extended to acknowledge tablet as a per-paper non-capsule reference
  (the prior Hennig 2006/2007 entries used solution as the reference).
