utilities: Print and Summary Method Utilities for Mixed Effects

In lme4: Linear Mixed-Effects Models using 'Eigen' and S4

Print and Summary Method Utilities for Mixed Effects

Description

The print, summary methods (including the print for the summary() result) in lme4 are modular, using about ten small utility functions. Other packages, building on lme4 can use the same utilities for ease of programming and consistency of output.

Notably see the Examples.

llikAIC() extracts the log likelihood, AIC, and related statics from a Fitted LMM.

formatVC() “format()”s the VarCorr matrix of the random effects – for print()ing and show()ing; it is also the “workhorse” of .prt.VC(), and returns a character matrix.

.prt.*() all use cat and print to produce output.

Usage

llikAIC(object, devianceFUN = devCrit, chkREML = TRUE,
        devcomp = object@devcomp)

methTitle(dims)

.prt.methTit(mtit, class)
.prt.family (famL)
.prt.resids (resids, digits, title = "Scaled residuals:", ...)
.prt.call   (call, long = TRUE)
.prt.aictab (aictab, digits = 1)
.prt.grps   (ngrps, nobs)
.prt.warn   (optinfo, summary = FALSE, ...)

.prt.VC (varcor, digits, comp = "Std.Dev.", corr = any(comp == "Std.Dev."),
         formatter = format, ...)
formatVC(varcor, digits = max(3, getOption("digits") - 2),
         comp = "Std.Dev.", corr = any(comp == "Std.Dev."),
         formatter = format,
         useScale = attr(varcor, "useSc"), ...)

Arguments

object	a LMM model fit
devianceFUN	the function to be used for computing the deviance; should not be changed for lme4 created objects.
chkREML	optional logical indicating if object maybe a REML fit.
devcomp	for lme4 always the equivalent of object@devcomp; here a list
dims	for lme4 always the equivalent of object@devcomp$dims, a named vector or list with components "GLMM", "NLMM", "REML", and "nAGQ" of which the first two are logical scalars, and the latter two typically are FALSE or numeric.
mtit	the result of methTitle(object)
class	typically class(object).
famL	a list with components family and link, each a character string; note that standard R family objects can be used directly, as well.
resids	numeric vector of model residuals.
digits	non-negative integer of (significant) digits to print minimally.
title	character string.
...	optional arguments passed on, e.g., to residuals().
call	the call of the model fit; e.g., available via (generic) function getCall().
long	logical indicating if the output may be long, e.g., printing the control part of the call if there is one.
aictab	typically the AICtab component of the result of llikAIC().
varcor	typically the result of VarCorr().
comp	optional character vector of length 1 or 2, containing "Std.Dev." and/or "Variance", indicating the columns to use.
corr	logical indicating if correlations or covariances should be used for vector random effects.
formatter	a function used for formatting the numbers.
ngrps	integer (vector), typically the result of ngrps(object).
nobs	integer; the number of observations, e.g., the result of nobs.
optinfo	typically object @ optinfo, the optimization infos, including warnings if there were.
summary	logical
useScale	(logical) whether the parent model estimates a scale parameter.

Value

llikAIC() returns a list with components

logLik	which is logLik(object), and
AICtab	a “table” of AIC, BIC, logLik, deviance and df.residual() values.

Examples

## Create a few "lme4 standard" models ------------------------------
fm1 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
fmM <- update(fm1, REML=FALSE) # -> Maximum Likelihood
fmQ <- update(fm1, . ~ Days + (Days | Subject))

gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
             data = cbpp, family = binomial)
gmA <- update(gm1, nAGQ = 5)


(lA1 <- llikAIC(fm1))
(lAM <- llikAIC(fmM))
(lAg <- llikAIC(gmA))

(m1 <- methTitle(fm1 @ devcomp $ dims))
(mM <- methTitle(fmM @ devcomp $ dims))
(mG <- methTitle(gm1 @ devcomp $ dims))
(mA <- methTitle(gmA @ devcomp $ dims))

.prt.methTit(m1, class(fm1))
.prt.methTit(mA, class(gmA))

.prt.family(gaussian())
.prt.family(binomial())
.prt.family( poisson())

.prt.resids(residuals(fm1), digits = 4)
.prt.resids(residuals(fmM), digits = 2)

.prt.call(getCall(fm1))
.prt.call(getCall(gm1))

.prt.aictab ( lA1 $ AICtab ) # REML
.prt.aictab ( lAM $ AICtab ) # ML --> AIC, BIC, ...

V1 <- VarCorr(fm1)
m <- formatVC(V1)
stopifnot(is.matrix(m), is.character(m), ncol(m) == 4)
print(m, quote = FALSE) ## prints all but the first line of .prt.VC() below:
.prt.VC( V1, digits = 4)
## Random effects:
##  Groups   Name        Std.Dev. Corr
##  Subject  (Intercept) 24.740
##           Days         5.922   0.07
##  Residual             25.592
p1 <- capture.output(V1)
p2 <- capture.output( print(m, quote=FALSE) )
pX <- capture.output( .prt.VC(V1, digits = max(3, getOption("digits")-2)) )
stopifnot(identical(p1, p2),
          identical(p1, pX[-1])) # [-1] : dropping 1st line

(Vq <- VarCorr(fmQ)) # default print()
print(Vq, comp = c("Std.Dev.", "Variance"))
print(Vq, comp = c("Std.Dev.", "Variance"), corr=FALSE)
print(Vq, comp = "Variance")

.prt.grps(ngrps = ngrps(fm1),
          nobs  = nobs (fm1))
## --> Number of obs: 180, groups:  Subject, 18

.prt.warn(fm1 @ optinfo) # nothing .. had no warnings
.prt.warn(fmQ @ optinfo) # (ditto)

lme4 documentation built on July 3, 2024, 5:11 p.m.