Overview

Description

Fit several structural PK models and select the best one based on a corrected Bayesian Information Criterion for mixed effects models.

Models to compare can be defined by rate constants and/or clearances and can include or not nonlinear elimination models.

Usage

pkbuild <- function(data=NULL, project=NULL, stat=FALSE, param="clearance", new.dir=".", 
                    MM=FALSE, level=NULL, settings.stat=NULL) 

Arguments

data
a list with fields
project
a Monolix project
stat
({FALSE}, TRUE): the statistical model is also built (using buildmlx)
param
parameterization ({“clearance”}, “rate”, “both)
new.dir
name of the directory where the created files are stored (default is the current working directory)
MM
({FALSE}, TRUE): tested models include or not Michaelis Menten elimination models
level
an integer between 1 and 9 (used by setSettings)
settings.stat
list of settings used by buildmlx (only if stat=TRUE)


Examples

Reading a data file

Let us use the warfarin PK data in this example:

head(read.csv('data/warfarinPK.csv'))
##   id time   y amount   wt age sex
## 1  1    0   .    100 66.7  50   1
## 2  1   24 9.2      . 66.7  50   1
## 3  1   36 8.5      . 66.7  50   1
## 4  1   48 6.4      . 66.7  50   1
## 5  1   72 4.8      . 66.7  50   1
## 6  1   96 3.1      . 66.7  50   1

warfarin is administrated orally:

library(Rsmlx)
warfarinPK <- list(
  dataFile = "data/warfarinPK.csv",
  headerTypes = c("id", "time", "observation", "amount", "contcov", "contcov", "catcov"),
  administration = "oral")

By default, PK models parameterized with clearance(s) are fitted.

warf.pk1 <- pkbuild(data=warfarinPK,  new.dir="warfarinPK")
## 
## [1] "warfarinPK/pk_kaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0V1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaV1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood...

Comparison of these models is based on a corrected Bayesian Information Criterion for mixed effects models:

\[{\rm BIC}_{\rm cor} = -2 LL + \log(N)\times d_N + \log(n) \times d_n\] where \(N\) is the number of individuals and \(n\) the total number of observations. Then, \(d_N\) is the number of parameters associated to the variability of the individual parameters and \(d_n\) the number of parameters associated to the variability of the observations.

According to this criterion, the best PK model for this data is oral0_1cpt_TlagTk0VCl.txt.

print(warf.pk1)
## $pop.ini
##      Tlag       Tk0         V        Cl 
## 0.4975827 2.6767010 8.1342567 0.1317972 
## 
## $project
## [1] "warfarinPK/pk_TlagTk0VCl.mlxtran"
## 
## $model
## [1] "C:/ProgramData/Lixoft/MonolixSuite2019R1/factory/library/pk/oral0_1cpt_TlagTk0VCl.txt"
## 
## $data
## $data$dataFile
## [1] "data/warfarinPK.csv"
## 
## $data$headerTypes
## [1] "id"          "time"        "observation" "amount"      "contcov"    
## [6] "contcov"     "catcov"     
## 
## $data$administration
## [1] "oral"
## 
## 
## $bic
##                           model  bic.cor
## 1     oral0_1cpt_TlagTk0VCl.txt 713.0319
## 2      oral1_1cpt_TlagkaVCl.txt 716.8131
## 3  oral1_2cpt_TlagkaClV1QV2.txt 718.1846
## 4 oral0_2cpt_TlagTk0ClV1QV2.txt 720.6063
## 5         oral0_1cpt_Tk0VCl.txt 840.9993
## 6          oral1_1cpt_kaVCl.txt 919.0924
## 
## $pop.est
##   Tlag_pop    Tk0_pop      V_pop     Cl_pop omega_Tlag  omega_Tk0 
## 0.83981762 1.33397176 8.01427029 0.13202203 0.59192544 0.62607546 
##    omega_V   omega_Cl          a          b 
## 0.22381333 0.29282128 0.33080886 0.07114829

Models parameterized with rate constants can be used instead:

