---
title: "Adalimumab (Nader 2017)"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Adalimumab (Nader 2017)}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

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

# Model and source

* Citation: Nader A, Beck D, Noertersheuser P, Williams D, Mostafa N. Population Pharmacokinetics and Immunogenicity of Adalimumab in Adult Patients with Moderate-to-Severe Hidradenitis Suppurativa. *Clin Pharmacokinet* 2017;56(9):1091-1102.
* Article: <https://doi.org/10.1007/s40262-016-0502-4>
* PubMed: <https://pubmed.ncbi.nlm.nih.gov/28066879/>

A one-compartment population PK model with first-order subcutaneous absorption
and linear elimination, parameterized in apparent clearance (CL/F) and apparent
volume of distribution (V/F). Baseline body weight enters as a power effect on
both CL/F and V/F (Nader 2017 Eqs. 5 and 6); baseline C-reactive protein enters
as a power effect on CL/F; and anti-adalimumab-antibody (AAA) positivity enters
as a multiplicative power effect on CL/F that activates 21 days after the first
adalimumab dose. Model time origin (t = 0) coincides with the first SC dose.
Time unit is days, dose is mg, and concentration is ug/mL (equivalent to mg/L
for adalimumab).

# Population

The estimation dataset pooled one phase II study (NCT00918255, n = 143) and two
phase III studies -- PIONEER-I (NCT01468207, n = 295) and PIONEER-II
(NCT01468233, n = 162) -- for a total of 600 adult patients with moderate-to-
severe hidradenitis suppurativa (HS) who had inadequate response to a trial of
oral antibiotics and were naive to anti-tumor-necrosis-factor-alpha treatment.
Median age was 35 years (range 18-67), 34% male, median body weight 94 kg
(range 43-221), and median baseline CRP 7.9 mg/L (range 0.1-189). Hurley stage
II 53%, stage III 43%. Of 787 randomized patients 187 were excluded from the
population PK dataset (175 randomized to placebo or discontinued before receipt
of adalimumab; 12 without measurable adalimumab concentrations above the assay
LLOQ). See Nader 2017 Table 1 for the full demographics.

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

# Source trace

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

| Equation / parameter | Value | Source location |
|----------------------|-------|-----------------|
| `lka` = log(ka)      | log(0.195 /day)     | Table 2, ka row (0.195, RSE 9.90%) |
| `lcl` = log(CL/F)    | log(0.667 L/day)    | Table 2, CL/F row (0.667, RSE 2.65%) |
| `lvc` = log(V/F)     | log(13.5 L)         | Table 2, Vd/F row (13.5, RSE 3.40%) |
| `e_ada_cl`           | 6.76 (unitless)     | Table 2, AAA on CL/F row (theta_4 of Eq. 5) |
| `e_crp_cl`           | 0.173 (unitless)    | Table 2, CRP on CL/F row (theta_5 of Eq. 5) |
| `e_wt_cl`            | 0.888 (unitless)    | Table 2, Weight on CL/F row (theta_6 of Eq. 5) |
| `e_wt_vc`            | 0.707 (unitless)    | Table 2, Weight on V/F row (theta_7 of Eq. 6) |
| `etalcl` variance    | 0.346 (58.8% CV)    | Table 2, IIV on CL/F row |
| `etalvc` variance    | 0.091 (30.2% CV)    | Table 2, IIV on V/F row |
| `propSd`             | sqrt(0.052) = 0.228 | Table 2, proportional residual variance |
| `addSd`              | sqrt(2.06) = 1.435  | Table 2, additive residual variance (ug/mL)^2 |
| CRP reference        | 7.9 mg/L            | Eq. (5) COV_2 explanatory text; Table 1 median |
| WT reference         | 94 kg               | Eq. (5)/(6) COV_3/COV_4 explanatory text; Table 1 median |
| AAA-effect delay     | 21 days             | Eq. (5) COV_1 definition; Section 3.3 ("3 weeks after first adalimumab dose") |
| Equation: CL/F       | theta_1 * COV_1 * COV_2 * COV_3 * exp(eta_1) | Eq. (5) |
| Equation: V/F        | theta_3 * COV_4 * exp(eta_2) | Eq. (6) |
| ODE: depot           | d/dt(depot) = -ka * depot | Section 2.4 ("one-compartment model with first-order absorption and elimination") |
| ODE: central         | d/dt(central) = ka * depot - (CL/V) * central | Section 2.4 |
| Observation          | Cc = central / V/F  | Standard 1-cpt observation |
| Residual error       | combined additive + proportional | Section 2.4 Eq. (2) |

