1
$\begingroup$

Why do the standard errors for gamlss are much lower than the standard errors from lmer?

library(lme4)
library(lmerTest)
library(gamlss)

fit <- lmer(Reaction ~ Days + (1|Subject), data = sleepstudy)
summary(fit)
# Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
# Formula: Reaction ~ Days + (1 | Subject)
#    Data: sleepstudy
# 
# REML criterion at convergence: 1786.5
# 
# Scaled residuals: 
#     Min      1Q  Median      3Q     Max 
# -3.2257 -0.5529  0.0109  0.5188  4.2506 
# 
# Random effects:
#  Groups   Name        Variance Std.Dev.
#  Subject  (Intercept) 1378.2   37.12   
#  Residual              960.5   30.99   
# Number of obs: 180, groups:  Subject, 18
# 
# Fixed effects:
#             Estimate Std. Error       df t value Pr(>|t|)    
# (Intercept) 251.4051     9.7467  22.8102   25.79   <2e-16 ***
# Days         10.4673     0.8042 161.0000   13.02   <2e-16 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 
# Correlation of Fixed Effects:
#      (Intr)
# Days -0.371

fit <- gamlss(Reaction ~ Days + re(random = ~ 1|Subject), 
              data = sleepstudy, 
              family = "NO")
summary(fit)
# ******************************************************************
#   Family:  c("NO", "Normal") 
# 
# Call:  gamlss(formula = Reaction ~ Days + re(random = ~1 |  
#                                                Subject), family = "NO", data = sleepstudy, opt = "optim") 
# 
# Fitting method: RS() 
# 
# ------------------------------------------------------------------
#   Mu link function:  identity
# Mu Coefficients:
#   Estimate Std. Error t value Pr(>|t|)    
# (Intercept) 251.4051     4.0759   61.68   <2e-16 ***
#   Days         10.4673     0.7635   13.71   <2e-16 ***
#   ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 
# ------------------------------------------------------------------
#   Sigma link function:  log
# Sigma Coefficients:
#   Estimate Std. Error t value Pr(>|t|)    
# (Intercept)   3.3817     0.0527   64.16   <2e-16 ***
#   ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 
# ------------------------------------------------------------------
#   NOTE: Additive smoothing terms exist in the formulas: 
#   i) Std. Error for smoothers are for the linear effect only. 
# ii) Std. Error for the linear terms maybe are not accurate. 
# ------------------------------------------------------------------
#   No. of observations in the fit:  180 
# Degrees of Freedom for the fit:  18.76598
# Residual Deg. of Freedom:  161.234 
# at cycle:  2 
# 
# Global Deviance:     1728.238 
# AIC:     1765.77 
# SBC:     1825.689 
# ******************************************************************
$\endgroup$

1 Answer 1

1
$\begingroup$

The standard errors of the gamlss are conditioned on the values of the random component. The lmer standard errors integrate out the random component. As stated in Section 5.2 "Functions to obtain standard errors" of Flexible Regression and Smoothing: Using GAMLSS in R

When smoothing terms are fitted, gen.likelihood() considers them as fixed at their fitted values, so the Hessian in this case does not take into account the variability in the fitting of the smoothers.

The gen.likelihood()-function is used to make the vcov matrix from where the standard errors originate. In gamlss random effects are carefully constructed as a special type of smoothing. Hence, random effects are kept at their fixed value and which makes the standard errors for gamlss small.

In far most use cases, one should use the lmer-standard errors.

$\endgroup$
2
  • $\begingroup$ Did you mean that the standard errors for gamlss are underestimated instead of inflated since they do not take into account the variability fitting of the smoothers? $\endgroup$ Jun 7, 2023 at 17:45
  • $\begingroup$ Ah yes, you are correct, thank you! I have edited the post $\endgroup$
    – svendvn
    Jun 8, 2023 at 21:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.