RxODE (version 0.9.0-8)

rxControl: Solving \& Simulation of a ODE/solved system (and solving options) equation

Description

This uses RxODE family of objects, file, or model specification to solve a ODE system.

Usage

rxControl(scale = NULL, method = c("liblsoda", "lsoda", "dop853"),
  transitAbs = NULL, atol = 1e-08, rtol = 1e-06, maxsteps = 70000L,
  hmin = 0L, hmax = NA, hmaxSd = 0, hini = 0, maxordn = 12L,
  maxords = 5L, ..., cores, covsInterpolation = c("locf", "linear",
  "nocb", "midpoint"), addCov = FALSE, matrix = FALSE, sigma = NULL,
  sigmaDf = NULL, sigmaLower = -Inf, sigmaUpper = Inf,
  nCoresRV = 1L, sigmaIsChol = FALSE, nDisplayProgress = 10000L,
  amountUnits = NA_character_, timeUnits = "hours", stiff,
  theta = NULL, thetaLower = -Inf, thetaUpper = Inf, eta = NULL,
  addDosing = FALSE, stateTrim = Inf, updateObject = FALSE,
  omega = NULL, omegaDf = NULL, omegaIsChol = FALSE,
  omegaLower = -Inf, omegaUpper = Inf, nSub = 1L, thetaMat = NULL,
  thetaDf = NULL, thetaIsChol = FALSE, nStud = 1L, dfSub = 0,
  dfObs = 0, returnType = c("rxSolve", "matrix", "data.frame",
  "data.frame.TBS"), seed = NULL, nsim = NULL, minSS = 7,
  maxSS = 70000, atolSS = 1e-09, rtolSS = 1e-09, params = NULL,
  events = NULL, istateReset = TRUE, subsetNonmem = TRUE,
  linLog = FALSE, maxAtolRtolFactor = 0.1, from = NULL, to = NULL,
  by = NULL, length.out = NULL, iCov = NULL, keep = NULL)

rxSolve(object, ...)

# S3 method for default rxSolve(object, params = NULL, events = NULL, inits = NULL, ...)

# S3 method for rxSolve update(object, ...)

# S3 method for RxODE predict(object, ...)

# S3 method for rxSolve predict(object, ...)

# S3 method for rxEt predict(object, ...)

# S3 method for rxParams predict(object, ...)

# S3 method for RxODE simulate(object, nsim = 1L, seed = NULL, ...)

# S3 method for rxSolve simulate(object, nsim = 1L, seed = NULL, ...)

# S3 method for rxParams simulate(object, nsim = 1L, seed = NULL, ...)

# S3 method for rxSolve solve(a, b, ...)

# S3 method for RxODE solve(a, b, ...)

# S3 method for rxParams solve(a, b, ...)

# S3 method for rxEt solve(a, b, ...)

Arguments

scale

a numeric named vector with scaling for ode parameters of the system. The names must correstond to the parameter identifiers in the ODE specification. Each of the ODE variables will be divided by the scaling factor. For example scale=c(center=2) will divide the center ODE variable by 2.

method

The method for solving ODEs. Currently this supports:

  • "liblsoda" thread safe lsoda. This supports parallel thread-based solving, and ignores user Jacobian specification.

  • "lsoda" -- LSODA solver. Does not support parallel thread-based solving, but allows user Jacobian specification.

  • "dop853" -- DOP853 solver. Does not support parallel thread-based solving nor user Jacobain specification

transitAbs

boolean indicating if this is a transit compartment absorption

atol

a numeric absolute tolerance (1e-8 by default) used by the ODE solver to determine if a good solution has been achieved; This is also used in the solved linear model to check if prior doses do not add anything to the solution.

rtol

a numeric relative tolerance (1e-6 by default) used by the ODE solver to determine if a good solution has been achieved. This is also used in the solved linear model to check if prior doses do not add anything to the solution.

maxsteps

maximum number of (internally defined) steps allowed during one call to the solver. (5000 by default)

hmin

The minimum absolute step size allowed. The default value is 0.

hmax

The maximum absolute step size allowed. When hmax=NA (default), uses the average difference (+hmaxSd*sd) in times and sampling events. When hmax=NULL RxODE uses the maximum difference in times in your sampling and events. The value 0 is equivalent to infinite maximum absolute step size.

hmaxSd

The number of standard deviations of the time difference to add to hmax. The default is 0

hini

The step size to be attempted on the first step. The default value is determined by the solver (when hini = 0)

maxordn

The maximum order to be allowed for the nonstiff (Adams) method. The default is 12. It can be between 1 and 12.

maxords

The maximum order to be allowed for the stiff (BDF) method. The default value is 5. This can be between 1 and 5.

...