# Errata

No erratum has been published for Nader 2017 as of 2026-07-11.

Nader 2017 Section 3.3 states that "The final population PK model included
correlated exponential terms for inter-individual variability on CL/F and V/F"
but Table 2 reports only the two diagonal variances (0.346 for CL/F, 0.091 for
V/F) and does not publish the off-diagonal covariance or correlation coefficient.
The nlmixr2lib implementation therefore encodes the two etas as **independent**
(diagonal only). Simulated between-subject variability on CL/F and V/F is
therefore uncorrelated in this vignette, which slightly under-represents the
tail spread that the correlated etas would produce.

# Virtual cohort

We simulate two cohorts approximating the phase III (PIONEER-I / -II) dosing
regimen and one alternative regimen from the phase II open-label period:

- **Phase III 40 mg weekly** -- loading 160 mg SC at week 0 and 80 mg SC at
  week 2, then 40 mg SC weekly starting at week 4 through week 12 (the primary
  efficacy window; Nader 2017 Fig. 2).
- **Phase II 40 mg eow** -- loading 80 mg SC at week 0, then 40 mg SC every
  other week starting at week 1 through week 16 (Nader 2017 Fig. 1).

Each cohort is 150 subjects (below the 200-per-arm skill cap). The population
covariate mixture is set to the pooled Table 1 medians and marginals.

```{r cohort}
set.seed(20170109)  # paper published online 2017-01-09
n_sub <- 150L

make_covariates <- function(n, id_offset) {
  # Body weight: log-normal to match median 94 kg (range 43-221) in Table 1.
  # sd = 0.30 on log scale gives ~ 5th-95th of ~ 58-153 kg, a reasonable
  # trial-cohort spread.
  wt <- exp(rnorm(n, mean = log(94), sd = 0.30))
  wt <- pmin(pmax(wt, 43), 221)  # clip to the reported range

  # Baseline CRP: log-normal to match median 7.9 mg/L (range 0.1-189) in
  # Table 1. sd = 1.20 on log scale reproduces the heavy right tail.
  crp <- exp(rnorm(n, mean = log(7.9), sd = 1.20))
  crp <- pmin(pmax(crp, 0.1), 189)

  # AAA-positive fraction 6.5% pooled across phase III (Section 3.2); apply the
  # same marginal in the phase II cohort for simulation simplicity (the paper
  # reports 10.4% in phase II, but with only 16/154 patients and no separate
  # per-treatment estimate).
  ada_pos <- as.integer(runif(n) < 0.065)

  tibble::tibble(
    id      = id_offset + seq_len(n),
    WT      = wt,
    CRP     = crp,
    ADA_POS = ada_pos
  )
}

# Simulation grid: daily observations through week 16 (day 112), plus a
# pre-dose sample immediately before each maintenance dose so steady-state
# troughs are exact.
weeks_per_day <- 1 / 7
hours_per_day <- 24

# Cohort A -- Phase III 40 mg weekly
covA <- make_covariates(n_sub, id_offset =   0L) |>
  mutate(treatment = "Phase III 40 mg weekly")

doses_A_times <- c(
  0,                      # 160 mg week 0
  14,                     # 80 mg week 2
  seq(28, 12 * 7, by = 7) # 40 mg weekly starting week 4 through week 12
)
doses_A_amts <- c(160, 80, rep(40, length(doses_A_times) - 2))

# Cohort B -- Phase II 40 mg eow
covB <- make_covariates(n_sub, id_offset = n_sub) |>
  mutate(treatment = "Phase II 40 mg eow")

doses_B_times <- c(
  0,                       # 80 mg week 0
  seq(7, 16 * 7, by = 14)  # 40 mg every other week from week 1 through week 15
)
doses_B_amts <- c(80, rep(40, length(doses_B_times) - 1))

# Observation grid (same for both cohorts): weekly + 1 h pre-dose troughs.
obs_grid <- sort(unique(c(
  seq(0, 16 * 7, by = 1)  # daily observations to catch peaks and troughs
)))

build_events <- function(cov_df, dose_times, dose_amts) {
  cov_df |>
    rowwise() |>
    do({
      cov <- .
      dplyr::bind_rows(
        tibble::tibble(
          id        = cov$id,
          time      = dose_times,
          evid      = 1L,
          amt       = dose_amts,
          cmt       = "depot",
          WT        = cov$WT,
          CRP       = cov$CRP,
          ADA_POS   = cov$ADA_POS,
          treatment = cov$treatment
        ),
        tibble::tibble(
          id        = cov$id,
          time      = obs_grid,
          evid      = 0L,
          amt       = 0,
          cmt       = "central",
          WT        = cov$WT,
          CRP       = cov$CRP,
          ADA_POS   = cov$ADA_POS,
          treatment = cov$treatment
        )
      )
    }) |>
    ungroup()
}

events <- dplyr::bind_rows(
  build_events(covA, doses_A_times, doses_A_amts),
  build_events(covB, doses_B_times, doses_B_amts)
) |>
  arrange(id, time, desc(evid))

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

# Sanity: cohort ID ranges disjoint and matching n_sub.
stopifnot(length(unique(covA$id)) == n_sub)
stopifnot(length(unique(covB$id)) == n_sub)
stopifnot(length(intersect(covA$id, covB$id)) == 0L)
```

