Learn R Programming

spaMM (version 4.7.0)

hatvalues.HLfit: Hat values and leverages extractor

Description

hatvalues.HLfit is a method for the generic function hatvalues:
\(\bullet\) hatvalues(.,type="projection") will return “hat values”, defined as diagonal elements of the projection matrix of the augmented linear model;
\(\bullet\) hatvalues(.,type="std") will return generalized “leverages”, which are standardizing valeus that can be used for estimation of dispersion parameters (see Details for the differences).

To get the full hat matrix, see get_matrix(., which="hat_matrix").

Usage

# S3 method for HLfit
hatvalues(model, type = "projection", which = "resid", force=FALSE, ...)

Value

A list with separate components resid (leverages of the observations) and ranef if which="both", and a vector otherwise.

Arguments

model

An object of class HLfit, as returned by the fitting functions in spaMM.

type

Character: "projection", "std", or more cryptic values not documented.

which

Character: "resid" for the traditional leverages of the observations, "ranef" for random-effect leverages, or "both" for both.

force

Boolean: to force recomputation of the leverages even if they are available in the object, for checking purposes.

...

For consistency with the generic.

Details

Leverages may have distinct meanings depending on context. The textbook version for linear models is that leverages \((q_i)\) are the diagonal elements of a projection matrix (“hat matrix”), and that they may be used to standardize (“studentize”) residuals as follows. If the residual variance \(\phi\) is known, then the variance of each fitted residual \(\hat{e}_i\) is \(\phi(1-q_i)\). Standardized residuals, all with variance 1, are then \(\hat{e}_i/\)\(\sqrt{}\)\((\phi(1-q_i))\). This standardization of variance no longer holds exactly with estimated \(\phi\), but if one uses here an unbiased (REML) estimator of \(\phi\), the studentized residuals may still practically have a unit expected variance.

In mixed models, a generalized definition of the projection matrix is available, in terms of the matrix of the augmented linear model [its value \(\mathbf{P}_\mathrm{a}\) can be retrieved by get_matrix(., which="hat_matrix")]. The “hat values”, defined as the diagonal elements of \(\mathbf{P}_\mathrm{a}\), are also defined for the random effects. However, in general the standardizing “leverages” \(\mathbf{q}^*\) that can again be used to obtain estimates of \(\phi\) and \(\lambda\), by standardizing residuals (for \(\phi\)), or by standardizing predicted random effects (for \(\lambda\)), are no longer simply the “hat values”.

For illustration, the following simple differences hold between "projection" hat values and "std" values:
\(\bullet\) "std" values are defined so that they are 0 for ML fits of simple linear model.
\(\bullet\) In LMMs fitted by REML, "std" values are the same as "projection" values. For ML fits, these two forms of leverages differ, because "std" values are diagonal values of another matrix than \(\mathbf{P}_\mathrm{a}\). Note that the "projection" hat values also differ between ML and REML fits because the fitted \(\mathbf{P}_\mathrm{a}\) itself depends on \(\hat{\phi}\) and \(\hat{\lambda}\), which differ between ML and REML fits.
\(\bullet\) In GLMMs (or HGLMs), the "std" values may differ from the "projection" ones even for REML, because they may include corrections for non-Gaussian response family, for non-Gaussian random effects, and for taking into account the variation of the GLM weights in the logdet(info.mat) derivatives. Which corrections are included depend on the precise method used to fit the model (e.g., EQL vs PQL vs REML).

References

Lee, Y., Nelder, J. A. and Pawitan, Y. (2006) Generalized linear models with random effects: unified analysis via h-likelihood. Chapman & Hall: London.

Examples

Run this code
if (spaMM.getOption("example_maxtime")>0.8) {
data("Orthodont",package = "nlme")
rnge <- (107:108)

# all different:
#
hatvalues(rlfit <- fitme(distance ~ age+(age|Subject), 
                         data = Orthodont, method="REML"))[rnge]
hatvalues(mlfit <- fitme(distance ~ age+(age|Subject), 
                         data = Orthodont))[rnge] 
hatvalues(mlfit,type="std")[rnge]
}

Run the code above in your browser using DataLab