# Linear mixed model in R; modelling fixed effects with multiple levels and interactions. Help!

I am new to R and to mixed linear modelling. I have a dataset with variables from a cross-sectional study looking at fractional anisotrophy (a property of the brains white matter) in 6 different white matter fibre tracts in the brain. For each fibre tract there are 2 measures (one from each hemisphere). There are 66 participants divided into two groups. We want to control the group comparisons for differences in age, the average FA volume across the brain (wholebrain FA) and tract volume.

I assume that Group (Patients/Controls), Tract (CB/SLF1/SLF2/SLF3/UF/OFST) and Hemisphere (Left/Right) are fixed effects and that Subject (n=66) is a random effect. I also assume that Age, Wholebrain FA and Tract volume should be modelled as fixed effects. For Age and Wholebrain FA there is one value for each participant, but for Volume there is one value for each observation. The attached picture presents the table in the long format. For each subject there are 12 observations.

We hypothesized that there would be group difference in each of the six fibre tracts. We had no a priori assumptions about the hemispheres but would like to explore this post-hoc. We would also like to explore associations between Age and FA in different tracts.

My suggested model look like this

mixed.lmer <- lmer(FA ~ Age + Wholebrain_FA + Volume + GroupTractHemisphere + (1|Subjects), data = DTI)


Question 1: Given that Tract and Hemisphere are assumed to be fixed variables but also are within-subject variables, are they correctly modelled? I am having a hard time understanding how the model "understand" that these variables have multiple levels from the way it is written above.

Question 2: The Volume variable is a within-subject variable whereas Age is a between subject variable. Should they then not be modelled differently?

Question 3: Whether or not to include a three-way interaction is a major debate in my research group. Some saying that for practical purposes its impossible to really make sense of it. Other say it can guide the decision of whether or not to test differences between groups for each hemisphere in each tract. Including a three way interaction to the model is likely to change the results significantly so it seems pretty important to get it right the first time. Any thought on this? Is it being a criminal to include it?

• It is customary not to include "thank you" or other forms or greeting on all stackexchange sites. Although this may seem unusual, it helps to reduce the questions to the main content and improves readability. Feb 19, 2020 at 8:58

Question 1: Given that Tract and Hemisphere are assumed to be fixed variables but also are within-subject variables, are they correctly modelled? I am having a hard time understanding how the model "understand" that these variables have multiple levels from the way its written above.

In most software, such as lme4 or GLMMadaptive it is not necessary to specify at which level a variable varies because, contrary to your understanding, the software really does "know". The level at which a variable varies is a property of the data and it is easy to demonstrate with cross-tabulations.

You may also want to allow a within-subject fixed effect to vary randomly across subjects in which case you can also specify it as a random slope. For example:

lmer(FA ~ Age + Wholebrain_FA + Volume + Hemisphere + (Hemisphere | Subjects)


will estimate a fixed effect for Hemisphere and also allow it to vary by subject. The software will estimate a variance for the "random slope" of Hemisphere.

The difference between the model without random slopes and with random slopes is that in the former, the "within-subject" variable is estimated to have a fixed effect which is the same for all subjects, whereas fitting random slopes allows each subject to have it's own effect of that variable (a global fixed effect and a random offset)

Question 2: The Volume variable is a within-subject variable whereas Age is a between subject variable. Should they then not be modelled differently?

Fixed effects are estimated in the same way regardless of whether they vary within levels of a grouping variable (Subject in your case). This means that the entries in the model matrix of fixed effects will be quite different for within vs. between variables, but this is not something you need to worry about. These kinds of concerns often arise when people come from a traditional ANOVA background.

Question 3: Whether or not to include a three-way interaction is a major debate in my research group. Some saying that for practical purposes its impossible to really make sense of it. Other say it can guide the decision of whether or not to test differences between groups for each hemisphere in each tract. Including a three way interaction to the model is likely to change the results significantly so it seems pretty important to get it right the first time. Any thought on this? Is it being a criminal to include it?

In general there is no problem in interpreting statistical interactions. They have a fairly simple interpretation. This question is too broad to answer. I would suggest posting a new question about this, and including as much detail as possible.

• Dear Robert Long. Thank you very much for taking your time to answer. I am impressed that lme4 understand the variables. I will follow through with your suggestions i relation to Q3. I understand that you cannot give an answer as to whether or not the three way interaction should be modelled or not (or replaced by a groupXhemisphere and hemisphereXtract interaction). However, as a follow-up question I would like to ask; given that we want to test group differences for each tract and hemisphere while controlling for variations in age, wholebrain_FA and volume, is the suggested model correct? Feb 19, 2020 at 6:50
• I am not sure, what you mean by "In most software, such as lme4 or GLMMadaptive it is not necessary to specify at which level a variable varies because, contrary to your understanding, the software really does "know". ". The way it is currently modeled, lmer will only include a random intercept (1|Subject), but no random slopes (e.g. (Hemisphere|Subjects)). If you don't tell lmer that you want random slopes, it won't just assume you need them. Feb 19, 2020 at 9:06
• @LiKao This is concerning fixed effects that either vary within subjects, or between subjects, not about random slopes. I am saying that lme4 understands what a between variable is and what a within variable is. I am not talking about random slopes. In the case of a between variable, random slopes would not make any sense anyway. Feb 19, 2020 at 9:13
• @RobertLong Agreed. For a between-subject variable, a random slope makes absolutely no sense. The reason I pointed this out is, that without a random slope there really is no need for the software to "know" which variables are within and which are between. The formula for within and between are exactly the same in the final regression. Thus, I was not sure what really was the point in the original question. As it might (not sure my guess was right) have been about random slopes, I just felt the need to clarify, that your answer should not be taken to imply random slopes are auto-detected. Feb 19, 2020 at 11:41

You are currently only including a random intercept, that is the base rate of $$FA$$ can differ between subjects. While that is a good first step for a multilevel model, you should usually also include random intercepts, in other words, the effect sizes should be allowed to differ between subjects. This closely relates to your question about within-subject variables.

1. Your model understands that the categorical variables have multiple levels, because these variables can't be interpreted numerically directly. What number would correspond to "left" or "right"? However, multilevel models (basically an advanced form of regression analysis) always require numerical predictors. So whenever R encounters a categorical variable, it automatically recodes it into numerical values using dummy coding. Dummy coding is used because in most cases it is not too wrong. If you need more control over the coding (e.g. choosing the reference class of the dummy coding, centering, ANOVA coding, etc), there are additional methods that you may use.

2. Categorial variables (with multiple levels) and within-subject variables are two quite different concepts you seem to confuse somewhat. If you had included handedness or gender, those would have been between-subject but categorical (with different levels instead of numerically). It's good to keep those two concepts separate.

3. The (1|Subjects) term in your model indicates that the intercept may vary between the subjects (random intercept). So different subjects can have a different base level of $$FA$$ but the effects that switching from right to left side of the brain must be the same for each subject. A simple example: Assume I want to estimate the average weight gain per year for members of my family, so I measure my son and my cat at birth (don't think to much about this, it's just an example for illustration), year 1, etc. So if I then use weight ~ age + (1|subject), this model captures the fact that my cat was much smaller than my son at the time of birth. However, this model does not capture the fact that my son also grew much faster, so it will try to fit the same growth rate for both my cat and my son.

So you should also include random intercepts. Random intercepts can only be used on within-subject variables (think about this for a while, if it isn't obvious why this is the case than re-read the previous until you understand what random-intercept mean).

So you could use a model like FA ~ Age + Wholebrain_FA + Volume + Group + Tract + Hemisphere + (Tract + Hemisphere|Subjects). In this model, the effects of Tract and Hemisphere are also allowed to vary between subjects. If I understand your data correctly, the tracts have been measured for both hemispheres. So you could even nest factors: FA ~ Age + Wholebrain_FA + Volume + Group + Tract + Hemisphere + (Tract|Hemisphere:Subjects) + (Hemisphere | Subjects), which indicates that the effect of tract could vary for each hemisphere within each subject, and the effect of hemisphere can vary between subjects.

Warning: If you don't have enough data, a complicated model often cannot be computed. If you include all random intercepts and random slopes, the model may not converge. Nevertheless, it is often a good idea to start with the most complete model (all random intercepts and slopes) and then work from there to reduce parts of the model until you get convergence. Have a look at the between-subject variances and correlations of between-subject differences to see what should be dropped (often it is advisable to drop the smallest between-subject effects first) until you reach convergence. You should also describe this process in your publication: "We started with a model that included all random intercepts and slopes. Because this model did not converge, we dropped the random slope on Factor XYZ (SD=...). After dropping this random slope, the final model converged."

I found this page to be a very good primer on multilevel models in R.

• Very stimulating answer. With (Tract|Hemisphere:Subjects) I get: "Error: number of observations (=792) <= number of random effects (=792) for term (tract | hemisphere:subjects); the random-effects parameters and the residual variance (or scale parameter) are probably unidentifiable" I suppose there are not enough observations to support the model? If I add tract|subject the model converge. Is it reasonable to add a random effect factor with only two levels (Hemisphere)? stats.stackexchange.com/questions/37647/what-is-the-minimum-recommended-number-of-groups-for-a-random-effects-factor Feb 20, 2020 at 21:59