利用平方残差的Levene试验

Question

如何使用平方残差而不是绝对值来进行Levene的测试？

我尝试过在levene.test包中使用lawstat包，在car包中使用leveneTest，这两种方法都使用绝对残差。

其目的是再现SAS输出，默认情况下使用平方残差。

Answer 1

iris.lm <- lm(Petal.Width ~ Species, data = iris)
anova(lm(residuals(iris.lm)^2 ~ iris$Species))
## Analysis of Variance Table
## 
## Response: residuals(iris.lm)^2
##               Df Sum Sq Mean Sq F value  Pr(>F)    
## iris$Species   2  0.100  0.0500    14.8 1.4e-06 ***
## Residuals    147  0.497  0.0034                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

也许它有助于理解这是如何工作的。

正如这里所指出的，Levene's检验只是每一个观测值与其组中心之间距离的一个方差。Levene测试的不同实施方式因其对“距离”和“中心”的定义而不同。

“距离”可以是绝对差异，也可以是平方差异。

“中心”可以指平均值或中间值。

默认情况下，SAS使用平方差和平均值。leveneTest在R的car包中只使用绝对差，默认情况下使用中间值。levene.test在lawstat包中也是如此。

所有四种可能的组合都可以手工完成，如下所示。

require(plyr)
x <- ddply(iris, .(Species), summarize
            , abs.mean = abs(Petal.Width - mean(Petal.Width))
            , abs.median = abs(Petal.Width - median(Petal.Width))
            , squared.mean = (Petal.Width - mean(Petal.Width))^2
            , squared.median = (Petal.Width - median(Petal.Width))^2)

anova(lm(abs.mean ~ Species, data = x)) # Levene's test
## Analysis of Variance Table
## 
## Response: abs.mean
##            Df Sum Sq Mean Sq F value  Pr(>F)    
## Species     2   0.53  0.2648    19.6 2.7e-08 ***
## Residuals 147   1.98  0.0135                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

anova(lm(abs.median ~ Species, data = x)) # Brown-Forsythe test
## Analysis of Variance Table
## 
## Response: abs.median
##            Df Sum Sq Mean Sq F value  Pr(>F)    
## Species     2  0.642   0.321    19.9 2.3e-08 ***
## Residuals 147  2.373   0.016                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

anova(lm(squared.mean ~ Species, data = x)) # default SAS Levene's Test
## Analysis of Variance Table
## 
## Response: squared.mean
##            Df Sum Sq Mean Sq F value  Pr(>F)    
## Species     2  0.100  0.0500    14.8 1.4e-06 ***
## Residuals 147  0.497  0.0034                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

anova(lm(squared.median ~ Species, data = x)) # Who-Knows-Whose Test
## Analysis of Variance Table
## 
## Response: squared.median
##            Df Sum Sq Mean Sq F value  Pr(>F)    
## Species     2  0.096  0.0478    13.6 3.7e-06 ***
## Residuals 147  0.515  0.0035                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

为了显示上面的前两个复制leveneTest：

require(car)
leveneTest(Petal.Width ~ Species, data = iris, center = mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##        Df F value  Pr(>F)    
## group   2    19.6 2.7e-08 ***
##       147                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

leveneTest(Petal.Width ~ Species, data = iris, center = median)
## Levene's Test for Homogeneity of Variance (center = median)
##        Df F value  Pr(>F)    
## group   2    19.9 2.3e-08 ***
##       147                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

由于一个人通常有一个线性模型来处理剩余的(平均值)，所以这样做更方便。

iris.lm <- lm(Petal.Width ~ Species, data = iris)
anova(lm(residuals(iris.lm)^2 ~ iris$Species))
## Analysis of Variance Table
## 
## Response: residuals(iris.lm)^2
##               Df Sum Sq Mean Sq F value  Pr(>F)    
## iris$Species   2  0.100  0.0500    14.8 1.4e-06 ***
## Residuals    147  0.497  0.0034                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

因此，答案在顶部。