Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Standardized coefficients for lmer model

I used to use the code below to calculate standardized coefficients of a lmer model. However, with the new version of lme the structure of the returned object has changed.

How to adapt the function stdCoef.lmer to make it work with the new lme4 version?

# Install old version of lme 4
install.packages("lme4.0", type="both",
                 repos=c("http://lme4.r-forge.r-project.org/repos",
                         getOption("repos")[["CRAN"]]))

# Load package
detach("package:lme4", unload=TRUE)
library(lme4.0)

# Define function to get standardized coefficients from an lmer
# See: https://github.com/jebyrnes/ext-meta/blob/master/r/lmerMetaPrep.R
stdCoef.lmer <- function(object) {
  sdy <- sd(attr(object, "y"))
  sdx <- apply(attr(object, "X"), 2, sd)
  sc <- fixef(object)*sdx/sdy
  #mimic se.ranef from pacakge "arm"
  se.fixef <- function(obj) attr(summary(obj), "coefs")[,2]
  se <- se.fixef(object)*sdx/sdy
  return(list(stdcoef=sc, stdse=se))
}

# Run model
fm0 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)

# Get standardized coefficients
stdCoef.lmer(fm0)

# Comparison model with prescaled variables
fm0.comparison <- lmer(scale(Reaction) ~ scale(Days) + (scale(Days) | Subject), sleepstudy)
like image 723
majom Avatar asked Dec 12 '22 04:12

majom


1 Answers

The answer by @LeonardoBergamini works, but this one is more compact and understandable and only uses standard accessors — less likely to break in the future if/when the structure of the summary() output, or the internal structure of the fitted model, changes.

stdCoef.merMod <- function(object) {
  sdy <- sd(getME(object,"y"))
  sdx <- apply(getME(object,"X"), 2, sd)
  sc <- fixef(object)*sdx/sdy
  se.fixef <- coef(summary(object))[,"Std. Error"]
  se <- se.fixef*sdx/sdy
  return(data.frame(stdcoef=sc, stdse=se))
}
library("lme4")
fm1 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
fixef(fm1)
## (Intercept)        Days 
##   251.40510    10.46729
stdCoef.merMod(fm1)
##               stdcoef      stdse
## (Intercept) 0.0000000 0.00000000
## Days        0.5352302 0.07904178

(This does give the same results as stdCoef.lmer in @LeonardoBergamini's answer ...)

You can get partially scaled coefficients — scaled by 2 times the SD of x, but not scaled by SD(y), and not centred — using broom.mixed::tidy + dotwhisker::by_2sd:

library(broom.mixed)
library(dotwhisker)
(fm1 
    |> tidy(effect="fixed") 
    |> by_2sd(data=sleepstudy) 
    |> dplyr::select(term, estimate, std.error)
)
like image 70
Ben Bolker Avatar answered Dec 27 '22 17:12

Ben Bolker