正态性 Shapiro-Wilk 检验#

假设我们希望从测量数据中推断一项医学研究中成年男性的人体重是否不是正态分布 [1]。体重(磅)记录在下面的数组 x 中。

import numpy as np
x = np.array([148, 154, 158, 160, 161, 162, 166, 170, 182, 195, 236])

正态性检验 scipy.stats.shapiro,见 [1][2],首先计算一个统计量,该统计量基于观测值与正态分布的期望次序统计量之间的关系。

from scipy import stats
res = stats.shapiro(x)
res.statistic
0.7888146948631716

对于从正态分布中抽取的样本,此统计量的值往往较高(接近 1)。

通过将观察到的统计量值与零分布进行比较来执行检验:零分布是在零假设(即体重是从正态分布中抽取的)下形成的统计量值的分布。对于此正态性检验,零分布不容易精确计算,因此通常通过蒙特卡罗方法进行近似,即从正态分布中抽取许多与 x 大小相同的样本,并计算每个样本的统计量值。

def statistic(x):
    # Get only the `shapiro` statistic; ignore its p-value
    return stats.shapiro(x).statistic
ref = stats.monte_carlo_test(x, stats.norm.rvs, statistic,
                             alternative='less')
import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=(8, 5))
bins = np.linspace(0.65, 1, 50)

def plot(ax):  # we'll reuse this
    ax.hist(ref.null_distribution, density=True, bins=bins)
    ax.set_title("Shapiro-Wilk Test Null Distribution \n"
                 "(Monte Carlo Approximation, 11 Observations)")
    ax.set_xlabel("statistic")
    ax.set_ylabel("probability density")

plot(ax)
plt.show()
../../_images/582d197e8a9cfa342066dae13235ae905cc3f741841016f9a7373f0927ad70a1.png

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

fig, ax = plt.subplots(figsize=(8, 5))
plot(ax)
annotation = (f'p-value={res.pvalue:.6f}\n(highlighted area)')
props = dict(facecolor='black', width=1, headwidth=5, headlength=8)
_ = ax.annotate(annotation, (0.75, 0.1), (0.68, 0.7), arrowprops=props)
i_extreme = np.where(bins <= res.statistic)[0]
for i in i_extreme:
    ax.patches[i].set_color('C1')
plt.xlim(0.65, 0.9)
plt.ylim(0, 4)
plt.show()
../../_images/cc01e2268ef1fedbd994e9c3b6f3384a47c4a4bab3824e8e03ad3b732c2411fe.png
res.pvalue
0.006703814061898823

如果 p 值“很小”,也就是说,如果从正态分布的总体中采样数据产生如此极端的统计量值的概率很低,则可以将其作为反对零假设的证据,而支持备择假设:体重不是从正态分布中抽取的。请注意

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

  • 将被视为“小”的值的阈值是应该在分析数据之前做出的选择 [3],同时考虑到假阳性(错误地拒绝零假设)和假阴性(未能拒绝错误的零假设)的风险。

参考文献#