How does the REML fitting procedure know when a predictor is nested

Consider an experiment with two categorical predictors, A and B, where the levels of B are nested within those of A, and thus described by the following design information:

A B
1 1
1 1
1 2
1 2
2 3
2 3
2 4
2 4

A = random, B = random

If I treat both A and B to be random effect predictors, then I can specify the correct random effects model using lme4/lmer notation:

R ~ (1|A) + (1|B)


Which results in the following fixed and random effect design matrices:

A = fixed, B = random

If I now treat A as a fixed effect predictor, and B as a random effect predictor, then I can specify the correct mixed effects model using lme4/lmer notation:

R ~ A + (1|B)


Which results in the following fixed and random effect design matrices:

My questions are:

(1) How does lmer's REML fitting procedure "know" that predictor B is nested within A, for either model?

(2.1) When A and B are both random effects, their design columns both appear in the random effect design matrix. How does the REML procedure know to estimate the effects of A first, and then model the offsets of the nested B levels within those of A?

(2.2) In the random effects model, where A and B are both in the random effects design matrix, is the REML fitting procedure effectively performing a partitioning of variance, similar to that of an ANOVA?

(3) In the mixed effects model, where A is fixed and B is random, how does the REML procedure "look across" from the random effects design matrix to the fixed effects design matrix, to know that B's levels are nested relative to A's levels?

The reason why I think this question is interesting to ask specifically for random/mixed effects models, is because the way a random effects design matrix encodes nested predictors seems (to me) to be very different from the way a fixed effect design matrix encodes nested predictors.

In a fixed effect design matrix, only the "offset" levels of the nested predictor are included as columns in the design matrix.

Using the predictor information given above, but now with both A and B as fixed effects, the fixed effects design matrix (for the fixed effect model) would be:

Where only levels 2 and 4 of nested predictor B are even included as columns in the design matrix.

Thus, there is no possible confusion as to how a ML/OLS fitting procedure "knows" that predictor B is nested, since only a selected number levels are included as columns in the fixed effects design matrix.

This is not the case with the random effect design matrices described above, in which all 4 levels of nested predictor B are included as separate columns in the random effects design matrix. Despite this, the REML fitting procedure is still able to "know" that B's levels are nested, and furthermore nested relative to predictor A.

• 1, 2, and 4: You tell it all that in the call to lmer. Commented Jul 29, 2023 at 22:32
• Right, I know that the lmer formula is what is telling the REML procedure what to do, but I am more curious how the REML procedure works under the hood. I added the fixed effect model to compare, and highlight why I think the random effects design matrix is fundamentally different from fixed effects design matrix. Thanks! Commented Jul 29, 2023 at 23:17

I would actually argue that lmer doesn't know they are nested, doesn't estimate the effects of $$A$$ first, and isn't really performing a partitioning of variance.
lmer uses a penalised least squares approach to computation. That is, given values of the variance components, it fits the fixed and random effects by minimising something like $$\| Y-X\beta-Zb\|_2^2+\| u\|_2^2$$ where $$\beta$$ are the fixed effects, $$b$$ are the random effects, and $$u$$ are the underlying 'sphered' random effects distributed as iid $$N(0,\sigma^2)$$. This alternates with estimating the variance components using a profile (restricted) likelihood.
Because lmer doesn't care about nesting, it can be used for problems that aren't nested, including observational clustered or longitudinal data. Conversely, it doesn't take advantage of nesting when it's present in some designed experiments and still does the full iterative estimation.
• If $b\sim N(0,\sigma^2V)$ the sphered random effects are $V^{-1/2}u$, so $Zb$ is $ZV^{1/2}u$ Commented Jul 30, 2023 at 1:17