Extract or Get Generalized Components from a Fitted Mixed Effects Model

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".

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 of lmer(), glmer() or nlmer().

name

a character vector specifying the name(s) of the “component”. If length(name) > 1 or if name = "ALL", a named list 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 since lme4.0; to get a backward-compatible structure, use do.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 to getME(.,"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 of cmp 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 in isSingular.)

"devfun":

deviance function (so far only available for LMMs)

"glmer.nb.theta":

negative binomial \(\theta\) parameter, only for glmer.nb.

"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))

Extract or Get Generalized Components from a Fitted Mixed Effects Model

Arguments

Value

Details

See also

Examples