当您的数据不正常且方差可能不同时,如何进行方差分析?

机器算法验证 方差分析 正态分布 方差 数据转换 计数数据
2022-03-24 19:38:25

我希望对 13 个组进行 ANOVA:每组都有不同的样本量。另外,我的组不正常,几乎没有通过 Levene 的测试p=.052. 我已经尝试过通过平方根和 ln 进行转换,但这只会让情况变得更糟。我尝试运行 KWH 测试,但我的数据不符合方差测试。还有什么我可以做的吗?

链接到数据

3个回答

我注意到您的回答被称为“有机体数量”,并且所有值都是非负整数。我怀疑这些是计数数据它们不应被视为正态分布并使用传统的 ANOVA 进行分析。相反,计数 GLM 是合适的。我们可以尝试泊松回归

anova(glm(N.Organisms~as.factor(Group), data=d, family=poisson), test="LRT")
# Analysis of Deviance Table
# Model: poisson, link: log
# Response: N.Organisms
# Terms added sequentially (first to last)
# 
# Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
# NULL                               723     619.31              
# as.factor(Group) 12   78.907       711     540.41 6.669e-12 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

aggregate(N.Organisms~Group, data=d, FUN=function(x){ c(mean=mean(x),variance=var(x)) })
# Group N.Organisms.mean N.Organisms.variance
#     1         4.476190             2.161905
#     2         4.902439             1.790244
#     3         4.523077             2.690865
#     4         4.585714             2.883851
#     5         4.260000             2.767755
#     6         4.387097             3.060814
#     7         3.610169             2.690240
#     8         3.717391             3.540580
#     9         3.692308             2.334842
#    10         3.290909             3.284175
#    11         3.058824             2.256471
#    12         3.064516             2.520360
#    13         3.022222             2.471411

奇怪的是,您的数据似乎分散不足这是相当不寻常的,我不知道该怎么做。为了稳健性,我们可以尝试一个简单的卡方检验,它应该是动力不足的(因为它将计数视为类别),但它也非常重要:

with(d, table(Group, N.Organisms))
#      N.Organisms
# Group  1  2  3  4  5  6  7  8
#     1   0  2  3  7  3  4  2  0
#     2   0  3  2  8 17  5  6  0
#     3   2  4 14 12 14 10  8  1
#     4   2  5 13 15 14  9 10  2
#     5   2  3 14 13  4  7  7  0
#     6   6  1 13 12 10 13  7  0
#     7   7  5 21  9  9  4  4  0
#     8   8  5 10  5  7  9  2  0
#     9   5  5 16  9 10  6  1  0
#     10 13  5 16  6  6  7  2  0
#     11 10  6 20  5  7  2  1  0
#     12 12 10 21  8  5  4  2  0
#     13 22  8 33  7 14  5  1  0

set.seed(1)  # makes the simulated p-value reproducible
chisq.test(with(d, table(Group, N.Organisms)), simulate.p.value=T)
#  Pearson's Chi-squared test with simulated p-value (based on 2000 replicates)
# 
# data:  with(d, table(Group, N.Organisms))
# X-squared = 172.02, df = NA, p-value = 0.0004998

最好的检验可能是Kruskal-Wallis 检验,它类似于非正态数据的单向方差分析。它不会受到欠分散的影响,但会考虑到以下事实:8生物是>7,并通过 Dunn 检验轻松提供成对的事后比较。(Dunn 的 R 测试实现是由我们自己的@Alexis 开发的;有关 Dunn 测试的更多信息,请参见此处此处。)请注意,Kruskal-Wallis 测试并不假设方差相等,尽管您无法做到如果分布的形状不同,则将显着结果解释为中位数的简单变化。

kruskal.test(N.Organisms~as.factor(Group), data=d)
#  Kruskal-Wallis rank sum test
# 
# data:  N.Organisms by as.factor(Group)
# Kruskal-Wallis chi-squared = 99.948, df = 12,
# p-value = 5.699e-16

library(dunn.test)
with(d, dunn.test(N.Organisms, as.factor(Group), method="Holm", kw=FALSE))
#                       Comparison of N.Organisms by group                       
#                                     (Holm)                                     
# Col Mean-|
# Row Mean |          1          2          3          4          5          6
# ---------+------------------------------------------------------------------
#        2 |  -0.955248
#          |     1.0000
#          |
#        3 |  -0.012069   1.270114
#          |     0.4952     1.0000
#          |
#        4 |  -0.112368   1.161274  -0.144722
#          |     1.0000     1.0000     1.0000
#          |
#        5 |   0.586382   1.940375   0.826708   0.974481
#          |     1.0000     0.9420     1.0000     1.0000
#          |
#        6 |   0.216494   1.544994   0.324981   0.473740  -0.514632
#          |     1.0000     1.0000     1.0000     1.0000     1.0000
#          |
#        7 |   2.045019   3.816537   2.906720   3.098461   1.910109   2.556624
#          |     0.8171    0.0043*     0.0950     0.0515     0.9260     0.2537
#          |
#        8 |   1.678281   3.251401   2.309693   2.475995   1.417073   1.990415
#          |     1.0000     0.0316     0.4704     0.3122     1.0000     0.8844
#          |
#        9 |   1.755800   3.400923   2.456288   2.632395   1.522142   2.123496
#          |     1.0000    0.0191*     0.3229     0.2077     1.0000     0.7080
#          |
#       10 |   2.692297   4.589523   3.786058   3.987944   2.754012   3.433308
#          |     0.1774    0.0002*    0.0047*    0.0022*     0.1501    0.0176*
#          |
#       11 |   3.223104   5.206152   4.483634   4.691145   3.432921   4.131515
#          |     0.0342    0.0000*    0.0002*    0.0001*    0.0173*    0.0012*
#          |
#       12 |   3.321963   5.440195   4.741819   4.969661   3.610449   4.365576
#          |     0.0250    0.0000*    0.0001*    0.0000*    0.0093*    0.0004*
#          |
#       13 |   3.487721   5.846371   5.211207   5.479181   3.927490   4.789963
#          |    0.0146*    0.0000*    0.0000*    0.0000*    0.0027*    0.0001*
# 
# Col Mean-|
# Row Mean |          7          8          9         10         11         12
# ---------+------------------------------------------------------------------
#        8 |  -0.394796
#          |     1.0000
#          |
#        9 |  -0.345287   0.059170
#          |     1.0000     1.0000
#          |
#       10 |   0.912191   1.244372   1.223490
#          |     1.0000     1.0000     1.0000
#          |
#       11 |   1.652971   1.936170   1.936942   0.746270
#          |     1.0000     0.8984     0.9232     1.0000
#          |
#       12 |   1.754493   2.038829   2.046213   0.799659   0.016138
#          |     1.0000     0.8086     0.8351     1.0000     0.9871
#          |
#       13 |   1.943622   2.224783   2.246161   0.903324   0.054395   0.039280
#          |     0.9609     0.5611     0.5433     1.0000     1.0000     1.0000

如链接到电子表格中的变量名称所示,数据与 13 组中的生物体数量有关,某种计数数据也是如此。如果我们对数据有更多了解,我们可以在这里做得更好!但我不同意@gung 的回答,即我们应该使用一些计数数据模型,比如泊松回归。我将用图形说明我为什么这么说。

泊松分布的方差等于均值,因此简单的第一个分析是针对均值绘制(经验)方差。将数据读入data.frame(在R中)后,我做了:

> summary(dat)
     Number          Group    
 Min.   :1.000   13     : 90  
 1st Qu.:3.000   1      : 76  
 Median :4.000   4      : 70  
 Mean   :3.826   3      : 65  
 3rd Qu.:5.000   6      : 62  
 Max.   :8.000   12     : 62  
                 (Other):299 

>  s2<- with(dat,tapply(Number,Group,FUN=var))
>  m <- with(dat,tapply(Number,Group,FUN=mean))
> plot(m,s2,ylim=c(1.5,5.5))
> abline(0,1,col="red")

在此处输入图像描述

红线表示方差等于平均值​​。我们看到所有的点都在这条线以下,并且没有太多证据表明方差也随着均值的增加而增加。所以这不是泊松数据。此外,虽然方差有所不同,但它们并没有太大的变化(都在两倍之内),所以通常的等方差 AnovaFtest 可能可以使用,因为它对方差的微小变化具有相当的鲁棒性。米勒的经典著作:“超越方差分析”给出了可接受的四倍。因此,让我们在 R 中尝试 Anova,不管有没有这个假设:

> oneway.test(Number ~ Group,data=dat,var.equal=FALSE)

        One-way analysis of means (not assuming equal variances)

data:  Number and Group
F = 10.265, num df = 11.00, denom df = 265.99, p-value = 1.29e-15

> oneway.test(Number ~ Group,data=dat,var.equal=TRUE)

        One-way analysis of means

data:  Number and Group
F = 9.4014, num df = 11, denom df = 712, p-value = 6.497e-16

两者在这里都给出了基本相同的结论:均值相等的原假设被拒绝。请注意,F比率非常相似。

让我们还显示一个图形分析,一个带有“缺口”的箱线图,显示中位数的置信区间:

> boxplot(Number ~ Group, data=dat,notch=TRUE,col="blue",varwidth=TRUE)
Warning message:
In bxp(list(stats = c(1, 2, 3, 5, 7, 2, 4, 5, 6, 7, 1, 3, 5, 6,  :
  some notches went outside hinges ('box'): maybe set notch=FALSE

我们看到一些置信区间根本不重叠,这再次与 Anova 结果一致。我不会过分担心基于这些数据的 Anova 假设存在的问题,因此结论非常可靠。如果需要,可以使用引导程序进一步分析,但我会将其留给其他人。如果需要进一步分析,我们将需要有关数据的更多信息。

在此处输入图像描述

这在很大程度上是@kjetil 答案的脚注,我倾向于同意这一点。但我有一个额外的图表,它说得更多,不适合评论。

这是一个组合图。对于每个组,都有一个分位数图显示每个组的详细信息,一个以标准形式显示中位数和四分位数的叠加框,最后一个显示平均值的菱形图,为了清晰起见,显示得很大。图片本质上是相当恒定的散布、接近对称且缺乏异常值或其他病态或异常结构的图片。在这些情况下,我很高兴看到方差分析。

在此处输入图像描述

也就是说,如果对方差均值结构进行了一些仔细的处理(例如在标准误差计算中),我不认为泊松回归会对此类数据表现不佳或误导。

大概有(或应该有)其他生物和/或环境预测因子会出现在完整的分析中。