CODEWITHSUNDEEP
rplotrandomlme4effects
JoeN
JoeN

Reputation: 57

How to plot some terms of a lmer model

I am trying to plot only a few terms of random effects from a lmer object. I will borrow the dataset posted by oshun here.

Make up data:

tempEf <- data.frame(
  N = rep(c("Nlow", "Nhigh"), each=300),
  Myc = rep(c("AM", "ECM"), each=150, times=2),
  TRTYEAR = runif(600, 1, 15),
  site = rep(c("A","B","C","D","E"), each=10, times=12)   #5 sites
)

Make up some response data

tempEf$r <- 2*tempEf$TRTYEAR +                   
  -8*as.numeric(tempEf$Myc=="ECM") +
  4*as.numeric(tempEf$N=="Nlow") +
  0.1*tempEf$TRTYEAR * as.numeric(tempEf$N=="Nlow") +
  0.2*tempEf$TRTYEAR*as.numeric(tempEf$Myc=="ECM") +
  -11*as.numeric(tempEf$Myc=="ECM")*as.numeric(tempEf$N=="Nlow")+ 
  0.5*tempEf$TRTYEAR*as.numeric(tempEf$Myc=="ECM")*as.numeric(tempEf$N=="Nlow")+ 
  as.numeric(tempEf$site) +  #Random intercepts; intercepts will increase by 1
  tempEf$TRTYEAR/10*rnorm(600, mean=0, sd=2)    #Add some noise

library(lme4)
model <- lmer(r ~ Myc * N * TRTYEAR + (1|site), data=tempEf)

I can plot random effects by using type = "re" as follows:

plot_model(model, type = "re")

I would like to show only A and E so I add the 'terms' argument as follows:

plot_model(model, type = "re", terms = c("A", "E"))

But this does not work. Any guidance on how I can show "A" and "E" only??

Upvotes: 2

Views: 482

Answers (1)

StupidWolf
StupidWolf

Reputation: 46978

For some reason, the terms option doesn't work. Should be reflected to the authors of sjplot. Here's two workarounds:

1) manually define the terms using ggplot2:

library(ggplot2)
library(sjPlot)

plot_model(model, type = "re") + scale_x_discrete(limits=c("A","D"))

enter image description here

You get the following warnings because you throw out data

Scale for 'x' is already present. Adding another scale for 'x', which will replace the existing scale. Warning messages: 1: Removed 3 rows containing missing values (geom_point). 2: Removed 3 rows containing missing values (geom_errorbar).

2) Plot it from get_model()

    df = get_model_data(model,type="re")
    df = subset(df,term %in% c("A","D"))
    ggplot(df,aes(x=term,y=estimate,col=group)) +
    geom_point(show.legend=FALSE) + 
    geom_segment(aes(xend=term,y=conf.low,yend=conf.high),show.legend=FALSE)+
    scale_colour_brewer(palette = "Set1")+
    coord_flip()+ggtitle("Random Effects")

enter image description here

Upvotes: 1

Related Questions