Time varying F with ALAG with mrgsolve > 1.0.5

[FelicienLL]

Hi Kyle,
I detected an other issue with time-varying bioavailability in presence of lag time. See the reprex below, that I derived from a more complex model that used to run in my validation settings in mapbayr, but now it fails. Here three patients simulated, ID1 is with the bioavailability on the administration record, ID2 is with the bioavailability "lagged" of one row so that the value taken into account is the one interpolated with locf (the way I see it is how NONMEM works) ; ID3 is like ID2, but the the value of F1 at TIME = 0 is filled up with 0.1. The trick here is not to ask for a simulation record between the two doses.

My guess is that the difference comes from the fact that the behavior of $MAIN calls in presence of lag-times was changed in mrgsolve 0.5.0.

Setup

cod <- '$PARAM F1 = 1, ALAG1 = 0.5
$CMT CENT
$MAIN
F_CENT = F1 ; 
ALAG_CENT = ALAG1 ;
'
mod <-  mrgsolve::mcode("mod", cod)
#> Building mod ... done.
dat <- data.frame(ID = 1, time = c(0,1,2), amt = c(100,100,0), cmt = 1, evid = c(1,1,0))
dat <- dplyr::bind_rows(
  dplyr::mutate(dat, ID = 1, F1 = c(0.1, 1, 1)), 
  dplyr::mutate(dat, ID = 2, F1 = c(0, 0.1, 1)), 
  dplyr::mutate(dat, ID = 3, F1 = c(0.1, 0.1, 1))
)
print(dat)
#>   ID time amt cmt evid  F1
#> 1  1    0 100   1    1 0.1
#> 2  1    1 100   1    1 1.0
#> 3  1    2   0   1    0 1.0
#> 4  2    0 100   1    1 0.0
#> 5  2    1 100   1    1 0.1
#> 6  2    2   0   1    0 1.0
#> 7  3    0 100   1    1 0.1
#> 8  3    1 100   1    1 0.1
#> 9  3    2   0   1    0 1.0

mrgsolve 1.0.6

packageVersion("mrgsolve")
#> [1] '1.0.6'
mrgsolve::mrgsim(mod, dat, end = -1, output = "df", carry_out = "F1")
#>   ID time  F1 CENT
#> 1  1    0 0.1    0
#> 2  1    1 1.0   10
#> 3  1    2 1.0  110
#> 4  2    0 0.0    0
#> 5  2    1 0.1    0
#> 6  2    2 1.0   10
#> 7  3    0 0.1    0
#> 8  3    1 0.1   10
#> 9  3    2 1.0   20

mrgsolve 1.0.4

packageVersion("mrgsolve")
#> [1] '1.0.4'
mrgsolve::mrgsim(mod, dat , end = -1, output = "df", carry_out = "F1")
#>   ID time  F1 CENT
#> 1  1    0 0.1    0
#> 2  1    1 1.0  100
#> 3  1    2 1.0  200
#> 4  2    0 0.0    0
#> 5  2    1 0.1    0
#> 6  2    2 1.0  100
#> 7  3    0 0.1    0
#> 8  3    1 0.1   10
#> 9  3    2 1.0  110

NONMEM 7.4.4

nonmem code

$PROBLEM mrgsolvelag
$INPUT ID TIME AMT CMT BIOAV EVID DV
$DATA data.csv IGNORE=@
$SUBROUTINES ADVAN13 TOL=8
$MODEL
COMP=(CENT) ; 1
$PK
BLA = ETA(1)
F1 = BIOAV
ALAG1 = 0.5
S1 = 1
$OMEGA 1
$SIGMA 1
$DES
DADT(1) = 0
$ERROR
BLU = EPS(1)
Y = F
CONC = Y 
$SIM (123) ONLYSIM
$TABLE ID TIME BIOAV CONC NOPRINT NOAPPEND FILE=tab001

nmtab
#>   ID TIME  F1 CENT
#> 1  1    0 0.1    0
#> 2  1   24 1.0  100
#> 3  1   25 1.0  200
#> 4  2    0 0.0    0
#> 5  2   24 0.1   10
#> 6  2   25 1.0  110
#> 7  3    0 0.1    0
#> 8  3   24 0.1   10
#> 9  3   25 1.0  110

I'm struggling a bit to switch from a package version to another, so I hope I did not fail my reprex and you will be able to reproduce it.

Best regards

Félicien

[kylebaron]

Thanks for the reprex, @FelicienLL. Lag time is really tricky and this is really helpful to see what the behavior should be. I'll try to get this addressed in the next release.

[kylebaron]

Establishing that the outputs are the same between mrgsolve and NONMEM when there is no lag time.

$PROBLEM mrgsolvelag
$INPUT ID TIME AMT CMT EVID BIOAV DV
$DATA data.csv IGNORE=@
$SUBROUTINES ADVAN13 TOL=8
$MODEL
COMP=(CENT) ; 1
$PK
BLA = ETA(1)
F1 = BIOAV
ALAG1 = 0
S1 = 1
$OMEGA 1
$SIGMA 1
$DES
DADT(1) = 0
$ERROR
BLU = EPS(1)
Y = F
CONC = Y 
$SIM (123) ONLYSIM
$TABLE ID TIME BIOAV CONC NOPRINT NOAPPEND FILE=tab001

ID,time,amt,cmt,evid,F1,DV
1,0,100,1,1,0.1,.
1,1,100,1,1,1,.
1,2,0,1,0,1,.
2,0,100,1,1,0,.
2,1,100,1,1,0.1,.
2,2,0,1,0,1,.
3,0,100,1,1,0.1,.
3,1,100,1,1,0.1,.
3,2,0,1,0,1,.

