nlmixrdevelopment / nlmixr

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

1.1.0-3 Installer version: nlmixr invents MDV=0 when EVID=1 and stops execution #152

Closed nholford closed 5 years ago

nholford commented 5 years ago

pk.fit.ka1t <- nlmixr(one.compartment.saem, data.pk, est="saem", control=saemControl(print=100000))

-- Fit Data (object pk.fit.ka1t is a modified tibble): ------------------------- Compiling VPC model...done done (2.67 sec) Error in RxODE::etTrans(data, paste(paste(.tmp, collapse = "\n"), "\n", : MDV cannot be 0 when EVID=1 Calls: source ... nlmixr.function -> do.call -> -> Execution halted

data.pk does not have an MDV column:

image

so it seems nlmixr has invented an MDV value of 0 in a data record with EVID=1.

mattfidler commented 5 years ago

If you look at the low level functions, invent a model, and use the data above I cannot reproduce this issue

library(RxODE)

df <- data.frame(ID=1, TIME=c(0, 0.5, 1), AMT=c(100, 0, 0), DVID=1, EVID=c(1, 0, 0), WT=66.7, AGE=50, SEX=1, CMT=c(1, 3, 3))

print(df)
#>   ID TIME AMT DVID EVID   WT AGE SEX CMT
#> 1  1  0.0 100    1    1 66.7  50   1   1
#> 2  1  0.5   0    1    0 66.7  50   1   3
#> 3  1  1.0   0    1    0 66.7  50   1   3

mod1 <-RxODE({
    C2 = centr/V2;
    C3 = peri/V3;
    d/dt(depot) =-KA*depot;
    d/dt(centr) = KA*depot - CL*C2 - Q*C2 + Q*C3;
    d/dt(peri)  =                    Q*C2 - Q*C3;
    d/dt(eff)  = Kin - Kout*(1-C2/(EC50+C2))*eff;
});

tmp <- RxODE::etTrans(df, mod1)

Created on 2019-04-11 by the reprex package (v0.2.1)

Can you provide the data for me to test.

mattfidler commented 5 years ago

Thank you for providing the data and script.

Here is the data causing the issue:

> head(data.pd)
   ID TIME AMT DV DVID EVID   WT AGE SEX      IKTR       IKA       ICL       IV CMT MDV
1   1    0 100  0    1    1 66.7  50   1 0.4897436 0.8194047 0.2762888 8.003758   1   0
10  1   24   0 44    2    0 66.7  50   1 0.4897436 0.8194047 0.2762888 8.003758   4   0
12  1   36   0 27    2    0 66.7  50   1 0.4897436 0.8194047 0.2762888 8.003758   4   0
14  1   48   0 28    2    0 66.7  50   1 0.4897436 0.8194047 0.2762888 8.003758   4   0
16  1   72   0 31    2    0 66.7  50   1 0.4897436 0.8194047 0.2762888 8.003758   4   0
17  1   96   0 60    2    0 66.7  50   1 0.4897436 0.8194047 0.2762888 8.003758   4   0

Note the first line MDV=0 and EVID=1, which violates the idea of a dose. Is this acceptable in NONMEM?

mattfidler commented 5 years ago

Also this MDV column is created by:

data.pd$MDV <- ifelse(is.na(data.pd$DV), 1, 0)
data.pd <- data.pd[data.pd$MDV==0,]
data.pd <- data.pd[data.pd$AMT!=0 | data.pd$DVID==2, ]

It really should be:


data.pd$MDV <- ifelse(is.na(data.pd$DV) | data.pd$DV == 0, 1, 0)
data.pd <- data.pd[data.pd$MDV==0,]
data.pd <- data.pd[data.pd$AMT!=0 | data.pd$DVID==2, ]
nholford commented 5 years ago

Thanks Matt for correctly identifying the data set triggering the error message was data.pd and not data.pk. I had wrongly identified the data set by looking at the R command line output when the output only showed standout without echoing commands.

I appreciate your suggestion to try to correct the code generating data.pd but this line data.pd <- data.pd[data.pd$MDV==0,] throws the baby out with the bathwater because the dose records are no longer included. The original code had worked with previous versions of nlmixr which did not use MDV in the way it is interpreted by NONMEM and now recognized by nlmixr.

I've changed the code as follows -- data.pd is then correctly structured with MDV consistent with EVID and is used without error with nlmixr.

#Create MDV=1 for dose record or when DV is recorded as zero
data.pd$MDV <- ifelse(data.pd$AMT>0 | data.pd$DV==0, 1, 0)
#select dosing records and PCA effect records
data.pd <- data.pd[data.pd$AMT!=0 | data.pd$DVID==2, ]
mattfidler commented 5 years ago

I'm glad it works!