问题描述
我想要做的是执行一个回归循环,该循环始终具有相同的预测变量但循环响应(此处:y1、y2 和 y3)。问题是我希望也为分组变量的每个类别完成它。在下面的示例数据中,我想对所有三个 y 变量进行回归 y_i=x,这将导致三个回归。但我希望这对 group=a、group=b 和 group=c 分别完成,从而产生 9 个不同的回归(最好存储为列表)。想不出怎么做!任何人都知道如何做到这一点?
到目前为止,我的想法是将 for 循环或 lapply 与 dplyr::group_by 结合使用,但无法使其工作。
set.seed(123)
dat <- data.frame(group=c(rep("a",10),rep("b",rep("c",10)),x=rnorm(30),y1=rnorm(30),y2=rnorm(30),y3=rnorm(30))
解决方法
1) 在 nlme 中使用 lmList
(它随 R 一起提供,因此您不必安装它)。
library(nlme)
regs <- lmList(cbind(y1,y2,y3) ~ x | group,dat)
给出一个 lmList
对象,每个对象都有一个组件。我们展示了组 a
的组件,其他组类似。
> regs$a
Call:
lm(formula = object,data = dat,na.action = na.action)
Coefficients:
y1 y2 y3
(Intercept) 0.2943 0.1395 0.4539
x 0.3721 -0.2206 -0.2255
2) 另一种方法是执行一个整体 lm
,给出具有与上述相同系数的 lm
对象。
lm(cbind(y1,y3) ~ group + x:group + 0,dat)
3) 我们也可以使用几种列表理解包中的一种。这给出了 9 个组件的列表。组件的名称标识了所使用的组合,就像每个主要组件中的调用组件(显示在输出的 Call: 行中)一样。请注意,当前的 CRAN 版本是 0.1.0,但下面的代码依赖于 listcompr 0.1.1,它可以从 github 获得,直到它被放到 CRAN 上。
# install.github("patrickroocks/listcompr")
library(listcompr)
packageVersion("listcompr") # need version 0.1.1 or later
regs <- gen.named.list("{y}.{g}",do.call("lm",list(reformulate("x",y),quote(dat),subset = bquote(dat$group == .(g)))
),y = c("y1","y2","y3"),g = unique(dat$group)
)
如果您不介意输出中的 Call: 行描述性较差,那么可以将其简化为:
gen.named.list("{y}.{g}",lm(reformulate("x",dat,subset = group == g),g = unique(dat$group))
注意
从有两个 y2 的问题中纠正的输入。
set.seed(123)
dat <- data.frame(group=c(rep("a",10),rep("b",rep("c",10)),x=rnorm(30),y1=rnorm(30),y2=rnorm(30),y3=rnorm(30))