partialResiduals: Compute partial residuals for linear and linear mixed models.

In bozenne/butils: Various useful functions

Compute partial residuals for linear and linear mixed models.

Description

Compute partial residuals for linear and linear mixed models.

Usage

partialResiduals(
  model,
  var,
  keep.intercept = FALSE,
  conditional = FALSE,
  interval = "confidence",
  level = 0.95,
  npoints = 100,
  quantile.range = c(0, 1),
  FUN.df,
  ...
)

Arguments

model	Model object (e.g. lm)
var	[character] Independent variable(s) whose effect should be kept in the partial residuals.
keep.intercept	[logical] should the partial residual be computed keeping the contribution of the reference level? Can also a vector of lenght 2, the second element indicating whether the uncertainty calculation should account for the estimation of the intercept.
conditional	[logical] should the partial residuals be computed without the random effects? (if any)
interval	[character] Type of interval calculation ("confidence" or "prediction").
level	[numeric, 0-1] Level of confidence limits (default 95%)
npoints	[integer] Length of the vector of unique values relative to each continuous variable.
quantile.range	[numeric vector, 0-1] the quantiles (for the continous covariates) between which the fitted values will be computed.
FUN.df	Optional function returning the residual degree of freedoms.
...	additional arguments to lower level functions

Details

In a linear model:

Y = \alpha + \beta X + \gamma Z + \varepsilon

The partial residuals are defined by:

\varepsilon_{X} = \beta X + \varepsilon

or

\varepsilon_{X} = \alpha + \beta X + \varepsilon

depending on the value of the argument keep.intercept. The X matrix is defined contains the variables defined by the var argument.

Confidence intervals are only valid in homoschedastic models. When using mixed models, the confidence and prediction intervals ignore the uncertainty of the covariance parameters/random effects.

Value

list with following members:

data	Original dataset with an additional column containing the partial residuals
partialFit	Fitted values according to the range of values of var
var	variable(s) for which the partial residuals have been computed
level	Level of confidence
name.Y	name of the response variable

Author(s)

Brice Ozenne

See Also

plot.partialResiduals for a graphical display and more examples.'

Examples

library(lava)
set.seed(10)
m.lvm <- lvm(Y~X1+X2+Id)
categorical(m.lvm, K = 5) <- ~Id
d <- lava::sim(n = 1e2, m.lvm)
d$Id <- as.factor(d$Id)
m <- lm(Y~X1+X2+Id, data = d)

#### 1- partial residuals relative to no variable matches the residuals
pres0 <- partialResiduals(m, var = NULL)
range(residuals(m) - pres0$pResiduals)

#### 2- partial residuals regarding the variable X1
pres1 <- partialResiduals(m, var = "X1")
fit1 <- coef(m)["(Intercept)"] + coef(m)["X2"] * d$X2 +  coef(m)["Id"] * d$Id
range((d$Y - fit1) - pres1$pResiduals)
cor(d$X1, pres1$pResiduals)
lava::partialcor(~X2+Id,d[,c("Y","X1","X2","Id")])

## Not run: 
## graphical display
autoplot(partialResiduals(m, var = "X1", keep.intercept = TRUE))
autoplot(partialResiduals(m, var = "X1", keep.intercept = FALSE))
autoplot(partialResiduals(m, var = "X1", keep.intercept = c(TRUE,FALSE)))
autoplot(partialResiduals(m, var = "Id", keep.intercept = TRUE))

## End(Not run)

#### 3-  partial residuals regarding both X1 and X2
pres2 <- partialResiduals(m, var = c("X1","X2"))
fit2 <- coef(m)["(Intercept)"] + d$Id * coef(m)["Id"]
range((d$Y - fit2) - pres2$pResiduals)

## Not run: 
## graphical display
autoplot(partialResiduals(m, var = c("X1","Id"), keep.intercept = TRUE))
autoplot(partialResiduals(m, var = c("X1","Id"), keep.intercept = FALSE))
autoplot(partialResiduals(m, var = c("X1","Id"), keep.intercept = c(TRUE,FALSE)))

## End(Not run)

## partial residuals binary variable
pres3 <- partialResiduals(m, var = "Id")

if(require(ggplot2)){
   pres3$Id.char <- factor(pres3$Id,
                                levels = 0:4,
                                labels = c("a","b","c","d","e"))
   gg <- ggplot(pres3, aes(y = pResiduals, group = Id.char, x = Id.char))
   gg <- gg + geom_boxplot()
   gg <- gg + geom_dotplot(binaxis = "y", stackdir = "center")
   gg
}

## partial residuals in presence of interactions
m.I <- lm(Y~X1*X2+Id, data = d)
pres.I <- partialResiduals(m.I, var = c("X1","X2"))
fit.I <- coef(m.I)["(Intercept)"] + d$Id * coef(m.I)["Id"]
range((d$Y - fit.I) - pres.I$pResiduals)

bozenne/butils documentation built on Oct. 14, 2023, 6:19 a.m.