Individual model to simualte other events

[Martin-Umpierrez]

Dear Cyril,

First of all, thank you for this amazing package.

I couldn't find a place to post Q&A, so I am reaching out here.

After using the function poso_estim_map(), an object with ETAs, the model, and the event is generated. We can then perform optimization for that ID using poso_dose_conc(), which I assume utilizes the Individual Model with the estimated ETAs. My question is whether there is a possibility to use or extract the Individual Model from the generated object to simulate not only the optimized dose but also other scenarios. I tried some manipulation but couldn't achieve this.

Maybe you can help me!

Thanks again.

[levenc]

Dear Martin,
Thank you for the kind words.
The model you want to use is directly accessible from the output of poso_estim_map() as you assumed. You can store the result of the estimation in an object and access the model:

myEstimate <- poso_estim_map(patient_data,mod_run001)

The full output is a list

myEstimate

$eta
    ETA_Cl     ETA_Vc     ETA_Ka 
 0.6019034 -0.4291750  0.1278449 

$model
── Solved rxode2 object ──
── Parameters ($params): ──
  THETA_Cl   THETA_Vc   THETA_Ka     ETA_Cl     ETA_Vc     ETA_Ka 
 4.0000000 70.0000000  1.0000000  0.6019034 -0.4291750  0.1278449 
── Initial Conditions ($inits): ──
depot centr   AUC 
    0     0     0 
── First part of data (object): ──
# A tibble: 151 × 12
   time  TVCl  TVVc  TVKa    Cl    Vc    Ka   K20     Cc depot centr    AUC
  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl>  <dbl>
1   0       4    70     1  7.30  45.6  1.14 0.160  0        0    0   0     
2   0.1     4    70     1  7.30  45.6  1.14 0.160  0.478  378.  21.8 0.0161
3   0.2     4    70     1  7.30  45.6  1.14 0.160  1.83   716.  83.5 0.125 
4   0.3     4    70     1  7.30  45.6  1.14 0.160  3.95  1017. 180.  0.408 
5   0.4     4    70     1  7.30  45.6  1.14 0.160  6.75  1286. 307.  0.938 
6   0.5     4    70     1  7.30  45.6  1.14 0.160 10.1   1526. 461.  1.78  
# ℹ 145 more rows
# ℹ Use `print(n = ...)` to see more rows

$event
        id  time   amt  evid   dur
     <int> <num> <num> <int> <num>
  1:     1   0.0    NA     0    NA
  2:     1   0.0  2000     1   0.5
  3:     1   0.1    NA     0    NA
  4:     1   0.2    NA     0    NA
  5:     1   0.3    NA     0    NA
 ---                              
148:     1  14.6    NA     0    NA
149:     1  14.7    NA     0    NA
150:     1  14.8    NA     0    NA
151:     1  14.9    NA     0    NA
152:     1  15.0    NA     0    NA

$model is an element of this list

myEstimate$model

── Solved rxode2 object ──
── Parameters ($params): ──
  THETA_Cl   THETA_Vc   THETA_Ka     ETA_Cl     ETA_Vc     ETA_Ka 
 4.0000000 70.0000000  1.0000000  0.6019034 -0.4291750  0.1278449 
── Initial Conditions ($inits): ──
depot centr   AUC 
    0     0     0 
── First part of data (object): ──
# A tibble: 151 × 12
   time  TVCl  TVVc  TVKa    Cl    Vc    Ka   K20     Cc depot centr    AUC
  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl>  <dbl>
1   0       4    70     1  7.30  45.6  1.14 0.160  0        0    0   0     
2   0.1     4    70     1  7.30  45.6  1.14 0.160  0.478  378.  21.8 0.0161
3   0.2     4    70     1  7.30  45.6  1.14 0.160  1.83   716.  83.5 0.125 
4   0.3     4    70     1  7.30  45.6  1.14 0.160  3.95  1017. 180.  0.408 
5   0.4     4    70     1  7.30  45.6  1.14 0.160  6.75  1286. 307.  0.938 
6   0.5     4    70     1  7.30  45.6  1.14 0.160 10.1   1526. 461.  1.78  
# ℹ 145 more rows
# ℹ Use `print(n = ...)` to see more rows

It's a model built by rxode2, with the values of the individual parameters. The flexibility of rxode2's event tables means that almost any scenario can be simulated. See for reference: https://nlmixr2.github.io/rxode2/articles/rxode2-event-table.html

Please don't hesitate to get in touch with me if you have a problem with a particular case.

[Martin-Umpierrez]

Dear Cyril :
I was trying some coding prior to your answer and it actually worked well, but i ll see if it is easier the way you are showing me here

What I did was the following:

TacrolimusIOV<- list(
  ppk_model   = rxode2::rxode({
    depot(0) = 0;
    centr(0) = 0;
    #### TV VALUES###
    #An RxODE model specification (this line is a comment).
    if (CYP3A5==1) { 
      beta_polimorf = 0.91;
    } else if (CYP3A5==2) {
      beta_polimorf = 0.63;
    } else {
      beta_polimorf = 0;
    }

    TVKa  = THETA_Ka ;
    TVCl  = THETA_Cl *exp(beta_polimorf) * (HCT/33.5)^(-1.0) * (LBW/60)^(0.75);
    TVVc  = THETA_Vc                                         * (LBW/60)       ;
    TVVp  = THETA_Vp                                         * (LBW/60)       ;
    TVQ   = THETA_Q                                          * (LBW/60)^(0.75);
    
    Ka    = TVKa *exp(ETA_Ka + KAPPA_Ka);
    Cl    = TVCl *exp(ETA_Cl + KAPPA_Cl);
    Vc    = TVVc *exp(ETA_Vc)           ;
    Q     = TVQ *exp(ETA_Q)             ;
    Vp    = TVVp *exp(ETA_Vp)           ;
    ke    = Cl/Vc;
    k12   = Q/Vc;
    k21   = Q/Vp;
    Cc    = centr/Vc;
    d/dt(depot)  = - Ka*depot                         ;
    d/dt(centr)  = + Ka*depot- ke*centr - k12*centr + k21*periph;
    d/dt(periph) =                      + k12*centr - k21*periph;
    d/dt(AUC)    =   Cc;
  }),
  error_model = function(f,sigma){
    g <- sigma[1] + sigma[2]*f
    return(g)
  },
  theta = c(THETA_Ka= 2.49, THETA_Cl=19.84, THETA_Vc=330.4, THETA_Vp=118.2,THETA_Q=35.44),
  omega = lotri::lotri({ETA_Ka + ETA_Cl + ETA_Vc + ETA_Q + ETA_Vp ~
      c(0.820,
        0     ,   0.360,
        0     ,       0,   1.190,
        0     ,       0,       0,    0.580,
        0     ,       0,       0,        0,   0)}),
  pi_matrix = lotri::lotri({KAPPA_Cl + KAPPA_Ka~
      c(0.260,
        0     ,   0.860)}),
  covariates  = c("CYP3A5","HCT","LBW"),
  sigma       = c(additive_a = 0, proportional_b = 0.227))


tdm_tacro <- data.frame(ID=1,
                        TIME=c(0,12.0,24.0,36.0,48.0, 60.0,72),
                      DV=c(NA,NA,NA,NA,NA,NA,8.0),
                      AMT=c(3000,3000,3000,3000,3000,3000,0),
                      OCC= 1 , 
                      EVID=c(1,1,1,1,1,1,0),
                      CYP3A5=3,LBW=60,HCT=40)

## MAP Estimation###
patA_map <-  poso_estim_map(dat=tdm_tacro,
                            prior_model=TacrolimusIOV)
covar <- utils::tail(tdm_tacro[,TacrolimusIOV$covariates],1)

### here extended_et comes from some time values that i need###
dfIOV= data.frame(KAPPA_Cl= extended_et$KAPPA_Cl[1],
                  KAPPA_Ka= extended_et$KAPPA_Ka[1])
indiv_param <- cbind(patA_map$model$params,
                     covar,
                     dfIOV)

new_scenario<- as.et(extended_et) %>%  add.dosing(dose=1000, 
                                                  start.time = from,
                                                  nbr.doses= 12, 
                                                  dosing.interval=12)

runscenario_ppk_model <- rxode2::rxSolve(object=TacrolimusIOV$ppk_model,
                                         params=indiv_param,
                                         event= new_scenario,
                                         nDisplayProgress=200)

Thx for your answer again!
Regards, Martin

[levenc]

Hi Martin,
IOV and covariates do complicate things a bit, that's true. Your approach seems right to me. If you want to retrieve the latest KAPPA and covariate values, you can also get them directly from the model produced after MAP-BE:

> tail(patA_map$model[,c(‘KAPPA_Cl’, ‘KAPPA_Ka’, ‘CYP3A5’, ‘HCT’, ‘LBW’)],1)
     KAPPA_Cl KAPPA_Ka CYP3A5 HCT LBW
731 0.1181386 0.006395748 3 40 60

Please note: in posologyr dose selection functions, IOV values are not used (KAPPA is set to 0 for all parameters). See #32

[levenc]

Hi Martin,

There is now an easier way to simulate other scenarios. The new function poso_replace_et() [1] enables updating a model with events from a new rxode2 event table, while accounting for and interpolating any covariates or inter-occasion variability.

It is now possible to do this:

# convert your dataset to an event table
et_tacro <- rxode2::as.et(tdm_tacro)

# change it to your needs
new_scenario<- et_tacro$add.dosing(dose=1000, 
                                   start.time = 72,
                                   nbr.doses= 12, 
                                   dosing.interval=12)

# update the model from poso_estim_map() using this new scenario
poso_replace_et(target_model=patA_map$model,
                prior_model=TacrolimusIOV,
                event_table = new_scenario)

[1] https://levenc.github.io/posologyr/reference/poso_replace_et.html

[Martin-Umpierrez]

Dear Cyril,
Sorry for the delay, I was busy with other projects. This function sounds much more convenient than what I was using in the past. I’ll try to implement it, and if I have any issues, I’ll come back to this post.

Thank you very much!