# singular fit issue for simple mixed model

I'm trying to run a simple MLM, but I'm bumping into singular fit warnings:

    y <- Reduce(c,replicate(10, scale(rnorm(700, 0, 1))))
x <- rep(c("A","B"), each=700, times=10)
g <- rep(c("g1", "g2", "g3", "g4", "g5", "g6", "g7", "g8",
"g9", "g10"), each=7000)

df <- data.frame(y=y, x=x, g=g)

m <- lmer(y ~ x + (1|g), data=df)
boundary (singular) fit: see ?isSingular


I believe the reason for this warning is that there seems to be no variation of the random intercept as ranef(m) shows.

I read several answers here in CV suggesting to reduce the complexity of the model, often by removing random slopes. However, this model is already as simple as it can get, and there are no random slopes. I also have over 200 groups in my actual data set so it cannot be that I have too few groups.

What are my options? Should I remove the random intercept entirely?

• This nicely reproducible example is completely balanced, with the same number of observations in each group and the same number of A and B cases in each group. Are your actual data similarly balanced?
– EdM
Nov 26, 2021 at 15:40
• @EdM, yes, my actual data are also balanced Nov 26, 2021 at 18:55
• This example strongly violates a basic assumption of lmer: namely, that of independence of errors. You use the same 1400 responses y in a dataframe of 70000 observations. Could you describe the real-world data you are attempting to model in such an unusual fashion?
– whuber
Nov 26, 2021 at 20:05
• @whuber, sorry, yes you're right. I made a typo and already edited the question. The problem however remains even with different y values Nov 26, 2021 at 23:58
• Your creation of y is so strange that I wonder whether it reflects your intentions. What would be the problem with y <- rnorm(70000) in place of the first line?
– whuber
Nov 27, 2021 at 17:41

Your example model samples fine, if a little slowly (~45 minutes on my laptop), using rstanarm. So one solution is to use a Bayesian MLM/LMM:

library(rstanarm)
sm <- stan_lmer(y ~ x + (1|g), data=df,
chains= 4,
cores = 4)


which recovers the expected values:

                                   mean   sd   10%   50%   90%
(Intercept)                       0.0    0.0  0.0   0.0   0.0
xB                                0.0    0.0 -0.1   0.0   0.0
b[(Intercept) g:g1]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g10]              0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g2]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g3]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g4]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g5]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g6]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g7]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g8]               0.0    0.0  0.0   0.0   0.0
b[(Intercept) g:g9]               0.0    0.0  0.0   0.0   0.0
sigma                             1.0    0.0  1.0   1.0   1.0


This is just using the default priors, which you should tune to your specific question. My experience has been that mildly regularizing priors allow these otherwise singular models in lmer() to sample adequately.