How to plot random intercept and slope in a mixed model with multiple predictors?

lmer, mixed-models, plot, r, random-effects

Solution

## generate one response, two predictors and one factor (random effect)
set.seed(101)
resp <- runif(100,1,100)
pred1<- rnorm(100, 
           mean=rep(resp[1:50],2)+rep(c(-10,20),each=50),
           sd=rep(c(10,5),each=50))
pred2<- rnorm(100, resp-10, 10)

NOTE that you should probably not be trying to fit a random effect for an grouping variable with only two levels -- this will almost invariably result in an estimated random-effect variance of zero, which will in turn put your predicted lines right on top of each other -- I'm switching from `gl(2,50)` to `gl(10,10)` ...

RF1<-gl(10,10)
d <- data.frame(resp,pred1,pred2,RF1)

#lmer
library(lme4)
mod <- lmer(resp ~ pred1 + pred2 + (1|RF1),data=d)

The development version of `lme4` has a `predict()` function that makes this a little easier ...

- Predict for a range of `pred1` with `pred2` equal to its mean, and vice versa. This is all a little bit cleverer than it needs to be, since it generates all the values for both focal predictors and plots them with ggplot in one go ...

()

nd <- with(d,
           rbind(data.frame(expand.grid(RF1=levels(RF1),
                      pred1=seq(min(pred1),max(pred1),length=51)),
                      pred2=mean(pred2),focus="pred1"),
                 data.frame(expand.grid(RF1=levels(RF1),
                      pred2=seq(min(pred2),max(pred2),length=51)),
                      pred1=mean(pred1),focus="pred2")))
nd$val <- with(nd,pred1[focus=="pred1"],pred2[focus=="pred2"])
pframe <- data.frame(nd,resp=predict(mod,newdata=nd))
library(ggplot2)
ggplot(pframe,aes(x=val,y=resp,colour=RF1))+geom_line()+
         facet_wrap(~focus,scale="free")

- Alternatively, focusing just on `pred1` and generating predictions for a (small/discrete) range of `pred2` values ...

()

nd <- with(d,
           data.frame(expand.grid(RF1=levels(RF1),
                      pred1=seq(min(pred1),max(pred1),length=51),
                      pred2=seq(-20,100,by=40))))
pframe <- data.frame(nd,resp=predict(mod,newdata=nd))
ggplot(pframe,aes(x=pred1,y=resp,colour=RF1))+geom_line()+
         facet_wrap(~pred2,nrow=1)

You might want to set `scale="free"` in the last `facet_wrap()` ... or use `facet_grid(~pred2,labeller=label_both)`

For presentation you might want to replace the `colour` aesthetic, with `group`, if all you want to do is distinguish among groups (i.e. plot separate lines) rather than identify them ...

Problem

Is it possible to plot the random intercept or slope of a mixed model when it has more than one predictor? With one predictor I would do like this: ``` #generate one response, two predictors and one factor (random effect) resp<-runif(100,1, 100) pred1<-c(resp[1:50]+rnorm(50, -10, 10),resp[1:50]+rnorm(50, 20, 5)) pred2<-resp+rnorm(100, -10, 10) RF1<-gl(2, 50) #gamm library(mgcv) mod<-gamm(resp ~ pred1, random=list(RF1=~1)) plot(pred1, resp, type="n") for (i in ranef(mod$lme)[[1]]) { abline(fixef(mod$lme)[1]+i, fixef(mod$lme)[2]) } #lmer library(lme4) mod<-lmer(resp ~ pred1 + (1|RF1)) plot(pred1, resp, type="n") for (i in ranef(mod)[[1]][,1]) { abline(fixef(mod)[1]+i, fixef(mod)[2]) } ``` But what if I have a model like this instead?: ``` mod<-gamm(resp ~ pred1 + pred2, random=list(RF1=~1)) ``` Or with lmer ``` mod<-lmer(resp ~ pred1 + pred2 + (1|RF1)) ``` Should I consider all the coefficients or only the ones of the variable that I'm plotting? Thanks

Original source