【计算机程序设计艺术 TAOCP】卷2·3.3 统计检验：频数正常，能说明多少问题

对应原书：第2卷《半数值算法》第3版，§3.3，书页35至104。本节包含一般检验方法、经验检验、理论检验与谱检验。本文以卡方、KS和相邻组合为主线，最后解释谱检验的几何含义。

拿到一万项输出后，最容易做的是画直方图。如果每根柱子高度差不多，我们通常会觉得它“挺随机”。

然而，直方图把顺序丢掉了。把同一批数字从小到大排序，直方图一根都不会改变。§3.3 的出发点正是：要检验一种性质，必须先写清楚它对应的概率模型和统计量。

1. 检验给出的结论有边界

假设我们关心独立、一致的十进制数字。一个频数检验可以问：在 $N$ 次观察中，各数字的出现次数是否与每种概率 $1/10$ 相容？

如果结果明显异常，它为“数据不符合这一模型”提供证据。原因可能是实现错误、参数选择、调用方式，也可能是模型下偶然出现的罕见事件。

如果结果不异常，只能说这项检查没有揭示足够明显的偏差。它没有覆盖的模式仍然可能存在。增加样本量能提高对某些偏差的敏感度，但无法把一个频数检查变成全部随机性检查。

尤其要避免反复更换种子，直到出现一次好看的结果。这样的筛选改变了被报告样本的来源。比较多种方法时，应事先固定样本长度、检验规则和报告方式。

2. 卡方统计量为什么要除以期望次数

把结果分成 $k$ 类，第 $i$ 类的理论概率为 $p_i>0$，观察次数为 $O_i$。总样本数为 $N$，期望次数为 $E_i=N{p}_i$。

Pearson 卡方统计量为

$$V=\sum _{i=1}^k\frac{(O_i-E_i{)}^2}{E_i}.$$

分母让“偏差有多大”相对于各类的规模来衡量。同样多出10次，对期望1000次的一类和期望5次的一类，意义显然不同。

展开平方，并利用 $\sum O_i=N$、$\sum p_i=1$，得到计算上的另一种形式：

$$V=\frac{1}{N}\sum _{i=1}^k\frac{O_i^2}{p_i}-N.$$

在观察独立、类别概率事先给定、样本量足以支持近似等条件下，$V$ 的分布可用自由度为 $k-1$ 的卡方分布近似。少掉的一个自由度来自总次数固定。

如果参数是用同一批数据估计出来的，自由度等校准也要随之调整。若某些期望次数很小，直接套用近似表可能不合适；“每类期望至少5次”是常见经验规则，不是所有情况下都足够的定理。

普通拟合优度检验常把过大的 $V$ 视为偏离证据。随机数诊断还有另一个角度：多次实验中的统计量若总是小得异常，也可能表明输出被强行平衡。单次得到0可以发生，不能仅凭这一项就宣布程序错误。

3. 相邻数字要另做检查

考虑400位的串 0101...。两种数字各出现200次，一维卡方统计量恰好为0。

现在把数据分成不重叠的二元组：

$$(X_0,X_1),(X_2,X_3),\dots ,(X_{398},X_{399}).$$

如果原序列真是独立公平比特，这200个二元组也是独立观察，00、01、10、11 的理论概率各为 $1/4$。

实际却只有01，计数为 $(0,200,0,0)$。每类期望50次，因此

$$V=\frac{50^2+150^2+50^2+50^2}{50}=600.$$

一个一维检查完全看不到的缺陷，换一个观察对象就暴露了。

这里特意使用不重叠二元组。如果改成 $(X_0,X_1),(X_1,X_2),\dots$，相邻观察共享变量，不能再不加说明地使用同一个独立多项模型和自由度。重叠计数本身可以研究，但统计量的分布需要另外推导。

原书的经验检验还包括间隔、扑克式分组、集券、排列、运行、极值等检查。它们分别保留了不同的信息：

检验对象	典型问题
落入某一区间的间隔	事件之间的等待长度是否合理
小组内不同值的个数	重复模式是否异常
一组数的相对顺序	某些排列是否过多
连续上升或下降的段	次序中是否有持续趋势
一组数的最大值	多项组合后的分布是否正确

这些检查彼此可能相关。把它们的“通过数”简单相加，并不会自动得到一个有概率意义的质量总分。

4. KS检验比较整条分布函数

直方图需要先决定分箱。对于连续分布，还可以直接比较经验分布函数与理论分布函数。

设排序后的样本为 $x_{(1)}\le \cdots \le x_{(N)}$，理论分布函数为连续的 $F$。定义

$$D_N^{+}=\underset{1\le i\le N}{\max }\left(\frac{i}{N}-F(x_{(i)})\right),$$

$$D_N^{-}=\underset{1\le i\le N}{\max }\left(F(x_{(i)})-\frac{i-1}{N}\right).$$

