Simulate responses from an exposure-response model

Generates simulated response datasets from a fitted erglm model, propagating uncertainty in the parameter estimates. Useful for simulation-based confidence bands, predictive checks, or bootstrapping downstream analyses. Implements the standard stats::simulate() generic, so it is called as simulate(object, ...) rather than through an erglm-specific function name.

Usage

# S3 method for class 'erglm_model'
simulate(object, nsim = 1, seed = NULL, ...)

Arguments

object

An erglm model, as returned by erglm_model()

nsim

Number of replicates. Must be a single positive whole number. Defaults to 1.

seed

Used to set the RNG seed. If NULL, a random seed is chosen and reported.

...

Ignored

Value

A data frame with one row per observation per simulated replicate, containing:

• dat_id, sim_id: identifiers for the original observation and the simulation replicate

• mu: the expected response (response scale) at the sampled parameter vector

• val: the simulated response value (mu plus family-appropriate noise)

• one coef_* column per model coefficient (e.g. coef_`(Intercept)`, coef_aucss), giving the sampled parameter values used for that replicate – prefixed to avoid colliding with predictor columns of the same name

• the model's predictor columns (not including the response)

Details

Samples new parameter values from the multivariate normal distribution implied by the model's variance-covariance matrix (via mvtnorm::rmvnorm()), evaluates the expected response at each sampled parameter vector using erglm_fun(), then draws a simulated response at each prediction using family-appropriate residual noise (the same .erglm_draw_response() mechanism used elsewhere in the package: Bernoulli draws for binomial, Poisson draws for poisson, normal draws for gaussian, gamma draws for Gamma). The dispersion parameter used for that noise is a single point estimate (summary(object)$dispersion), not resampled per replicate. Other glm() families are not currently supported and will raise an informative error.

For a VPC-style plot comparing observed and simulated response rates, see the companion erplots package's er_vpc() mini-grammar – specifically erplots::er_vpc_add_simulated(), which can build its simulated replicates directly from a fitted model (its model = argument) via the same .erglm_draw_response() noise mechanism this method uses, without needing to call simulate() yourself.

Examples

mod <- erglm_model(ae1 ~ aucss + sex, erglm_data, family = binomial())
sim <- simulate(mod, nsim = 5, seed = 963)
head(sim)
#>   dat_id sim_id        mu val coef_(Intercept)  coef_aucss coef_sexMale   aucss
#> 1      1      1 0.8961156   1        -1.768383 0.006181968   -0.2378474  673.09
#> 2      2      1 0.9999998   1        -1.768383 0.006181968   -0.2378474 2806.12
#> 3      3      1 0.1457435   0        -1.768383 0.006181968   -0.2378474    0.00
#> 4      4      1 0.9957585   1        -1.768383 0.006181968   -0.2378474 1169.04
#> 5      5      1 0.5808257   1        -1.768383 0.006181968   -0.2378474  377.29
#> 6      6      1 0.5630660   1        -1.768383 0.006181968   -0.2378474  327.08
#>      sex
#> 1   Male
#> 2 Female
#> 3 Female
#> 4 Female
#> 5   Male
#> 6 Female