FelicienLL / mapbayr

Easy Maximum A Posteriori Bayesian Estimation of PK parameters in R.
19 stars 2 forks source link

Allow >2 types of DV #51

Open FelicienLL opened 3 years ago

FelicienLL commented 3 years ago

Currently, only two types of DV are admitted for postprocessing purpose, internally they are refered as PAR and MET.

What to do ?

Impacted functions:

Not-impacted functions (among others):

mhismail commented 3 years ago

Hi @FelicienLL, nice package.

For the situation when estimating more than one DV, I wanted to share how we approached this in a similar package for general parameter optimization with mrgsolve (rather than map bayes estimation) because I think it may be helpful here too. In the R function we allow users to specify the output column (only one output column allowed), and in the mrgsolve code itself, users are responsible for including a switch-like statement if there are more than one compartments/outputs. See example below, where there are both PK and PD measurements in the NONMEM formatted dataset:

mrgsolve model

code<-"
$PROB 2 cmt PK Model, Emax PD model

$PARAM
CL=10
VC = 20
VP = 20
Q=20
Emax = 60
BL = 50
EC50 = 10
gamma =1
sigma1 = 0.1
sigma2 = 0.1

$CMT X1 X2 

$ODE
dxdt_X1 = -(Q+CL)/VC*X1+Q/VP*X2;
dxdt_X2 = Q/VC*X1-Q/VP*X2;

$TABLE
capture PK = X1/VC;
capture varPK = (PK*sigma1)*(PK*sigma1);

capture PD = BL-(pow(PK,gamma)*Emax)/(pow(PK,gamma)+pow(EC50,gamma));
capture varPD = (PD*sigma2)*(PD*sigma2);

capture ipred = NAN;
capture var = NAN; 

if(self.cmt == 1) {
   ipred = PK;
   var = varPK;
}

if(self.cmt == 2) {
   ipred = PD;
   var = varPD;
}"

mod <- mcode("2cmtPK-Emax", code)

Optimization function

fit <- mod %>%
  data_set(data) %>%
  mrgoptim(output = "ipred",
           var = "var",
           prms = c("CL",
                    "VC",
                    "VP",
                    "Q",
                    "Emax",
                    "BL",
                    "EC50",
                    "gamma"),
           v_prms = c("sigma1", "sigma2"),
           method = "newuoa")
FelicienLL commented 3 years ago

Hi @mhismail,

Thanks for sharing, this is interesting. I'll have a look on your work too.

Note that mapbayr already handles optimization of multiple outcomes, also with a switch-like statement in order to change the values of DV (unique outcome) as function of the compartment set in data. My problem is more with "post-processing", where we need to capture predictions for each outcome, i.e. PK and PD in your example. This is required because methods like plot() represents predictions vs. time for each outcome with a time grid different than the one defined in the data. Currently, this is managed by forcing them to be named PAR and MET, meaning maximum 2 outcomes allowed, and it is also not very convenient if working on anything else than a "parent + 1 metabolite" model.

I need to find a way to access these variables programmatically. Probably either by forcing the naming in model code (PRED1, PRED2 ...), or by adding a "prediction outcome(s)" argument in the estimation function.

Félicien