07-solving.Rmd

# Solving and solving options {#solve}

In general, ODEs are solved using a combination of:

- A compiled model specification from `rxode2()`, specified with
  `object=`

- Input parameters, specified with `params=` (and could be blank)

- Input data or event table, specified with `events=`

- Initial conditions, specified by `inits=` (and possibly in the model
  itself by `state(0)=`)


The solving options are given in the sections below:


## General Solving Options

### object
`object` is a either a rxode2 family of objects, or a file-name
with a rxode2 model specification, or a string with a rxode2
model specification.

### params
`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
`events` an `eventTable` object describing the input
(e.g., doses) to the dynamic system and observation sampling
time points (see [eventTable()]);

### inits
`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);

### sigdig
`sigdig` Specifies the "significant digits" that the ode
solving requests.  When specified this controls the relative and
absolute tolerances of the ODE solvers.  By default the tolerance
is `0.5*10^(-sigdig-2)` for regular ODEs. For the
sensitivity equations the default is `0.5*10\^(-sigdig-1.5)`
(sensitivity changes only applicable for liblsoda).  This also
## lsoda/dop solving options

controls the `atol`/`rtol` of the steady state solutions. The
`ssAtol`/`ssRtol` is `0.5*10\^(-sigdig)` and for the sensitivities
`0.5*10\^(-sigdig+0.625)`.  By default
this is unspecified (`NULL`) and uses the standard `atol`/`rtol`.

### atol
`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
`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.

### atolSens
`atolSens` Sensitivity atol, can be different than atol with
liblsoda.  This allows a less accurate solve for gradients (if desired)

### rtolSens
`rtolSens` Sensitivity rtol, can be different than rtol with
liblsoda.  This allows a less accurate solve for gradients (if desired)

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

### hmin
`hmin` The minimum absolute step size allowed. The default
value is 0.

### hmax
`hmax` The maximum absolute step size allowed.  When
`hmax=NA` (default), uses the average difference +
hmaxSd*sd in times and sampling events. The `hmaxSd` is a user
specified parameter and which defaults to zero.  When
`hmax=NULL` rxode2 uses the maximum difference in times in
your sampling and events.  The value 0 is equivalent to infinite
maximum absolute step size.

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

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

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

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

### mxhnil
`mxhnil` maximum number of messages printed (per problem)
warning that `T + H = T` on a step (`H` = step size).  This must
be positive to result in a non-default value.  The default
value is 0 (or infinite).

### hmxi
`hmxi` inverse of the maximum absolute value of `H` to are used.
hmxi = 0.0 is allowed and corresponds to an infinite `hmax1
(default).  `hmin` and `hmxi` may be changed at any time, but will
not take effect until the next change of `H` is considered.
This option is only considered with `method="liblsoda"`.

### istateReset
`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.

## Inductive Linerization Options

### indLinMatExpType
`indLinMatExpType` This is them matrix exponential type that
is use for rxode2.  Currently the following are supported:

* `Al-Mohy` Uses the exponential matrix method of Al-Mohy Higham (2009)

* `arma` Use the exponential matrix from RcppArmadillo

* `expokit` Use the exponential matrix from Roger B. Sidje (1998)


### indLinMatExpOrder
`indLinMatExpOrder` an integer, the order of approximation to
be used, for the `Al-Mohy` and `expokit` values.
The best value for this depends on machine precision (and
slightly on the matrix). We use `6` as a default.

### indLinPhiTol
`indLinPhiTol` the requested accuracy tolerance on
exponential matrix.

### indLinPhiM
`indLinPhiM`  the maximum size for the Krylov basis

## Steady State Solving Options

### minSS
`minSS` Minimum number of iterations for a steady-state dose

### maxSS
`maxSS` Maximum number of iterations for a steady-state dose

### strictSS
`strictSS` Boolean indicating if a strict steady-state is
required. If a strict steady-state is (`TRUE`) required
then at least `minSS` doses are administered and the
total number of steady states doses will continue until
`maxSS` is reached, or `atol` and `rtol` for
every compartment have been reached.  However, if ODE solving
problems occur after the `minSS` has been reached the
whole subject is considered an invalid solve. If
`strictSS` is `FALSE` then as long as `minSS`
has been reached the last good solve before ODE solving
problems occur is considered the steady state, even though
either `atol`, `rtol` or `maxSS` have not
been achieved.

