In MikeKSmith/rspeaksnonmem: Read, Modify and run NONMEM from within R

rspeaksnonmem is designed to allow the user to craft workflows based on a given NONMEM model.

For more information on how to use rspeaksnonmem to run NONMEM and PsN please review the vignette "Introduction to rspeaksnonmem".

After importing and parsing a control stream to an R object using the RNMImport package, rspeaksnonmem allows the user to easily change initial estimates, change data attributes or change task properties (estimation or simulation settings) without having to change the model. This saves a lot of textual searching by the user. The resulting data, parameter, model and task information is then written back out to file using an existing (possibly the original) model as a template. This ensures that the written NONMEM control stream mirrors as closely as possible the structure of the original file.

Installation

Dependencies

rspeaksnonmem relies on the package RNMImport function importNmMod which reads and parses the NONMEM control stream. rspeaksnonmem then works with the data, parameter values, and task information separately from the model.

Install RNMIMport from Github using the devtools package. rspeaksnonmem requires the version of RNMImport on MikeKSmith's Github repository which includes handling of OMEGA and SIGMA blocks and identifying which elements are fixed for these parameters.

  #devtools::install_git("MikeKSmith/RNMImport")

Install rspeaksnonmem

Eventually, rspeaksnonmem will be released to CRAN, but while still in development rspeaksnonmem can most easily be installed from GitHub using the devtools package:

    #devtools::install_github("MikeKSmith/rspeaksnonmem")

 devtools::load_all(pkg = ".")
 devtools::load_all(pkg = "C:\\Users\\smith_mk\\Documents\\Working documents\\RNMImport")

Load the rspeaksnonmem package

#library(rspeaksnonmem)

Copy an example dataset and model to a directory of your choice

getwd()

file.copy(from = file.path(system.file("exdata", package="rspeaksnonmem"),"warfarin_conc_pca.csv"),
            to = getwd(), overwrite=T )
file.copy(from = file.path(system.file("exdata", package="rspeaksnonmem"),"warfarin.ctl"),
            to = getwd(), overwrite=T )

The initial model within the workflow should act as a "template" for modifications. The best practice is to use a model where all possible parameters are defined (including OMEGAs and covariances / correlations) but where it is possible to fix parameters to zero or some null value. We can then run and test a wide variety of models simply by updating the $THETA, $OMEGA and $SIGMA parameters to allow estimation.

Read the control stream using RNMImport

First, we need to read the NONMEM control stream into R using the importNmMod function of RNMImport.

warfModel <- importNmMod("warfarin.ctl")
class(warfModel)
names(warfModel)

This creates a list of R objects containing the Raw NONMEM control stream as a vector of characters:

warfBaseRaw <- warfModel$Raw
cat(warfModel$Raw, sep="\n")

It also creates a parsed list of control statements in the controlFile slot:

warfTemplateModel <- warfModel$problemContents[[1]]
warfModel$problemContents[[1]]

Using rspeaksnonmem the user can then update elements of the model using this parsed set of commands. rspeaksnonmem provides some additional functions to extract certain elements of the parsed control stream.

"Object" | rspeaksnonmem Function | NMTRAN blocks ---------|---------------------------|----------------------------- Data | getNMDataObjects | $DATA, $INPUT Parameter| getNMParameterObjects | $THETA, $OMEGA, $SIGMA Task | getNMTaskPropertiesObjects| $EST,$COV, $TAB Model | getNMModelObjects | everything else

For example: getNMParameterObjects returns the $THETA, $OMEGA and $SIGMA records. Since the THETA parameters are represented as a data frame, we can easily update the initial values, lower or upper bounds, fix or estimate (unfix) THETA parameters. With each object returned by the various getNM<...> functions we get not only the parsed objects but the raw control stream as well.

NOTE: Generally speaking, rspeaksnonmem is set up to allow the user to quickly and easily update Data, Parameter or Task items, while leaving the Model unchanged.

The rspeaksnonmem function getNMDataObjects returns the RAW $DATA and $INPUT lines from the NONMEM control stream and the associated parsed objects.

getNMDataObjects(warfModel)

Similarly, getNMParameterObjects returns the $THETA, $OMEGA and $SIGMA records.

getNMParameterObjects(warfModel)

Finally, getNMTasProperties returns the $EST, $COV and $TABLES records.

getNMTaskPropertiesObjects(warfModel)

Together the Data, Parameters, Model and Task information forms a "Modelling Object Group" or MOG which is used for a specific estimation task. We can update the MOG using the rspeaksnonmem function updateModel.

newTheta <- warfTemplateModel$Theta
newTheta["POP_TLAG",] <- list(Lower=0, Est=0.75, Upper=1.5, FIX=FALSE, comments="POP_TLAG")
newWarfModel <- updateModel(parsedObject = warfTemplateModel, theta=newTheta)
newWarfModel

Note how the THETA settings for POP_TLAG have been updated, but all other elements are unchanged.

It is also possible to change the initial estimates while retaining the specified lower and upper bounds:

newTheta2 <- warfTemplateModel$Theta
newInits <- c(0.25, 10, 4, 0.25, 0, 0.1)
newTheta2[,"Est"] <- newInits
newTheta2

We can finally write out the updated control stream using the function writeNMControlStream. This function uses an existing NONMEM control stream (or raw records) as a template for writing the updated, parsed statements. This ensures that the order of statements corresponds between the original and updated files, which allows easier checking of changes.

writeNMControlStream(templateModel = warfModel$Raw, 
                     parsedControl = list(newWarfModel), 
                     outputFile = "warfarin2", 
                     outFileExtension = "ctl")
cat(readLines("warfarin2.ctl", warn = F), sep="\n")