Other arguments including scaling factors for each compartment. This includes S# = numeric will scale a compartment # by a dividing the compartment amount by the scale factor, like NONMEM.

cores

Number of cores used in parallel ODE solving. This defaults to the number or system cores determined by rxCores for methods that support parallel solving (ie thread-safe methods like "liblsoda").

covsInterpolation

specifies the interpolation method for time-varying covariates. When solving ODEs it often samples times outside the sampling time specified in events. When this happens, the time varying covariates are interpolated. Currently this can be:

  • "linear" interpolation (the default), which interpolates the covariate by solving the line between the observed covariates and extrapolating the new covariate value.

  • "constant" -- Last observation carried forward.

  • "NOCB" -- Next Observation Carried Backward. This is the same method that NONMEM uses.

  • "midpoint" Last observation carried forward to midpoint; Next observation carried backward to midpoint.

addCov

A boolean indicating if covariates should be added to the output matrix or data frame. By default this is disabled.

matrix

A boolean inticating if a matrix should be returned instead of the RxODE's solved object.

sigma

Named sigma covariance or Cholesky decomposition of a covariance matrix. The names of the columns indicate parameters that are simulated. These are simulated for every observation in the solved system.

sigmaDf

Degrees of freedom of the sigma t-distribution. By default it is equivalent to Inf, or a normal distribution.

sigmaLower

Lower bounds for simulated unexplained variability (by default -Inf)

sigmaUpper

Upper bounds for simulated unexplained variability (by default Inf)

nCoresRV

Number of cores used for the simulation of the sigma variables. By default this is 1. This uses the package rmvn and rmvt. To reproduce the results you need to run on the same platform with the same number of cores. This is the reason this is set to be one, regardless of what the number of cores are used in threaded ODE solving.

sigmaIsChol

Boolean indicating if the sigma is in the Cholesky decomposition instead of a symmetric covariance

nDisplayProgress

An integer indicating the minimum number of c-based solves before a progress bar is shown. By default this is 10,000.

amountUnits

This supplies the dose units of a data frame supplied instead of an event table. This is for importing the data as an RxODE event table.

timeUnits

This supplies the time units of a data frame supplied instead of an event table. This is for importing the data as an RxODE event table.

stiff

a logical (TRUE by default) indicating whether the ODE system is stiff or not.

For stiff ODE sytems (stiff = TRUE), RxODE uses the LSODA (Livermore Solver for Ordinary Differential Equations) Fortran package, which implements an automatic method switching for stiff and non-stiff problems along the integration interval, authored by Hindmarsh and Petzold (2003).

For non-stiff systems (stiff = FALSE), RxODE uses DOP853, an explicit Runge-Kutta method of order 8(5, 3) of Dormand and Prince as implemented in C by Hairer and Wanner (1993).

theta

