在 R 中运行统计模型(例如,
glm
、lm
、lme4::lmer
等)后,我使用 summary()
运行 corr=TRUE
命令来获取系数相关性表。它具有给定模型中截距和自变量的相关矩阵。我想做的是提取使这些相关性成为可能的列。例如,如果模型是 x ~ a + b
,并且表格看起来像这样......
(Intercept) a
a 0.##
b 0.## 0.##
我想提取允许
a
和 b
之间相关的列。
我已经找到了各种命令来提取拟合值和残差,但还没有任何命令可以为我提供与该表中的相关性相对应的
a
和 b
值。
您想要的是
summary(.)$correlation
。
如果您想从
object
中提取某些内容,请查看结构 str(object)
以找到所需的元素。 (如果您在 RStudio 中工作,也许可以从 str(object, max.levels=1)
开始,因为如果输出太大,它往往会挂起或崩溃。)
查看
fitted()
等 stats:::fitted.lm
方法的源代码发现,它们本质上是提取相应对象的 $fitted.values
元素(类似适用于 resid()
)。
我们可以用同样的方式编写函数
ex_corr_lm
。我们可以使用公式元素上的 all.vars()
自动获取 RHS 变量。
> ex_corr_lm <- \(x, vars) {
+ if (missing(vars)) {
+ vars <- all.vars(x$call$formula)[-1]
+ }
+ if (inherits(x, c("summary.lm", "summary.glm"))) {
+ x$correlation[vars, vars]
+ } else if (inherits(x, "summary.merMod")) {
+ as.matrix(x$vcov@factors$correlation)[vars, vars]
+ } else {
+ stop('not implemented.')
+ }
+ }
>
> ex_corr_lm(s_lm)
a b
a 1.0000000 -0.4798395
b -0.4798395 1.0000000
>
> ex_corr_lm(s_lm, vars=c('(Intercept)', 'a', 'b')) ## explicitly specifying vars
(Intercept) a b
(Intercept) 1.00000000 -0.7514856 -0.08217174
a -0.75148557 1.0000000 -0.47983950
b -0.08217174 -0.4798395 1.00000000
>
> ex_corr_lm(s_glm)
a b
a 1.0000000 -0.4798395
b -0.4798395 1.0000000
>
> ex_corr_lm(s_lmer)
a b
a 1.000000 -0.288444
b -0.288444 1.000000
为了获得特定值,我们可以进行子集化。
> ex_corr_lm(s_lm)['a', 'b']
[1] -0.4798395
数据:
> set.seed(42)
> d <- data.frame(a=1:10, b=(1:10)^2);d$x <- .5 + 2*d$a + .1*d$b + rnorm(10,,.5)
> s_lm <- summary(f <- lm(x ~ a + b, d,), corr=TRUE)
> s_glm <- summary(f <- glm(x ~ a + b, d, fam=gaussian()), corr=TRUE)
> s_lmer <- summary(lme4::lmer(x ~ a + b + (1|a), d), corr=TRUE)