A Brief hands-on Introduction to PMX simulation with RxODE"
In RxODE: Facilities for Simulating from ODE-Based Models

Addendum based on current version

Some things have changed based on the new package;

Notes: - RxODE's event tables can be NONMEM input files directly - NONMEM EVIDs do not need to be translated to RxODE evids by the user; They are automatically recognized - Note: as.et is not required and you can directly solve from theo_sd; It supports NONMEM-style events and data frames. Also as.et drops covariate information - Note: This lecture uses R to simulate varaibility, you can also simulate by specifying the omega and sigma matrices as well. - as.et will drop any co-variate information. - RxODE can be used with PopED (See https://andrewhooker.github.io/PopED/articles/model_def_other_pkgs.html#ode-solution-using-rxode-1)

Day 1 Exercises

library(RxODE)
m1 <- 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
   eff(0) = Kin/Kout
})

# system parameters
theta <- 
   c(KA=2.94E-01, 
     CL=1.86E+01, V2=4.02E+01,    # central 
     Q=1.05E+01, V3=2.97E+02,     # peripheral
     Kin=1, Kout=1, EC50=200)     # effects

# simulation & plot
simu_plot = function()
{
    x <- solve(m1, theta, ev)
    plot(x, C2)
}


#---------------------------
# single dose
ev <- eventTable()
ev$add.dosing(dose=10000, nbr.doses=1)
ev$add.sampling(0:24)
simu_plot()

# multiple doses
ev <- eventTable()
ev$add.dosing(dose=10000, nbr.doses=5, dosing.interval = 24)
ev$add.sampling(0:120)
simu_plot()

# bid for 5 days, followedby qd for 5 days
ev <- eventTable()
ev$add.dosing(dose=10000, nbr.doses=10, dosing.interval=12)
ev$add.dosing(dose=20000, nbr.doses=5, start.time=120, dosing.interval = 24)
ev$add.sampling(0:240)
simu_plot()

# infusion for 5 days, followedby oral for 5 days
ev <- eventTable()
ev$add.dosing(
    dose=10000, 
    nbr.doses=5, 
    dosing.to=2,
    rate=5000,
    dosing.interval = 24
)
ev$add.dosing(dose=10000, nbr.doses=5, start.time=120, dosing.interval = 24)
ev$add.sampling(0:240)
simu_plot()


# 2wk-on, 1wk-off
ncyc = 4
ev <- eventTable()
for (i in 1:ncyc)
    ev$add.dosing(dose=10000, nbr.doses=14, start.time=(i-1)*21*24, dosing.interval = 24)
ev$add.sampling(0:(ncyc*21*24))
simu_plot()

library(RxODE)
m1 <- 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
   eff(0) = Kin/Kout
})

# system parameters
nsub=100
theta.all <- 
   cbind(KA=2.94E-01, CL=1.86E+01*exp(rnorm(nsub,0,.1)),
     V2=4.02E+01, Q=1.05E+01, V3=2.97E+02,
     Kin=1, Kout=1, EC50=200)

# dosing & sampling
ev <- eventTable()
ev$add.dosing(dose=10000, nbr.doses=5, dosing.interval=24)
ev$add.sampling(0:120)

# simu
x <- solve(m1, theta.all, ev)
plot(x, C2)

#plot quantiles
y = aggregate(C2~time, x, quantile, prob=c(.05, .5, .95))
y = as.data.frame(as.matrix(y))
names(y)[2:4] = c("C.05", "C.5", "C.95")
require(ggplot2)
ggplot(y, aes(time, C.5)) +
  geom_line(col="red") +
  geom_line(aes(time, C.05)) +
  geom_line(aes(time, C.95))

Day 1 Lecture

Day 2 exercises

Weight based dosing using f(depot)

nsub=100
theta.all <- 
   cbind(KA=2.94E-01, CL=1.86E+01*exp(rnorm(nsub,0,.3)),
     V2=4.02E+01, Q=1.05E+01, V3=2.97E+02,
     Kin=1, Kout=1, EC50=200, 
     WT=rnorm(nsub, 70, 10)
)

x = list()
for (i in 1:nsub) {
    theta = theta.all[i,]
    ev <- eventTable()
    ev$add.dosing(dose=1000*theta["WT"], nbr.doses=5, dosing.interval=24)
    ev$add.sampling(0:120)
    x[[i]] = cbind(id=i, solve(m1, theta, ev))
}

y = do.call("rbind", x)
require(ggplot2)
ggplot(y, aes(time, C2, group=id)) + geom_line()


#------------------------------
m2 <- 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
   eff(0) = Kin/Kout
   f(depot) = WT
})

i=1; theta = theta.all[i,]

ev <- eventTable()
ev$add.dosing(dose=1000*theta["WT"], nbr.doses=5, dosing.interval=24)
ev$add.sampling(0:120)
x1 = solve(m1, theta, ev)

ev <- eventTable()
ev$add.dosing(dose=1000, nbr.doses=5, dosing.interval=24)
ev$add.sampling(0:120)
x2 = solve(m2, theta, ev)

x1$C2 - x2$C2

Simultaneous 1st order and 0 order absorption

## ex6.R
library(RxODE)
m3 <- 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
  f(depot) = 1-F2 
  f(centr) = F2 
  dur(centr) = D2
})

# system parameters
theta <- 
  c(KA=2.94E-01, 
    CL=1.86E+01, V2=4.02E+01,    # central 
    Q=1.05E+01, V3=2.97E+02,     # peripheral
    F2=.6, D2=2)

# simulation & plot
simu_plot = function()
{
  x <- solve(m3, theta, ev)
  plot(x, C2)
}


#---------------------------
# simutaneous 0-order & 1-order absorption
ev <- eventTable()
ev$add.dosing(
  dose=10000, 
  nbr.doses=5, 
  dosing.to=2,
  rate=-2,
  dosing.interval = 24
)
ev$add.dosing(
  dose=10000, 
  nbr.doses=5, 
  dosing.interval = 24
)
ev$add.sampling(0:120)
simu_plot()