# Simulation

```{r simulate}
mod <- readModelDb("Nader_2017_adalimumab")

sim <- rxode2::rxSolve(
  mod,
  events = events,
  keep   = c("treatment", "WT", "CRP", "ADA_POS")
) |>
  as.data.frame()
```

# Replicate published figures

## Figure 2 -- Phase III 40 mg weekly steady-state concentrations

Nader 2017 Fig. 2 shows mean (SD) serum adalimumab concentrations for the two
PIONEER trials under 40 mg weekly dosing (starting at week 4 with 160 mg loading
at week 0 and 80 mg at week 2). The reported concentrations reach steady state
of 7-9 ug/mL by week 2 and are maintained through week 12.

```{r figure-2}
sim |>
  dplyr::filter(treatment == "Phase III 40 mg weekly", time > 0) |>
  dplyr::mutate(week = time / 7) |>
  dplyr::group_by(week) |>
  dplyr::summarise(
    mean_Cc = mean(Cc, na.rm = TRUE),
    sd_Cc   = sd(Cc,   na.rm = TRUE),
    .groups = "drop"
  ) |>
  ggplot(aes(week, mean_Cc)) +
  geom_ribbon(aes(ymin = pmax(mean_Cc - sd_Cc, 0), ymax = mean_Cc + sd_Cc),
              alpha = 0.2) +
  geom_line() +
  geom_hline(yintercept = c(7, 9), linetype = "dashed", alpha = 0.5) +
  labs(
    x       = "Time (weeks)",
    y       = "Mean serum adalimumab (ug/mL)",
    title   = "Phase III 40 mg weekly (replicates Nader 2017 Fig. 2)",
    subtitle = "Dashed lines mark the paper's 7-9 ug/mL steady-state band."
  )
```

## Figure 1 -- Phase II 40 mg weekly vs 40 mg eow (weeks 4, 8, 16)

Nader 2017 Fig. 1 reports phase II double-blind mean concentrations of
12.4 +/- 9.16 ug/mL (40 mg weekly, week 16) and 5.9 +/- 5.2 ug/mL (40 mg eow). Our
Cohort B ("Phase II 40 mg eow") reproduces the eow arm; the weekly arm of the
Phase II study used the same 160/80 loading + 40 mg weekly maintenance as the
Phase III studies, so the Phase III weekly cohort simulated above is the
appropriate comparator.

