Problem of transit model for dual absorption.

[Lee21A]

Hi!
I am testing a dual absorption scenario consisting of transit and first-order absorption.

library(rxode2)

mod <- function() {
  ini({
    ## Table 3 from Savic 2007
    cl  <- 17.2 # (L/hr)
    vc  <- 45.1 # L
    ka  <- 0.38 # 1/hr
    mtt <- 1.37 # hr
    
    f2 <-0.5    # Fraction of 1st Order portion
    n   <- 20.1
  })
  model({
    k           <- cl/vc
    bio <- 1-f2
    #ktr         <- (n+1)/mtt
    d/dt(depot1) <- transit(n,mtt,bio)-ka*depot1
    d/dt(depot2) <- -ka*depot2
    f(depot2) <-f2
    d/dt(cen)   <- ka*depot1 + ka*depot2-k*cen
  })
}

ev1 <- et(0, 7, length.out=200) %>%
  et(amt=20, cmt='depot1', evid=7) %>%
    et(amt=20, cmt='depot2', evid=1)

ev2 <- et(0, 7, length.out=200) %>%
  et(amt=20, cmt='depot1', time=0.001, evid=7) %>%
    et(amt=20, cmt='depot2', time=0, evid=1)

case1 <- rxSolve(mod, ev1) 
case2 <- rxSolve(mod, ev2)

plot(case1, depot1, depot2, cen, ylab="Central Concentration")
plot(case2, depot1, depot2, cen, ylab="Central Concentration")

But the answer is

Case1, ev1

However, interestingly, if I artificially adjust the absorption start time from zero to 0.001(Case2, ev2), it works well.

Case2, ev2

Anyway, if I adjust the absorption start time like this, the transit model can be used for the double absorption scenario, but I am worried that this may also cause a problem during the fitting using nlmixr2.

I will wait for your valuable comments.

[mattfidler]

Hi @Lee21A still looking into this.

[mattfidler]

Hi @Lee21A

Thank you for your patience. I have tracked down the issue. With the fix in #819, this gives the following:

library(rxode2)
#> rxode2 3.0.2.9000 using 8 threads (see ?getRxThreads)
#>   no cache: create with `rxCreateCache()`

mod <- function() {
  ini({
    ## Table 3 from Savic 2007
    cl  <- 17.2 # (L/hr)
    vc  <- 45.1 # L
    ka  <- 0.38 # 1/hr
    mtt <- 1.37 # hr

    f2 <-0.5    # Fraction of 1st Order portion
    n   <- 20.1
  })
  model({
    k           <- cl/vc
    bio <- 1-f2
    #ktr         <- (n+1)/mtt
    d/dt(depot1) <- transit(n,mtt,bio)-ka*depot1
    d/dt(depot2) <- -ka*depot2
    f(depot2) <-f2
    d/dt(cen)   <- ka*depot1 + ka*depot2-k*cen
  })
}

ev1 <- et(0, 7, length.out=200) %>%
  et(amt=20, cmt='depot1', evid=7) %>%
  et(amt=20, cmt='depot2', evid=1)

case1 <- rxSolve(mod, ev1)
#> ℹ parameter labels from comments are typically ignored in non-interactive mode
#> ℹ Need to run with the source intact to parse comments
ev2 <- et(0, 7, length.out=200) %>%
  et(amt=20, cmt='depot1', time=0.001, evid=7) %>%
  et(amt=20, cmt='depot2', time=0, evid=1)

library(ggplot2)
library(patchwork)

case2 <- rxSolve(mod, ev2)
#> ℹ parameter labels from comments are typically ignored in non-interactive mode
#> ℹ Need to run with the source intact to parse comments

p1 <- plot(case1, depot1, depot2, cen, ylab="Central Concentration")  +
  ggtitle("Case 1, start at zero")

p2 <- plot(case2, depot1, depot2, cen, ylab="Central Concentration")  +
  ggtitle("Case 2, start at 0.001")


p1/p2

^{Created on 2024-12-02 with reprex v2.1.1}