ISoP ISoP on Twitter | ACoP ACoP on Twitter | PAGE | WCoP | PAGANZ | PAGJa

Mixed effects modeling alternatives in R


#1

Is anybody able suggest a mixed effects modeling package in R? I’ve tried {nlmeODE}, I’m looking for something with a bit more documentation.

Any help is appreciated!
Justin


#2

Hi Justin -

I did some non-linear mixed effects modeling in R with plain old nlme + mrgsolve. The poster is here http://metrumrg.com/assets/pubs/asbmroct13.pdf but it has limited methodological details (presented to a more clinical audience). I was able to fit the model and get what I needed.

If I were to do it over today, I’d probably take the time to look at saemix (https://github.com/ecomets/saemix or https://cran.r-project.org/web/packages/saemix/index.html) first.

Best Regards,
Kyle


#3

Thanks, Kyle! I’ll give that a shot. Look forward to hearing your mrgsolve presentation in the study group, too.

Cheers
Justin


#4

Dear Justin,

As of 23 October, nlmixr is available for testing on GitHub at https://github.com/nlmixrdevelopment/nlmixr

nlmixr builds on Wenping Wang’s RxODE package for simulation of nonlinear mixed effect models using ordinary differential equations, by implementing parameter estimation algorithms like nlme, gnlmm and SAEM. nlmixr greatly expands the utility of existing packages (like nlme) by providing an efficient and versatile way to specify pharmacometric models and dosing scenarios, with rapid execution due to compilation to C.

Give it a try and tells us how you like it!

Rik Schoemaker
on behalf of the nlmixr team


#5

Thanks Rik! Already forked it.


#6

Dear Rik,
I have the nlmixr installed successfully. I go thru some of the examples in the introduction document and found this can be very valuable new tool for pkpd modeler. Is there a formal forum for nlmixr users to discuss issues/questions on the use of nlmixr?

I am in the process of trying one of my pk model in nlmixr and I do have some problem to be resolved. How do you code the nlme_ode to allow log-transformed drug PK data?


#7

Dear Sam,

There is no forum yet for discussions: something we should think about. Does GitHub allow discussions of this sort?

I checked and I believe I may have an implementation of additive on log scale error. You do not need to specify weight (because there is none) and you can define the log transformation inside the ode.

You can find the data file in the following zip:
https://github.com/nlmixrdevelopment/demos/blob/master/GitHubnlmixrFiles.zip

Kind regards,

Rik

library(nlmixr)
source(“print.summary.lme.R”)
library(data.table)

datr <- read.csv(“BOLUS_1CPT.csv”, header=TRUE,stringsAsFactors=F)
datr$EVID<-ifelse(datr$EVID==1,101,datr$EVID)
datr<-data.table(datr)
datr<-datr[EVID!=2]

ode1 <- “
d/dt(centr) = -(CL/V)*centr;
logC1=log(centr/V);

mypar1 <- function(lCL, lV )
{
CL = exp(lCL)
V = exp(lV)
}

specs1<-list(fixed=lCL+lV~1, random = pdDiag(lCL+lV~1), start=c(lCL=1.6,lV=4.5))

################################################################################################

1 compartment IV Bolus: single dose, ODE with additive on log-scale error

################################################################################################

dat<-datr[SD==1]
dat[,DV:=log(DV)]

fitODElog <- nlme_ode(dat, model=ode1, par_model=specs1, par_trans=mypar1, response=“logC1”,
verbose=TRUE,control = nlmeControl(pnlsTol = .01, msVerbose = TRUE))

summary(fitODElog)
intervals(fitODElog)

################################################################################################

1 compartment IV Bolus: single dose, ODE with proportional error on linear scale

################################################################################################

dat2<-datr[SD==1]

fitODE <- nlme_ode(dat2, model=ode1, par_model=specs1, par_trans=mypar1, response=“centr”, response.scaler=“V”,
verbose=TRUE,weight=varPower(fixed=c(1)),control = nlmeControl(pnlsTol = .01, msVerbose = TRUE))

summary(fitODE)
intervals(fitODE)


#8

Dear @rikschoemaker - github issues is where all the discussion regarding usage should happen, although you can link those discussions here. e.g. checkout the issues page in the mrgsolve development here


#9

Dear Rik,
Thanks a lot. This example models in “GitHubnlmixrFiles” can be very useful for new nlmixr users.


#10

Dear Rik,
is there a documentation on the nlmixr input data format? A brief summary of any deviation from NONMEM data format would be very helpful.
Is there a way to show the prediction from the model to debug the new PK model?


#11

Dear Sam,
I agree the documentation can be improved big time :-). In the meantime, the nlmixr vignette gives instructions on how to specify the EVID record (in the appendix of the document). Examples of most data modifications to NONMEM datasets are provided in the https://github.com/nlmixrdevelopment/nlmixr/tree/master/inst/examples/ACoP7poster section that gives original NONMEM-style input data (the .csv files) while the nlmixr_examples_161020.R file gives syntax both to modify the data and to run the nlmixr models associated with the ACoP7 poster.
The main changes are in the more advanced EVID structure:

• EVID=0 denotes an observation event
• EVID>0 denotes an dosing event. In general, an EVID in nlmixr has five digits:
a. The right-most two digits are reserved for defining different events;
b. The next two digits point to which state variable (or compartment) this event is applied to;
c. The fifth digit from the right takes a value 1 if the dosing is an infusion and 0 if the dosing is a bolus.

So as an example, EVID is 101 if a bolus is entered into compartment 1 and EVID is 10101 if an infusion is started in compartment 1.

A big difference is that infusions need to be stopped with a separate record, and where AMT is set to the negative AMT of the start record, and steady state dosing is not implemented. In the code SS-dosing is simulated by running 7 doses in succession.

From the vignete:

• AMT is the drug amount with a bolus dosing. In case of an infusion, a positive number denotes the start of an infusion at a particular time with such an infusion rate; a negative number denotes the end of a previously started infusion

Hope that helps!

Kind regards,

Rik