Metamodels in RsNLME • Certara.RsNLME

Introduction

Metamodel is a file with .mmdl extension that makes integration of RsNLME into Pirana seamlessly. It contains all the information needed to run a PML (Pharmacometric Modeling Language behind Phoenix NLME and RsNLME) model. In such sense, it is similar to the control file in NONMEM. The purpose of this vignette is to give an overview of metamodels and to demonstrate how they can be useful for easy storing, fitting, and exploring PML models. To do this, we need to understand how PML models are run in Phoenix NLME and RsNLME.

Here we assume that both input data and PML model files are in the current working directory. In addition, the directory with NLME Engine is assumed to be given as an environment variable (INSTALLDIR), which can be checked with the following commands:

  # loading the package
  library(Certara.RsNLME)
  Sys.getenv("INSTALLDIR")
[1] "C:\\Program Files\\Certara\\NLME_Engine"

These assumptions will be applied to the rest of this document.

PML models in Phoenix NLME

Phoenix NLME users may know that the PML model used in Phoenix infrastructure has an .mdl extension with statements/equations/assignments. However, it does not contain any information on the dataset to be used, columns to be mapped to model variables, or tables to be generated (e.g., see below for an example of PML models whose syntax can be found here).


    test(){
        cfMicro(A1, Cl / V)
        dosepoint(A1)
        C = A1 / V
        error(CEps = 0.1)
        observe(CObs = C * (1 + CEps))
        stparm(V = tvV * exp(nV))
        stparm(Cl = tvCl * exp(nCl))
        fixef(tvV = c(, 5, ))
        fixef(tvCl = c(, 1, ))
        ranef(diag(nV, nCl) = c(1, 1))
    }

The reason behind this is that Phoenix project keeps all the information (e.g., model, data, and mapping), which is loaded by NLME Exectutables only during execution. However, if the user wants to run PML models in command line, data and column definition (relationship between data columns and model variables, tables to be generated, additional input options such as ADDL and SS) files should be provided along with engine arguments:

%INSTALLDIR%\runNLME.bat 5 100 OneCpt_IVInfusion.mdl columnMapping.txt OneCpt_IVInfusionData.csv

PML models in RsNLME

Alternatively, the user can run PML models in RsNLME command-line interface. To do this, the user needs to provide input data and PML model files to create a model object. Then, mappings between data columns and model variables need to be provided if they are not automatically mapped. In addition, if any additional input options (such as ADDL, SS, Reset, MDV, and infusion) is desired, it should be provided too. An example to run a PML model is given as follows:

PMLModelCodeOutput  <-" 
    test(){
        cfMicro(A1, Cl / V)
        dosepoint(A1)
        C = A1 / V
        error(CEps = 0.1)
        observe(CObs = C * (1 + CEps))
        stparm(V = tvV * exp(nV))
        stparm(Cl = tvCl * exp(nCl))
        fixef(tvV = c(, 5, ))
        fixef(tvCl = c(, 1, ))
        ranef(diag(nV, nCl) = c(1, 1))
    } "

  PMLModelCodeFile <- file.path(tempdir(TRUE), "OneCpt_IVInfusion.mdl")
  writeLines(PMLModelCodeOutput, PMLModelCodeFile)
  # for the input data description please refer to ?OneCpt_IVInfusionData
  # Create the model object
  model <- textualmodel(modelName = "OneCpt_IVInfusion",
                        mdl = PMLModelCodeFile,
                        data = OneCpt_IVInfusionData)
  # View the model and its associated column mappings 
  print(model)

 Model Overview 
 ------------------------------------------- 
Model Name        :  OneCpt_IVInfusion
Working Directory :  C:/Users/jcraig/Documents/GitHub/R-RsNLME/vignettes/OneCpt_IVInfusion
Model Type        :  Textual

 PML 
 ------------------------------------------- 

    test(){
        cfMicro(A1, Cl / V)
        dosepoint(A1)
        C = A1 / V
        error(CEps = 0.1)
        observe(CObs = C * (1 + CEps))
        stparm(V = tvV * exp(nV))
        stparm(Cl = tvCl * exp(nCl))
        fixef(tvV = c(, 5, ))
        fixef(tvCl = c(, 1, ))
        ranef(diag(nV, nCl) = c(1, 1))
    }

 Structural Parameters 
 ------------------------------------------- 
 V Cl
 ------------------------------------------- 
 Column Mappings 
 ------------------------------------------- 
