R中的回归与循环

Question

我需要在 R 中使用 Lm() 运行一个简单的回归。它很简单，因为我只有一个自变量。然而，问题是我需要为数据框中的列的许多依赖项测试这个自变量。

所以基本上我有一个常见的 X 和许多 Y，我需要提取截距和斜率并将它们全部存储在一个数据框中。

在 excel 中，这可以使用截距和斜率函数，然后跨列拖动。我需要在 R 中做一些基本相同的事情，我当然可以运行单独的回归，但要求是我需要在一个循环中运行所有这些，并将截距和斜率的估计值存储在一起。

我还在学习 R，任何关于这方面的帮助都会很棒。谢谢:)

Answer 1

nlme 包中的 lmList 函数就是为此而设计的。

我们以iris 数据集为例：

DF <- iris[, 1:4]
#  Sepal.Length Sepal.Width Petal.Length Petal.Width
#1          5.1         3.5          1.4         0.2
#2          4.9         3.0          1.4         0.2
#3          4.7         3.2          1.3         0.2
#4          4.6         3.1          1.5         0.2
#5          5.0         3.6          1.4         0.2
#6          5.4         3.9          1.7         0.4
#...

首先我们必须重塑它。在此示例中，我们希望 Sepal.Length 作为依赖项，其他列作为预测变量。

library(reshape2)
DF <- melt(DF, id.vars = "Sepal.Length")
#  Sepal.Length    variable value
#1          5.1 Sepal.Width   3.5
#2          4.9 Sepal.Width   3.0
#3          4.7 Sepal.Width   3.2
#4          4.6 Sepal.Width   3.1
#5          5.0 Sepal.Width   3.6
#6          5.4 Sepal.Width   3.9
#...

现在我们可以适应了。

library(nlme)
mods <- lmList(Sepal.Length ~ value | variable, 
               data = DF, pool = FALSE)

我们现在可以提取每个模型的截距和斜率。

coef(mods)
#             (Intercept)      value
#Sepal.Width     6.526223 -0.2233611
#Petal.Length    4.306603  0.4089223
#Petal.Width     4.777629  0.8885803

并获得通常的 t 表：

summary(mods)
# Call:
#   Model: Sepal.Length ~ value | variable 
# Data: DF 
# 
# Coefficients:
#   (Intercept) 
#              Estimate Std. Error  t value      Pr(>|t|)
# Sepal.Width  6.526223 0.47889634 13.62763  6.469702e-28
# Petal.Length 4.306603 0.07838896 54.93890 2.426713e-100
# Petal.Width  4.777629 0.07293476 65.50552 3.340431e-111
#   value 
#                Estimate Std. Error   t value     Pr(>|t|)
# Sepal.Width  -0.2233611 0.15508093 -1.440287 1.518983e-01
# Petal.Length  0.4089223 0.01889134 21.646019 1.038667e-47
# Petal.Width   0.8885803 0.05137355 17.296454 2.325498e-37 

或 R 平方值：

summary(mods)$r.squared
#[1] 0.01382265 0.75995465 0.66902769

但是，如果您需要更高效的东西，您可以将包 data.table 与lm 的主力lm.fit 一起使用：

library(data.table)
setDT(DF)
DF[, setNames(as.list(lm.fit(cbind(1, value), 
                             Sepal.Length)[["coefficients"]]), 
              c("intercept", "slope")), by = variable]
#       variable intercept      slope
#1:  Sepal.Width  6.526223 -0.2233611
#2: Petal.Length  4.306603  0.4089223
#3:  Petal.Width  4.777629  0.8885803

当然，这些模型的 R.squared 值只是 Pearson 相关系数的平方：

DF[, .(r.sq = cor(Sepal.Length, value)^2), by = variable]
#       variable       r.sq
#1:  Sepal.Width 0.01382265
#2: Petal.Length 0.75995465
#3:  Petal.Width 0.66902769