warf.pk2 <- pkbuild(data=warfarinPK,  new.dir="warfarinPK", param="rate")
## 
## [1] "warfarinPK/pk_kaVk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0Vk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0Vk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0Vk12k21k.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVk12k21k.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVk12k21k13k31k.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood...
print(warf.pk2$bic)
##                                 model  bic.cor
## 1       oral1_2cpt_TlagkaVkk12k21.txt 702.7145
## 2      oral0_2cpt_TlagTk0Vkk12k21.txt 706.4141
## 3            oral0_1cpt_TlagTk0Vk.txt 709.6980
## 4             oral1_1cpt_TlagkaVk.txt 715.5151
## 5 oral1_3cpt_TlagkaVkk12k21k13k31.txt 729.4074
## 6                oral0_1cpt_Tk0Vk.txt 836.8413
## 7                 oral1_1cpt_kaVk.txt 918.5904
print(warf.pk2$pop.est)
##   Tlag_pop     ka_pop      V_pop      k_pop    k12_pop    k21_pop 
## 0.84167886 1.25322089 7.35773190 0.01526672 0.00846664 0.04796776 
## omega_Tlag   omega_ka    omega_V    omega_k  omega_k12  omega_k21 
## 0.68084643 0.67159549 0.19227345 0.06875926 0.84023440 2.75131335 
##          a          b 
## 0.30642230 0.06354331

or both:

warf.pk3 <- pkbuild(data=warfarinPK,  new.dir="warfarinPK", param="both")
## 
## [1] "warfarinPK/pk_kaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0V1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaV1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_kaVk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0Vk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0Vk.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0Vk12k21k.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVk12k21k.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVk12k21k13k31k.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood...
print(warf.pk3$bic)
##                                  model  bic.cor
## 1        oral1_2cpt_TlagkaVkk12k21.txt 702.7145
## 2       oral0_2cpt_TlagTk0Vkk12k21.txt 706.4141
## 3             oral0_1cpt_TlagTk0Vk.txt 709.6980
## 4            oral0_1cpt_TlagTk0VCl.txt 713.0319
## 5              oral1_1cpt_TlagkaVk.txt 715.5151
## 6             oral1_1cpt_TlagkaVCl.txt 716.8131
## 7         oral1_2cpt_TlagkaClV1QV2.txt 718.1846
## 8        oral0_2cpt_TlagTk0ClV1QV2.txt 720.6063
## 9  oral1_3cpt_TlagkaVkk12k21k13k31.txt 729.4074
## 10                oral0_1cpt_Tk0Vk.txt 836.8413
## 11               oral0_1cpt_Tk0VCl.txt 840.9993
## 12                 oral1_1cpt_kaVk.txt 918.5904
## 13                oral1_1cpt_kaVCl.txt 919.0924


Using a Monolix project

If a Monolix project has been already created using the data file, it can be used for providing data information

warf.pk4  <- pkbuild(project="projects/warfarinPK.mlxtran", new.dir="warfarinPK")
## 
## [1] "warfarinPK/pk_kaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0V1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaV1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood...

Of course, the results are the same as those obtained using the data file

print(warf.pk4$bic)
##                           model  bic.cor
## 1     oral0_1cpt_TlagTk0VCl.txt 703.7249
## 2  oral1_2cpt_TlagkaClV1QV2.txt 705.1475
## 3      oral1_1cpt_TlagkaVCl.txt 706.1825
## 4 oral0_2cpt_TlagTk0ClV1QV2.txt 711.8775
## 5         oral0_1cpt_Tk0VCl.txt 816.8908
## 6          oral1_1cpt_kaVCl.txt 895.0166


Building the statistical model also

