# Diagnose singular fit error in linear mixed effects model

Below is a snippet of real data that I am analyzing for a project. I typically use Bayesian models and am not that familiar with LMEMs but I've been asked to fit them here. I know there is no association between the response and the predictor, but I want to understand if any variation in the predictor is explained by a group factor--or at least develop an intuition for why LMEMs won't work here. Here is a reproducible example.

library(lme4)
response <- c(0.787710494836211, -2.26848461101917, 0.335538079423925, 1.01509294449791,
-0.936751069402168, 0.507819877762778, 0.704504970975397, -0.0943414852337048,
1.95710731281598, -0.627030992426569, -0.28854682169693, -0.852092774084596,
-0.189045885168638, 0.113786773723001, -0.495804931259161, 1.43325791198439,
-1.03114568704162, -0.0715741086870386, 1.29835268746477, 0.242053455499142,
-0.862187705028092, -0.137175939400688, -0.706729051540225, 0.679704540427017,
-0.429808555375643, 1.35848150548549, 0.0350388950514984, -0.750208363378522,
-1.88703297079487, 0.575026753749559, 0.645356897860597, -0.109038367699859,
-0.638219281675638, 1.44008221529131, -0.474565415915845, -0.240759479672387,
0.565007054418043, 0.527573249579339, -1.19209881470121, 1.34022284508866,
0.125736715160556, -0.633348268842764, 0.215920047701383, -0.928636406436253,
0.818675674446728, -0.76956505520075, 0.356773783283443, 1.70417265999518,
-1.57985011976655, -0.256089253779126, 0.742572425688287, 0.894849443061542,
0.728009178885092, -2.9823525423482, 0.284055563418989)
predictor <- c(-0.501621272662627, 0.84601796732652, 0.441726195329776, -0.726227812660818,
-0.321936040664074, 2.95731944330952, 0.08235573133267, -0.366857348663712,
-0.726227812660818, -0.00748688466660642, -0.501621272662627,
-0.726227812660818, -0.726227812660818, 0.306962271330861, -0.726227812660818,
-0.366857348663712, -0.726227812660818, 1.78936543531892, 0.694728391628287,
-0.0510632086690036, 1.04276447176702, 0.943325591727384, -0.697415928926656,
-0.697415928926656, -0.697415928926656, 0.147814551410274, 1.09248391178684,
-0.697415928926656, -0.399099288807739, -0.200221528728462, -0.697415928926656,
-0.697415928926656, -0.697415928926656, -0.697415928926656, -0.697415928926656,
-0.150502088708642, -0.697415928926656, 0.694728391628287, -0.697415928926656,
-0.697415928926656, -0.697415928926656, -0.697415928926656, 0.44613119152919,
-0.448818728827559, -0.34937984878792, 0.794167271667926, 0.147814551410274,
3.82705311287691, -0.697415928926656, -0.697415928926656, 0.296972871469732,
2.68350599242106, -0.697415928926656, 0.197533991430093, 0.44613119152919
)
group <- c("A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A",
"A", "A", "A", "A", "A", "A", "B", "B", "B", "B", "B", "B", "B",
"B", "B", "B", "B", "B", "B", "B", "B", "B", "B", "B", "B", "B",
"B", "B", "B", "B", "B", "B", "B", "B", "B", "B", "B", "B", "B",
"B", "B", "B", "B")

testdat<-data.frame(group=factor(group), predictor=predictor, response=response)
lmer(response ~ predictor + (1|group), data=testdat)


This produces the following result:

boundary (singular) fit: see ?isSingular
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ predictor + (1 | group)
Data: testdat
REML criterion at convergence: 156.7178
Random effects:
Groups   Name        Std.Dev.
group    (Intercept) 1.255e-18
Residual             9.844e-01
Number of obs: 55, groups:  group, 2
Fixed Effects:
(Intercept)    predictor
-8.244e-18    1.760e-01
optimizer (nloptwrap) convergence code: 0 (OK) ; 0 optimizer warnings; 1 lme4 warnings


I'm not worried about alternative ways to analyze the data, but I would appreciate help developing an intuition for why the LMEM failed to fit these data. The real dataset has many more groups (~20) but about the same mean number of samples per group (~23).

The warning occurs because one of variance components (In this case, group) is approximately 0. We can do a quick analysis to identify the root cause.

> testdat%>%group_by(group)%>%summarise(mean=mean(response),sd=sd(response))
# A tibble: 2 x 3
group      mean    sd
<fct>     <dbl> <dbl>
1 A     -1.87e-16    1.
2 B      9.30e-17    1.


The mean of group A and the mean of group B are extremely close (in the magnitude of 10^-17). Therefore, there is not much variability in the group-level. In addition, the mean and sd of group A and B are essentially the same. There is no point to use mixed model since it is impossible to differentiate the two groups. This indicates the result from mixed model would be similar to the result from linear model.

> lm(response ~ predictor, data=testdat)

Call:
lm(formula = response ~ predictor, data = testdat)

Coefficients:
(Intercept)          predictor
9.67631665817e-17  1.75978567043e-01
$$$$
`
• Also, I assume the response is normalized, which is not an issue. The main problem is group A and B are (nearly) identical.
– one
Jul 7, 2022 at 16:27