How do I add a population parameter using model piping?

[john-harrold]

I know that based on patterns, population parameters can be added when a model portion is added through piping. I find this a little confusing and I would rather add the parameter explicitly. So say I want to add an end point Cp but I want to add the error parameter (myerror) to the ini portion first then add the endpoint. I would think I could do something like this but it doesn't seem to work. I've tried lots of variations on this, but I'm kind of at a loss.

library(rxode2)
library(nlmixr2)

PK_2cmt_des <- function() {
  description <- "Two compartment PK model with linear clearance using differential equations"
  ini({
    lka <- 0.45 ; label("Absorption rate (Ka)")
    lcl <- 1 ; label("Clearance (CL)")
    lvc  <- 3 ; label("Central volume of distribution (V)")
    lvp  <- 5 ; label("Peripheral volume of distribution (Vp)")
    lq  <- 0.1 ; label("Intercompartmental clearance (Q)")
    propSd <- 0.5 ; label("Proportional residual error (fraction)")
  })
  model({
    ka <- exp(lka)
    cl <- exp(lcl)
    vc <- exp(lvc)
    vp <- exp(lvp)
    q  <- exp(lq)

    kel <- cl/vc
    k12 <- q/vc
    k21 <- q/vp

    d/dt(depot) <- -ka*depot
    d/dt(central) <-  ka*depot - kel*central - k12*central + k21*peripheral1
    d/dt(peripheral1) <- k12*central - k21*peripheral1
    Cc <- central / vc

    Cc ~ prop(propSd)
  })
}

PK_2cmt_des |>
  ini({myerror <- .1}) |>
  model({Cp     <- peripheral1/vp
         Cp ~ add(myerror)}) 

[billdenney]

You have to add it to the model first, and since you're not modifying a line, you need to add an append argument to model():

library(rxode2)
#> rxode2 2.1.3.9000 using 8 threads (see ?getRxThreads)
#>   no cache: create with `rxCreateCache()`
library(nlmixr2)
#> Loading required package: nlmixr2data

PK_2cmt_des <- function() {
  description <- "Two compartment PK model with linear clearance using differential equations"
  ini({
    lka <- 0.45 ; label("Absorption rate (Ka)")
    lcl <- 1 ; label("Clearance (CL)")
    lvc  <- 3 ; label("Central volume of distribution (V)")
    lvp  <- 5 ; label("Peripheral volume of distribution (Vp)")
    lq  <- 0.1 ; label("Intercompartmental clearance (Q)")
    propSd <- 0.5 ; label("Proportional residual error (fraction)")
  })
  model({
    ka <- exp(lka)
    cl <- exp(lcl)
    vc <- exp(lvc)
    vp <- exp(lvp)
    q  <- exp(lq)
    
    kel <- cl/vc
    k12 <- q/vc
    k21 <- q/vp
    
    d/dt(depot) <- -ka*depot
    d/dt(central) <-  ka*depot - kel*central - k12*central + k21*peripheral1
    d/dt(peripheral1) <- k12*central - k21*peripheral1
    Cc <- central / vc
    
    Cc ~ prop(propSd)
  })
}

PK_2cmt_des |>
  model({
    Cp <- peripheral1/vp
    Cp ~ add(myerror)
  },
  append = TRUE
  ) |>
  ini(myerror <- .1)
#> ℹ add residual parameter `myerror` and set estimate to 1
#> ℹ change initial estimate of `myerror` to `0.1`
#>  ── rxode2-based free-form 3-cmt ODE model ────────────────────────────────────── 
#>  ── Initalization: ──  
#> Fixed Effects ($theta): 
#>     lka     lcl     lvc     lvp      lq  propSd myerror 
#>    0.45    1.00    3.00    5.00    0.10    0.50    0.10 
#> 
#> States ($state or $stateDf): 
#>   Compartment Number Compartment Name
#> 1                  1            depot
#> 2                  2          central
#> 3                  3      peripheral1
#>  ── Multiple Endpoint Model ($multipleEndpoint): ──  
#>   variable               cmt               dvid*
#> 1   Cc ~ … cmt='Cc' or cmt=4 dvid='Cc' or dvid=1
#> 2   Cp ~ … cmt='Cp' or cmt=5 dvid='Cp' or dvid=2
#>   * If dvids are outside this range, all dvids are re-numered sequentially, ie 1,7, 10 becomes 1,2,3 etc
#> 
#>  ── Model (Normalized Syntax): ── 
#> function() {
#>     description <- "Two compartment PK model with linear clearance using differential equations"
#>     ini({
#>         lka <- 0.45
#>         label("Absorption rate (Ka)")
#>         lcl <- 1
#>         label("Clearance (CL)")
#>         lvc <- 3
#>         label("Central volume of distribution (V)")
#>         lvp <- 5
#>         label("Peripheral volume of distribution (Vp)")
#>         lq <- 0.1
#>         label("Intercompartmental clearance (Q)")
#>         propSd <- c(0, 0.5)
#>         label("Proportional residual error (fraction)")
#>         myerror <- c(0, 0.1)
#>     })
#>     model({
#>         ka <- exp(lka)
#>         cl <- exp(lcl)
#>         vc <- exp(lvc)
#>         vp <- exp(lvp)
#>         q <- exp(lq)
#>         kel <- cl/vc
#>         k12 <- q/vc
#>         k21 <- q/vp
#>         d/dt(depot) <- -ka * depot
#>         d/dt(central) <- ka * depot - kel * central - k12 * central + 
#>             k21 * peripheral1
#>         d/dt(peripheral1) <- k12 * central - k21 * peripheral1
#>         Cc <- central/vc
#>         Cc ~ prop(propSd)
#>         Cp <- peripheral1/vp
#>         Cp ~ add(myerror)
#>     })
#> }

^{Created on 2024-08-05 with reprex v2.1.1}

[john-harrold]

Thanks. This was breaking my brain.