### infSSstep
`infSSstep` Step size for determining if a constant infusion
has reached steady state.  By default this is large value,
12.

### ssAtol
`ssAtol` Steady state atol convergence factor.  Can be
a vector based on each state.

### ssRtol
`ssRtol` Steady state rtol convergence factor.  Can be a
vector based on each state.

### ssAtolSens
`ssAtolSens` Sensitivity absolute tolerance (atol) for
calculating if steady state has been achieved for sensitivity compartments.

### ssRtolSens
`ssRtolSens` Sensitivity relative tolerance (rtol) for
calculating if steady state has been achieved for sensitivity compartments.

## rxode2 numeric stability options

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

### stateTrim
`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 large negative
amounts/concentrations.  If you want to trim between a range
say `c(0, 2000000)` you may specify 2 values with a lower and
upper range to make sure all state values are in the
reasonable range.

### safeZero
`safeZero` Use safe zero divide and log routines.  By default
this is turned on but you may turn it off if you wish.

### sumType
`sumType` Sum type to use for `sum()` in
rxode2 code blocks.

`pairwise` uses the pairwise sum (fast, default)

`fsum` uses the PreciseSum package's fsum function (most accurate)

`kahan` uses Kahan correction

`neumaier` uses Neumaier correction

`c` uses no correction: default/native summing

### prodType
`prodType` Product to use for `prod()` in rxode2 blocks

`long double` converts to long double, performs the
multiplication and then converts back.

`double` uses the standard double scale for multiplication.

### maxwhile
`maxwhile` represents the maximum times a while loop is
evaluated before exiting.  By default this is 100000

## Linear compartment model sensitivity options

### sensType
`sensType` Sensitivity type for `linCmt()` model:

`advan` Use the direct advan solutions

`autodiff` Use the autodiff advan solutions

`forward` Use forward difference solutions

`central` Use central differences

### linDiff
`linDiff` This gives the linear difference amount for all the
types of linear compartment model parameters where sensitivities
are not calculated. The named components of this numeric vector are:

* `"lag"` Central compartment lag
* `"f"` Central compartment bioavailability
* `"rate"` Central compartment modeled rate
* `"dur"` Central compartment modeled duration
* `"lag2"` Depot compartment lag
* `"f2"` Depot compartment bioavailability
* `"rate2"` Depot compartment modeled rate
* `"dur2"` Depot compartment modeled duration

### linDiffCentral
`linDiffCentral` This gives the which parameters use central
differences for the linear compartment model parameters.  The
are the same components as `linDiff`

## Covariate Solving Options

### iCov
`iCov` A data frame of individual non-time varying covariates
to combine with the `events` dataset by merge.

### covsInterpolation
`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, which interpolates the covariate
by solving the line between the observed covariates and extrapolating the new
covariate value.

* `"constant"` -- Last observation carried forward (the default).

* `"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
`addCov` A boolean indicating if covariates should be added
to the output matrix or data frame. By default this is
disabled.

## Simulation options

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

### nsim
`nsim` represents the number of simulations.  For rxode2, if
you supply single subject event tables (created with
`[eventTable()]`)

### thetaMat
`thetaMat` Named theta matrix.

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

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

### thetaDf
`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
`thetaIsChol` Indicates if the `theta` supplied is a
Cholesky decomposed matrix instead of the traditional
symmetric matrix.

### nStud
`nStud` Number virtual studies to characterize uncertainty in estimated
parameters.

