Simulate from a fitted TAM

Generates possible catch and survey observations from a fitted model, with options to include parameter uncertainty and new population histories. Use these draws to explore uncertainty or test whether simulated data resemble the observations. Results are tidy data frames with a sim = 1..n column.

Usage

sim_tam(
  fit,
  n = 10,
  par_uncertainty = c("joint", "fixed", "none"),
  redraw_random = FALSE,
  progress = TRUE,
  globals = NULL,
  seed = TRUE
)

Arguments

fit

A fitted TAM object returned by fit_tam().

n

Integer; number of simulations (default 10).

par_uncertainty

Character; one of "joint", "fixed", "none". Controls how the parameter list is sampled before each simulation (see Details).

redraw_random

Logical; if TRUE, re-draw random-effect fields from their process models over the whole modeled period on each run. If FALSE, keep random effects and simulate observations only.

progress

Logical; show a progress bar using progressr::with_progress() (default TRUE).

globals

Optional character vector of global objects for parallel workers (forwarded to furrr::furrr_options()).

seed

Seed passed to furrr::furrr_options() for reproducibility.

Value

A named list of data frames (e.g., catch, index, N, F, M, Z, ssb, …), each stacked across n simulations with a sim column.

Details

Missing observations, including zeros treated as missing by the likelihood, remain NA in the simulated catch and index tables. Population states and predictions are still calculated throughout the modeled period, including projection years. Missing observations are not estimated as extra parameters.

For mono() effects, negative dq or dmu_f values from Gaussian parameter draws are clipped to zero before generating predictions or observations. This projects the normal approximation onto the feasible increments and adds probability at zero; it is not a boundary-corrected posterior sampler. This can put probability at the asymptotic or Gaussian endpoint; normal uncertainty approximations near either boundary require particular care.

The simulation has two orthogonal controls:

• Parameter uncertainty via par_uncertainty:

  • "none" — use point estimates (û, θ̂).

  • "fixed" — sample fixed effects θ ~ MVN(sdrep$par.fixed, sdrep$cov.fixed).

  • "joint" — sample (random + fixed) jointly from the Laplace Gaussian approximation using the joint precision (sparse Cholesky). When no random effects are integrated, this is the same as "fixed".

• Random-effect handling via redraw_random:

  • FALSE — keep the sampled/fitted random effects and simulate observations only (an approximate predictive draw when par_uncertainty = "joint").

  • TRUE — generate new process fields for the random effects and re-simulate across the entire modeled history, including historical years. Recruitment follows rec_settings. Stock-recruit draws use parent SSB from the same simulated population; fixed early recruitment boundary states are retained. See recruitment_formulas. Initial-age states are retained: latent initial states have no proper distribution to redraw. Joint parameter draws can vary them through their fitted uncertainty. Fixed r0 is retained in process redraws. Temporal RW fields retain their supplied first process row and draw subsequent deviations conditionally; see nll_fun() for the N cohort construction. No distribution is imposed on the RW starting row. Fixed anchors may still vary through the chosen par_uncertainty draw. This is not a future-only forecast conditional on historical states. See tinyAM-model for the exact simulation and projection assumptions.

Parallel execution is supported via furrr::future_map(). Call future::plan() beforehand if you want parallel workers.

See also

fit_tam(), tidy_rep(), stack_nested()

Examples

fit <- fit_tam(
  cod_obs,
  years = 1983:2024,
  ages = 2:14,
  N_settings = list(process = "iid", init = "exp"),
  F_settings = list(process = "rw", mu_form = NULL),
  M_settings = list(process = "off", mu_supplied = ~ I(0.3)),
  catch_settings = list(sd_form = ~ 1),
  index_settings = list(sd_form = ~ 1, q_form = ~ q_block),
  silent = TRUE
)
# \donttest{
if (interactive()) {
  # Set-up parallel workers and fit model
  future::plan(future::multisession, workers = 4)

  # Draw fixed-effects uncertainty and redraw random process fields
  sims1 <- sim_tam(fit, n = 10, par_uncertainty = "fixed", redraw_random = TRUE)
  plot_trend(sims1$ssb, split = ~sim, line = list(width = 0.5))

  # Joint draw of random and fixed effects, keep random effects, simulate obs only
  sims2 <- sim_tam(fit, n = 10, par_uncertainty = "joint", redraw_random = FALSE)
  plot_trend(sims2$ssb, split = ~sim, line = list(width = 0.5))

  # Conditional/fast: point estimates, obs only
  sims3 <- sim_tam(fit, n = 10, par_uncertainty = "none", redraw_random = FALSE)
  plot_trend(sims3$index, y = ~obs, frame = ~age, split = ~sim, line = list(width = 0.5))
}
# }