```{r figure-1}
phase2 <- sim |>
  dplyr::filter(time %in% (c(4, 8, 16) * 7)) |>
  dplyr::mutate(week = time / 7) |>
  dplyr::group_by(treatment, week) |>
  dplyr::summarise(
    mean_Cc = mean(Cc, na.rm = TRUE),
    sd_Cc   = sd(Cc,   na.rm = TRUE),
    .groups = "drop"
  )

phase2 |>
  dplyr::rename(
    "Treatment"                       = treatment,
    "Week"                            = week,
    "Mean serum adalimumab (ug/mL)"  = mean_Cc,
    "SD"                              = sd_Cc
  ) |>
  knitr::kable(
    digits  = 2,
    caption = "Simulated mean concentrations at weeks 4, 8, 16 (compare with Nader 2017 Fig. 1 phase II mean 12.4 +/- 9.16 ug/mL weekly, 5.9 +/- 5.2 ug/mL eow at week 16)."
  )
```

## AAA effect on serum adalimumab

Nader 2017 Section 3.2 and Fig. 5-6 show that AAA-positive patients have
substantially lower serum adalimumab concentrations from ~ week 4 onward (the
AAA effect activates 21 days after the first dose). We replicate this by
splitting the Phase III weekly cohort by AAA status.

```{r ada-effect}
sim |>
  dplyr::filter(treatment == "Phase III 40 mg weekly", time > 0) |>
  dplyr::mutate(
    week       = time / 7,
    ada_status = ifelse(ADA_POS == 1L, "AAA-positive", "AAA-negative")
  ) |>
  dplyr::group_by(week, ada_status) |>
  dplyr::summarise(
    mean_Cc = mean(Cc, na.rm = TRUE),
    sd_Cc   = sd(Cc,   na.rm = TRUE),
    n       = dplyr::n(),
    .groups = "drop"
  ) |>
  ggplot(aes(week, mean_Cc, colour = ada_status, fill = ada_status)) +
  geom_ribbon(aes(ymin = pmax(mean_Cc - sd_Cc, 0), ymax = mean_Cc + sd_Cc),
              alpha = 0.2, colour = NA) +
  geom_line() +
  geom_vline(xintercept = 3, linetype = "dashed", alpha = 0.5) +
  labs(
    x        = "Time (weeks)",
    y        = "Mean serum adalimumab (ug/mL)",
    colour   = "AAA status",
    fill     = "AAA status",
    title    = "Phase III 40 mg weekly by AAA status (replicates Nader 2017 Fig. 6)",
    subtitle = "Dashed line marks the AAA effect activation at day 21 (3 weeks)."
  )
```

# PKNCA validation

The paper does not report a tabulated NCA, but we can still compute the
steady-state Cmax/Cmin/AUC over the last complete weekly interval to summarize
the simulated exposure per treatment arm.

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

# Guarantee a t=0 row per (id, treatment); for extravascular pre-dose Cc = 0.
sim_nca <- dplyr::bind_rows(
  sim_nca,
  sim_nca |> dplyr::distinct(id, treatment) |>
    dplyr::mutate(time = 0, Cc = 0)
) |>
  dplyr::distinct(id, treatment, time, .keep_all = TRUE) |>
  dplyr::arrange(id, treatment, time)

dose_df <- events |>
  dplyr::filter(evid == 1L) |>
  dplyr::select(id, time, amt, treatment)
```

```{r pknca}
conc_obj <- PKNCA::PKNCAconc(sim_nca, Cc ~ time | treatment + id,
                             concu = "ug/mL", timeu = "day")
dose_obj <- PKNCA::PKNCAdose(dose_df, amt ~ time | treatment + id,
                             doseu = "mg")

# One steady-state interval per cohort.
#   Phase III weekly: last full weekly interval (days 77 -> 84) after week 4
#                     maintenance start.
#   Phase II eow:     last full 2-weekly interval (days 98 -> 112).
intervals <- data.frame(
  treatment = c("Phase III 40 mg weekly", "Phase II 40 mg eow"),
  start     = c(77, 98),
  end       = c(84, 112),
  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)
```

## Steady-state NCA summary vs. published mean concentrations

```{r nca-summary}
nca_tbl <- as.data.frame(nca_res$result)

nca_summary <- nca_tbl |>
  dplyr::group_by(treatment, PPTESTCD) |>
  dplyr::summarise(
    median = median(PPORRES, na.rm = TRUE),
    q05    = quantile(PPORRES, 0.05, na.rm = TRUE),
    q95    = quantile(PPORRES, 0.95, na.rm = TRUE),
    .groups = "drop"
  ) |>
  dplyr::arrange(treatment, PPTESTCD)