warf.pk5  <- pkbuild(data=warfarinPK, stat=TRUE, new.dir="warfarinPK")
## 
## [1] "warfarinPK/pk_kaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0V1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaV1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## 
## direction = "full" will be used for the covariate search
## 
## ____________________________________________
## Initialization:
## 
## Covariate model:
##      sex age wt
## Tlag   0   0  0
## Tk0    0   0  0
## V      0   0  0
## Cl     0   0  0
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## Sampling of the conditional distribution using the initial model ... 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 657.94   0.11 677.94 704.86 692.60 
## ____________________________________________
## ____________________________________________
## Iteration 1:
## 
## Covariate model:
##      sex age wt
## Tlag   0   0  0
## Tk0    0   0  0
## V      0   0  1
## Cl     0   0  1
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## Run scenario for model 1 ... 
## Estimation of the population parameters... 
## Sampling from the conditional distribution... 
## Estimation of the log-likelihood... 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 627.90   0.11 651.90 681.75 669.49 
## ____________________________________________
## Iteration 2:
## 
## Covariate model:
##      sex age wt
## Tlag   0   0  0
## Tk0    0   0  0
## V      0   0  1
## Cl     0   1  1
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## Run scenario for model 2 ... 
## Estimation of the population parameters... 
## Sampling from the conditional distribution... 
## Estimation of the log-likelihood... 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 622.97   0.04 648.97 680.29 668.03 
## ____________________________________________
## Iteration 3:
## No difference between two successive iterations
## ____________________________________________
## Final model:
## 
## Covariate model:
##      sex age wt
## Tlag   0   0  0
## Tk0    0   0  0
## V      0   0  1
## Cl     0   1  1
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 622.97   0.04 648.97 680.29 668.03 
## ____________________________________________
## total time: 22.9s
## ____________________________________________
print(warf.pk5)
## $pop.ini
##      Tlag       Tk0         V        Cl 
## 0.4975827 2.6767010 8.1342567 0.1317972 
## 
## $project
## [1] "warfarinPK/pk_TlagTk0VCl_built.mlxtran"
## 
## $model
## [1] "C:/ProgramData/Lixoft/MonolixSuite2019R1/factory/library/pk/oral0_1cpt_TlagTk0VCl.txt"
## 
## $data
## $data$dataFile
## [1] "data/warfarinPK.csv"
## 
## $data$headerTypes
## [1] "id"          "time"        "observation" "amount"      "contcov"    
## [6] "contcov"     "catcov"     
## 
## $data$administration
## [1] "oral"
## 
## 
## $bic
##                           model  bic.cor
## 1     oral0_1cpt_TlagTk0VCl.txt 713.0319
## 2      oral1_1cpt_TlagkaVCl.txt 716.8131
## 3  oral1_2cpt_TlagkaClV1QV2.txt 718.1846
## 4 oral0_2cpt_TlagTk0ClV1QV2.txt 720.6063
## 5         oral0_1cpt_Tk0VCl.txt 840.9993
## 6          oral1_1cpt_kaVCl.txt 919.0924
## 
## $pop.est
##     Tlag_pop      Tk0_pop        V_pop   beta_V_cwt       Cl_pop 
##  0.800605379  1.523056962  7.978402734  0.013493997  0.132219081 
## beta_Cl_cage  beta_Cl_cwt   omega_Tlag    omega_Tk0      omega_V 
##  0.009411299  0.009340206  0.552894204  0.537567272  0.141679307 
##     omega_Cl            a            b 
##  0.253304510  0.329102249  0.070179168 
## 
## $covariate.model
## $covariate.model$Tlag
##   sex   age    wt 
## FALSE FALSE FALSE 
## 
## $covariate.model$Tk0
##   sex   age    wt 
## FALSE FALSE FALSE 
## 
## $covariate.model$V
##   sex   age    wt 
## FALSE FALSE  TRUE 
## 
## $covariate.model$Cl
##   sex   age    wt 
## FALSE  TRUE  TRUE 
## 
## 
## $correlation.model
## list()
## 
## $error.model
##           y 
## "combined2"

Define the settings of buildmlx in a list

warf.pk6  <- pkbuild(data=warfarinPK, stat=TRUE, new.dir="warfarinPK",
                     settings.stat =list(covToTransform = c("wt"),
                                         criterion ="AIC"))
## 
## [1] "warfarinPK/pk_kaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_Tk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaVCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0VCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagTk0V1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## [1] "warfarinPK/pk_TlagkaV1V2QCl.mlxtran"
## Estimation of the population parameters...
## Estimation of the log-likelihood... 
## 
## 
## direction = "full" will be used for the covariate search
## 
## ____________________________________________
## Initialization:
## 
## Covariate model:
##      sex age wt
## Tlag   0   0  0
## Tk0    0   0  0
## V      0   0  0
## Cl     0   0  0
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## Sampling of the conditional distribution using the initial model ... 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 657.94   0.11 677.94 704.86 692.60 
## ____________________________________________
## ____________________________________________
## Iteration 1:
## 
## Covariate model:
##      sex age wt lwt
## Tlag   0   0  0   0
## Tk0    0   0  0   0
## V      0   0  0   1
## Cl     0   1  0   1
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## Run scenario for model 1 ... 
## Estimation of the population parameters... 
## Sampling from the conditional distribution... 
## Estimation of the log-likelihood... 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 621.35   0.08 647.35 678.67 666.40 
## ____________________________________________
## Iteration 2:
## No difference between two successive iterations
## ____________________________________________
## Final model:
## 
## Covariate model:
##      sex age wt lwt
## Tlag   0   0  0   0
## Tk0    0   0  0   0
## V      0   0  0   1
## Cl     0   1  0   1
## 
## Correlation model:
## list()
## 
## Residual error model:
##           y 
## "combined2" 
## 
## Estimated criteria (importanceSampling):
##   -2LL   s.e.    AIC   BICc    BIC 
## 621.35   0.08 647.35 678.67 666.40 
## ____________________________________________
## total time: 18.6s
## ____________________________________________