Skip to contents

Simulate from a fitted sdmTMB model

Source: R/tmb-sim.R
simulate.sdmTMB.Rd

simulate.sdmTMB is an S3 method for producing a matrix of simulations from a fitted model. This is similar to lme4::simulate.merMod() and glmmTMB::simulate.glmmTMB(). It can be used with the DHARMa package among other uses.

Usage

# S3 method for sdmTMB
simulate(
  object,
  nsim = 1L,
  seed = sample.int(1e+06, 1L),
  type = c("mle-eb", "mle-mvn"),
  model = c(NA, 1, 2),
  re_form = NULL,
  mcmc_samples = NULL,
  ...
)

Arguments

object

sdmTMB model

nsim

Number of response lists to simulate. Defaults to 1.

seed

Random number seed

type

How parameters should be treated. "mle-eb": fixed effects are at their maximum likelihood (MLE) estimates and random effects are at their empirical Bayes (EB) estimates. "mle-mvn": fixed effects are at their MLEs but random effects are taken from a single approximate sample. This latter option is a suggested approach if these simulations will be used for goodness of fit testing (e.g., with the DHARMa package).

model

If a delta/hurdle model, which model to simulate from? NA = combined, 1 = first model, 2 = second mdoel.

re_form

NULL to specify a simulation conditional on fitted random effects (this only simulates observation error). ~0 or NA to simulate new random affects (smoothers, which internally are random effects, will not be simulated as new).

mcmc_samples

An optional matrix of MCMC samples. See extract_mcmc() in the sdmTMBextra package.

...

Extra arguments (not used)

Value

Returns a matrix; number of columns is nsim.

See also

sdmTMB_simulate()

Examples


# start with some data simulated from scratch:
set.seed(1)
predictor_dat <- data.frame(X = runif(300), Y = runif(300), a1 = rnorm(300))
mesh <- make_mesh(predictor_dat, xy_cols = c("X", "Y"), cutoff = 0.1)
dat <- sdmTMB_simulate(
  formula = ~ 1 + a1,
  data = predictor_dat,
  mesh = mesh,
  family = poisson(),
  range = 0.5,
  sigma_O = 0.2,
  seed = 42,
  B = c(0.2, -0.4) # B0 = intercept, B1 = a1 slope
)
fit <- sdmTMB(observed ~ 1 + a1, data = dat, family = poisson(), mesh = mesh)

# simulate from the model:
s1 <- simulate(fit, nsim = 300)
dim(s1)
#> [1] 300 300

# test whether fitted models are consistent with the observed number of zeros:
sum(s1 == 0)/length(s1)
#> [1] 0.3297667
sum(dat$observed == 0) / length(dat$observed)
#> [1] 0.3466667

# simulate with random effects sampled from their approximate posterior
s2 <- simulate(fit, nsim = 1, params = "mle-mvn")
# these may be useful in conjunction with DHARMa simulation-based residuals

# simulate with new random fields:
s3 <- simulate(fit, nsim = 1, re_form = ~ 0)