nca_summary |>
  dplyr::rename(
    "Treatment"           = treatment,
    "NCA parameter"       = PPTESTCD,
    "Median"              = median,
    "5th percentile"      = q05,
    "95th percentile"     = q95
  ) |>
  knitr::kable(
    digits  = 2,
    caption = "Simulated steady-state NCA over the last dosing interval by treatment (n = 150 virtual subjects per cohort). Concentration units ug/mL, time units days."
  )
```

The paper does not tabulate NCA, so the direct cross-check is against
Nader 2017 Section 3.1 and Figures 1-2:

- **Phase III 40 mg weekly steady state** -- Nader 2017 reports mean
  serum concentration 7-9 ug/mL from week 2 through week 12 (Fig. 2 and
  Section 3.1). The simulation median Cmax and Cmin over the last weekly
  interval should sit within, or just above, this band.
- **Phase II 40 mg eow week 16** -- Nader 2017 reports mean 5.9 +/- 5.2 ug/mL
  (Section 3.1). The simulation median Cav over the last eow interval should
  fall inside this range.

```{r cross-check}
phaseIII_ss <- sim |>
  dplyr::filter(treatment == "Phase III 40 mg weekly",
                time >= 14, time <= 84) |>
  dplyr::summarise(mean_Cc = mean(Cc, na.rm = TRUE))

phaseII_eow_wk16 <- sim |>
  dplyr::filter(treatment == "Phase II 40 mg eow", time == 16 * 7) |>
  dplyr::summarise(
    mean_Cc = mean(Cc, na.rm = TRUE),
    sd_Cc   = sd(Cc,   na.rm = TRUE)
  )

data.frame(
  Comparison = c(
    "Phase III 40 mg weekly, mean Cc (weeks 2-12)",
    "Phase II 40 mg eow, mean (SD) Cc at week 16"
  ),
  Published  = c("7-9 ug/mL (Fig. 2)", "5.9 +/- 5.2 ug/mL"),
  Simulated  = c(
    sprintf("%.1f ug/mL",           phaseIII_ss$mean_Cc),
    sprintf("%.1f +/- %.1f ug/mL", phaseII_eow_wk16$mean_Cc, phaseII_eow_wk16$sd_Cc)
  )
) |>
  knitr::kable(caption = "Simulated vs. published mean concentrations in the two dosing regimens.")
```

# Assumptions and deviations

- **Correlated IIV covariance not published.** Nader 2017 Section 3.3 states
  that the final model included correlated exponential etas on CL/F and V/F,
  but Table 2 reports only the two diagonal variances. The model file encodes
  the etas as independent (diagonal only); simulated between-subject variability
  on CL/F and V/F is therefore uncorrelated. See the Errata section.
- **Body weight distribution.** Table 1 reports only the median (94 kg) and
  range (43-221 kg). The virtual cohort samples WT from a log-normal with
  median 94 kg and log-sd 0.30, clipped to the reported range.
- **Baseline CRP distribution.** Table 1 reports median 7.9 mg/L and
  range 0.1-189. The virtual cohort samples CRP from a log-normal with
  median 7.9 mg/L and log-sd 1.20 to reproduce the heavy right tail.
- **AAA prevalence.** Set at 6.5% to match the pooled phase III rate reported
  in Section 3.2. Phase II reported 10.4% (16/154), but per-cohort simulation
  of that difference is not required for the pooled model, so the same 6.5%
  marginal is used in Cohort B.
- **AAA effect timing.** The AAA multiplier of 6.76 activates once the model's
  ODE time exceeds 21 days, per Nader 2017 Eq. (5). This assumes model time
  origin (t = 0) coincides with the first adalimumab dose (SC loading), which
  matches the paper's Fig. 2 / Fig. 6 x-axis convention.
- **Race / ethnicity.** Not reported in the population PK Table 1; not used
  as a covariate.
- **Sex.** Not retained as a covariate in the final model; not simulated.
- **Bioavailability.** F is not estimated. The published CL/F and V/F are
  apparent values that already absorb F. The depot compartment is dosed with
  implicit F = 1.
