knitr::opts_chunk$set(comment='.', message=FALSE, 
                      fig.path="inst/maintenance/img/README-")

mrgsolvetk

A toolkit to be used with mrgsolve

Installation

library(devtools)
install_github("mrgsolve/mrgsolvetk")

To install branch with mrgoptim function, until merged with mrgsolve/master:

library(devtools)
install_github("mhismail/mrgsolvetk",ref="mrgoptim")

Examples

library(ggplot2)
library(dplyr)
library(mrgsolve)
library(mrgsolvetk)

theme_set(theme_bw())

mod <- mread_cache("pk1cmt",modlib())

mod <- ev(mod, amt=100) %>% Req(CP) %>% update(end = 48, delta = 0.25)

param(mod)

Sensitivity analyses

sens_unif

out <- 
  mod %>% 
  select(CL,VC,KA1) %>%
  sens_unif(.n=10, lower=0.2, upper=3)

out

sens_plot(out, CP)

We can also make a univariate version of this

mod %>% 
  select(CL,VC,KA1) %>%
  sens_unif(.n=10, lower=0.2, upper=3, univariate = TRUE) %>%
  sens_plot(CP, split = TRUE)

sens_norm

mod %>% 
  select(CL,VC) %>%
  sens_norm(.n=10, cv=30) %>%
  sens_plot(CP)

sens_seq

mod %>% sens_seq(CL = seq(2,12,2), VC = seq(30,100,10)) %>% sens_plot(CP)

sens_range

mod %>%
  select(CL,VC) %>%
  sens_range(.n = 5, .factor = 4) %>%
  sens_plot(CP, split = TRUE)

or

mod %>%
  sens_range(CL = c(0.5, 1.5), VC = c(10,40), .n = 5) %>%
  sens_plot(CP)

sens_grid

mod %>%  sens_grid(CL = seq(1,10,1), VC = seq(20,40,5)) %>% sens_plot(CP)

sens_covset

cov1 <- dmutate::covset(CL ~ runif(1,3.5), VC[0,] ~ rnorm(50,25))

cov1
out <- mod %>% sens_covset(cov1) 
out

distinct(out,ID,CL,VC)

Maximum Likelihood Parameter Optimization

mrgoptim

This example shows a simultaneous fit of PK and PD data from five dose levels.

Data structure

The data to be fit is an mrgsolve dataset. Required columns for fitting are:

data <-  read.csv("inst/maintenance/data/optim-example.csv")

head(data)

Plot the data to get an idea of the profiles to be fit. cmt 1 is plasma concentration data and cmt 2 is PD data

ggplot(data,aes(x=time,y=dv,color=as.factor(ID)))+
  geom_point()+
  geom_line()+
  facet_wrap("cmt")+
  guides(color=F)

The following model will be fit to these data:

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);"

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

Here, the plasma concentrations, response, and variances were captured in the PK, PD, varPK, and varPD outputs, respectively.

Let's check how the initial parameter values fit the data.

out <- mod%>%
  data_set(data)%>%
  obsonly()%>%
  mrgsim()%>%
  as.data.frame()

ggplot(out,aes(x=time,y=PK,color=as.factor(ID)))+
  geom_line()+
  geom_point(data=filter(data,cmt==1),aes(y=dv))+
  guides(color=F)

ggplot(out,aes(x=time,y=PD,color=as.factor(ID)))+
  geom_line()+
  geom_point(data=filter(data,cmt==2),aes(y=dv))+
  guides(color=F)

Not terrible, should be good enough for initial estimates.

Now let's use mrgoptim to optimize the parameters and return parameter values and precision. The input, output,and var arguments map the observed values to the output compartments and variances of the outputs. Since compartment 1 in the input dataset corresponds to plasma concentration data, which is the PK column in the output, they both need to be the first value in their respective vectors. Similarly, since varPK corresponds to the variance of the PK data, it is the first value in the var argument vector.

fit<- mod%>%
  data_set(data)%>%
  mrgoptim(input=c(1,2),
           output=c("PK","PD"),
           var= c("varPK","varPD"),
           prms=c("CL",
             "VC",
             "VP",
             "Q",
             "Emax",
             "BL",
             "EC50",
             "gamma"),
           v_prms=c("sigma1","sigma2"),
           method="newuoa")

The function returns a list with some information about the optimization, the final objective function value (-2LL), final parameter estimates, covariance and correlation matrices, CV percent, and output dataset.

print(fit)

Lets check how the optimized parameters fit the data.

out_fit <- mod%>%
  param(fit$par)%>%
  data_set(data)%>%
  obsonly()%>%
  mrgsim()%>%
  as.data.frame()


ggplot(out_fit,aes(x=time,y=PK,color=as.factor(ID)))+
  geom_line()+
  geom_point(data=filter(data,cmt==1),aes(y=dv))+
  guides(color=F)

ggplot(out_fit,aes(x=time,y=PD,color=as.factor(ID)))+
  geom_line()+
  geom_point(data=filter(data,cmt==2),aes(y=dv))+
  guides(color=F)


mrgsolve/mrgsolvetk documentation built on May 11, 2019, 4:19 p.m.