Extract or Get Generalized Components from a Fitted Mixed Effects Model
getME.Rd
Extract (or “get”) “components” – in a
generalized sense – from a fitted mixed-effects model,
i.e., (in this version of the package) from an object of
class "merMod"
.
Usage
getME(object, name, ...)
# S3 method for class 'merMod'
getME(object,
name = c("X", "Z", "Zt", "Ztlist", "mmList", "y", "mu", "u", "b",
"Gp", "Tp", "L", "Lambda", "Lambdat", "Lind", "Tlist",
"A", "RX", "RZX", "sigma", "flist",
"fixef", "beta", "theta", "ST", "REML", "is_REML",
"n_rtrms", "n_rfacs", "N", "n", "p", "q",
"p_i", "l_i", "q_i", "k", "m_i", "m",
"cnms", "devcomp", "offset", "lower", "devfun", "glmer.nb.theta"),
...)
Arguments
- object
a fitted mixed-effects model of class
"merMod"
, i.e., typically the result oflmer()
,glmer()
ornlmer()
.- name
a character vector specifying the name(s) of the “component”. If
length(name) > 1
or ifname = "ALL"
, a namedlist
of components will be returned. Possible values are:"X"
:fixed-effects model matrix
"Z"
:random-effects model matrix
"Zt"
:transpose of random-effects model matrix. Note that the structure of
Zt
has changed sincelme4.0
; to get a backward-compatible structure, usedo.call(Matrix::rBind,getME(.,"Ztlist"))
"Ztlist"
:list of components of the transpose of the random-effects model matrix, separated by individual variance component
"mmList"
:list of raw model matrices associated with random effects terms
"y"
:response vector
"mu"
:conditional mean of the response
"u"
:conditional mode of the “spherical” random effects variable
"b"
:conditional mode of the random effects variable
"Gp"
:groups pointer vector. A pointer to the beginning of each group of random effects corresponding to the random-effects terms, beginning with 0 and including a final element giving the total number of random effects
"Tp"
:theta pointer vector. A pointer to the beginning of the theta sub-vectors corresponding to the random-effects terms, beginning with 0 and including a final element giving the number of thetas.
"L"
:sparse Cholesky factor of the penalized random-effects model.
"Lambda"
:relative covariance factor \(\Lambda\) of the random effects.
"Lambdat"
:transpose \(\Lambda'\) of \(\Lambda\) above.
"Lind"
:index vector for inserting elements of \(\theta\) into the nonzeros of \(\Lambda\).
"Tlist"
:vector of template matrices from which the blocks of \(\Lambda\) are generated.
"A"
:Scaled sparse model matrix (class
"dgCMatrix"
) for the unit, orthogonal random effects, \(U\), equal togetME(.,"Zt") %*% getME(.,"Lambdat")
"RX"
:Cholesky factor for the fixed-effects parameters
"RZX"
:cross-term in the full Cholesky factor
"sigma"
:residual standard error; note that
sigma(object)
is preferred."flist"
:a list of the grouping variables (factors) involved in the random effect terms
"fixef"
:fixed-effects parameter estimates
"beta"
:fixed-effects parameter estimates (identical to the result of
fixef
, but without names)"theta"
:random-effects parameter estimates: these are parameterized as the relative Cholesky factors of each random effect term
"ST"
:A list of S and T factors in the TSST' Cholesky factorization of the relative variance matrices of the random effects associated with each random-effects term. The unit lower triangular matrix, \(T\), and the diagonal matrix, \(S\), for each term are stored as a single matrix with diagonal elements from \(S\) and off-diagonal elements from \(T\).
"n_rtrms"
:number of random-effects terms
"n_rfacs"
:number of distinct random-effects grouping factors
"N"
:number of rows of
X
"n"
:length of the response vector,
y
"p"
:number of columns of the fixed effects model matrix,
X
"q"
:number of columns of the random effects model matrix,
Z
"p_i"
:numbers of columns of the raw model matrices,
mmList
"l_i"
:numbers of levels of the grouping factors
"q_i"
:numbers of columns of the term-wise model matrices,
ZtList
"k"
:number of random effects terms
"m_i"
:numbers of covariance parameters in each term
"m"
:total number of covariance parameters, i.e., the same as
dims@nth
below."cnms"
:the “component names”, a
list
."REML"
:0
indicates the model was fitted by maximum likelihood, any other positive integer indicates fitting by restricted maximum likelihood"is_REML"
:same as the result of
isREML(.)
"devcomp"
:a list consisting of a named numeric vector,
cmp
, and a named integer vector,dims
, describing the fitted model. The elements ofcmp
are:- ldL2
twice the log determinant of
L
- ldRX2
twice the log determinant of
RX
- wrss
weighted residual sum of squares
- ussq
squared length of
u
- pwrss
penalized weighted residual sum of squares, “wrss + ussq”
- drsum
sum of residual deviance (GLMMs only)
- REML
REML criterion at optimum (LMMs fit by REML only)
- dev
deviance criterion at optimum (models fit by ML only)
- sigmaML
ML estimate of residual standard deviation
- sigmaREML
REML estimate of residual standard deviation
- tolPwrss
tolerance for declaring convergence in the penalized iteratively weighted residual sum-of-squares (GLMMs only)
The elements of
dims
are:- N
number of rows of
X
- n
length of
y
- p
number of columns of
X
- nmp
n-p
- nth
length of
theta
- q
number of columns of
Z
- nAGQ
see
glmer
- compDev
see
glmerControl
- useSc
TRUE
if model has a scale parameter- reTrms
number of random effects terms
- REML
0
indicates the model was fitted by maximum likelihood, any other positive integer indicates fitting by restricted maximum likelihood- GLMM
TRUE
if a GLMM- NLMM
TRUE
if an NLMM
"offset"
:model offset
"lower"
:lower bounds on random-effects model parameters (i.e, "theta" parameters). In order to constrain random effects covariance matrices to be semi-positive-definite, this vector is equal to 0 for elements of the
theta
vector corresponding to diagonal elements of the Cholesky factor,-Inf
otherwise. (getME(.,"lower")==0
can be used as a test to identify diagonal elements, as inisSingular
.)"devfun"
:deviance function (so far only available for LMMs)
"glmer.nb.theta"
:negative binomial \(\theta\) parameter, only for
glmer.nb
.
%% -- keep at the very end:
"ALL"
:get all of the above as a
list
.
- ...
currently unused in lme4, potentially further arguments in methods.
Value
Unspecified, as very much depending on the name
.
Details
The goal is to provide “everything a user may
want” from a fitted "merMod"
object as far
as it is not available by methods, such as
fixef
, ranef
,
vcov
, etc.
Examples
## shows many methods you should consider *before* using getME():
methods(class = "merMod")
#> [1] PBmodcomp PBrefdist VarCorr
#> [4] anova as.function coef
#> [7] confint cooks.distance deviance
#> [10] df.residual drop1 extractAIC
#> [13] family fitted fixef
#> [16] formula getData getL
#> [19] getME hatvalues influence
#> [22] isGLMM isLMM isNLMM
#> [25] isREML logLik model.frame
#> [28] model.matrix model2restriction_matrix ngrps
#> [31] nobs plot predict
#> [34] print profile qqmath
#> [37] ranef rePCA refit
#> [40] refitML residuals rstudent
#> [43] show sigma simulate
#> [46] summary terms update
#> [49] vcov weights
#> see '?methods' for accessing help and source code
(fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
#> Linear mixed model fit by REML ['lmerMod']
#> Formula: Reaction ~ Days + (Days | Subject)
#> Data: sleepstudy
#> REML criterion at convergence: 1743.628
#> Random effects:
#> Groups Name Std.Dev. Corr
#> Subject (Intercept) 24.741
#> Days 5.922 0.07
#> Residual 25.592
#> Number of obs: 180, groups: Subject, 18
#> Fixed Effects:
#> (Intercept) Days
#> 251.41 10.47
Z <- getME(fm1, "Z")
stopifnot(is(Z, "CsparseMatrix"),
c(180,36) == dim(Z),
all.equal(fixef(fm1), b1 <- getME(fm1, "beta"),
check.attributes=FALSE, tolerance = 0))
## A way to get *all* getME()s :
## internal consistency check ensuring that all work:
parts <- getME(fm1, "ALL")
str(parts, max=2)
#> List of 45
#> $ X : num [1:180, 1:2] 1 1 1 1 1 1 1 1 1 1 ...
#> ..- attr(*, "dimnames")=List of 2
#> ..- attr(*, "assign")= int [1:2] 0 1
#> ..- attr(*, "msgScaleX")= chr(0)
#> $ Z :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> $ Zt :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> $ Ztlist :List of 2
#> ..$ Subject.(Intercept):Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> ..$ Subject.Days :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> $ mmList :List of 1
#> ..$ Days | Subject: num [1:180, 1:2] 1 1 1 1 1 1 1 1 1 1 ...
#> .. ..- attr(*, "dimnames")=List of 2
#> .. ..- attr(*, "assign")= int [1:2] 0 1
#> $ y : num [1:180] 250 259 251 321 357 ...
#> $ mu : num [1:180] 254 273 293 313 332 ...
#> $ u : num [1:36] 2.34 39.68 -41.79 -34.58 -40.3 ...
#> $ b :Formal class 'dgeMatrix' [package "Matrix"] with 4 slots
#> $ Gp : int [1:2] 0 36
#> $ Tp : Named num [1:2] 0 3
#> ..- attr(*, "names")= chr [1:2] "beg__" "Subject"
#> $ L :Formal class 'dCHMsimpl' [package "Matrix"] with 11 slots
#> $ Lambda :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> $ Lambdat :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> $ Lind : int [1:54] 1 2 3 1 2 3 1 2 3 1 ...
#> $ Tlist :List of 1
#> ..$ Subject: num [1:2, 1:2] 0.9667 0.0152 0 0.2309
#> $ A :Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
#> $ RX : num [1:2, 1:2] 3.79 0 2.3 16.56
#> ..- attr(*, "dimnames")=List of 2
#> $ RZX : num [1:36, 1:2] 3.022 0.269 3.022 0.269 3.022 ...
#> ..- attr(*, "dimnames")=List of 2
#> $ sigma : num 25.6
#> $ flist :List of 1
#> ..$ Subject: Factor w/ 18 levels "308","309","310",..: 1 1 1 1 1 1 1 1 1 1 ...
#> ..- attr(*, "assign")= int 1
#> $ fixef : Named num [1:2] 251.4 10.5
#> ..- attr(*, "names")= chr [1:2] "(Intercept)" "Days"
#> $ beta : num [1:2] 251.4 10.5
#> $ theta : Named num [1:3] 0.9667 0.0152 0.2309
#> ..- attr(*, "names")= chr [1:3] "Subject.(Intercept)" "Subject.Days.(Intercept)" "Subject.Days"
#> $ ST :List of 1
#> ..$ Subject: num [1:2, 1:2] 0.9667 0.0157 0 0.2309
#> $ REML : int 2
#> $ is_REML : logi TRUE
#> $ n_rtrms : int 1
#> $ n_rfacs : int 1
#> $ N : int 180
#> $ n : int 180
#> $ p : int 2
#> $ q : int 36
#> $ p_i : Named int 2
#> ..- attr(*, "names")= chr "Days | Subject"
#> $ l_i : Named int 18
#> ..- attr(*, "names")= chr "Subject"
#> $ q_i : Named int 36
#> ..- attr(*, "names")= chr "Days | Subject"
#> $ k : int 1
#> $ m_i : Named num 3
#> ..- attr(*, "names")= chr "Days | Subject"
#> $ m : int 3
#> $ cnms :List of 1
#> ..$ Subject: chr [1:2] "(Intercept)" "Days"
#> $ devcomp :List of 2
#> ..$ cmp : Named num [1:10] 7.60e+01 8.28 9.89e+04 1.77e+04 1.17e+05 ...
#> .. ..- attr(*, "names")= chr [1:10] "ldL2" "ldRX2" "wrss" "ussq" ...
#> ..$ dims: Named int [1:12] 180 180 2 178 36 3 1 1 0 2 ...
#> .. ..- attr(*, "names")= chr [1:12] "N" "n" "p" "nmp" ...
#> $ offset : num [1:180] 0 0 0 0 0 0 0 0 0 0 ...
#> $ lower : Named num [1:3] 0 -Inf 0
#> ..- attr(*, "names")= chr [1:3] "Subject.(Intercept)" "Subject.Days.(Intercept)" "Subject.Days"
#> $ devfun :function (theta)
#> $ glmer.nb.theta: logi NA
stopifnot(identical(Z, parts $ Z),
identical(b1, parts $ beta))