In nlmixrdevelopment/nlmixr.examples: Nonlinear Mixed Effects Models in Population Pharmacokinetics and Pharmacodynamics

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  message = FALSE,
  warning = FALSE,
  out.width = "100%"
  )

## These options cache the models and the model simulations in R
## To run the actual models on your system, take the save options off.
options(nlmixr.save=TRUE,
        nlmixr.save.dir=system.file(package="nlmixr.examples"));

Adding Covariances between random effects

You can simply add co-variances between two random effects by adding the effects together in the model specification block, that is eta.cl+eta.v ~. After that statement, you specify the lower triangular matrix of the fit with c().

An example of this is the Phenobarbital data:

## Load Phenobarb data
library(nlmixr)

Model Specification

pheno <- function() {
  ini({
    tcl <- log(0.008) # typical value of clearance
    tv <-  log(0.6)   # typical value of volume
    ## var(eta.cl)
    eta.cl + eta.v ~ c(1, 
                       0.01, 1) ## cov(eta.cl, eta.v), var(eta.v)
                      # interindividual variability on clearance and volume
    add.err <- 0.1    # residual variability
  })
  model({
    cl <- exp(tcl + eta.cl) # individual value of clearance
    v <- exp(tv + eta.v)    # individual value of volume
    ke <- cl / v            # elimination rate constant
    d/dt(A1) = - ke * A1    # model differential equation
    cp = A1 / v             # concentration in plasma
    cp ~ add(add.err)       # define error model
  })
}

Fit with SAEM

fit <- nlmixr(pheno, pheno_sd, "saem", table=list(cwres=TRUE, npde=TRUE))

print(fit)

Basic Goodness of Fit Plots

plot(fit)

Those individual plots are not that great, it would be better to see the actual curves; You can with augPred

plot(augPred(fit))

Two types of VPCs

library(ggplot2)
p1 <- nlmixr::vpc(fit, show=list(obs_dv=TRUE));
p1 <- p1+ ylab("Concentrations")

## A prediction-corrected VPC
p2 <- nlmixr::vpc(fit, pred_corr = TRUE, show=list(obs_dv=TRUE))
p2 <- p2+ ylab("Prediction-Corrected Concentrations")

library(gridExtra)
grid.arrange(p1,p2)

nlmixrdevelopment/nlmixr.examples documentation built on Nov. 4, 2019, 10:08 p.m.