### omega
`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.  When `omega` is `NA` and you are using it with a
`rxode2` ui model, the between subject variability described by
the `omega` matrix are set to zero.

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

### omegaSeparation
`omegaSeparation` Omega separation strategy

Tells the type of separation strategy when
simulating covariance with parameter uncertainty with standard
deviations modeled in the `thetaMat` matrix.

* `"lkj"` simulates the correlation matrix from the
`rLKJ1` matrix with the distribution parameter `eta`
equal to the degrees of freedom `nu` by `(nu-1)/2`

*  `"separation"` simulates from the identity inverse Wishart
covariance matrix with `nu` degrees of freedom.  This is then
converted to a covariance matrix and augmented with the modeled
standard deviations.  While computationally more complex than the
`"lkj"` prior, it performs better when the covariance matrix
size is greater or equal to 10

* `"auto"` chooses `"lkj"` when the dimension of the
matrix is less than 10 and `"separation"` when greater
than equal to 10.

### omegaXform
`omegaXform` When taking `omega` values from the `thetaMat`
simulations (using the separation strategy for covariance
simulation), how should the `thetaMat` values be turned int
standard deviation values:

- `identity` This is when standard deviation values are
directly modeled by the `params` and `thetaMat` matrix

- `variance` This is when the `params` and `thetaMat`
simulates the variance that are directly modeled by the
`thetaMat` matrix

- `log` This is when the `params` and `thetaMat`
simulates `log(sd)`

- `nlmixrSqrt` This is when the `params` and
`thetaMat` simulates the inverse cholesky decomposed matrix
with the `x\^2` modeled along the diagonal.  This only works
with a diagonal matrix.

- `nlmixrLog` This is when the `params` and
`thetaMat` simulates the inverse cholesky decomposed matrix
with the `exp(x\^2)` along the diagonal.  This only works
with a diagonal matrix.

- `nlmixrIdentity` This is when the `params` and
`thetaMat` simulates the inverse cholesky decomposed matrix.
This only works with a diagonal matrix.

### omegaLower
`omegaLower` Lower bounds for simulated ETAs (by default -Inf)

### omegaUpper
`omegaUpper` Upper bounds for simulated ETAs (by default Inf)

### omegaDf
`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.

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

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

### sigma
`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. When `sigma` is `NA` and you are using it with a
`rxode2` ui model, the unexplained variability described by
the `sigma` matrix are set to zero.

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

### sigmaUpper
`sigmaUpper` Upper bounds for simulated unexplained variability (by default Inf)

### sigmaXform
`sigmaXform` When taking `sigma` values from the `thetaMat`
simulations (using the separation strategy for covariance
simulation), how should the `thetaMat` values be turned int
standard deviation values:

- `identity` This is when standard deviation values are
directly modeled by the `params` and `thetaMat` matrix

- `variance` This is when the `params` and `thetaMat`
simulates the variance that are directly modeled by the
`thetaMat` matrix

- `log` This is when the `params` and `thetaMat`
simulates `log(sd)`

- `nlmixrSqrt` This is when the `params` and
`thetaMat` simulates the inverse cholesky decomposed matrix
with the `x\^2` modeled along the diagonal.  This only works
with a diagonal matrix.

- `nlmixrLog` This is when the `params` and
`thetaMat` simulates the inverse cholesky decomposed matrix
with the `exp(x\^2)` along the diagonal.  This only works
with a diagonal matrix.

- `nlmixrIdentity` This is when the `params` and
`thetaMat` simulates the inverse cholesky decomposed matrix.
This only works with a diagonal matrix.


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

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

### sigmaSeparation
`sigmaSeparation` separation strategy for sigma;

Tells the type of separation strategy when
simulating covariance with parameter uncertainty with standard
deviations modeled in the `thetaMat` matrix.

* `"lkj"` simulates the correlation matrix from the
`rLKJ1` matrix with the distribution parameter `eta`
equal to the degrees of freedom `nu` by `(nu-1)/2`

*  `"separation"` simulates from the identity inverse Wishart
covariance matrix with `nu` degrees of freedom.  This is then
converted to a covariance matrix and augmented with the modeled
standard deviations.  While computationally more complex than the
`"lkj"` prior, it performs better when the covariance matrix
size is greater or equal to 10

*  `"auto"` chooses `"lkj"` when the dimension of the
matrix is less than 10 and `"separation"` when greater
than equal to 10.

### dfObs
`dfObs` Degrees of freedom to sample the unexplained variability matrix from the
inverse Wishart distribution (scaled) or scaled inverse chi squared distribution.

### resample
`resample` A character vector of model variables to resample
from the input dataset; This sampling is done with replacement.
When `NULL` or `FALSE` no resampling is done.  When
`TRUE` resampling is done on all covariates in the input
dataset

### resampleID
`resampleID` boolean representing if the resampling should be
done on an individual basis `TRUE` (ie. a whole patient is
selected) or each covariate is resampled independent of the
subject identifier `FALSE`.  When `resampleID=TRUE`
correlations of parameters are retained, where as when
`resampleID=FALSE` ignores patient covariate correaltions.
Hence the default is `resampleID=TRUE`.

## rxode2 output options

### returnType
`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 rxode2,  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`

* `"data.table"` -- returns a `data.table`; The `data.table` is
created by reference (ie `setDt()`), which should be fast.

* `"tbl"` or `"tibble"` returns a tibble format.

### addDosing
`addDosing` Boolean indicating if the solve should add rxode2
EVID and related columns.  This will also include dosing
information and estimates at the doses.  Be default, rxode2
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
`rxode2` EVID events are returned. When `addDosing` is `TRUE`
add the event information in NONMEM-style format; If
`subsetNonmem=FALSE` rxode2 will also include 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.

### keep
`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.

