Getting singular fit error on lmer model after standardizing the response variable

Question

I'm running a mixed model with the lmerfunction in R, and am running into an issue with singular fits. My dataset is comprised of 48,538 observations of sleep duration values (in hours) from 220 participants (identifier is snapshot_id) collected over the course of a year. If I run a baseline model using the raw sleep duration metrics, my model fits fine:

mod_RQ1_base <- lmer(adjusted_sleep_duration ~ (1|snapshot_id), REML = F, data=sleep_daily_summary_complete_year_meq_week, na.action = na.exclude)
jtools::summ(mod_RQ1_base, scale = TRUE, digits = 5)

However, in my case, I want to build a model to predict the changes in sleep duration based on an individuals baseline (i.e. their mean sleep). To do so, I thought I would standardize the sleep duration metric by subtracting the mean sleep duration and dividing by the standard deviation of sleep duration for that participant. However, when I try to run the same baseline model with the standardized sleep duration, I get a singular fit error:

mod_RQ1_base <- lmer(standardized_adjusted_sleep_duration ~ (1|snapshot_id), REML = F, data=sleep_daily_summary_complete_year_meq_week, na.action = na.exclude)
jtools::summ(mod_RQ1_base, scale = TRUE, digits = 5)

I get the same error when I add other fixed-effect variables I am planning to test to the models as well, but as I'm running into this error already in the baseline model, I think there is a problem here I am not understanding.

My questions are as follows:

1. Am I getting this error because I normalized sleep duration by participant identifier (snapshot_id) and am adding a random intercept effect on snapshot_id as well?

2. Since I want to model what predicts variation sleep duration that are relative to each person's individual baseline, how can I do this realistically and also avoid the singular fit error using mixed models?

Any insight into why I'm getting this error just when I standardize the response (and why I'm not when I don't standardize the response), and how to model the data properly would be greatly appreciated.

Thanks in advance.

Answer 1

Normalizing the dependent variable as you have does not make sense. It would be one thing to take a regular z-score of it (based on the sample mean and standard deviation). Sometimes people do that to their predictors and outcome to get their coefficients into an effect size metric - 1 standard deviation increase in predictor is associated with XX standard deviation change in outcome.

The multilevel model was designed to analyze outcomes that, because of nesting or other non-nested grouping, have both within- and between-group variance. The proportion of variance that is within/between is quantified by the ICC as reported in your output. But the normalization you carried out removes all the individual variation at level 1 (person's standard deviation) and at between person variation at level 2 (their mean). So what is left over for the model to estimate?

Just to show that I can replicate your results using the sleepstudy data in lme4, I ran models where I transformed the outcome (Reaction) in the ways mentioned and ran so-called empty models.

First the Reaction variable in its raw metric:

library(lme4)
library(dplyr)
data("sleepstudy")
raw <- lmer(Reaction ~ 1 + (1|Subject), sleepstudy)
jtools::summ(raw)

MODEL INFO:
Observations: 180
Dependent Variable: Reaction
Type: Mixed effects linear regression 

MODEL FIT:
AIC = 1910.33, BIC = 1919.91
Pseudo-R² (fixed effects) = 0.00
Pseudo-R² (total) = 0.39 

FIXED EFFECTS:
---------------------------------------------------------
                      Est.   S.E.   t val.    d.f.      p
----------------- -------- ------ -------- ------- ------
(Intercept)         298.51   9.05    32.98   17.00   0.00
---------------------------------------------------------

p values calculated using Satterthwaite d.f.

RANDOM EFFECTS:
------------------------------------
  Group      Parameter    Std. Dev. 
---------- ------------- -----------
 Subject    (Intercept)     35.75   
 Residual                   44.26   
------------------------------------

Grouping variables:
---------------------------
  Group    # groups   ICC  
--------- ---------- ------
 Subject      18      0.39 
---------------------------

Next, the traditional z-score, using the sample mean and standard deviation:

sleepstudy <- sleepstudy %>% mutate(s_mn_reaction=mean(Reaction), s_sd_reaction=sd(Reaction), z_reaction=(Reaction-s_mn_reaction)/s_sd_reaction)

z <- lmer(z_reaction ~ 1 + (1|Subject), sleepstudy)

jtools::summ(z)
MODEL INFO:
Observations: 180
Dependent Variable: z_reaction
Type: Mixed effects linear regression 

MODEL FIT:
AIC = 467.16, BIC = 476.73
Pseudo-R² (fixed effects) = 0.00
Pseudo-R² (total) = 0.39 

FIXED EFFECTS:
-------------------------------------------------------
                    Est.   S.E.   t val.    d.f.      p
----------------- ------ ------ -------- ------- ------
(Intercept)         0.00   0.16     0.00   17.00   1.00
-------------------------------------------------------

p values calculated using Satterthwaite d.f.

RANDOM EFFECTS:
------------------------------------
  Group      Parameter    Std. Dev. 
---------- ------------- -----------
 Subject    (Intercept)     0.63    
 Residual                   0.79    
------------------------------------

Grouping variables:
---------------------------
  Group    # groups   ICC  
--------- ---------- ------
 Subject      18      0.39 
---------------------------

Note that the ICC is the same as in the raw metric model (good). The standard deviation of the random intercepts has changed because of the z-scoring, which could be good or bad depending on one's perspective. Personally, I like to keep dependent variables in their original metric.

Now for the results when I normalize Reaction in the same manner as you proposed:

sleepstudy <- sleepstudy %>% group_by(Subject) %>% mutate(mn_reaction = mean(Reaction), sd_reaction = sd(Reaction)) %>% 
  ungroup() %>% mutate(p_z_reaction = (Reaction-mn_reaction)/sd_reaction) 

p_z <- lmer(p_z_reaction ~ 1 + (1|Subject), sleepstudy)
jtools::summ(p_z)

MODEL INFO:
Observations: 180
Dependent Variable: p_z_reaction
Type: Mixed effects linear regression 

MODEL FIT:
AIC = 501.31, BIC = 510.89
Pseudo-R² (fixed effects) = 0.00
Pseudo-R² (total) = 0.00 

FIXED EFFECTS:
---------------------------------------------------------
                     Est.   S.E.   t val.     d.f.      p
----------------- ------- ------ -------- -------- ------
(Intercept)         -0.00   0.07    -0.00   179.00   1.00
---------------------------------------------------------

p values calculated using Satterthwaite d.f.

RANDOM EFFECTS:
------------------------------------
  Group      Parameter    Std. Dev. 
---------- ------------- -----------
 Subject    (Intercept)     0.00    
 Residual                   0.95    
------------------------------------

Grouping variables:
---------------------------
  Group    # groups   ICC  
--------- ---------- ------
 Subject      18      0.00 
---------------------------

As you can see, all variation between groups is gone, or alternatively, through the by-person standardization, we removed all within-person variation in Reaction. Bottom line is that the normalization you proposed removes the very thing the multilevel model is designed to analyze.