方差分析和回归在 R 中给出相反的结果

机器算法验证 r 回归 方差分析 相互作用
2022-03-30 02:05:15

我在因子设计中进行了一项实验：我测量了三种食草动物处理和六种营养处理的光 (PAR)。实验被阻止了。

我已经按如下方式运行线性模型（您可以从我的网站下载数据进行复制）

dat <- read.csv('http://www.natelemoine.com/testDat.csv')
mod1 <- lm(light ~ Nutrient*Herbivore + BlockID, dat)

残差图看起来不错

par(mfrow=c(2,2))
plot(mod1)

当我查看 ANOVA 表时，我看到了 Nutrient 和 Herbivore 的主要影响。

anova(mod1)

Analysis of Variance Table 

Response: light 
                    Df  Sum Sq Mean Sq F value    Pr(>F)     
Nutrient             5  4.5603 0.91206  7.1198 5.152e-06 *** 
Herbivore            2  2.1358 1.06791  8.3364 0.0003661 *** 
BlockID              9  5.6186 0.62429  4.8734 9.663e-06 *** 
Nutrient:Herbivore  10  1.7372 0.17372  1.3561 0.2058882     
Residuals          153 19.5996 0.12810                       
--- 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

然而，回归表显示了不显着的主效应和显着的交互作用。

summary(mod1)

Call: 
lm(formula = light ~ Nutrient * Herbivore + BlockID, data = dat) 

Residuals: 
     Min       1Q   Median       3Q      Max  
-0.96084 -0.19573  0.01328  0.24176  0.74200  

Coefficients: 
                           Estimate Std. Error t value Pr(>|t|)     
(Intercept)                1.351669   0.138619   9.751  < 2e-16 *** 
Nutrientb                  0.170548   0.160064   1.066  0.28833     
Nutrientc                 -0.002172   0.160064  -0.014  0.98919     
Nutrientd                 -0.163537   0.160064  -1.022  0.30854     
Nutriente                 -0.392894   0.160064  -2.455  0.01522 *   
Nutrientf                  0.137610   0.160064   0.860  0.39129     
HerbivorePaired           -0.074901   0.160064  -0.468  0.64049     
HerbivoreZebra            -0.036931   0.160064  -0.231  0.81784     
... 
Nutrientb:HerbivorePaired  0.040539   0.226364   0.179  0.85811     
Nutrientc:HerbivorePaired  0.323127   0.226364   1.427  0.15548     
Nutrientd:HerbivorePaired  0.642734   0.226364   2.839  0.00513 **  
Nutriente:HerbivorePaired  0.454013   0.226364   2.006  0.04665 *   
Nutrientf:HerbivorePaired  0.384195   0.226364   1.697  0.09168 .   
Nutrientb:HerbivoreZebra   0.064540   0.226364   0.285  0.77594     
Nutrientc:HerbivoreZebra   0.279311   0.226364   1.234  0.21913     
Nutrientd:HerbivoreZebra   0.536160   0.226364   2.369  0.01911 *   
Nutriente:HerbivoreZebra   0.394504   0.226364   1.743  0.08338 .   
Nutrientf:HerbivoreZebra   0.324598   0.226364   1.434  0.15362     
--- 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.3579 on 153 degrees of freedom 
Multiple R-squared:  0.4176,    Adjusted R-squared:  0.3186  
F-statistic: 4.219 on 26 and 153 DF,  p-value: 8.643e-09

我知道这个问题之前已经在多个帖子中被问过和回答。在之前的帖子中，问题围绕 anova() 和 lm() 中使用的不同类型的 SS。但是，我认为这不是这里的问题。首先，设计是平衡的：

with(dat, tapply(light, list(Nutrient, Herbivore), length))

其次，使用 Anova() 选项不会更改 anova 表。这并不奇怪，因为设计是平衡的。

Anova(mod1, type=2)
Anova(mod1, type=3)

改变对比度不会改变结果（定性地）。我仍然从 anova() 与 summary() 中得到了相当多的反向解释。

options(contrasts=c("contr.sum","contr.poly"))
mod2 <- lm(light ~ Nutrient*Herbivore + BlockID, dat)
anova(mod2)
summary(mod2)

我很困惑，因为我在回归上阅读的所有内容都与 ANOVA 不一致，这暗示了 R 将 SS 用于 summary() 和 anova() 函数的方式存在差异。但是，在平衡设计中，SS 类型是等价的，这里的结果不会改变。根据我使用的输出，我怎样才能有完全相反的解释？

2个回答

本质上，问题是，为什么线性模型中的一个系数与 0 显着不同，但 ANOVA 显示没有显着影响，反之亦然。

为此，让我们考虑一个更简单的例子。

set.seed( 123 )
data <- data.frame( x= rnorm( 100 ), g= rep( letters[1:10], each= 10 ) )
data$x[ data$g == "d" ] <- data$x[ data$g == "d" ] + 0.5
boxplot( x ~ g, data )
l <- lm( x ~ 0 + g, data )
summary( l )
anova( l )

您可以看到只有一组 (d) 脱颖而出（系数显着不同于零）。但是，鉴于其他九个组没有显示效果，方差分析返回。但是，让我们删除一些组：p>0.1

data2 <- data[ data$g %in% c( "a", "d" ), ]
anova( lm( x ~ 0 + g, data2 )

返回

          Df  Sum Sq Mean Sq F value  Pr(>F)  
g          2  6.8133  3.4066  5.7363 0.01182 *
Residuals 18 10.6898  0.5939

ANOVA 考虑组内和组间的总体方差。在第一种情况（10 组）中，组之间的差异较小，因为许多组没有效果。第二，只有两组，所有组间方差都来自这两组之间的差异。

反过来呢？这更容易：想象三个组的均值等于 -1、0、1。总平均值为 0。每个组单独与 0 不一定有显着差异，但组 1 和组 3 之间有足够的差异来解释显着性组间总方差。

这里发生的是多重比较问题。您有 10 个 df 用于交互，因此您可以查看 10 个独立的交互效果（尽管我怀疑您实际查看的 10 个，即 10 个回归效果，不是独立的）。

当且仅当 Scheffe 的多重比较检验可以找到在该水平上显着的交互作用效应时，10 df 交互作用检验才会在某个水平上显着。因此，使用 Scheffe 方法，您将无法找到显着的交互回归系数。所报告的回归系数的 P 值等同于查看 Fisher 的 LSD 多重比较方法，众所周知，这种方法更容易声明显着性。所以基本上你有一种方法声明没有效果，另一种方法发现一些效果，但由于它们是不同的方法，这并不奇怪。您需要决定要使用的标准。（除非整体测试显着，否则更复杂的 LSD 使用不会查看单个系数。）

另一种思考方式是 10df 交互测试是 10 个 1df 测试的平均值，如果交互效果不是很显着，它们可能会在平均它们的过程中丢失。但是，如果您单独查看它们，您可以看到它们的效果。

我不会讨论主要影响问题。但我认为 R 最强烈地告诉你关于相互作用的内容 (P=.00513) 是使用营养素 a 和 d 的不同效果取决于你使用的是未命名的草食动物还是配对的草食动物。如果 a 和 d 的微分效应可以改变，那么对 a 和 d 一定有一些影响，但是 Nut 的回归系数。d （它确实看到了它们的差异）似乎是在说没有——但它似乎只是在说，因为存在交互作用的主要影响变得如此复杂，以至于它们不值得试图弄清楚。

其它你可能感兴趣的问题
上一篇如何在预测中使用 rpart 的结果 下一篇每个受试者需要多少次观察才能在混合模型中拟合随机斜率？