allEffects() returns effects in incorrect order

Question

A few months ago, I computed a logistic mixed-effects model in R using the lme4 package:

mymodel = glmer(cbind(nr_corr,maximum-nr_corr) ~ (condition|Participant) + condition + CEFR_level + Other variable 1 + Other variable 2 + Other variable 3 + Other variable 4 + condition:CEFR_level, mydata, control=glmerControl(optimizer="bobyqa", optCtrl=list(maxfun=100000)), family = binomial)

I saved the final model as .rda object, as well as the effects of this model which I had compiled with the help of the allEffects() function from the effects package, i.e. simply:

myeffects <- allEffects(mymodel)
save(myeffects, file="myeffects.rda")

Since the model contains an interaction of two categorical variables (Condition x CEFR level), allEffects() also gave me the following table in the output:

condition*CEFR_level effect
                  CEFR_level
condition                 C2         A1         A2        B1        B2        C1
  Correct_items_C1 0.3193710 0.07661543 0.09713965 0.2116102 0.2502621 0.2928581
  Correct_items_C2 0.7539475 0.26099725 0.39254774 0.4964012 0.6800637 0.6915149
  Correct_items_C3 0.9664311 0.42636993 0.58915857 0.7566142 0.9076332 0.9324940
  Correct_items_C4 0.9255429 0.34131708 0.48417129 0.5688259 0.7985903 0.8905785

(C2 and Correct_items_C1 were set as the reference categories for the two variables when computing the model).

When playing around with some plotting options, I now ran the code allEffects(mymodel) again (I loaded the mymodel.rda file into Rstudio), and suddenly, I get a different output for the interaction:

condition*CEFR_level effect
                  CEFR_level
condition                 A1         A2         B1        B2        C1        C2
  Correct_items_C4 0.3193710 0.07661543 0.09713965 0.2116102 0.2502621 0.2928581
  Correct_items_C3 0.7539475 0.26099725 0.39254774 0.4964012 0.6800637 0.6915149
  Correct_items_C2 0.9664311 0.42636993 0.58915857 0.7566142 0.9076332 0.9324940
  Correct_items_C1 0.9255429 0.34131708 0.48417129 0.5688259 0.7985903 0.8905785

Upon closer inspection, I noticed that it's actually the same values as when I did it the first time, but the labeling is completely off - the column that should be 'C2' (the reference category in the model object I loaded) is 'A1', and the order for the variable 'condition' has been reversed. I suspect this has something to do with the the version of the effects package I'm using (back then 4.1-1, now 4.2-1). The original version displayed the accurate effects, which is obvious when looking at the model summary.

I also noticed that I'm receiving slightly different values for the other effects, e.g. for a continuous variable 'Other variable 1': Back then:

 Other variable 1 effect
Other variable 1
       -5        -4        -2      -0.5         1 
0.5273364 0.5728093 0.6594978 0.7184645 0.7707691 

Now:

 Other variable 1 effect
Other variable 1
       -5        -4        -2      -0.5         1 
0.4182078 0.4634982 0.5551416 0.6218171 0.6841844 

The same is true for 'Other variable 2-4' in the model. When plotting the effect, it looks fairly similar with the values from 'back then' and 'now', just the starting point on the y-axis is different:

I suspect something is different in the more recent R package containing the allEffects() function and that maybe the way it extracts the higher order terms of the model has changed. I'm afraid it now gives me inaccurate output. I've already tried to uninstall and re-install the effects package, but to no avail. How can I get allEffects () to display the output correctly? Maybe I'm overlooking something?

Answer 1

I was able to fix the issue by making sure that I was working with the newest version of R, Rtools, Rstudio and the newest version of the lme4 package. I then computed the model again and of course made sure to set the reference categories correctly beforehand:

#setting the reference categories:    
levels(mydata$CEFR_level)  #checking the reference level
mydata$CEFR_level=as.factor(mydata$CEFR_level)
mydata$CEFR_level <- relevel(mydata$CEFR_level, ref="C2")
levels(mydata$CEFR_level)   #checking whether it worked

levels(mydata$condition) #checking the reference level
mydata$condition <- relevel(as.factor(mydata$condition), ref="Correct_items_C1")
levels(mydata$condition)  #checking whether it worked

####computing the model:
mymodel = glmer(cbind(nr_corr,maximum-nr_corr) ~ (condition|Participant) + 
   condition + 
   CEFR_level + Other variable 1 + Other variable 2 +
   Other variable 3 + Other variable 4 + 
   condition:CEFR_level, mydata, 
       control=glmerControl(optimizer="bobyqa",
                  optCtrl=list(maxfun=100000)), 
       family = binomial)

Then directly afterwards, I computed the effects:

allEffects(mymodel)

I then got the accurate results, as I did when computing everything a few months ago. It works with both the earlier (4.1-1) and the current version of effects (4.2-1). So it seems that for allEffects() to recognize the reference categories accurately, it's not enough to reload the model into R and apply allEffects() to it (even when the reference categories have been set when computing the model), but it is necessary to use relevel(), then compute the model, then use allEffects().