A vector of parameters that will be named THETA[#] and added to parameters

thetaLower

Lower bounds for simulated population parameter variability (by default -Inf)

thetaUpper

Upper bounds for simulated population unexplained variability (by default Inf)

eta

A vector of parameters that will be named ETA[#] and added to parameters

addDosing

Boolean indicating if the solve should add RxODE EVID and related columns. This will also include dosing information and estimates at the doses. Be default, RxODE only includes estimates at the observations. (default FALSE). When addDosing is NULL, only include EVID=0 on solve and exclude any model-times or EVID=2. If addDosing is NA the classic RxODE EVID events. When addDosing is TRUE add the event information in NONMEM-style format; If subsetNonmem=FALSE RxODE will also extra event types (EVID) for ending infusion and modeled times:

  • EVID=-1 when the modeled rate infusions are turned off (matches rate=-1)

  • EVID=-2 When the modeled duration infusions are turned off (matches rate=-2)

  • EVID=-10 When the specified rate infusions are turned off (matches rate>0)

  • EVID=-20 When the specified dur infusions are turned off (matches dur>0)

  • EVID=101,102,103,... Modeled time where 101 is the first model time, 102 is the second etc.

stateTrim

When amounts/concentrations in one of the states are above this value, trim them to be this value. By default Inf. Also trims to -stateTrim for lage negative amounts/concentrations

updateObject

This is an internally used flag to update the RxODE solved object (when supplying an RxODE solved object) as well as returning a new object. You probably should not modify it's FALSE default unless you are willing to have unexpected results.

omega

Estimate of Covariance matrix. When omega is a list, assume it is a block matrix and convert it to a full matrix for simulations.

omegaDf

The degrees of freedom of a t-distribution for simulation. By default this is NULL which is equivalent to Inf degrees, or to simulate from a normal distribution instead of a t-distribution.

omegaIsChol

Indicates if the omega supplied is a Cholesky decomposed matrix instead of the traditional symmetric matrix.

omegaLower

Lower bounds for simulated ETAs (by default -Inf)

omegaUpper

Upper bounds for simulated ETAs (by default Inf)

nSub

Number between subject variabilities (ETAs) simulated for every realization of the parameters.

thetaMat

Named theta matrix.

thetaDf

The degrees of freedom of a t-distribution for simulation. By default this is NULL which is equivalent to Inf degrees, or to simulate from a normal distribution instead of a t-distribution.

thetaIsChol

Indicates if the theta supplied is a Cholesky decomposed matrix instead of the traditional symmetric matrix.

nStud

Number virtual studies to characterize uncertainty in estimated parameters.

dfSub

Degrees of freedom to sample the between subject variaiblity matrix from the inverse Wishart distribution (scaled) or scaled inverse chi squared distribution.

dfObs

Degrees of freedom to sample the unexplained variaiblity matrix from the inverse Wishart distribution (scaled) or scaled inverse chi squared distribution.

returnType

This tells what type of object is returned. The currently supported types are:

  • "rxSolve" (default) will return a reactive data frame that can change easily change different pieces of the solve and update the data frame. This is the currently standard solving method in RxODE, is used for rxSolve(object, ...), solve(object,...),

  • "data.frame" -- returns a plain, non-reactive data frame; Currently very slightly Faster than returnType="matrix"

  • "matrix" -- returns a plain matrix with column names attached to the solved object. This is what is used object$run as well as object$solve

seed

an object specifying if and how the random number generator should be initialized

nsim

represents the number of simulations. For RxODE, if you supply single subject event tables (created with eventTable)

minSS

Minimum number of iterations for a steady-state dose

maxSS

Maximum number of iterations for a steady-state dose

atolSS

Absolute tolerance to check if a solution arrived at steady state.

rtolSS

Relative tolerance to check if a solution arrived at steady state.

params

a numeric named vector with values for every parameter in the ODE system; the names must correspond to the parameter identifiers used in the ODE specification;

events

an eventTable object describing the input (e.g., doses) to the dynamic system and observation sampling time points (see eventTable);

istateReset

When TRUE, reset the ISTATE variable to 1 for lsoda and liblsoda with doses, like deSolve; When FALSE, do not reset the ISTATE variable with doses.

subsetNonmem

subset to NONMEM compatible EVIDs only. By default TRUE.

linLog

Boolean indicating if linear compartment models be calculated more accurately in the log-space (slower) By default this is off (FALSE)

maxAtolRtolFactor

The maximum atol/rtol that FOCEi and other routines may adjust to. By default 0.1

from

When there is no observations in the event table, start observations at this value. By default this is zero.

to

When there is no observations in the event table, end observations at this value. By default this is 24 + maximum dose time.

by

When there are no observations in the event table, this is the amount to increment for the observations between `from` and `to`.

length.out

The number of observations to create if there isn't any observations in the event table. By default this is 200.

iCov

A data frame of individual non-time varying covariates to combine with the params to form a parameter data.frame.

keep

Columns to keep from either the input dataset or the iCov dataset. With the iCov dataset, the column is kept once per line. For the input dataset, if any records are added to the data LOCF (Last Observation Carried forward) imputation is performed.

object

is a either a RxODE family of objects, or a file-name with a RxODE model specification, or a string with a RxODE model specification.

inits

a vector of initial values of the state variables (e.g., amounts in each compartment), and the order in this vector must be the same as the state variables (e.g., PK/PD compartments);

a

when using solve, this is equivalent to the object argument. If you specify object later in the argument list it overwrites this parameter.

b

when using solve, this is equivalent to the params argument. If you specify params as a named argument, this overwrites the output

Value

An “rxSolve” solve object that stores the solved value in a matrix with as many rows as there are sampled time points and as many columns as system variables (as defined by the ODEs and additional assignments in the RxODE model code). It also stores information about the call to allow dynamic updating of the solved object.

The operations for the object are simialar to a data-frame, but expand the $ and [[""]] access operators and assignment operators to resolve based on different parameter values, initial conditions, solver parameters, or events (by updaing the time variable).

You can call the eventTable methods on the solved object to update the event table and resolve the system of equations.

References

Hindmarsh, A. C. ODEPACK, A Systematized Collection of ODE Solvers. Scientific Computing, R. S. Stepleman et al. (Eds.), North-Holland, Amsterdam, 1983, pp. 55-64.

Petzold, L. R. Automatic Selection of Methods for Solving Stiff and Nonstiff Systems of Ordinary Differential Equations. Siam J. Sci. Stat. Comput. 4 (1983), pp. 136-148.

Hairer, E., Norsett, S. P., and Wanner, G. Solving ordinary differential equations I, nonstiff problems. 2nd edition, Springer Series in Computational Mathematics, Springer-Verlag (1993).

See Also

RxODE