Learn R Programming

spaMM (version 4.1.20)

residuals.HLfit: Extract model residuals

Description

Extracts several types of residuals from an object of class HLfit. Note that the default type ("deviance") of returned residuals differs from the default (response residuals) of equivalent functions in base R.

Usage

# S3 method for HLfit
residuals(object, 
  type = c("deviance", "pearson", "response", "working", "std_dev_res"), force=FALSE, ...)

Value

A vector of residuals

Arguments

object

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

type

The type of residuals which should be returned. See Details for additional information.

force

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

...

For consistency with the generic.

Details

The four types "deviance" (default), "pearson", "response" are "working" are, for GLM families, the same that are returned by residuals.glm. "working" residuals may be returned only for fixed-effect models. "deviance" residuals are the signed square root of those returned by dev_resids.

Following Lee et al. (2006, p.52), the standardized deviance residuals returned for type="std_dev_res" are defined as the deviance residuals divided by \(\phi\sqrt(1-q)\), where the deviance residuals are defined as for a GLM, \(\phi\) is the dispersion parameter of the response family (a vector of values, for heteroscedastic cases), and \(q\) is a vector of leverages given by hatvalues(., type="std") (see hatvalues for details about these specific standardizing leverages).

These definitions are special cases of more general ones holding for non-GLM response families. In the latter case, the deviance residuals are as defined in Details of llm.fit, and "std_dev_res" residuals are defined from them as above for GLM response families, with the additional convention that \(\phi=1\). Pearson residuals and response residuals are defined as in stats:::residuals.glm. The "working" residuals are defined for each response as \(- [d \log(clik)/d \eta]/[d^2 \log(clik)/d \eta^2]\) where clik is the conditional likelihood.

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
data("wafers")
fit <- fitme(y ~X1+(1|batch) ,data=wafers, init=list(phi=NaN))  # : this 'init' 
#                 implies that standardized deviance residuals are saved in the 
#                 fit result, allowing the following comparison: 

r1 <- residuals(fit, type="std_dev_res") # gets stored value
r2 <- residuals(fit, type="std_dev_res", force=TRUE) # forced recomputation
if (diff(range(r1-r2))>1e-14) stop()

Run the code above in your browser using DataLab