Edit: Just adding a relevant blog post that discusses checking if a random effect should be included or not, but my question is more specifically based on deciding if the intercepts and slopes of random effects should be included or not (particularly in the case of a crossed design with 2+ random effects).

Note: I'm using the same sample data as my previous (different) question here. I am aware of previous related questions such as this and this but I don't think the question was asked/answered in a general sense.

Consider the following data (R code):

sample <- read.table(text = " Session ID Response Predictor
                     1: 1 1 60.47 0.012
                     2: 2 1 42.53 -0.50
                     3: 3 1 64.45  0.01
                     4: 1 2 67.64 0.01
                     5: 2 2 51.77 0.04
                     6: 3 2 68.84 0.09
                     7: 1 3 79.80 -0.05
                     8: 2 3 46.95 0.43
                     9: 3 3 83.3 -0.05 ", h = T)
sample$Session <- factor(sample$Session)
sample$ID <- factor(sample$ID)

Where Session and ID are crossed random effects and we are examining the relationship between Predictor and Response. Now considering the possible combinations of random slopes and intercepts I believe (based on this thread) that there are a number of ways to model this situation.

model.1 <- lmer(Response ~ Predictor + (1 + Predictor| ID) + (1 + Predictor| Session), data=sample, REML = FALSE)
model.2 <- lmer(Response ~ Predictor + (1 | ID) + (1 + Predictor| Session), data=sample, REML = FALSE)
model.3 <- lmer(Response ~ Predictor + (1 + Predictor| ID) + (1 | Session), data=sample, REML = FALSE)
model.4 <- lmer(Response ~ Predictor + (0 + Predictor| ID) + (0 + Predictor| Session), data=sample, REML = FALSE)
model.5 <- lmer(Response ~ Predictor + (0 + Predictor | ID) + (1 + Predictor| Session), data=sample, REML = FALSE)
model.6 <- lmer(Response ~ Predictor + (1 + Predictor| ID) + (0 + Predictor | Session), data=sample, REML = FALSE)

This doesn't even include all of the possibilities where we find one of the random effects to not be significant e.g.

model.7 <- lmer(Response ~ Predictor + (1 + Predictor| ID) , data=sample, REML = FALSE)

So the question is: What is a systematic way of deciding which of these models should be used and thus should different random intercepts and slopes be included?

Some possibilities I've considered are:

  • Using the anova function and choosing the lowest AIC
  • Comparing to the null model e.g. model.null <- lmer(Response ~ 1 + (1 | ID) + (1 | Session), data=sample, REML = FALSE)
  • Some version of backwards/forwards stepwise elimination based on likelihood ratio tests
  • Rejecting models that have a standard deviation of 0 for any random effect
  • $\begingroup$ If a random effect has a standard deviation that is actually zero, it is likely that all of its covariance with the outcome is modeled by another term in the same model. Model.3 is the most extensive specification you have. If the grouping variables are crossed, it's possible you could additionally have (1 + Predictor | ID:Session). If all these specifications make sense and you have sufficient information in the data to estimate them, you can attempt the more flexible specifications. In your example above, you can't do this as the interaction of ID and Session is a single data point. $\endgroup$ – Heteroskedastic Jim Jul 10 '18 at 12:35
  • $\begingroup$ Thanks for the response @user162986. What exactly do you mean by extensive specification? From what I read I got the impression that using ID:Session would make sense if they were nested rather than crossed? How is it determined if these specifications make sense exactly, is it just a matter of stating some justification for them? $\endgroup$ – Seraf Fej Jul 10 '18 at 13:45
  • $\begingroup$ When you add more terms to a model, the model becomes more complicated, so extensive specification. By makes sense, I mean justify them in some way. The claim is that a coefficient is different depending on the group the case belongs to. And you can have the interaction with crossed factors too. $\endgroup$ – Heteroskedastic Jim Jul 10 '18 at 13:58
  • $\begingroup$ @user162986 In that case would model.1 not be the most extensive specification in that case as air has slopes and intercepts for both of the random effects? The linked blog post seems to indicate (see "Our Second Mixed Model") that interactions don't make sense for a crossed design (although it could well be my interpretation that is wrong). ourcodingclub.github.io/2017/03/15/mixed-models.html $\endgroup$ – Seraf Fej Jul 10 '18 at 14:35
  • $\begingroup$ You're right about model 1. You can interact crossed factors. Think standard analysis of variance. $\endgroup$ – Heteroskedastic Jim Jul 10 '18 at 15:28

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.