Julie Charbonneau
Julie Charbonneau

Reputation: 35

Extracting partial residuals from effects package for a linear mixed effects model with an interaction term

I’m running a linear mixed effects model with an interaction term.

library(lme4)
library(effects)
data<-mtcars
model<-lme4::lmer(mpg~cyl + hp*gear+ disp + (1|carb), REML=T, data=data)

# The partial residuals appear on the default effect() plot
plot(effect("hp:gear",partial.residuals=TRUE, model))

# But partial residual values are not included when saving as a dataframe
residuals <- as.data.frame(effect("hp:gear",partial.residuals=TRUE, model))

Any idea how to extract the partial residuals for a LMEM with interaction term from effects()?

Upvotes: 1

Views: 1820

Answers (2)

John Fox
John Fox

Reputation: 384

The partial residuals are computed by the plot() method, not by Effect(), because it's necessary to know which points are in each panel, information that's available to a panel function in a lattice plot, before computing the partial residuals.

The method used is described in J. Fox and S. Weisberg (2018), "Visualizing Fit and Lack of Fit in Complex Regression Models with Predictor Effect Plots with Partial Residuals," Journal of Statistical Software 87:9, 1–27 (which is referenced in ?Effect).

Upvotes: 2

Rui Barradas
Rui Barradas

Reputation: 76621

First, run the modeling code in the question. Then, save the output of function effects.

eff <- effect("hp:gear", partial.residuals = TRUE, model)

Now, resid(eff) gives the residuals. Since the question asks for the partial residuals, go look for them:

grep("resid", names(eff))
#[1] 12 13

names(eff)[grep("resid", names(eff))]
#[1] "residuals"               "partial.residuals.range"

And use the second name to extract the list member "partial.residuals.range".

eff$partial.residuals.range
#[1] 10.4 33.9

For all the partial residuals, they are defined as

predict(model) + resid(model)

range(predict(model) + resid(model))
#[1] 10.4 33.9

Upvotes: 0

Related Questions