Model Variable Name : Data Column name
id                  : ?
time                : Time
A1                  : ?
CObs                : CObs
  # Manually map the un-mapped model variables
  model <- colMapping(model, mappings = c(id = "Subject", A1 = "Dose"))
  
  # Add infusion information for dosing compartment A1
  model <- addInfusion(model, "A1", isDuration = TRUE, colName = "Duration")

Overall, we see that, to run PML models in command line, the user should provide at least two files: a .mdl file for PML model, and a file for column definitions/mappings. Metamodel concept allows users to integrate these parts in one file.

Metamodel Overview

Metamodels consist of different blocks, which start with doubled number sign (##) and the name of the block. The syntax for comments in metamodels is the same as the one in PML models: single line can be commented by ‘#’ or ‘//’, and multiple lines can be commented by /*…*/.

Author:

The author of the metamodel. An example of usage:
```
##Author: User
```
This block is not used for estimation, and hence can be omitted.
Description:

Used by Pirana for the description column presentation with the purpose to understand the goal of the metamodel created.

An example of usage:
```
##Description: User
```
This block is not used for estimation, and hence can be omitted.
Based on

Provide the name of the reference metamodel, which is used by Pirana for the reference trees building.

An example of usage:
```
##Based on mmdl1
```
This block is not used for estimation, and hence can be omitted.
DATA

A mandatory block specifying the path of the input data. Both absolute and relative paths are acceptable. The data.table::fread() function with default settings is used to load the data.

An example of usage:
```
##DATA ./data.csv
```
MAP

This block is used to specify mappings between model variables and data columns through “=” sign: “variableName=columnName”. If the mapping for a model variable is not specified, then the name of the model variable is assumed to be the same as the one for its corresponding data column: CObs=CObs is equal to CObs.

There are some special variables that are not presented in the model but are needed to be mapped to their corresponding data columns can also be specified in this block. These variables are listed below:
- id
  
  It is required to be mapped for population models to identify individual data profiles for different subjects data distinction. The user can map up to five data columns to it with column names separated by comma.
- time
  
  It is required to be mapped for time-based models.
- dosingCompartmentName_Rate
  
  Mapping this variable to its corresponding data column will indicate that the associated dosing compartment involves infusion with rate information provided.
- dosingCompartmentName_Duration
  
  Mapping this variable to its corresponding data column will indicate that the associated dosing compartment involves infusion with duration information provided.
- SS
  
  Mapping this variable to its corresponding data column if the input data contains a column indicating steady state. This specification will be translated to sscol statement in the column definition file.
- SSOffSet
  
  Mapping this variable to its corresponding data column if the input data contains a column indicating SSOffSet. Only applicable to the case where SS appears in the MAP block. This specification will be translated to ssoffcol statement in the column definition file.
- ADDL
  
  Mapping this variable to its corresponding data column if the input data contains an additional identical dose column. This specification will be translated to addlcol statement in the column definition file.
- II
  
  It denotes the inter-dose interval, and must be mapped to its corresponding data column when either SS or ADDL appears in the MAP block. This specification will be translated to iicol statement in the column definition file.
- MDV
  
  Mapping this variable to its corresponding data column if the input data contains a MDV column. If mapped, those observations corresponding to the rows in the MDV column having non-zero numeric values will be ignored.
- Reset
  
  Mapping this variable to its corresponding data column if the input data contains a column indicting reset with lowest and highest values being 4.
An example in which id is mapped to two data columns, Study and Subject:
```
## MAP id = Study, Subject time = Time Aa = Dose CObs
```
An example of usage where dosing compartment, A1, involves infusion with rate information provided in the data column Rate:
```
##MAP id = Subject time = Time A1 = Dose A1_Rate = Rate CObs
```
For a categorical covariate (specified through either fcovariate or covariate with the covariate name followed by an empty parenthesis), if its corresponding data column is of class character, then the user can also use the ##MAP block to define the label name for each of its category values. Specifically, label names are given in round brackets right after the column name (Note: label names with spaces are not supported by NLME): both with and without quotation marks around the label name are supported. See below for an example in which categorical covariate, Sex, is mapped to data column “Gender” that has values “Male” and “Female”.
```
##MAP id = Subject time = Time Aa = Dose CObs 

Sex = Gender(Male = 0, Female = 1)
```

DOSING CYCLE

This block provides an alternative way to define ADDL or steady state (SS) by specifying dosing cycle information for each of the involved dosing compartments (more information can be found here for ADDL and here for SS on this topic). The syntax for defining SS and ADDL are respectively given as follows:

SS = [COL]  Dosepoint = [CMT] Amount = [NUM/COL] Delta = [NUM/COL] Duration = [NUM/COL] Rate = [NUM/COL] 

ADDL = [COL] Delta = [NUM/COL] Dosepoint = [CMT] Amount = [NUM/COL] Duration = [NUM/COL] Rate = [NUM/COL]

Here [COL] is a column name with ADDL/SS flag,

[CMT] is a dosing compartment name from the model,

[NUM/COL] is a column name or a numeric value.

and argument Delta denotes the inter-dose interval.

An example of usage:

##DOSING CYCLE SS=SS Dosepoint = Aa Delta = 10 Amount = Dose Duration = 5
# A single line comment could be started with a simple number sign “#”
# Here it is assumed that the input dataset contains 'SS' and 'Dose' columns
# and the model has 'dosepoint(Aa)' statement.
# Each dosing cycle lasts 10 time units with the amount and
# the duration of the infusion set to ‘Dos’ and 5 time units, respectively.

Another example:

##DOSING CYCLE ADDL = ADDL DELTA = II AMOUNT = Dose Dosepoint = A1
# It is assumed that the input data contains 'ADDL', 'II' and 'Dose' columns
# and the model has 'dosepoint(A1)' statement. 
# The dose is administered to ‘A1’ as a bolus with amount set to be ‘Dose’, 
# and the inter-dose interval is set to ‘II’.

COLDEF

This block is used to define column definitions, which can be provided through the syntax given here. It is particularly useful for the case where column definitions cannot be defined through neither ##MAP nor ##DOSING CYCLE blocks. For example, see below for an example on using ##COLDEF block to define the column definition for the case where the input data contains a Reset column with lowest and highest values being 1 and 4, respectively.
```
##COLDEF
reset("Reset", c(1, 4))
```
Experienced users may want to omit ##MAP block and provide all definitions in this block. It is also possible to use both ##MAP and ##COLDEF blocks, in such case, they are appended. An example of usage where all the column definitions are defined in this block:
```
##COLDEF id("id")
         time("time")
         dose(A1<-"dose")
         covr(sex<-"sex"("male" = 0, "female" = 1))
         covr(wt<-"wt")
         obs(CObs<-"conc")
```

MODEL

This block is mandatory, and is used to specify the PML model. See “Modeling Syntax” for details on the available statements.

An example of usage:

##MODEL 
    test(){
        cfMicro(A1, Cl / V)
        dosepoint(A1)
        C = A1 / V
        error(CEps = 0.1)
        observe(CObs = C * (1 + CEps))
        stparm(V = tvV * exp(nV))
        stparm(Cl = tvCl * exp(nCl))
        fixef(tvV = c(, 5, ))
        fixef(tvCl = c(, 1, ))
        ranef(diag(nV, nCl) = c(1, 1))
    }

ESTARGS

This block is used to specify the values of engine arguments through the syntax provided by the R function, engineParams(), in the Certara.RsNLME package with arguments separated by either comma or space. If not provided, then the default values will be used.

An example of usage:
```
##ESTARGS
method = "QRPEM",
iSample = 1200,
maxStepsODE = 50000000,
mapAssist = 1,
ODE = "AutoDetect",
numIterations = 0
```
The user can specify multiple ESTARGS blocks, in such case they will be executed sequentially, the model for the next run will use the final estimates from the previous one
SIMARGS

This block is used to specify arguments for model simulation. The following arguments are supported:
- numReplicates - Number of replicates to simulate (default is 100),
- seed - Random number generator seed (default is 1234),
- sort, ODE, rtolODE, atolODE, maxStepsODE - see engineParams().
Multiple ESTARGS and SIMARGS are supported, but applied sequentially: all ESTARGS runs first, all SIMARGS the second.
TABLES

This block is used to define additional tables to be outputted (Note: posthoc.csv table is created by default with structural parameters and covariates outputted at each data row). The syntax for table statement can be found here.

An example of usage:
```
##TABLES
table(file=“table01.csv”, time(0,10,seq(2,8,0.1)), 
      dose(A1), covr(BW), obs(Conc), BW, C, cObs, V, Ke)
```
It is worth pointing out that tables are part of column definition, and hence they can be defined in the ##COLDEF block instead of the ##TABLES block. If they are defined in both blocks, then all of them will be added to the column definition file.

Metamodel Example

Now we create a simple metamodel and name it as OneCpt_IVInfusion.mmdl.

## Description: A one-compartment model with IV infusion

# The model is fitted using the default FOCE-ELS engine.
# Note: the default values for the relevant NLME engine arguments are chosen based on the model, type ?engineParams for details.

## DATA OneCpt_IVInfusion.csv
## MAP id = Subject time = Time A1 = Dose A1_Rate = Rate CObs

## MODEL
test(){
    cfMicro(A1, Cl / V)
    dosepoint(A1)
    C = A1 / V
    error(CEps = 0.1)
    observe(CObs = C * (1 + CEps))
    stparm(V = tvV * exp(nV))
    stparm(Cl = tvCl * exp(nCl))
    fixef(tvV = c(, 5, ))
    fixef(tvCl = c(, 1, ))
    ranef(diag(nV, nCl) = c(1, 1))
}

## ESTARGS
numIterations = 1 # one iteration only
stdErr = "None" # no standard error estimation requested

There are 5 blocks used in this model:

Description

This block provides descriptions of the model.
DATA

This block says that data file OneCpt_IVInfusion.csv (in the current working directory) will be used as the input data.
MAP

This block maps five variables to their corresponding data columns:
- id = Subject says that id is mapped to data column Subject.
- time = Time says that time is mapped to data column Time.
- A1 = Dose says that model variable A1 is mapped to data column Dose (which denoting the amount of drug administered to dosing compartment A1).
- A1_Rate = Rate says that data column Rate specifies the infusion rate associated with the dose administered to dosing compartment A1
  
  (Note: If rate is zero or unspecified, then it is assumed to be a bolus).
- CObs says that model variable CObs is mapped to data column CObs
(Note: Here the name of the model observable is the same as the column name, and hence the equal sign is omitted. It can be also rewritten as CObs = CObs).
MODEL

This block says that an one-compartment population model with clearance parameterization is used as PML model.
ESTARGS

This block specifies the values for engine arguments “numIterations” and “stdErr”. Specifically,
- numIterations = 1 sets the number of optimization iterations to 1.
- stdErr = “None” means that standard error is not requested.

Metamodel run

To run metamodels, one can use the function run_metamodel(). Here we demonstrate how to run the metamodel created above with host set to be a local one without any parallelization (Note: For the case where host is not specified, if MPI is found, then it will be automatically parallelized over 4 threads; otherwise, a simple one core execution is performed).

Note: We will additionally use the argument directoryToRun, in which case the model output files will be created in a new sub-folder of your working directory, named OneCpt_IVInfusion.

host <- hostParams(parallelMethod = "None",
                   hostName = "local",
                   numCores = 1)

OneCpt_IVInfusionFile <-
  system.file("vignettesdata/OneCpt_IVInfusion.mmdl",
              package = "Certara.RsNLME",
              mustWork = TRUE)

OneCpt_IVInfusionResults <-
  run_metamodel(mmdlfile = OneCpt_IVInfusionFile,
                directoryToRun = "OneCpt_IVInfusion",
                host = host)
#> 
#> NLME Job
#> 
#> Compiling 1 of 1 NLME models
#> TDL5 version: 23.10.1.283
#> Servername: no-net
#> Sentinel License file: C:\PROGRA~1\Certara\NLME_E~1\lservrc
#> 
#> License Name: NU_PMLCommandLine_N
#> License Type: Commercial
#> Expiration Date: 31 Jan 2024
#> Current Date: 27 Nov 2023
#> Days until program expires: 65
#> The model compiled
#> 
#> 
#>  Iteration    -2LL    tvV     tvCl nSubj nObs
#>          1 2171.54 4.7378 0.847792   100  700
#> 
#> Trying to generate job results...
#> 
#> Generating Overall.csv
#> Generating EtaEta.csv
#> Generating Eta.csv
#>  Generating EtaStacked.csv
#>  Generating bluptable.dat
#> Generating ConvergenceData.csv
#> Generating initest.csv
#> Generating doses.csv
#> Generating omega.csv
#> Generating theta.csv
#>  Generating thetaCorrelation.csv
#>  Generating thetaCovariance.csv
#>  Generating Covariance.csv
#> Generating Residuals.csv
#> Generating posthoc.csv
#> 
#> Finished summarizing results. Transferring data and loading the results...
#> Done generating job results.

Loading metamodel for other run modes

If the user wants to execute a run mode other than simple fit, he/she should load the metamodel into R and then provide the resulted NlmePmlModel object and engine parameters (if specified in the metamodel) to the desired model execution function as values for model and params arguments, respectively. To load the metamodel, the user should run create_model_from_metamodel() function. The resulted list will have two objects: the model object and engine parameters. For example, the codes below are used to perform bootstrap run for the metamodel created above.

  ModelParamsList <- 
    create_model_from_metamodel(mmdlfile = OneCpt_IVInfusionFile)
  bootParams <- BootstrapParams(numReplicates = 5,
                                randomNumSeed = 1234)
  bootstrap(model = ModelParamsList$model,
            params = ModelParamsList$params,
            bootParams = bootParams)
#> Boot Job
#> 
#> NLME Job
#> 
#> Preparing files for Bootstrap run
#> TDL5 version: 23.10.1.283
#> Servername: no-net
#> Sentinel License file: C:\PROGRA~1\Certara\NLME_E~1\lservrc
#> 
#> License Name: NU_PMLCommandLine_N
#> License Type: Commercial
#> Expiration Date: 31 Jan 2024
#> Current Date: 27 Nov 2023
#> Days until program expires: 65
#> The model compiled
#> 
#> Processing replicates
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/0/0
#> 
#> Num Jobs/Completed/Failed:5/2/0
#> 
#> Num Jobs/Completed/Failed:5/2/0
#> 
#> Num Jobs/Completed/Failed:5/3/0
#> 
#> Num Jobs/Completed/Failed:5/3/0
#> 
#> Num Jobs/Completed/Failed:5/3/0
#> 
#> Num Jobs/Completed/Failed:5/4/0
#> 
#> Num Jobs/Completed/Failed:5/4/0
#> 
#> Num Jobs/Completed/Failed:5/4/0
#> 
#> Num Jobs/Completed/Failed:5/4/0
#> 
#> Num Jobs/Completed/Failed:5/5/0
#> 
#> 
#> Summarizing bootstrap results
#> $BootOverall
#>    Scenario Replicate ReturnCode       LL
#> 1:      (B)         1          4 -1089.74
#> 2:      (B)         2          4 -1101.49
#> 3:      (B)         3          4 -1107.09
#> 4:      (B)         4          4 -1058.07
#> 5:      (B)         5          4 -1116.01
#> 
#> $BootTheta
#>    Scenario Parameter      Mean      Stderr      CV%    Median      2.5%
#> 1:      (B)       tvV 4.7780358 0.010690711 0.223747 4.7720495 4.7682849
#> 2:      (B)      tvCl 0.8549129 0.010439461 1.221114 0.8508297 0.8473296
#> 3:      (B)      CEps 0.0941178 0.007253549 7.706883 0.0912145 0.0849065
#>        97.5%
#> 1: 4.7949070
#> 2: 0.8731682
#> 3: 0.1038233
#> 
#> $BootOmega
#>    Scenario Label        nV      nCl
#> 1:      (B)    nV 0.4137216 0.000000
#> 2:      (B)   nCl 0.0000000 0.491738
#> 
#> $BootOmegaCorrelation
#>    Scenario Label nV nCl
#> 1:      (B)    nV  1   0
#> 2:      (B)   nCl  0   1
#> 
#> $BootOmegaStderr
#>    Scenario Label          nV         nCl
#> 1:      (B)    nV 0.009838646 0.000000000
#> 2:      (B)   nCl 0.000000000 0.007960785
#> 
#> $BootVarCoVar
#>    Scenario Var Name           tvV          tvCl          CEps
#> 1:      (B)      tvV  1.142913e-04 -3.831222e-05  1.534593e-05
#> 2:      (B)     tvCl -3.831222e-05  1.089823e-04 -4.863693e-05
#> 3:      (B)     CEps  1.534593e-05 -4.863693e-05  5.261397e-05