多独立样本秩检验:Kruskal-Wallis检验的理论与实践
一、引言
在统计学中,当数据不满足正态分布或方差齐性假设时,传统的参数检验(如方差分析ANOVA)可能失效。此时,非参数检验方法(如秩检验)成为更可靠的选择。本文将详细介绍多独立样本秩检验的核心方法——Kruskal-Wallis检验,包括其理论基础、公式推导、案例分析及Python实现。
二、理论基础
1. 问题定义
假设我们有 k k k 个独立样本 X 1 , X 2 , … , X k X_1, X_2, \dots, X_k X1,X2,…,Xk,每个样本包含 n i n_i ni 个观测值。我们需要检验这些样本是否来自相同的总体分布。
原假设:所有样本的总体分布相同
备择假设:至少有一个样本的总体分布不同
2. Kruskal-Wallis检验原理
Kruskal-Wallis检验通过以下步骤实现:
- 数据排序:将所有样本数据合并后统一排序,赋予秩次(相同数值取平均秩)
- 计算秩和:分别计算每个样本的秩和 R i R_i Ri
- 构造统计量:计算H统计量并进行显著性检验
H统计量公式:
H
=
12
N
(
N
+
1
)
∑
i
=
1
k
R
i
2
n
i
−
3
(
N
+
1
)
(
1
)
H = \frac{12}{N(N+1)} \sum_{i=1}^k \frac{R_i^2}{n_i} - 3(N+1) \quad (1)
H=N(N+1)12i=1∑kniRi2−3(N+1)(1)其中
N
=
∑
i
=
1
k
n
i
N = \sum_{i=1}^k n_i
N=∑i=1kni
打结处理:当存在相同数值时,使用修正公式:
H
a
d
j
=
H
1
−
∑
j
=
1
g
(
t
j
3
−
t
j
)
N
3
−
N
(
2
)
H_{adj} = \frac{H}{1 - \frac{\sum_{j=1}^g (t_j^3 - t_j)}{N^3 - N}} \quad (2)
Hadj=1−N3−N∑j=1g(tj3−tj)H(2)其中
t
j
t_j
tj 是第
j
j
j 个打结组的大小,
g
g
g 是打结组数
3. 渐近分布
当样本量较大时,H统计量近似服从自由度为 k − 1 k-1 k−1 的卡方分布。对于小样本,可采用精确检验或蒙特卡洛模拟。
三、数学推导
1. 秩和的期望与方差
假设所有样本来自同一总体,则每个样本的秩和期望为:
E
(
R
i
)
=
n
i
(
N
+
1
)
2
(
4
)
E(R_i) = \frac{n_i(N+1)}{2} \quad (4)
E(Ri)=2ni(N+1)(4)方差为:
Var
(
R
i
)
=
n
i
(
N
+
1
)
(
N
−
n
i
)
12
(
5
)
\text{Var}(R_i) = \frac{n_i(N+1)(N - n_i)}{12} \quad (5)
Var(Ri)=12ni(N+1)(N−ni)(5)
2. 统计量的构造
通过标准化秩和差,构造服从卡方分布的统计量:
H
=
12
N
(
N
+
1
)
∑
i
=
1
k
(
R
i
−
E
(
R
i
)
)
2
Var
(
R
i
)
(
6
)
H = \frac{12}{N(N+1)} \sum_{i=1}^k \frac{(R_i - E(R_i))^2}{\text{Var}(R_i)} \quad (6)
H=N(N+1)12i=1∑kVar(Ri)(Ri−E(Ri))2(6)展开后简化为式(1)。
3. 打结修正的推导
当存在打结时,总体方差会减小,修正因子为:
Correction
=
1
−
∑
j
=
1
g
(
t
j
3
−
t
j
)
N
3
−
N
(
7
)
\text{Correction} = 1 - \frac{\sum_{j=1}^g (t_j^3 - t_j)}{N^3 - N} \quad (7)
Correction=1−N3−N∑j=1g(tj3−tj)(7)因此,修正后的统计量为:
H
a
d
j
=
H
Correction
(
8
)
H_{adj} = \frac{H}{\text{Correction}} \quad (8)
Hadj=CorrectionH(8)
四、案例分析
1. 数据背景
某医院测试三种药物对患者的疗效评分(满分10分),数据如下:
药物A | 药物B | 药物C |
---|---|---|
7 | 5 | 3 |
8 | 6 | 4 |
9 | 7 | 5 |
6 | 8 | 2 |
2. 计算步骤
- 合并排序:
原始数据:2,3,4,5,5,6,6,7,7,8,8,9 秩次: 1,2,3,4.5,4.5,6.5,6.5,8.5,8.5,10.5,10.5,12
- 计算秩和:
R A = 8.5 + 10.5 + 12 + 6.5 = 37.5 ( 9 ) R B = 4.5 + 6.5 + 8.5 + 10.5 = 30 ( 10 ) R C = 1 + 2 + 3 + 4.5 = 10.5 ( 11 ) R_A = 8.5 + 10.5 + 12 + 6.5 = 37.5 \quad (9) \\ R_B = 4.5 + 6.5 + 8.5 + 10.5 = 30 \quad (10) \\ R_C = 1 + 2 + 3 + 4.5 = 10.5 \quad (11) RA=8.5+10.5+12+6.5=37.5(9)RB=4.5+6.5+8.5+10.5=30(10)RC=1+2+3+4.5=10.5(11)3. 计算H统计量:
H = 12 12 × 13 ( 37. 5 2 4 + 3 0 2 4 + 10. 5 2 4 ) − 3 × 13 = 6.857 ( 12 ) H = \frac{12}{12 \times 13} \left( \frac{37.5^2}{4} + \frac{30^2}{4} + \frac{10.5^2}{4} \right) - 3 \times 13 = 6.857 \quad (12) H=12×1312(437.52+4302+410.52)−3×13=6.857(12)4. 打结修正:
∑ ( t j 3 − t j ) = ( 2 3 − 2 ) + ( 2 3 − 2 ) + ( 2 3 − 2 ) + ( 2 3 − 2 ) = 6 + 6 + 6 + 6 = 24 ( 13 ) H a d j = 6.857 1 − 24 1 2 3 − 12 = 6.857 0.996 = 6.884 ( 14 ) \sum (t_j^3 - t_j) = (2^3-2)+(2^3-2)+(2^3-2)+(2^3-2) = 6+6+6+6=24 \quad (13) \\ H_{adj} = \frac{6.857}{1 - \frac{24}{12^3 - 12}} = \frac{6.857}{0.996} = 6.884 \quad (14) ∑(tj3−tj)=(23−2)+(23−2)+(23−2)+(23−2)=6+6+6+6=24(13)Hadj=1−123−12246.857=0.9966.857=6.884(14)5. 显著性检验:自由度 d f = 2 df=2 df=2,临界值 χ 0.05 2 ( 2 ) = 5.991 \chi^2_{0.05}(2)=5.991 χ0.052(2)=5.991,拒绝原假设( p = 0.032 p=0.032 p=0.032)
3. Python实现
import numpy as np
from scipy.stats import kruskal
data = [
[7,8,9,6],
[5,6,7,8],
[3,4,5,2]
]
# 手动计算秩次
combined = np.concatenate(data)
sorted_indices = np.argsort(combined)
ranks = np.empty_like(sorted_indices)
ranks[sorted_indices] = np.arange(1, len(combined)+1)
# 处理打结
unique_values = np.unique(combined)
for val in unique_values:
indices = np.where(combined == val)[0]
if len(indices) > 1:
avg_rank = np.mean(ranks[indices])
ranks[indices] = avg_rank
# 计算秩和
n = [len(d) for d in data]
ranks_split = np.split(ranks, np.cumsum(n)[:-1])
R = [np.sum(r) for r in ranks_split]
# 计算H统计量
N = sum(n)
H = (12 / (N*(N+1))) * np.sum([(ri**2)/ni for ri, ni in zip(R, n)]) - 3*(N+1)
# 打结修正
def tie_correction(ranks):
unique, counts = np.unique(ranks, return_counts=True)
return 1 - np.sum(counts**3 - counts) / (N**3 - N)
correction = tie_correction(ranks)
H_adj = H / correction
print(f"修正前H统计量: {H:.3f}, 修正后H统计量: {H_adj:.3f}")
五、方法比较
方法 | 适用场景 | 优点 | 缺点 |
---|---|---|---|
Kruskal-Wallis | 多独立样本非参数检验 | 无需分布假设 | 小样本效率较低 |
ANOVA | 正态分布且方差齐性的参数检验 | 统计效力高 | 对分布假设敏感 |
Friedman | 相关样本检验 | 处理重复测量数据 | 依赖秩次转换 |
六、注意事项
- 样本量要求:每个样本至少包含5个观测值
- 多重比较:显著结果后需进一步进行两两比较(如Dunn检验)
z i j = R ˉ i − R ˉ j N ( N + 1 ) 12 ( 1 n i + 1 n j ) ( 15 ) z_{ij} = \frac{\bar{R}_i - \bar{R}_j}{\sqrt{\frac{N(N+1)}{12} \left( \frac{1}{n_i} + \frac{1}{n_j} \right)}} \quad (15) zij=12N(N+1)(ni1+nj1)Rˉi−Rˉj(15)3. 数据预处理:异常值可能影响秩次分配,需谨慎处理
七、参考文献
- Hollander, M., & Wolfe, D. A. (1999). Nonparametric Statistical Methods. Wiley.
- Conover, W. J. (1999). Practical Nonparametric Statistics. Wiley.
- Kruskal, W. H., & Wallis, W. A. (1952). Use of ranks in one-criterion variance analysis. Journal of the American Statistical Association, 47(260), 583-621.
- Lehmann, E. L. (1998). Nonparametrics: Statistical Methods Based on Ranks. Springer.
八、总结
Kruskal-Wallis检验是处理多独立样本非参数检验的重要工具。通过秩次转换,它能够有效克服数据分布的限制,在医学、社会科学等领域有广泛应用。实际应用中需注意数据预处理和多重比较问题,结合专业背景进行结果解读。