Kruskal-Wallis 内置 R 函数与手动计算略有不一致

机器算法验证 r 克鲁斯卡尔-沃利斯测试 领带
2022-03-14 10:35:16

我对以下内容感到困惑,我无法在其他地方找到答案。

我正在尝试在做一些统计数据的同时学习 R,并且,作为练习,我尝试通过在 R 中“手动”执行这些操作来仔细检查内置 R 函数的结果。但是,对于 Kruskal-Wallis 测试,我不断得到不同的结果,我不知道为什么。

例如,我正在查看在练习中分发的以下数据

activity <- c(2, 4, 3, 2, 3, 3, 4, 0, 4, 3, 4, 0, 0, 1, 3, 1, 2, 0, 3, 1, 0, 3, 4, 0, 1, 2, 2, 2, 3, 2) 
group <- c(rep("A", 11), rep("B", 10), rep("C", 9))
group <- factor(group)
data.raw <- data.frame(activity, group)

我想按组分析活动。首先,我使用内置的 R 函数运行 Kruskal-Wallis 测试

kruskal.test(activity ~ group, data = data.raw)

返回H=8.9056

为了仔细检查,我尝试在 R 中“手动”执行相同的操作,并使用以下(无疑是无助的)代码

rank <- rank(activity)
data.rank <- data.frame(rank, group)
rank.sum <- aggregate(rank ~ group, data = data.rank, sum)

x <- rank.sum[1,2]^2 / 11 + rank.sum[2,2]^2 / 10 + rank.sum[3,2]^2 / 9
H <- (12 / (length(activity) * (length(activity) + 1))) * x - 3 * (length(activity) + 1)
H

这旨在反映以下公式:

H=12N(N+1)i=1g(Ri2ni)3(N+1)

其中是观察总数,是组数,是第组中的观察数,组的秩和NgniiRii

现在我得到,这让我更加困惑,这也是为相关练习提供的答案。我已经对几个不同的数据集进行了尝试,并且我倾向于使用内置函数获得稍高的值。H=8.499H

我试图找出我做错了什么或未能理解,但无济于事。谁能帮我理解为什么内置kruskal.test函数返回的值与我通过拼写得到的值不同?

1个回答

kruskal.test应用此 Wikipedia 文章(第 4 点)中描述的关系更正

如果使用前一点中描述的快捷公式,可以通过将 H 除以 , ...1i=1G(ti3ti)N3N

继续您的代码:

TIES <- table(activity)
H / (1 - sum(TIES^3 - TIES)/(length(activity)^3 - length(activity)))
#[1] 8.9056

您可以通过仔细研究代码来了解 R 函数的作用,您可以使用getAnywhere(kruskal.test.default).