Simulation: Simulation

In RTMB: 'R' Bindings for 'TMB'

Simulation

Description

An RTMB objective function can be run in 'simulation mode' where standard likelihood evaluation is replaced by corresponding random number generation. This facilitates automatic simulation under some restrictions, notably regarding the order of likelihood accumulation in the user template - see details. Simulations can be obtained directly from the model object by obj$simulate() or used indirectly via checkConsistency. In both cases a simulation is carried out from the full model i.e. both random effects and observations are simulated. The OBS function is used to mark which data items are observations - see examples.

Usage

simref(n)

## S3 replacement method for class 'simref'
dim(x) <- value

## S3 method for class 'simref'
length(x)

## S3 method for class 'simref'
dim(x)

## S3 method for class 'simref'
is.array(x)

## S3 method for class 'simref'
is.matrix(x)

## S3 method for class 'simref'
as.array(x, ...)

## S3 method for class 'simref'
is.na(x)

## S3 method for class 'simref'
x[...]

## S3 replacement method for class 'simref'
x[...] <- value

## S3 method for class 'simref'
Ops(e1, e2)

## S3 method for class 'simref'
Math(x, ...)

## S3 method for class 'simref'
t(x)

## S3 method for class 'simref'
diff(x, lag = 1L, differences = 1L, ...)

## S3 method for class 'simref'
Summary(..., na.rm = FALSE)

Arguments

n	Length
x	Object of class 'simref'
value	Replacement (numeric)
...	Extra arguments
e1	First argument
e2	Second argument
lag	As diff
differences	As diff
na.rm	Ignored

Details

In simulation mode all log density evaluation, involving either random effects or observations, is interpreted as probability assignment.

direct vs indirect Assignments can be 'direct' as for example

dnorm(u, log=TRUE) ## u ~ N(0, 1)

or 'indirect' as in

dnorm(2*(u+1), log=TRUE) ## u ~ N(-1, .25)

Indirect assignment works for a limited set of easily invertible functions - see methods(class="simref").

Simulation order Note that probability assignments are sequential: All information required to draw a new variable must already be simulated. It follows that, for the simulation to work, one cannot assume likelihood accumulation is commutative!

Vectorized assignment implicitly occurs elementwise from left to right. For example the assignment

dnorm(diff(u), log=TRUE)

is not valid without a prior assignment of u[1], e.g.

dnorm(u[1], log=TRUE)

Supported distributions Assignment must use supported density functions. I.e.

dpois(N, exp(u), log=TRUE)

cannot be replaced by

N * u - exp(u)

The latter will have no effect in simulation mode (the simulation will be NA).

Return value Note that when in simulation mode, the density functions all return zero. The actual simulation is written to the input argument by reference. This is very unlike standard R semantics.

Value

An object with write access to store the simulation.

Functions

• simref(): Construct simref

• dim(simref) <- value: Equivalent of dim<-

• length(simref): Equivalent of length

• dim(simref): Equivalent of dim

• is.array(simref): Equivalent of is.array

• is.matrix(simref): Equivalent of is.matrix

• as.array(simref): Equivalent of as.array

• is.na(simref): Equivalent of is.na

• [: Equivalent of [

• `[`(simref) <- value: Equivalent of [<-

• Ops(simref): Equivalent of Ops

• Math(simref): Equivalent of Math

• t(simref): Equivalent of t

• diff(simref): Equivalent of diff

• Summary(simref): Summary operations are not invertible and will throw an error.

Examples

## Basic example simulating response marked by OBS()
func <- function(par) {
  getAll(par, tmbdata)
  y <- OBS(y)
  ans <- -sum(dbinom(y, prob = plogis(a+b*x), size = 1, log = TRUE))
}
set.seed(101)
tmbdata <- list(x = seq(-3, 3, length = 25))
tmbdata$y <- rbinom(25, size = 1, prob = plogis(2 - 0.1*tmbdata$x))
obj <- MakeADFun(func, list(a = 0, b = 0), silent=TRUE)
with(obj, nlminb(par, fn, gr))
obj$simulate()
## Basic example simulating random effects
func <- function(p) {
  u <- p$u
  ans <- -dnorm(u[1], log=TRUE) ## u[1] ~ N(0,1)
  ans <- ans - sum(dnorm(diff(u), log=TRUE)) ## u[i]-u[i-1] ~ N(0,1)
}
obj <- MakeADFun(func, list(u=numeric(20)), random="u")
obj$simulate()
## Demonstrate how a 'simref' object works
s <- simref(4)
s2 <- 2 * s[1:2] + 1
s2[] <- 7
s ## 3 3 NA NA

RTMB documentation built on Nov. 5, 2025, 6:19 p.m.