Simulate from covariate/metadata in the absence of a real data set (EXPERIMENTAL)
Source:R/utils.R
simulate_new.RdSee vignette("sim", package = "glmmTMB") for more details and examples,
and vignette("covstruct", package = "glmmTMB")
for more information on the parameterization of different covariance structures.
Usage
simulate_new(
object,
nsim = 1,
seed = NULL,
family = gaussian,
newdata,
newparams = NULL,
...,
return_val = c("sim", "pars", "object")
)Arguments
- object
a one-sided model formula (e.g.
~ a + b + c(peculiar naming is for consistency with the generic function, which typically takes a fitted model object)- nsim
number of simulations
- seed
random-number seed
- family
a family function, a character string naming a family function, or the result of a call to a family function (variance/link function) information. See
familyfor a generic discussion of families orfamily_glmmTMBfor details ofglmmTMB-specific families.- newdata
a data frame containing all variables listed in the formula, including the response variable (which needs to fall within the domain of the conditional distribution, and should probably not be all zeros, but whose value is otherwise irrelevant)
- newparams
a list of parameters containing sub-vectors (
beta,betazi,betadisp,theta, etc.) to be used in the model. Ifbis specified in this list, then the conditional modes/BLUPs will be set to these values; otherwise they will be drawn from the appropriate Normal distribution. Seevignette("covstruct", package = "glmmTMB")for details on the parameterizations used for various random-effects models (i.e.,theta).- ...
other arguments to
glmmTMB(e.g.family)- return_val
what information to return: "sim" (the default) returns a list of vectors of simulated outcomes; "pars" returns the default parameter vector (this variant does not require
newparamsto be specified, and is useful for figuring out the appropriate dimensions of the different parameter vectors); "object" returns a fakeglmmTMBobject (useful, e.g., for retrieving the Z matrix (getME(simulate_new(...), "Z")) or covariance matrices (VarCorr(simulate_new(...))) implied by a particular set of input data and parameter values)
Details
Use the weights argument to set the size/number of trials per observation for binomial-type models; the default is 1 for every observation (i.e., Bernoulli trials)
See also
glmmTMB, family_glmmTMB (for conditional distribution parameterizations [betadisp]), put_cor (for correlation matrix parameterizations)
Examples
## use Salamanders data for observational design and covariate values
## parameters used here are sensible, but do not fit the original data
params <- list(beta = c(2, 1),
betazi = c(-0.5, 0.5), ## logit-linear model for zi
betadisp = log(2), ## log(NB dispersion)
theta = log(1)) ## log(among-site SD)
sim_count <- simulate_new(~ mined + (1|site),
newdata = Salamanders,
zi = ~ mined,
family = nbinom2,
seed = 101,
newparams = params
)
## simulate_new with return="sim" always returns a list of response vectors
Salamanders$sim_count <- sim_count[[1]]
summary(glmmTMB(sim_count ~ mined + (1|site), data=Salamanders, ziformula=~mined, family=nbinom2))
#> Family: nbinom2 ( log )
#> Formula: sim_count ~ mined + (1 | site)
#> Zero inflation: ~mined
#> Data: Salamanders
#>
#> AIC BIC logLik -2*log(L) df.resid
#> 3324.1 3350.9 -1656.0 3312.1 638
#>
#> Random effects:
#>
#> Conditional model:
#> Groups Name Variance Std.Dev.
#> site (Intercept) 0.669 0.8179
#> Number of obs: 644, groups: site, 23
#>
#> Dispersion parameter for nbinom2 family (): 1.69
#>
#> Conditional model:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.6481 0.2572 6.407 1.48e-10 ***
#> minedno 1.3861 0.3550 3.904 9.44e-05 ***
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Zero-inflation model:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) -0.7164 0.1510 -4.745 2.08e-06 ***
#> minedno 0.6331 0.1864 3.397 0.000682 ***
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
## return a glmmTMB object
sim_obj <- simulate_new(~ mined + (1|site),
return_val = "object",
newdata = Salamanders,
zi = ~ mined,
family = nbinom2,
newparams = params)
## simulate Gaussian data, multivariate random effect
data("sleepstudy", package = "lme4")
sim_obj <- simulate_new(~ 1 + (1|Subject) + ar1(0 + factor(Days)|Subject),
return_val = "pars",
newdata = sleepstudy,
family = gaussian,
newparams = list(beta = c(280, 1),
betad = log(2), ## log(residual std err)
theta = c(log(2), ## log(SD(subject))
log(2), ## log(SD(slope))
## AR1 correlation = 0.2
put_cor(0.2, input_val = "vec"))
)
)
#> Warning: unmatched parameter names: betad
#> Warning: length mismatch in component beta (1 != 2); not setting
#> Warning: unmatched parameter names: betad
#> Warning: length mismatch in component beta (1 != 2); not setting
#> Warning: length mismatch in component betad (0 != 1); not setting