从 R 中的同一线性模型中提取特定的交互系数

Question

我对线性模型的特定效果（或参数）感兴趣。例如，以 iris 数据集为例，我会感兴趣：

1. Sepal.Width 如何调节 versicolor 和 virginica 的效果（与参考水平 setosa 相比）。
2. Sepal.Width 和Sepal.Length 在species 的三个级别中的每个级别中的线性关系。

我发现这样做的唯一方法是拟合两个模型，一个带有交互，另一个带有嵌套交互：

fit <- lm(Sepal.Length ~ Species * Sepal.Width, data=iris)
summary(fit)  # The two last lines

fit <- lm(Sepal.Length ~ Species / Sepal.Width, data=iris)
summary(fit)  # The three last lines

尽管如此，将 相同 模型拟合两次似乎相对低效（尤其是在拟合需要很长时间计算的模型时）。

有没有办法从一个模型中获得第一个模型的最后两行和第二个模型的最后三行的等价物？谢谢

Answer 1

如果您对摘要中的系数估计感兴趣，那么从单个模型中，可以得到第一个模型的最后两行和第二个模型的最后三行的等效值通过加法。让我解释一下。

实际上，两个摘要都以不同的方式为您提供相同的信息。对于第一个摘要，这是我的输出：

                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     2.6390     0.5715   4.618 8.53e-06 ***
Speciesversicolor               0.9007     0.7988   1.128    0.261    
Speciesvirginica                1.2678     0.8162   1.553    0.123    
Sepal.Width                     0.6905     0.1657   4.166 5.31e-05 ***
Speciesversicolor:Sepal.Width   0.1746     0.2599   0.672    0.503    
Speciesvirginica:Sepal.Width    0.2110     0.2558   0.825    0.411 

这是从模型中删除 setosa 时的拟合结果。这会在 R 中自动发生，因为当您将一个因子传递给线性模型时，您只需要 n-1 部分来描述它，因为基础级部分（在本例中为 setosa）与其他两个因素。

在第二个摘要拟合中，这就是我所看到的。

                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     2.6390     0.5715   4.618 8.53e-06 ***
Speciesversicolor               0.9007     0.7988   1.128    0.261    
Speciesvirginica                1.2678     0.8162   1.553    0.123    
Speciessetosa:Sepal.Width       0.6905     0.1657   4.166 5.31e-05 ***
Speciesversicolor:Sepal.Width   0.8651     0.2002   4.321 2.88e-05 ***
Speciesvirginica:Sepal.Width    0.9015     0.1948   4.628 8.16e-06 ***

在此模型中，尽管是完全共线的，但基本级别的 setosa 仍保留在模型中。通常，如果要在模型中保留一个完全共线的变量，R 会返回一个错误，因为它会抱怨它无法反转矩阵以找到所需公式的最佳拟合（这是因为矩阵不是全等级）。

但是 R 对模型没问题，因为它成功返回了它。这意味着矩阵是可逆的，因为尽管在交互项中具有因子的所有三个水平，它仍然是一个满秩矩阵。使用这些知识，我认为第二个模型中“缺失”的任何信息都被添加到了第一个模型中。

这正是我所看到的。看下面的方程，方程左边是模型 1 的系数，方程右边是模型 2 的系数。

Sepal.Width = Speciessetosa:Sepal.Width
Sepal.Width + Speciesversicolo:Sepal.Width = Speciesversicolor:Sepal.Width
Sepal.Width + Speciesvirginica:Sepal.Width = Speciesvirginica:Sepal.Width

因此，要从第一个模型的最后两行和第二个模型的最后三行获得系数估计，只需简单地对第一个模型进行汇总，然后执行我在上面写的三个方程来获得存储在最后三行中的“额外”信息。