Compute standard error for a variable, for all variables of a data frame, for joint random and fixed effects coefficients of (non-/linear) mixed models, the adjusted standard errors for generalized linear (mixed) models, or for intraclass correlation coefficients (ICC).
se(x, nsim = 100, type = c("fe", "re"))
x | (Numeric) vector, a data frame, an |
---|---|
nsim | Numeric, the number of simulations for calculating the
standard error for intraclass correlation coefficients, as
obtained by the |
type | Type of standard errors for generalized linear mixed models.
|
The standard error of x
.
For linear mixed models, and generalized linear mixed models with
type = "re"
, this function computes the standard errors
for joint (sums of) random and fixed effects coefficients (unlike
se.coef
, which returns the standard error
for fixed and random effects separately). Hence, se()
returns the appropriate standard errors for coef.merMod
.
For generalized linear models or generalized linear mixed models,
approximated standard errors, using the delta method for transformed
regression parameters are returned (Oehlert 1992). For generalized
linear mixed models, by default, the standard errors refer to the
fixed effects only. Use type = "re"
to compute standard errors
for joint random and fixed effects coefficients. However,
computation for the latter is not based on the delta method,
so standard errors from type = "re"
are on the scale of the
link-function (and not back transformed).
The standard error for the icc
is based on bootstrapping,
thus, the nsim
-argument is required. See 'Examples'.
se()
also returns the standard error of an estimate (regression
coefficient) and p-value, assuming a normal distribution to compute
the z-score from the p-value (formula in short: b / qnorm(p / 2)
).
See 'Examples'.
Computation of standard errors for coefficients of mixed models
is based on this code.
Standard errors for generalized linear (mixed) models, if
type = "re"
, are approximations based on the delta
method (Oehlert 1992).
A remark on standard errors:
“Standard error represents variation in the point estimate, but
confidence interval has usual Bayesian interpretation only with flat prior.”
(Gelman 2017)
Oehlert GW. 1992. A note on the delta method. American Statistician 46(1).
Gelman A 2017. How to interpret confidence intervals? http://andrewgelman.com/2017/03/04/interpret-confidence-intervals/
# compute standard error for vector se(rnorm(n = 100, mean = 3))#> [1] 0.1074613# compute standard error for each variable in a data frame data(efc) se(efc[, 1:3])#> c12hour e15relat e16sex #> 1.69162290 0.06942207 0.01565588# compute standard error for merMod-coefficients library(lme4) fit <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy) se(fit)#> $Subject #> (Intercept) Days #> 1 13.86649 2.7752 #> 2 13.86649 2.7752 #> 3 13.86649 2.7752 #> 4 13.86649 2.7752 #> 5 13.86649 2.7752 #> 6 13.86649 2.7752 #> 7 13.86649 2.7752 #> 8 13.86649 2.7752 #> 9 13.86649 2.7752 #> 10 13.86649 2.7752 #> 11 13.86649 2.7752 #> 12 13.86649 2.7752 #> 13 13.86649 2.7752 #> 14 13.86649 2.7752 #> 15 13.86649 2.7752 #> 16 13.86649 2.7752 #> 17 13.86649 2.7752 #> 18 13.86649 2.7752 #># compute odds-ratio adjusted standard errors, based on delta method # with first-order Taylor approximation. data(efc) efc$services <- sjmisc::dicho(efc$tot_sc_e, dich.by = 0) fit <- glm(services ~ neg_c_7 + c161sex + e42dep, data = efc, family = binomial(link = "logit")) se(fit)#> term estimate std.error #> 1 (Intercept) 0.5946382 0.22373756 #> 2 neg_c_7 1.0427929 0.02035406 #> 3 c161sex 0.8033751 0.12974572 #> 4 e42dep 1.2379767 0.09724651# compute odds-ratio adjusted standard errors for generalized # linear mixed model, also based on delta method library(lme4) library(sjmisc) # create binary response sleepstudy$Reaction.dicho <- dicho(sleepstudy$Reaction, dich.by = "median") fit <- glmer(Reaction.dicho ~ Days + (Days | Subject), data = sleepstudy, family = binomial("logit")) se(fit)#> term estimate std.error #> 1 (Intercept) 0.02201714 0.02582094 #> 2 Days 2.43719045 0.57200380# compute standard error from regression coefficient and p-value se(list(estimate = .3, p.value = .002))#> [1] 0.09708008# NOT RUN { # compute standard error of ICC for the linear mixed model icc(fit) se(icc(fit)) # the standard error for the ICC can be computed manually in this way, # taking the fitted model example from above library(dplyr) library(purrr) dummy <- sleepstudy %>% # generate 100 bootstrap replicates of dataset bootstrap(100) %>% # run mixed effects regression on each bootstrap replicate # and compute ICC for each "bootstrapped" regression mutate( models = map(strap, ~lmer(Reaction ~ Days + (Days | Subject), data = .x)), icc = map_dbl(models, ~icc(.x)) ) # now compute SE and p-values for the bootstrapped ICC, values # may differ from above example due to random seed boot_se(dummy, icc) boot_p(dummy, icc) # }