library(mrgsolve)
#> 
#> Attaching package: 'mrgsolve'
#> The following object is masked from 'package:stats':
#> 
#>     filter
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union



cod <- '$PARAM F1 = 1, ALAG1 = 0
$CMT CENT
$MAIN
F_CENT = F1 ; 
ALAG_CENT = ALAG1 ;
'
mod <-  mrgsolve::mcode("mod", cod)
#> Building mod ...
#> done.
#> Building mod ... done.
dat <- data.frame(ID = 1, time = c(0,1,2), amt = c(100,100,0), cmt = 1, evid = c(1,1,0))
dat <- dplyr::bind_rows(
  dplyr::mutate(dat, ID = 1, F1 = c(0.1, 1, 1)), 
  dplyr::mutate(dat, ID = 2, F1 = c(0, 0.1, 1)), 
  dplyr::mutate(dat, ID = 3, F1 = c(0.1, 0.1, 1)), 
) 

mrgsim(mod, dat, carry_out = "F1") %>% as.data.frame()
#>   ID time  F1 CENT
#> 1  1    0 0.1   10
#> 2  1    1 1.0  110
#> 3  1    2 1.0  110
#> 4  2    0 0.0    0
#> 5  2    1 0.1   10
#> 6  2    2 1.0   10
#> 7  3    0 0.1   10
#> 8  3    1 0.1   20
#> 9  3    2 1.0   20


# > fread("tab001")
# ID TIME BIOAV CONC
# 1:  1    0   0.1   10
# 2:  1    1   1.0  110
# 3:  1    2   1.0  110
# 4:  2    0   0.0    0
# 5:  2    1   0.1   10
# 6:  2    2   1.0   10
# 7:  3    0   0.1   10
# 8:  3    1   0.1   20
# 9:  3    2   1.0   20

^{Created on 2023-02-02 with reprex v2.0.2}

[kylebaron]

mrgsolve output with and without lagtime

Still processing this but I want to say that the output with the lag time is sort of what I'd expect based on the output with no lagtime. Sort of like the answer is the same with and without lag time, but the result is bumped forward in time.

library(mrgsolve)
#> 
#> Attaching package: 'mrgsolve'
#> The following object is masked from 'package:stats':
#> 
#>     filter
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union



cod <- '$PARAM F1 = 1, ALAG1 = 0.5
$CMT CENT
$MAIN
F_CENT = F1 ; 
ALAG_CENT = ALAG1 ;
'
mod <-  mrgsolve::mcode("mod", cod)
#> Building mod ...
#> done.
#> Building mod ... done.
dat <- data.frame(ID = 1, time = c(0,1,2), amt = c(100,100,0), cmt = 1, evid = c(1,1,0))
dat <- dplyr::bind_rows(
  dplyr::mutate(dat, ID = 1, F1 = c(0.1, 1, 1)), 
  dplyr::mutate(dat, ID = 2, F1 = c(0, 0.1, 1)), 
  dplyr::mutate(dat, ID = 3, F1 = c(0.1, 0.1, 1)), 
) 

mrgsim(mod, dat, carry_out = "F1") %>% as.data.frame()
#>   ID time  F1 CENT
#> 1  1    0 0.1    0
#> 2  1    1 1.0   10
#> 3  1    2 1.0  110
#> 4  2    0 0.0    0
#> 5  2    1 0.1    0
#> 6  2    2 1.0   10
#> 7  3    0 0.1    0
#> 8  3    1 0.1   10
#> 9  3    2 1.0   20

cod2 <- '$PARAM F1 = 1, ALAG1 = 0
$CMT CENT
$MAIN
F_CENT = F1 ; 
ALAG_CENT = ALAG1 ;
'
mod2 <-  mrgsolve::mcode("mod2", cod2)
#> Building mod2 ... done.
mrgsim(mod2, dat, carry_out = "F1") %>% as.data.frame()
#>   ID time  F1 CENT
#> 1  1    0 0.1   10
#> 2  1    1 1.0  110
#> 3  1    2 1.0  110
#> 4  2    0 0.0    0
#> 5  2    1 0.1   10
#> 6  2    2 1.0   10
#> 7  3    0 0.1   10
#> 8  3    1 0.1   20
#> 9  3    2 1.0   20

^{Created on 2023-02-02 with reprex v2.0.2}

[kylebaron]

From the post here:

When we look at the NONMEM run with no lag time, we'd say that, when the data set sets the bioavailability for a dose, we use the BIOAV on the same record as the dose. That was this run here:

# > fread("tab001")
# ID TIME BIOAV CONC
# 1:  1    0   0.1   10
# 2:  1    1   1.0  110
# 3:  1    2   1.0  110
# 4:  2    0   0.0    0
# 5:  2    1   0.1   10
# 6:  2    2   1.0   10
# 7:  3    0   0.1   10
# 8:  3    1   0.1   20
# 9:  3    2   1.0   20

At time zero with no lag time, BIOAV was 0.1 and the dose was 100 and we got 10 units in the the compartment.

I'm not sure why this pattern should change when there is a lag time: with the lagtime, we don't see the dose right at time 0, but we know what happened with that dose when we look at time 1. NONMEM is already looking ahead at the next record in the data set for the BIOAV for the next dose. The way mrgsolve is handling lagged doses is they behave as regular doses, but they don't actually take effect until the apparent dose time + whatever the lag time is at the apparent dose time.