How to do partial regression plots with linear mixed-effect models?

Question

I have a linear mixed-effect model in R with two continuous fixed-effects and one random effect, like this:

                           model<-lmer(y~x1+x2+(1|r),data)

To graphically display the independent effect of x1 on y, while controlling the effects of x2 (fixed effect) and r (random effect), is it appropriate to do a partial regression plot using the same logic used for multiple linear regression models? I.e.:

                          #removing the effect of x2 and r on y

                         res.y<-residuals(lmer(y~x2+(1|r),data)) 

                          #removing the effect of x2 and r on x1

                         res.x1<-residuals(lmer(x1~x2+(1|r),data)) 

              #partial regression plot to display the pure effect of x1 on y

                                     plot(res.x1,res.y)

Also, I used the "plotLMER.fnc function" from the "LMERConvenienceFunctions" R package to plot the partial effect size of each fixed effect as follows:

                                      plotLMER.fnc(model) 

However, I am not sure what this package means by "effect size". Is it β1 and β2?

I will be very grateful for any help in this issue.

Answer 1

You can use the remef package to prepare data for the visualization of partial effects. The package can be installed with devtools.

install.packages("devtools")
devtools::install_github("hohenstein/remef")

To remove the influece of x2 and the random effects from the dependent variable, you can use

library(remef)
y_partial <- remef(model, fix = "x2", ran = "all")

This will create a modified version of y based on the partial effect while the residuals are still present. Hence, you can still visualize the deviations from the predictions.

With the adjusted data y_partial you can, for example, create a plot of y_partial as a function of x1 together with a linear regression line.

---

Disclaimer: I am the author of the remef package

Answer 2

I used the Effect function from the effect package for this.

model1<-lmer(y~x1+x2+(1|r),data)

est<-Effect("x1", partial.residuals=T, model1)
plot(est)