In john-harrold/ubiquity: PKPD, PBPK, and Systems Pharmacology Modeling Tools

knitr::opts_chunk$set(echo = TRUE, message=FALSE, eval=FALSE)
require(ubiquity)
require(deSolve)
require(ggplot2)
require(foreach)
require(doParallel)
require(rhandsontable)

Working with rxode2, nlmixr2, NONMEM, and Monolix

Say you have a model in ubiquity and you want to use it with nlmixr2. Or perhaps you want to try out the simulation engine of rxode2, or you want to hand the model off to someone using Monolix. This should provide you with a way to easily convert your model into different formats in an automated fashion.

Getting your model in rxode2 format

As an example consider the system file for the two compartment model (system_2cmt.txt below). When you build the system an rxode2 output target will be created. You can use the system_fetch_template() function to access this file. Below I'm placing the script into the temporary directory. This is a script that will generate an rxode2 model function and store the function in the object my_model. This is the default input for nlmixr2. You can then use this function and the piping methodology to do things like add IIV terms to parameters. So you can use ubiquity to build your structural model and then piping, rxode2, and nlmixr2 for analysis.

sf_2cmt_txt = "
# Two compartment model with absorption compartment (At). The central
# compartment (Cc) had a volume (Vc), and the tissue compartment (Cp) has a
# volume (Vt). The system is parameterized in terms of macro constants (Q and
# CL)
#
# While the the default dosing is a 1 mg IV dose into the central compartment,
# the model is written to accept dosing into an absorption compartment and
# through continuous IV infusion.
#       _________ 
#      |         |
#      |   At    |
#      |         |
#      |_________|
#           |
#           | ka, fb
#           |
#           V
#       _________         _________  
#      |         |   Q   |         | 
#      |   Cc    |------>|   Cp    |
#      |   Vc    |<------|   Vt    |
#      |_________|       |_________|
#           |
#           | CL
#           |
#           V  
#
# System Units:
#                
#   mass          [=] mg 
#   volume        [=] ml
#   concentration [=] mg/ml
#   time          [=] hr


# #-------------#
# | Parameters  |
# #-------------#
#
# System parameters
#    name              value     lower      upper   units editable    grouping
#                                bound      bound

<P>  Vc                  1.0       eps        Inf   ml         yes     System
<P>  Vt                  1.0       eps        Inf   ml         yes     System
<P>  CL                  1.0       eps        Inf   ml/hr      yes     System
<P>  Q                   1.0       eps        Inf   ml/hr      yes     System
<P>  ka                  1.0       eps        Inf    1/hr      yes     System
<P>  fb                  1.0       eps        Inf   --         yes     System


# Bolus Events
# ------------
# times/events state   values        scale      units
<B:times>;             [  0  ];      1;          hours
<B:events>;      Cc;   [1.0  ];      1/Vc;       mg     
<B:events>;      At;   [0.0  ];      1;          mg     

# Infusion Rates 
# ------------
#  name     time/levels  values  scale    units
<R:Dinf>;     times;       [0];     1;    hours
<R:Dinf>;     levels;      [0];     1;    mg/hour

<ODE:At> -ka*At           
<ODE:Cc>  ka*At*fb/Vc - CL/Vc*Cc -Q*(Cc-Cp)/Vc  + Dinf/Vc
<ODE:Cp>                         +Q*(Cc-Cp)/Vt

<O> Cc_mg_ml    = Cc
<VP> prop_err   0.1            eps    inf      --     yes      Variance
<VP> add_err    0.1            eps    inf      ng/ml  yes      Variance
<OE:Cc_mg_ml> add=add_err; prop=prop_err

<TS:hours> 1.0
<TS:days>  1.0/24"
sf_2cmt = file.path("system_2cmt.txt")
sf_2cmt_full = file.path(tempdir(), sf_2cmt)

fileConn = file(sf_2cmt_full)
writeLines(sf_2cmt_txt, fileConn)
close(fileConn)

cfg = build_system(system_file         = sf_2cmt_full, 
                   temporary_directory = file.path(tempdir(), "transient"))

fr = system_fetch_template(cfg, 
                           template         = "nlmixr2", 
                           output_directory = tempdir(),
                           overwrite       =  TRUE)

library(rxode2)
source(file.path(tempdir(), "system_nlmixr2.R"))

You can see the contents of system_nlmixr2.R and what the my_model object looks like below. Note that to include the output with ubiquity you will need to have some error model. So along with the <O> Descriptor you will need to define some variance parameters (<VP>) and also an output error model (<OE:?>).

Getting your model in NONMEM format

Once you have your model in rxode2 format you can convert it further into NONMEM format with the babelmixr2 package (via the ruminate package). To do this you need at least one inter-individual variability term. You can do this in the system file with the <IIV:?> delimiter. In this example we will do that by using model piping. Here we specify that the typical value of Vc, the parameter TV_Vc, should have an error that is lognormally distributed. This will add the IIV term eta_Vc.

my_model = my_model |> 
  model({Vc = exp(TV_Vc + eta_Vc)}) 

Next we just need to load the ruminate package (>=0.2.4) and use the rx2other() function to create a NONMEM control stream:

library(ruminate)
my_model_nm = rx2other(my_model, out_type="nonmem")

Getting your model in Monolix format

You can also get the model in Monolix by specifying monolix as the output type:

my_model_mlx = rx2other(my_model, out_type="monolix")

Further reading

For more discussion on this topic see the following:

• https://github.com/nlmixr2/babelmixr2/issues/107

Models in different formats {.tabset}

ubiquity

Contents of system_2cmt.txt file:

cat(readLines(sf_2cmt_full), sep="\n")

rxode2 script

Contents of system_nlmixr2.R file produced by system_fetch_template():

cat(readLines( file.path(tempdir(), "system_nlmixr2.R")) , sep="\n")

rxode2 object

The my_model object:

my_model

NONMEM

The file containing the control stream is found here:

my_model_nm$files$ctl$fn_full.

cat(readLines(  my_model_nm$files$ctl$fn_full ) , sep="\n")

Monolix

The file containing the mlxtran output is found here:

my_model_mlx$files$mlxtran$fn_full

cat(readLines(  my_model_mlx$files$mlxtran$fn_full ) , sep="\n")

john-harrold/ubiquity documentation built on March 1, 2025, 12:35 p.m.