Methods have been written that allow glmmTMB objects to be used with
several downstream packages that enable different forms of inference.
For some methods (Anova and emmeans, but not effects at present),
set the component argument
to "cond" (conditional, the default), "zi" (zero-inflation) or "disp" (dispersion) in order to produce results
for the corresponding part of a glmmTMB model.
Support for emmeans also allows additional options
component = "response" (response means taking both the cond and
zi components into account), and component = "cmean" (mean of the
[possibly truncated] conditional distribution).
In particular,
car::Anova constructs type-II and type-III Anova tables
for the fixed effect parameters of any component
the emmeans package computes estimated marginal means (previously known as least-squares means)
for the fixed effects of any component, or predictions with type = "response" or
type = "component". Note: In hurdle models,
component = "cmean" produces means
of the truncated conditional distribution, while
component = "cond", type = "response" produces means of the untruncated
conditional distribution.
the effects package computes graphical tabular effect displays
(only for the fixed effects of the conditional component)
Anova.glmmTMB(
mod,
type = c("II", "III", 2, 3),
test.statistic = c("Chisq", "F"),
component = "cond",
vcov. = vcov(mod)[[component]],
singular.ok,
include.rankdef.cols = FALSE,
ddf = c("asymptotic", "kenward-roger", "satterthwaite"),
...
)Effect.glmmTMB(focal.predictors, mod, ...)
a glmmTMB model
type of test, "II", "III", 2, or 3. Roman numerals are equivalent to the corresponding Arabic numerals. See Anova for details.
"Chisq" (default; a Wald chi-squared test) or "F"
(only available together with ddf != "asymptotic", see ddf below). An explicit
test.statistic = "Chisq" combined with ddf != "asymptotic" is an error, since
Kenward-Roger/Satterthwaite always produce an F table
which component of the model to test/analyze ("cond", "zi", or "disp") or, in emmeans only, "response" or "cmean" as described in Details.
variance-covariance matrix (usually extracted automatically); not
currently supported together with ddf != "asymptotic"
OK to do ANOVA with singular models (unused) ?
include all columns of a rank-deficient model matrix?
denominator degrees-of-freedom calculation, as in summary.glmmTMB
and anova.glmmTMB. The default "asymptotic" gives the classical Wald
chi-squared table; "kenward-roger" or "satterthwaite" instead give an F-ratio
table, with each term's denominator df computed via the Kenward-Roger or Satterthwaite
approximation (see dof_KR, dof_satt). "kenward-roger"
additionally requires a family with an estimated dispersion parameter, and throws an error
for families such as binomial or poisson that lack one ("satterthwaite"
has no such restriction). Not currently supported together with a user-supplied vcov.,
component != "cond", or models with aliased/rank-deficient or map-fixed
conditional coefficients.
Additional parameters that may be supported by the method.
a character vector of one or more predictors in the model in any order.
For models fitted with the ordinal family, emmeans() accepts
the same mode and rescale arguments as for MASS::polr fits
(ordinal::clm's method additionally offers "scale", which
does not apply here; see the clm/polr entries in
vignette("models", package = "emmeans")): "latent" (the default; means on the latent scale,
centered on the average threshold, with no back-transformation, so
type = "response" has no effect), "linear.predictor"
(\(\theta_j - x'\beta\) for each threshold \(j\), with a grid
variable cut), "cum.prob" (cumulative probabilities
\(P(Y \le j)\)), "exc.prob" (exceedance probabilities
\(P(Y > j)\)), "prob" (probabilities of each response
category, indexed by the response variable) and "mean.class"
(the expected category index). Standard errors combine the
fixed-effect covariance with the delta-method covariance of the
thresholds (as in summary()), and denominator degrees of freedom
are always asymptotic (a ddf request for
"satterthwaite" or "kenward-roger" warns and is
ignored). In "latent" mode, rescale = c(a, b) reports
\(a + b \mu\) in place of the latent mean \(\mu\), as for
MASS::polr. A user-supplied vcov. must be the joint
covariance matrix of the conditional fixed effects (in the order of
fixef(.)$cond, omitting coefficients dropped for rank
deficiency) followed by the \(K-1\) thresholds on the threshold
scale, as in the thresholds element of summary(.).
For models with random effects, the ddf argument to emmeans()
(default taken from getOption("glmmTMB.df", "asymptotic")) additionally accepts
"satterthwaite" and "kenward-roger" (see dof_KR and
dof_satt for the underlying calculations), matching the same argument
to summary.glmmTMB and anova.glmmTMB. ddf = "kenward-roger"
requires a model fitted with REML = TRUE: for an ML fit (glmmTMB's default),
it throws an error rather than silently substituting another method; it also requires a
family with an estimated dispersion parameter, and throws an error for families such as
binomial or poisson that lack one. For families other than gaussian,
"kenward-roger" and "satterthwaite" are allowed but emit a warning, because
their performance (and theoretical justification) for GLMMs is poorly understood.
For Gaussian models without random effects, emmeans() defaults
to the residual degrees of freedom for a plain fit, i.e. one with
dispformula = ~1 and an estimated dispersion parameter. Any other such
fit defaults to "asymptotic" (infinite df): a non-trivial
dispformula, dispformula = ~0, or a dispersion parameter held
fixed via the map argument to glmmTMB. In the last case
there is no variance parameter left to estimate, so the residual variance is
known and the Wald statistics are exactly standard normal. (As elsewhere in
glmmTMB, the residual degrees of freedom count the dispersion
parameter, so they are one lower than lm() reports for the same
fixed-effect model.) Those are defaults, as is a value taken from
getOption("glmmTMB.df"); a ddf passed in the call itself is
respected where possible. "kenward-roger" and "satterthwaite"
need random effects, so for models without them they fall back to the
residual degrees of freedom with a message, exactly as in
summary.glmmTMB. That fallback also applies to a model whose
dispersion parameter is fixed, where it overrides the infinite-df default
described above, so that emmeans() and summary() give the same
answer to the same request.
While the examples below are disabled for earlier versions of
R, they may still work; it may be necessary to refer to private
versions of methods, e.g. glmmTMB:::Anova.glmmTMB(model, ...).
warp.lm <- glmmTMB(breaks ~ wool * tension, data = warpbreaks)
salamander1 <- up2date(readRDS(system.file("example_files","salamander1.rds",package="glmmTMB")))
if (require(emmeans)) withAutoprint({
emmeans(warp.lm, poly ~ tension | wool)
emmeans(salamander1, ~ mined, type="response") # conditional means
emmeans(salamander1, ~ mined, component="cmean") # same as above, but re-gridded
emmeans(salamander1, ~ mined, component="zi", type="response") # zero probabilities
emmeans(salamander1, ~ mined, component="response") # response means including both components
})
if (getRversion() >= "3.6.0") {
if (require(car)) withAutoprint({
Anova(warp.lm,type="III")
Anova(salamander1)
Anova(salamander1, component="zi")
})
if (require(effects)) withAutoprint({
plot(allEffects(warp.lm))
plot(allEffects(salamander1))
})
}
Run the code above in your browser using DataLab