4

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

如何调整功能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)
4

2 回答 2

9

@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

stdCoef.lmer（这确实给出了与@LeonardoBergamini 的答案相同的结果......）

于 2014-10-05T18:59:42.257 回答
4

这应该工作：

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))
}
于 2014-10-02T14:09:01.540 回答