lme4: adding covariates and interpreting output

I am new to mixed models and lme4 package. I think my data is too complicated for me to understand things using it... what I want is for someone to please explain this to me in layman's terms using a simple example (if one exists).

My data: I have a lot of data (multiple experiments), but generally speaking... response variables = continuous, normally distributed (some are not, I will get to those later), explanatory variables = mostly discrete. Random effects = group/individual (I "nest" these).

First, in model setup, I need to add mass as a covariate for one of my models. I don't know how to do that, so a basic example (with code) of how one adds a covariate to a lmm (with lme4 package) would be great. I read their help manual on this but it wasn't very helpful to me.

Second, I get the output from lme4 and I have no idea what it is telling me. When I rerun lme4 things with the afex package for the oh-so-controversial p-values, it tells me everything is significant but I still don't know what that means.

Any help or any direction toward a YouTube video, website tutorial, etc. that can explain this to a complete idiot (which is what I feel like I am at the moment) would be greatly appreciated.

Thank you!

• "response variables = continuous, normally distributed (some are not, I will get to those later)" This quote is already concerning. Before you proceed to mixed-effects models, you need to understand that the distribution of your repsonse variable is generally irrelevant for regression analysis. What is important is the conditional distribution (conditional on the predictors). For a model with a Gaussian distribution family this means you need to check the distribution of the residuals because that is assumed to be a normal distribution. Oct 8, 2021 at 7:02
• Adding a covariate just means adding another independent variable. Oct 8, 2021 at 7:02

I need to add mass as a covariate for one of my models. I don't know how to do that

If your model looks something like

y ~ x1 + (1 | G)


where x1 is a fixed effect and G is a grouping variable for which you want to fit random intercepts, then if mass is the name of the mass variable you would fit:

y ~ x1 + mass + (1 | G)


When I rerun lme4 things with the afex package for the oh-so-controversial p-values, it tells me everything is significant but I still don't know what that means.

Each p-value refers to a statistical test of a null hypothesis - for example, a test of whether a regression coefficient is zero or not. The p-values tell you the probability of observing your data, or data more extreme, if (and only if) the null hypothesis is true. If this doesn't make sense to you, then you are not alone, and that is part of the reason for the controversy surrounding p-values. Another issue regarding mixed models in particular is that the p-values are actually just approximations, because the number of degrees of freedom for the relevant test is not known and can only be approximated (eg Satterthwaite's or Kenward-Roger's method), so the idea of using p-values with respect to some arbitrary cut-off value for "significance", is folly. I suggest erring on the side of caution and don't use p-values. Or if you really want a probability that makes sense, then you can obtain probabilities that the null hypotheses are true (which is typically what most people would like to know, rather than what a p-value tells them), then you can adopt a fully Bayesian approach using a package like brms.

Once final point on your model. You mention that you fit nested random effects for individuals within groups. Depending on how many groups you have, this probably does not make sense. Eg if it's a treatment group where participants are nested within 1 of 2 groups, then the group variable should be a fixed effect, not random. 2 is too few levels to fit random intercepts.

I hope this helps:

Question 1:

This is code assuming you have 3 predictors (pred1, pred2, pred3). A'covariate' is just treated the same as other predictors in the code. I use 'Response' as your dependent (outcome) variable. I use 'X' for your grouping variable. 'data.df' is your data frame, with a column for each of the predictors, outcome variable, and random effect variable.

#Load the packages needed for the analysis

library(lme4) # for running the mixed model
library (lmerTest) # for obtaining p-values in the output tables

# Run the first linear mixed effects model

model1 <- lmer(Response ~ pred1+(1|X), data = data.df ,control=lmerControl(optimizer="bobyqa"))

summary(model1) #prints the model output table

# Run the second model, with added covariates

model2 <- lmer(Response ~ pred1+pred2+pred3+(1|X), data = data.df ,control=lmerControl(optimizer="bobyqa"))

summary(model2)


• You can add in interactions between variables using * (e.g., Response ~ pred1*pred2+pred3)

• 1|X is a random intercept. Where theoretically appropriate, you can add in random slopes. For example, if individual participants in a study complete a spelling test under 2 conditions, then you could run the model

lmer(Spelling_performance ~ condition + (condition|participant),data = data.df ,control=lmerControl(optimizer="bobyqa"))

This would include a random intercept for participants (assuming participants each have a different 'baseline' performance, i.e., spelling ability). It would also include a by-participant random slope for condition, meaning that the effect of condition varies between participants: some participants might perform similarly in both conditions, while others might be more susceptible to changes in condition.

• If you have convergence problems, you may need to consider scaling the variable or centering (particularly if you have interactions as this can reduce multicollinearity between main effects and interactions). See references 1 and 2 below. You can also look to simplify the random effects structure.

Question 2: I would need to see your code and output table in order to interpret the results. However, I have made some points below that might help. I have also attached a general tutorial for mixed modelling in lme4.

• If you have a predictor with more than 2 levels (e.g., 'Colours' with levels green, blue, pink) then dummy coding is automatically used for contrasts in R (e.g., blue vs pink, blue vs green in the output table). If you want to compare different groups, there are many in-built functions in R that you can use to get different contrasts. See reference 3 below.

References: