Lmer模型的标准化系数

Question

我曾经使用下面的代码来计算lmer模型的标准化系数。然而，随着新版lme4的发布，返回对象的结构已经发生了变化。

如何调整stdCoef.lmer函数，使其与新版lme4兼容？

# 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)

- majom

Answer 1

@LeonardoBergamini的回答是可行的，但这个回答更加简洁易懂，并且只使用了标准访问器，如果/当summary()输出的结构或已拟合模型的内部结构发生改变时，它更不容易出错。

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

(这会得到与@LeonardoBergamini答案中的stdCoef.lmer相同的结果...)

你可以使用broom.mixed::tidy+dotwhisker::by_2sd获得部分缩放系数——被x的2倍标准差缩放，但未缩放y的标准差，并且未居中：

library(broom.mixed)
library(dotwhisker)
(fm1 
    |> tidy(effect="fixed") 
    |> by_2sd(data=sleepstudy) 
    |> dplyr::select(term, estimate, std.error)
)

- Ben Bolker

Answer 2

这应该可以工作：

stdCoef.lmer <- function(object) {
  sdy <- sd(attr(object, "resp")$y) # the y values are now in the 'y' slot 
  ###                                 of the resp attribute
  sdx <- apply(attr(object, "pp")$X, 2, sd) # And the X matriz is in the 'X' slot of the pp attr
  sc <- fixef(object)*sdx/sdy
  #mimic se.ranef from pacakge "arm"
  se.fixef <- function(obj) as.data.frame(summary(obj)[10])[,2] # last change - extracting 
  ##             the standard errors from the summary
  se <- se.fixef(object)*sdx/sdy
  return(data.frame(stdcoef=sc, stdse=se))
}

- Leonardo Bergamini