Bartlett’s 等方差检验#

[1]中,研究了维生素 C 对豚鼠牙齿生长的影响。在一项对照研究中,60只受试豚鼠被分为小剂量组、中剂量组和大剂量组,每天分别接受0.5、1.0和2.0毫克维生素C的剂量。42天后,测量了牙齿生长情况。

下面的 small_dosemedium_doselarge_dose 数组记录了三组牙齿生长的微米测量值。

import numpy as np
small_dose = np.array([
    4.2, 11.5, 7.3, 5.8, 6.4, 10, 11.2, 11.2, 5.2, 7,
    15.2, 21.5, 17.6, 9.7, 14.5, 10, 8.2, 9.4, 16.5, 9.7
])
medium_dose = np.array([
    16.5, 16.5, 15.2, 17.3, 22.5, 17.3, 13.6, 14.5, 18.8, 15.5,
    19.7, 23.3, 23.6, 26.4, 20, 25.2, 25.8, 21.2, 14.5, 27.3
])
large_dose = np.array([
    23.6, 18.5, 33.9, 25.5, 26.4, 32.5, 26.7, 21.5, 23.3, 29.5,
    25.5, 26.4, 22.4, 24.5, 24.8, 30.9, 26.4, 27.3, 29.4, 23
])

scipy.stats.bartlett 统计量对样本间方差的差异很敏感。

from scipy import stats
res = stats.bartlett(small_dose, medium_dose, large_dose)
res.statistic
np.float64(0.6654670663030519)

当方差差异较大时,统计量的值往往较高。

我们可以通过将观察到的统计量值与零分布(在三组总体方差相等的零假设下导出的统计量值分布)进行比较,来检验各组方差的不等性。

对于此检验,零分布遵循卡方分布,如下所示。

import matplotlib.pyplot as plt
k = 3  # number of samples
dist = dist = stats.chi2(df=k-1)
val = np.linspace(0, 5, 100)
pdf = dist.pdf(val)
fig, ax = plt.subplots(figsize=(8, 5))

def plot(ax):  # we'll reuse this
    ax.plot(val, pdf, color='C0')
    ax.set_title("Bartlett Test Null Distribution")
    ax.set_xlabel("statistic")
    ax.set_ylabel("probability density")
    ax.set_xlim(0, 5)
    ax.set_ylim(0, 1)

plot(ax)
plt.show()
../../_images/af67e6b59f8b7f4a269711bc3b78463ffe99f1f3e17ae93116d5716e31a37b2a.png

这种比较通过 p 值进行量化:零分布中大于或等于观察到的统计量值的比例。

fig, ax = plt.subplots(figsize=(8, 5))
plot(ax)
pvalue = dist.sf(res.statistic)
annotation = (f'p-value={pvalue:.3f}\n(shaded area)')
props = dict(facecolor='black', width=1, headwidth=5, headlength=8)
_ = ax.annotate(annotation, (1.5, 0.22), (2.25, 0.3), arrowprops=props)
i = val >= res.statistic
ax.fill_between(val[i], y1=0, y2=pdf[i], color='C0')
plt.show()
../../_images/bc151ad760d24a8618e5daa91344709dfe8ec51085d905a57d4157bd5218a6b0.png
res.pvalue
np.float64(0.71696121509966)

如果 p 值“小”——也就是说,从具有相同方差的分布中采样数据,产生如此极端的统计量值的概率很低——这可以作为反对零假设(支持备择假设:各组方差不相等)的证据。请注意

  • 反之则不成立;也就是说,此检验不能用于提供支持零假设的证据。

  • 在分析数据[2]之前,应选择“小”值的阈值,并考虑假阳性(错误地拒绝零假设)和假阴性(未能拒绝错误的零假设)的风险。

  • 小的 p 值并非是效应的证据;相反,它们只能提供“显著”效应的证据,这意味着它们在零假设下不太可能发生。

请注意,当观测值呈正态分布时,卡方分布提供了零分布。对于从非正态总体中抽取的小样本,进行置换检验可能更合适:在零假设下,即所有三个样本都来自同一总体,每个测量值在三个样本中的任何一个中被观测到的可能性是相等的。因此,我们可以通过在许多随机生成的观测值分区中计算统计量来形成随机的零分布,将观测值分配到三个样本中。

def statistic(*samples):
    return stats.bartlett(*samples).statistic
ref = stats.permutation_test(
    (small_dose, medium_dose, large_dose), statistic,
    permutation_type='independent', alternative='greater'
)
fig, ax = plt.subplots(figsize=(8, 5))
plot(ax)
bins = np.linspace(0, 5, 25)
ax.hist(
    ref.null_distribution, bins=bins, density=True, facecolor="C1"
)
ax.legend(['asymptotic approximation\n(many observations)',
           'randomized null distribution'])
plot(ax)
plt.show()
../../_images/f5a7b3161101d954a66af030a3514183f71bea3ac9662aaacc934e5ec3af02c3.png
ref.pvalue  # randomized test p-value
np.float64(0.5362)

请注意,此处计算的 p 值与上方由 scipy.stats.bartlett 返回的渐近近似值之间存在显著差异。尽管可以从置换检验中严格推断出的统计推论是有限的,但在许多情况下,它们仍可能是首选方法[3]

参考文献#