nlmixrdevelopment / nlmixr

nlmixr: an R package for population PKPD modeling
https://nlmixrdevelopment.github.io/nlmixr/
GNU General Public License v2.0
114 stars 45 forks source link

Reproduce simulation from NONMEM. Error in configsaem: covariate(s) not found: A1, A2, q, centr, Vc #604

Closed AprilCCC closed 2 years ago

AprilCCC commented 2 years ago

Hi nlmixr team,

I just tried to reproduce a NONMEM popPK dose simulation from observed real data with nlmixr, and when I tried to simulate it, occurred 2 errors, and under the same dose, however, the NCA result from nlmixr model-simulated data was far from the original observed data NCA result. Can anyone help me to check what is wrong with the nlmixr model? The observed real data is in the attachment. Many thanks! observed data.xlsx

Error in .nlmixrUIModel(fun, ini, bigmodel) : There must be at least one prediction in the model({}) block. Use ~ for predictions Error in configsaem(model = model, data = dat, inits = inits, mcmc = .mcmc, : covariate(s) not found: A1, A2, q, centr, Vc

The NONMEN model code showed below: $SUBROUTINE ADVAN6 TOL=3 $MODEL NCOMP=2 COMP=(CENTRAL,DEFOBS,DEFDOSE) COMP=(PER) $PK

TVCLS = THETA(1)*(WT/70)*THETA(7) CLS = TVCLS EXP(ETA(1))

TVVC = THETA(2)*(WT/70)*THETA(8) VC = TVVC EXP(ETA(2))

TVCLP = THETA(3)*(WT/70)*THETA(7) CLP = TVCLP EXP(ETA(3))

TVVP = THETA(4)*(WT/70)*THETA(8) VP = TVVP EXP(ETA(4))

TVVMAX = THETA(5) VMAX = TVVMAX* EXP(ETA(5))

TVKM = THETA(6) KM = TVKM* EXP(ETA(6))

;S2=V2

$DES

C1 = A(1)/VC C2 = A(2)/VP

DADT(1) = -VMAXC1/(KM+C1)-(CLSC1)-(CLPC1)+(CLPC2) DADT(2) = (CLPC1)-(CLPC2)

$ERROR

R=A(1)/VC IPRED = R IRES = DV-IPRED IWRES = IRES/DV Y = R*(1+EPS(1))+EPS(2)

$THETA 0.0223 FIX ; 1CLS,L/h 3.26 FIX ; 2VC,L 0.0433 FIX ; 3CLP,L/h 2.72 FIX ; 4VP,L 295.4 FIX ; 5VMAX,ug/h 115 FIX ; 6KM,ng/ml 0.607 FIX ; 7WT on CL 0.371 FIX ; 8WT on V ;0 FIX ;WT on VMAX $OMEGA 0.093 FIX ; ETA_CLS 0 FIX ; ETA_VC 0 FIX ; ETA_CLP 0.87 FIX ; ETA_VP 0.161 FIX ; ETA_VAMX 0 FIX ; ETA_KM $SIGMA 0.284 FIX ; EPS(1) 14.4 FIX ; EPS(2) ng/ml $SIMULATION (889216570) ONLYSIM SUBPROBLEM=10

The nlmixr model code I created showed below: cmt2sad5 <- function(){ ini({ lv <- 3.26 lcl <- 0.0223 lQ <- 0.0433 # Q inter-compartmental CL (Q) lvp <- 2.72 # Vp peripheral volume (Vp) lVM <- 295.4 lKM <- 115 WTcl <- 0.607 WTv <-0.371 wt.est <- 0.0 eta.v ~ 0 eta.cl ~ 0.093 eta.Q ~ 0 eta.vp ~ 0.87 eta.VM ~ 0.161 eta.KM ~ 0 EPS1 ~ 0.284 EPS2 ~ 14.1 }) model({ cl <- exp(lcl + eta.cl + wt.est (WTcl/70)) vc <- exp(lv + eta.v + wt.est (WTv/70)) Q <- exp(lQ + eta.Q + wt.est (WTcl/70)) vp <- exp(lvp + eta.vp + wt.est (WTv/70)) VM <- exp(lVM + eta.VM) KM <- exp(lKM + eta.KM) c1 = A1 / vc c2 = A2 / vp d/dt(center) = -VM c1/(KM+c1) - (cl c1) -(q c1) + (q c2) d/dt(periph) = (q c1) + (q c2)

linCmt() ~ lnorm(logn.sd)

R=A1 / vc
IPRED = R
IRES = DV-IPRED
IWRES = IRES/DV
Y = R*(1+EPS1)+EPS2

}) } cmt2msad5 <- nlmixr(cmt2sad5)

print(cmt2msad5)

cmt2fitsad5.logn <- nlmixr(cmt2msad5, sad5, "saem", control=list(print=0), table=tableControl(cwres=TRUE, npde=TRUE)) plot(augPred(cmt2fitsad5.logn))

simulation with same dose, 5mg/kg in human with 61.5kg weight in avarage

ev <- eventTable(amount.units="mg", time.units="h") %>% et(amt=5 61.5 1000) %>% et(seq(0, 504, by=1)) %>% et(id=1:60)

sim1 <- simulate(cmt2fitsad5.logn, events=ev)

mattfidler commented 2 years ago

With nlmixr you cannot simulate from the un-fitted model. The upcoming nlmixr2 allows you to simulate.

However, regardless of the method you use, you need to specify residual error with the ~ syntax. This would be something like R ~ prop(e1.sd)+add(e2.sd). You also do not have to specify IPRED, IRES, IWRES in the nlmixr estimation. They are added automatically by nlmixr/nlmixr2

AprilCCC commented 2 years ago

Hi Matthew, thanks for the suggesion. I still have several questions regarding reproduce model from NONMEM code, Can you help me to check with this? Really appreciate it.

  1. NONMEM code: TVCLS = THETA(1)*(WT/70)*THETA(7) CLS = TVCLS EXP(ETA(1)) $THETA 0.607 FIX ; 7WT on CL Is it correct way to write in the nlmirx? cl <- exp(lcl + eta.cl + log(WT/70) * WTcl)
  2. NONMEM code:
    $DES
    C1 = A(1)/VC
    C2 = A(2)/VP
    DADT(1) = -VMAXC1/(KM+C1)-(CLSC1)-(CLPC1)+(CLPC2)
    DADT(2) = (CLPC1)-(CLPC2)

    Is it correct way to write in the nlmixr as below?

    c1 = center / vc
    c2 = periph / vp
    d/dt(center) = -VM * c1/(KM+c1) - (cl *c1) -(Q * c1) + (Q * c2)
    d/dt(periph) = (Q * c1) - (Q * c2)
  3. How to specify $SIGMA 0.284 FIX ; EPS(1). 14.4 FIX ; EPS(2) ng/ml. in the nlmixr estimation?