### drop
`drop` Columns to drop from the output

### idFactor
`idFactor` This boolean indicates if original ID values
should be maintained. This changes the default sequentially
ordered ID to a factor with the original ID values in the
original dataset.  By default this is enabled.

### subsetNonmem
`subsetNonmem` subset to NONMEM compatible EVIDs only.  By
default `TRUE`.

### scale
`scale` a numeric named vector with scaling for ode
parameters of the system.  The names must correspond 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.

### amountUnits
`amountUnits` This supplies the dose units of a data frame
supplied instead of an event table.  This is for importing the
data as an rxode2 event table.

### timeUnits
`timeUnits` This supplies the time units of a data frame
supplied instead of an event table.  This is for importing the
data as an rxode2 event table.

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

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

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

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

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

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

### warnIdSort
`warnIdSort` Warn if the ID is not present and rxode2 assumes
the order of the parameters/iCov are the same as the order of
the parameters in the input dataset.

### warnDrop
`warnDrop` Warn if column(s) were supposed to be dropped, but
were not present.

## Internal rxode2 options

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

### simVariability
`simVariability` determines if the variability is simulated.
When `NA` (default) this is determined by the solver.

### ...
`...` 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.

### a
`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
`b` when using `solve()`, this is equivalent to the
`params` argument.  If you specify `params` as a
named argument, this overwrites the output

### updateObject
`updateObject` This is an internally used flag to update the
rxode2 solved object (when supplying an rxode2 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.

## Parallel/Threaded Solve

### cores
`cores` Number of cores used in parallel ODE solving.  This
is equivalent to calling [setRxThreads()]

### nCoresRV
`nCoresRV` Number of cores used for the simulation of the
sigma variables.  By default this is 1. 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.

### nLlikAlloc
`nLlikAlloc` The number of log likelihood endpoints that are
used in the model.  This allows independent log likelihood per
endpoint in focei for nlmixr2.  It likely shouldn't be set,
though it won't hurt anything if you do (just may take up more
memory for larger allocations).


### useStdPow
`useStdPow` This uses C's `pow` for exponentiation instead of
R's `R_pow` or `R_pow_di`.  By default this is `FALSE`

### ss2cancelAllPending
`ss2cancelAllPending` When `TRUE` the `SS=2` event type
cancels all pending doses like `SS=1`.  When `FALSE` the pending
doses not canceled with `SS=2` (the infusions started before
`SS=2` occurred are canceled, though).

### addlKeepsCov
`addlKeepsCov` This determines if the additional dosing items
repeats the dose only (`FALSE`) or keeps the covariates at the
record of the dose (`TRUE`)

### addlDropSs
`addlDropSs` When there are steady state doses with an `addl`
specification the steady state flag is dropped with repeated
doses (when `TRUE`) or retained (when `FALSE`)

### ssAtDoseTime
`ssAtDoseTime` Boolean that when `TRUE` back calculates the
steady concentration at the actual time of dose, otherwise when
`FALSE` the doses are shifted

### naTimeHandle
`naTimeHandle` Determines what time of handling happens when
the time becomes `NA`: current options are:

- `ignore` this ignores the `NA` time input and passes it through.

- `warn` (default) this will produce a warning at the end of the
solve, but continues solving passing through the `NA` time

- `error` this will stop this solve if this is not a parallel
solved ODE (otherwise stopping can crash R)