一个方向看经验阶梯比理论曲线高多少，另一个方向看它低多少。原书使用的统计量是 $K_N^{\pm }=\sqrt{N}{D}_N^{\pm }$。常见的双侧统计量为 $D_N=\max (D_N^{+},D_N^{-})$。

以一致分布为目标时，区间内部有 $F(x)=x$。样本如果大量挤在低值区域，经验曲线会提前上升；挤在高值区域时则相反。

使用标准KS临界值时，独立性、连续性，以及分布是否预先指定都很重要。离散整数、有限网格上的大量重复值，或用样本估计参数的情形，需要相应的校准。计算出一个距离，并不等于已经获得适用的 $p$ 值。

5. 用精确算术复现两个检查

下面故意不输出“随机/不随机”标签，也不计算未经校准的 $p$ 值，只展示统计量捕捉了什么。Fraction 用于让教学例子的结果精确可核对。

from collections import Counter
from fractions import Fraction


def chi_uniform(counts):
    if not counts or any(c < 0 or int(c) != c for c in counts):
        raise ValueError("counts must be nonnegative integers")
    total = sum(counts)
    if total == 0:
        raise ValueError("empty sample")
    expected = Fraction(total, len(counts))
    return sum((c - expected) ** 2 / expected for c in counts)


def ks_uniform(values):
    xs = sorted(Fraction(x) for x in values)
    if not xs or any(x < 0 or x >= 1 for x in xs):
        raise ValueError("sample must lie in [0, 1)")
    n = len(xs)
    plus = max(Fraction(i, n) - x for i, x in enumerate(xs, 1))
    minus = max(x - Fraction(i - 1, n) for i, x in enumerate(xs, 1))
    return plus, minus


bits = [0, 1] * 200
one = Counter(bits)
pairs = Counter(zip(bits[::2], bits[1::2]))
pair_counts = [pairs[a, b] for a in (0, 1) for b in (0, 1)]
print(chi_uniform([one[0], one[1]]))
print(pair_counts, chi_uniform(pair_counts))

grid = [Fraction(2 * i + 1, 8) for i in range(4)]
print(tuple(str(x) for x in ks_uniform(grid)))
assert chi_uniform([10, 10, 10, 10]) == 0
assert chi_uniform([40, 0, 0, 0]) == 120

输出：

0
[0, 200, 0, 0] 600
('1/8', '1/8')

最后一组样本是 $1/8,3/8,5/8,7/8$，恰好落在四个等宽区间的中点。很小的KS距离说明它覆盖区间比较整齐；这四个按公式写出的点仍然不是一次独立随机抽样。

6. 谱检验看的是点落在哪些平面上

经验检验观察生成出来的一段数据。理论检验则尽量从递推结构预先推知某些统计行为。谱检验把二者联系起来：考察线性同余生成器一个周期中，连续 $t$ 项构成的点集。

在二维中画 $(U_n,U_{n+1})$，在三维中画 $(U_n,U_{n+1},U_{n+2})$。如果这些点只能落在少数几张间隔很大的平行平面上，某些区域就会长期空缺。

为什么会有这种结构？对

$$X_{n+1}\equiv a{X}_n+c(\mathrm{mod}m)$$

可写成 $X_{n+1}=a{X}_n+c-q_{n}m$。除以 $m$ 后得到

$$a{U}_n-U_{n+1}=q_n-\frac{c}{m}.$$

右边以整数为间隔变化，二维点于是落在一族平行直线上。谱检验要进一步寻找能够覆盖全部点、且相邻平面间距最大的方向；只看水平和垂直分箱还不够。

原书用平面最大间距的倒数描述相应维度的精度。间隔大意味着某个方向上的分辨能力粗；维度上升时，即使一维频数完美，联合点集也可能显得稀疏而有规则。

“谱”在这里不等于对时间序列做一次普通傅里叶变换。它针对的是同余递推所形成的格结构。不同维度和参数需要各自分析，不能用一个二维散点图替代完整结论。

7. 一份检验记录应保留什么

可复核的报告至少应说明：生成器与参数、种子和采样长度、是否跳步、观察如何分组、采用的理论分布、统计量及其校准方式。若做了许多检验，还应说明如何处理多重比较和重复试验。

应用层面还可以加入有已知答案的小问题。例如，若随机数被用于洗牌，除了看数值直方图，也要检查小规模排列的概率是否正确。某个发生器通过几项通用检查，不会替转换算法承担正确性。

下一篇§3.4讨论的正是这个接口：从一致随机数出发，如何生成整数、非一致分布、排列和样本。

---

阅读定位：一般方法见书页36至52；经验检验见53页起；理论检验从70页开始；谱检验从82页开始。本文程序只计算展示用统计量，未把小例子包装为完整检验套件。上一篇：卷2·3.2；下一篇：卷2·3.4。

专栏：计算机程序设计艺术 TAOCP。底本：Knuth 著，苏运霖译，《计算机程序设计艺术》第2卷《半数值算法》第3版，国防工业出版社，2002年版。