Zhang_2023_valproic_acid_nonsaturable <- function() {
  description <- "One-compartment population PK model with first-order absorption for valproic acid in Chinese children with epilepsy (Zhang 2023 Model IV, the linear non-saturable binding model - the strategy the authors conclude is best suited to valproate protein-binding non-linearity). The disposition is carried on the UNBOUND concentration Cu = central/(V/F) with linear unbound clearance; the observed TOTAL plasma concentration is reconstructed as Cc = Cu + Cb, where the bound concentration adds a linear non-saturable term to the Langmuir isotherm, Cb = Bm * Cu / (Kd + Cu) + NS * Cu, with Bm, Kd and NS FIXED to the values of Gu 2021. Formulation-specific absorption rate constants FIXED from the literature (syrup 2.64 1/h reference, conventional tablet 1.57 1/h, sustained-release tablet 0.46 1/h). Best simulation-based behaviour of the five strategies compared; see modellib('Zhang_2023_valproic_acid_langmuir') for the pure Langmuir form it extends."
  reference <- "Zhang L, Liu M, Qin W, Shi D, Mao J, Li Z. Modeling the protein binding non-linearity in population pharmacokinetic model of valproic acid in children with epilepsy: a systematic evaluation study. Front Pharmacol. 2023;14:1228641. doi:10.3389/fphar.2023.1228641. PMID 37860114. Model IV binding isotherm from Eq. 6 with Kd, Bm and NS quoted from Gu X, Zhu M, Sheng C, et al. Population pharmacokinetics of unbound valproic acid in pediatric epilepsy patients in China: a protein binding model. Eur J Clin Pharmacol. 2021;77(7):999-1009. doi:10.1007/s00228-020-03080-y; parameter estimates from Supplementary Table S3."
  vignette <- "Zhang_2023_valproic_acid_protein_binding"
  units <- list(time = "h", dosing = "mg", concentration = "mg/L")

  # Issue #482: what each ODE state holds, in what amount units, in what
  # biological matrix.
  compartmentData <- list(
    depot   = list(analyte = "valproic acid", units = "mg", specimen = "administration site", verified = TRUE),
    central = list(analyte = "valproic acid", units = "mg", specimen = "plasma", verified = TRUE,
                   notes = "Holds administered valproic acid amount. Dividing by V/F yields the UNBOUND plasma concentration Cu, not the total: V/F is the unbound-referenced apparent volume (492 L, about 18-fold the base model's 27.8 L, matching the 7.8% unbound fraction the isotherm produces). The observed total concentration is the algebraic sum Cu + Cb.")
  )

  covariateData <- list(
    FORM_TABLET = list(
      description        = "Conventional (immediate-release) oral tablet formulation indicator",
      units              = "(binary)",
      type               = "binary",
      reference_category = "0 (oral syrup, the reference formulation in this cohort)",
      notes              = "Selects the FIXED conventional-tablet absorption rate constant Ka = 1.57 1/h. Oral syrup (Ka = 2.64 1/h) is the reference formulation: FORM_TABLET = 0 and FORM_VPA_SR = 0. The Zhang 2023 evaluation cohort itself received only syrup (194 records) or sustained-release tablet (61 records); the conventional-tablet level is retained so the full literature Ka set the authors quote is reachable. Zhang 2023 Supplementary Table S3 footnote.",
      source_name        = "conventional tablet"
    ),
    FORM_VPA_SR = list(
      description        = "Sustained-release valproic acid tablet formulation indicator",
      units              = "(binary)",
      type               = "binary",
      reference_category = "0 (oral syrup, the reference formulation in this cohort)",
      notes              = "Selects the FIXED sustained-release-tablet absorption rate constant Ka = 0.46 1/h. 61 of 255 observation records in the Zhang 2023 cohort were on the sustained-release tablet (Table 2). Zhang 2023 Supplementary Table S3 footnote.",
      source_name        = "sustained release tablet"
    )
  )

  covariatesDataExcluded <- list(
    WT = list(
      description = "Body weight",
      units       = "kg",
      type        = "continuous",
      notes       = "Collected in the evaluation cohort (median 19.00 kg, range 4.00-70.00; Zhang 2023 Table 2) but NOT retained on CL/F or V/F in any of the authors' own six models. Zhang 2023 Supplementary Table S3 reports no weight term."
    ),
    ALB = list(
      description = "Serum albumin",
      units       = "g/L",
      type        = "continuous",
      notes       = "Collected in the evaluation cohort (mean 42.21 g/L, median 42.10, range 29.70-70.50; Zhang 2023 Table 2). Eq. 6 folds the saturable binding capacity into the constant Bm and the non-saturable component into NS, neither of which is scaled by albumin; albumin enters explicitly only in the Model I one-binding-site strategy (Eq. 3)."
    )
  )

  population <- list(
    species        = "human",
    n_subjects     = 202,
    n_studies      = 1,
    age_range      = "0.17-15.00 years",
    age_median     = "4.92 years",
    weight_range   = "4.00-70.00 kg",
    weight_median  = "19.00 kg",
    sex_female_pct = 31.2,
    race_ethnicity = "Chinese (single-centre Han-predominant cohort; sub-ethnicity not reported)",
    disease_state  = "Childhood epilepsy on maintenance valproic acid with routine therapeutic drug monitoring",
    dose_range     = "60-1250 mg/day (median 480); 8.70-57.69 mg/kg/day (median 23.44); oral syrup or sustained-release tablet given once, twice or three times daily",
    regions        = "China (Wuhan Children's Hospital, single centre, January 2016 - November 2018)",
    notes          = "255 total plasma valproic acid trough concentrations in 202 children (139 male / 63 female), measured by gas chromatography (LOQ 1 mg/L, calibration range 12.5-150 mg/L, CV < 10%). Observed concentrations 22.60-118.50 mg/L (median 50.40). All samples were troughs collected under steady-state conditions, so absorption and distribution parameters are only weakly identified. Lowest objective function value of the five protein-binding strategies bar the simple exponent model (1617.5 against a base-model 1752.1), MDPE -0.19%, MAPE 19.8%, F30 72.16% - the best F30 of any model in the paper. The authors' Conclusion names this the strategy most suitable for modelling valproate protein-binding saturation. Baseline demographics: Zhang 2023 Table 2."
  )

  ini({
    # ----------------------------------------------------------------
    # Absorption - FIXED to the literature Ka set the authors adopted.
    # Zhang 2023 Supplementary Table S3 footnote: "Ka was fixed to
    # 2.64, 1.57, 0.46 for syrup, conventional tablet and SR tablet,
    # respectively". Oral syrup is the reference formulation.
    # ----------------------------------------------------------------
    lka <- fixed(log(2.64)); label("Absorption rate constant, oral syrup reference (1/h)")            # Zhang 2023 Supplementary Table S3 footnote (syrup Ka = 2.64 1/h)
    e_form_tablet_ka <- fixed(log(1.57 / 2.64)); label("Log-ratio shift on Ka for conventional tablet vs syrup") # Zhang 2023 Supplementary Table S3 footnote (conventional tablet Ka = 1.57 1/h)
    e_form_vpa_sr_ka <- fixed(log(0.46 / 2.64)); label("Log-ratio shift on Ka for sustained-release tablet vs syrup") # Zhang 2023 Supplementary Table S3 footnote (SR tablet Ka = 0.46 1/h)

    # ----------------------------------------------------------------
    # Structural parameters (Zhang 2023 Supplementary Table S3, "Model
    # IV: Linear non-saturable binding equation" column; estimate
    # (% RSE)). Both are UNBOUND-referenced apparent (oral) parameters:
    # elimination acts on Cu, and Cu = central / (V/F). Their ratio to
    # the base model's total-concentration parameters (0.311 L/h and
    # 27.8 L) recovers the unbound fraction - 7.07% on CL and 5.65% on
    # V - which agrees with the 7.8% the isotherm itself gives at the
    # steady-state unbound concentration for the cohort mean daily dose.
    # That triple agreement is what identifies the unbound reading,
    # because the paper prints Eq. 6 without saying how it couples to
    # the PK.
    # ----------------------------------------------------------------
    lcl <- log(4.4);   label("Apparent unbound plasma clearance CLp/F (L/h)")    # Zhang 2023 Supplementary Table S3 Model IV (CLp/F = 4.4 L/h, 4.0% RSE)
    lvc <- log(492.0); label("Apparent unbound volume of distribution V/F (L)")  # Zhang 2023 Supplementary Table S3 Model IV (V/F = 492.0 L, 5.8% RSE)

    # ----------------------------------------------------------------
    # Linear non-saturable binding constants (Zhang 2023 Eq. 6), FIXED
    # to the values the authors quote from Gu 2021 and never estimated.
    # NS is the slope of the non-saturable component - the term that
    # distinguishes this model from the pure Langmuir Model II and, per
    # the authors' Discussion, the reason it outperforms it.
    # ----------------------------------------------------------------
    lbmax_pb <- fixed(log(67.3)); label("Maximum saturable binding-site concentration Bm (mg/L)")  # Zhang 2023 Eq. 6 text (Bm = 67.3 mg/L, fixed from Gu 2021)
    lkd_pb   <- fixed(log(2.12)); label("Dissociation constant Kd for the saturable site (mg/L)")  # Zhang 2023 Eq. 6 text (Kd = 2.12, fixed from Gu 2021)
    lns_pb   <- fixed(log(2.25)); label("Slope of the linear non-saturable binding component (unitless)")  # Zhang 2023 Eq. 6 text (NS = 2.25, fixed from Gu 2021)

    # ----------------------------------------------------------------
    # IIV - diagonal (no correlation reported). Zhang 2023 reports
    # between-subject variability as CV%, so the internal variance is
    # omega^2 = log(CV^2 + 1).
    #   BSV_CL/F = 55.5% -> omega^2 = log(0.555^2 + 1) = 0.26852
    #   BSV_V/F  = 77.7% -> omega^2 = log(0.777^2 + 1) = 0.47233
    # ----------------------------------------------------------------
    etalcl ~ 0.26852  # Zhang 2023 Supplementary Table S3 Model IV (BSV_CL/F = 55.5%, 10.7% RSE)
    etalvc ~ 0.47233  # Zhang 2023 Supplementary Table S3 Model IV (BSV_V/F = 77.7%, 8.9% RSE)

    # ----------------------------------------------------------------
    # Residual error - combined proportional + additive, on the TOTAL
    # plasma concentration (the measured quantity).
    # ----------------------------------------------------------------
    propSd <- 0.139; label("Proportional residual SD (fraction; 13.9%)")  # Zhang 2023 Supplementary Table S3 Model IV (residual error 13.9%, 11.0% RSE)
    addSd  <- 4.1;   label("Additive residual SD (mg/L)")                 # Zhang 2023 Supplementary Table S3 Model IV (residual error 4.1 mg/L, 25.6% RSE)
  })

  model({
    # 1. Formulation-specific absorption rate constant. Syrup is the
    #    reference (both indicators 0).
    ka <- exp(lka + e_form_tablet_ka * FORM_TABLET + e_form_vpa_sr_ka * FORM_VPA_SR)

    # 2. Individual PK parameters (unbound-referenced; no covariates)
    cl <- exp(lcl + etalcl)
    vc <- exp(lvc + etalvc)

    # 3. Fixed binding constants
    bmax_pb <- exp(lbmax_pb)
    kd_pb   <- exp(lkd_pb)
    ns_pb   <- exp(lns_pb)

    # 4. Micro-constant. Elimination acts on the unbound concentration,
    #    CL * Cu = (CL / V) * central, so the amount kinetics stay
    #    first-order; the non-linearity sits entirely in the observation
    #    equation below.
    kel <- cl / vc

    # 5. One-compartment ODE system with first-order oral absorption
    d/dt(depot)   <- -ka * depot
    d/dt(central) <-  ka * depot - kel * central

    # 6. Observation. Cu is the unbound plasma concentration; Cb adds a
    #    linear non-saturable component to the saturable Langmuir term
    #    (Zhang 2023 Eq. 6); the measured quantity is their sum.
    Cu <- central / vc
    Cb <- bmax_pb * Cu / (kd_pb + Cu) + ns_pb * Cu
    Cc <- Cu + Cb
    Cc ~ add(addSd) + prop(propSd)
  })
}
