Skip to contents

Simulate draws from a statistical model to return a data frame of estimates.

Usage

simulate_model(model, iterations = 1000, ...)

# Default S3 method
simulate_model(model, iterations = 1000, component = "all", ...)

Arguments

model

Statistical model (no Bayesian models).

iterations

The number of draws to simulate/bootstrap.

...

Arguments passed to insight::get_varcov(), e.g. to allow simulated draws to be based on heteroscedasticity consistent variance covariance matrices.

component

Should all parameters, parameters for the conditional model, for the zero-inflation part of the model, or the dispersion model be returned? Applies to models with zero-inflation and/or dispersion component. component may be one of "conditional", "zi", "zero-inflated", "dispersion" or "all" (default). May be abbreviated.

Value

A data frame.

Details

Technical Details

simulate_model() is a computationally faster alternative to bootstrap_model(). Simulated draws for coefficients are based on a multivariate normal distribution (MASS::mvrnorm()) with mean mu = coef(model) and variance Sigma = vcov(model).

Models with Zero-Inflation Component

For models from packages glmmTMB, pscl, GLMMadaptive and countreg, the component argument can be used to specify which parameters should be simulated. For all other models, parameters from the conditional component (fixed effects) are simulated. This may include smooth terms, but not random effects.

Model components

Possible values for the component argument depend on the model class. Following are valid options:

  • "all": returns all model components, applies to all models, but will only have an effect for models with more than just the conditional model component.

  • "conditional": only returns the conditional component, i.e. "fixed effects" terms from the model. Will only have an effect for models with more than just the conditional model component.

  • "smooth_terms": returns smooth terms, only applies to GAMs (or similar models that may contain smooth terms).

  • "zero_inflated" (or "zi"): returns the zero-inflation component.

  • "dispersion": returns the dispersion model component. This is common for models with zero-inflation or that can model the dispersion parameter.

  • "instruments": for instrumental-variable or some fixed effects regression, returns the instruments.

  • "nonlinear": for non-linear models (like models of class nlmerMod or nls), returns staring estimates for the nonlinear parameters.

  • "correlation": for models with correlation-component, like gls, the variables used to describe the correlation structure are returned.

Special models

Some model classes also allow rather uncommon options. These are:

  • mhurdle: "infrequent_purchase", "ip", and "auxiliary"

  • BGGM: "correlation" and "intercept"

  • BFBayesFactor, glmx: "extra"

  • averaging:"conditional" and "full"

  • mjoint: "survival"

  • mfx: "precision", "marginal"

  • betareg, DirichletRegModel: "precision"

  • mvord: "thresholds" and "correlation"

  • clm2: "scale"

  • selection: "selection", "outcome", and "auxiliary"

  • lavaan: One or more of "regression", "correlation", "loading", "variance", "defined", or "mean". Can also be "all" to include all components.

For models of class brmsfit (package brms), even more options are possible for the component argument, which are not all documented in detail here.

Examples

model <- lm(Sepal.Length ~ Species * Petal.Width + Petal.Length, data = iris)
head(simulate_model(model))
#>   (Intercept) Speciesversicolor Speciesvirginica Petal.Width Petal.Length
#> 1    3.436303        -0.7644174        -1.921772  0.81511491    0.8968868
#> 2    3.352833        -1.3050538        -2.057244  0.58321747    1.0497614
#> 3    3.724642        -0.7274907        -1.377786  0.65939731    0.7937357
#> 4    3.368219        -1.2476274        -2.599220  0.33148979    1.0075941
#> 5    3.383443        -1.5451483        -2.851981 -0.08413517    1.0958122
#> 6    3.591734        -1.2903453        -2.781200  0.32582409    0.9524036
#>   Speciesversicolor:Petal.Width Speciesvirginica:Petal.Width
#> 1                    -1.1498087                   -0.7648875
#> 2                    -1.0266680                   -0.8025617
#> 3                    -0.9864884                   -0.7338618
#> 4                    -0.6499459                   -0.2341848
#> 5                    -0.3796449                    0.1095252
#> 6                    -0.6257724                   -0.0805681
# \donttest{
if (require("glmmTMB", quietly = TRUE)) {
  model <- glmmTMB(
    count ~ spp + mined + (1 | site),
    ziformula = ~mined,
    family = poisson(),
    data = Salamanders
  )
  head(simulate_model(model))
  head(simulate_model(model, component = "zero_inflated"))
}
#>   (Intercept)   minedno
#> 1   1.0070148 -2.037876
#> 2   0.6898270 -1.774201
#> 3   0.7126957 -1.687886
#> 4   0.8666444 -1.847277
#> 5   0.9512715 -1.788466
#> 6   0.6543926 -1.515102
# }