主要先修:10.6 中心极限定理及其推广、16.1 统计模型、似然与可识别性、18.3 复合假设、单调似然比与干扰参数。
技术部分另用:17.2 矩估计与最大似然估计、17.4 Fisher信息与Cramér–Rao界、10.3 弱收敛、紧性与连续映射、19.1 随机阶与渐近展开,适用段落见下述说明。
本节分两次学习:可先读18.5.1–18.5.2的已知概率模型,以及18.5.5的Fisher精确检验;拟合参数的18.5.3–18.5.4另需17.2、17.4、10.3和19.1的随机阶。建议在19.3的统一似然理论之后返回这些证明;19.4供延伸比较。 查看全书依赖图与分支路线。
18.5.1 分类计数与Pearson统计量
分类数据记录每次观察属于哪一类。设 \(Y_1,\ldots,Y_n\) 独立同分布,取值于预先确定的 \(k\ge2\) 个类别,概率向量为 \(p=(p_1,\ldots,p_k)^\top\),\(\sum_i p_i=1\)。令 \(N_i\) 为第 \(i\) 类计数,则
本节首先检验完全指定的概率向量 \(p=p^0\),并要求每个 \(p_i^0>0\)。不同计数不独立,因为总数固定;单个类别的计数虽然是二项分布,不能把各类残差当作独立标准正态变量。
记经验比例 \(\widehat p_i=N_i/n\),原假设期望频数为 \(E_i=np_i^0\)。Pearson统计量为
每项把平方偏差按模型期望频数加权。较大的 \(Q_n\) 表示计数向量与原假设概率不相容的程度较高;它并非直接给出原假设为真的概率。
18.5.2 从多元CLT到 \(k-1\) 个自由度
若类别数 \(k\) 固定,原假设向量 \(p^0\) 固定且各分量严格为正,在上述iid抽样模型下,\(n\to\infty\) 时 \(Q_n\Rightarrow\chi^2_{k-1}\)。这是渐近分布,并非多项样本的有限样本精确分布。
证明。以下简写 \(p=p^0\),令 \(D=\operatorname{diag}(p_1,\ldots,p_k)\)、\(s=(\sqrt{p_1},\ldots,\sqrt{p_k})^\top\)。单次观察的类别指示向量均值为 \(p\),协方差为 \(D-pp^\top\)。由多元中心极限定理,
\(s^\top s=1\),故 \(P^\top=P\)、\(P^2=P\),它是到 \(s^\perp\) 的正交投影,秩为 \(k-1\)。总数约束还给出精确关系 \(s^\top Z_n=\sqrt n\sum_i(\widehat p_i-p_i)=0\)。取 \(k\times(k-1)\) 矩阵 \(B\),其列为 \(s^\perp\) 的一组标准正交基,则
正交坐标中的极限正态分量相互独立,其平方和才是相应卡方分布。减少的一个维度来自总数约束,原始 \(k\) 个标准化残差从未被假定独立。连续映射定理完成最后一步。证毕。
这个投影证明与Maryland Stat 701讲义的多项CLT推导对应。有限样本中仍有 \(\mathbb E_0 Q_n=\sum_i(1-p_i)=k-1\),但期望相同不能推出分布相同;习题18.5-1给出最小样本反例。
取假设数据 \(N=(30,20,50)\),样本量一百,检验已知概率 \(p^0=(0.2,0.3,0.5)\)。期望频数为 \(E=(20,30,50)\),于是 \(Q_n=10^2/20+(-10)^2/30=25/3\)。自由度为二,卡方参考分布的右尾为 \(e^{-Q_n/2}\),得到近似p值 \(e^{-25/6}\approx0.01550385\)。在5%水平的渐近检验中拒绝该已知概率向量。这里期望频数较大,但卡方p值仍是近似值;不能把近似有效性改称精确有限样本保证。
18.5.3 拟合参数如何消去切空间方向
下面用\(o_P(1)\)表示依概率趋于零的余项;系统定义与运算见19.1节。证明还使用MLE和Fisher信息。首次学习可先读本节的已知概率模型与Fisher精确检验,完成19.1并结合19.3节的一般似然理论后,再返回这里及下一小节。
现在检验 \(p=p(\theta)\) 是否属于一个正则模型,\(\theta\in\Theta\subset\mathbb R^r\),其中 \(\Theta\) 开放、\(1\le r<k-1\)。设真参数 \(\theta_0\) 为内点,所有 \(p_i(\theta_0)>0\),概率映射在其邻域二次连续可微,Jacobian \(J=Dp(\theta_0)\) 的列秩为 \(r\)。仍令 \(D=\operatorname{diag}(p(\theta_0))\),定义
\(A\) 的列空间是标准化概率的模型切空间,\(A^\top A\) 是这个分类模型的单次Fisher信息。\(H\) 把随机残差投影到模型能通过微小参数变化解释的方向。
在以上条件下,若 \(\widehat\theta\) 一致,且基于同一分类数据满足有效拟合展开
则使用拟合期望频数的Pearson统计量渐近服从 \(\chi^2_{k-1-r}\)。仅仅写出一个含 \(r\) 个参数的估计公式,或者仅有一致性,都不足以保证该展开。
证明与MLE充分条件。对同一组分类计数,令 \(\ell_n(\theta)=\sum_iN_i\log p_i(\theta)\)。若选取一致、内点、满足得分方程的最大似然估计,则 \(\sum_i\partial_\theta p_i=0\) 给出
第二式由有限类别的频率大数定律得到,概率二阶导数求和为零,剩下的一阶导数乘积给出信息矩阵。正概率与局部连续性也使沿 \(\theta_0\) 到 \(\widehat\theta\) 线段平均的Hessian有同一极限。在得分方程中使用积分形式Taylor展开,除以 \(\sqrt n\),并利用 \(A^\top A\) 正定,即得式(18.5-6)。这里一致性和内点解的存在是明示假设,没有把它们由“MLE”名称自动推出。
对任何满足该展开的拟合,Taylor展开 \(p(\widehat\theta)\),再用一致性替换标准化分母,有
因为 \(H\) 是秩 \(r\) 的正交投影,且 \(Hs=0\),投影后的极限协方差为 \(P-H\)。它仍满足对称和幂等,秩为 \(k-1-r\)。按定理18.5-1的正交坐标论证,
因此减少自由度来自有效拟合产生的正交投影,而不是单纯清点参数。证毕。
对连续观测先分组时,分类模型的MLE应使用分组计数似然;原始未分组观测的MLE通常对应不同信息,须重新验证式(18.5-6)。NIST拟合优度说明明确区分了这两种拟合。一致但低效的代入估计可能产生加权卡方极限,习题18.5-2将算出一个缩放因子为四的反例。
仍取 \(N=(30,20,50)\),改为检验 \(p(\theta)=(\theta,\theta,1-2\theta)\),\(0<\theta<1/2\)。计数对数似然为 \((N_1+N_2)\log\theta+N_3\log(1-2\theta)\),其内点最大值为 \(\widehat\theta=(N_1+N_2)/(2n)=0.25\)。拟合期望频数为 \((25,25,50)\),得到 \(Q_n^{\mathrm{fit}}=5^2/25+(-5)^2/25=2\)。模型切空间秩为一,剩余自由度也是一,近似p值约 \(0.15729921\),在5%水平不拒绝。检验这个一参数模型与检验例18.5-1的已知概率向量,是两个不同的原假设。
18.5.4 列联表的独立性模型
对 \(R\times C\) 表,设每位观察对象都具有两个分类变量,且对象之间iid抽样,只有总数 \(n\) 固定。单元格计数 \(N_{ij}\) 组成一个有 \(RC\) 类的多项样本。独立性假设要求 \(p_{ij}=a_i b_j\),其中两组边际概率均为正且分别求和为一;它指两个分类变量独立,并不指各单元格计数独立。
对数似然分解为行总数与列总数的两部分,故独立模型的MLE分别为 \(N_{i+}/n\)、\(N_{+j}/n\)。拟合期望和统计量为
在固定表尺寸、正概率和内点拟合条件下,模型维数为 \((R-1)+(C-1)\)。切空间变化写为 \(b_j\,\delta a_i+a_i\,\delta b_j\);若所有变化为零,分别按行与列求和,便得全部 \(\delta a_i,\delta b_j\) 为零,因此参数方向满秩。由前节定理,
行总数由设计固定时,需要改用各行独立的多项抽样模型分析齐性;计数公式相似不代表抽样模型相同。重复测量、聚类或调查权重也会改变计数协方差,不能直接沿用本节iid证明。
18.5.5 Fisher精确条件检验
小表可以不借助卡方近似。考虑上述多项抽样下的 \(2\times2\) 独立性假设,条件于第一行总数 \(m\) 和第一列总数 \(K\)。若左上角计数为 \(A=a\),整张表被确定为
在独立性原假设下,给定有正概率的边际,
证明。独立模型的多项概率中,关于边际参数的部分为 \(a_1^m a_2^{n-m}b_1^K b_2^{n-K}\),在所有具有同一边际的表中相同。条件化后这些未知参数全部约去,剩下的概率与四个单元格阶乘倒数的乘积成比例。将其写成 \(\binom Ka\binom{n-K}{m-a}\),再用Vandermonde恒等式求和,归一化常数为 \(\binom nm\),即得超几何分布。证毕。
边际可以由设计固定,也可以在当前多项模型中条件于观测值。后一种情况下,条件分布仍然精确;对每组边际控制错误概率后,再取全期望,便控制无条件错误概率。未知边际概率被条件化消去,不能将“条件精确”误解为仅在边际原先都固定时才有效。
若备择是正关联,即总体优势比 \(p_{11}p_{22}/(p_{12}p_{21})>1\),右尾p值取 \(\sum_{a\ge a_{\mathrm{obs}}}h(a)\)。双侧离散检验需要明确“更极端”的定义。本节采用概率排序:
它将原假设下不比观测表更常见的表全部计入,包括并列概率。不能普遍用右尾p值乘二来替代。按概率大小排列有限结果,式(18.5-14)满足条件超均匀性,因此其水平不超过标称值;离散性可能造成保守。关于条件超几何分布可核对Drton等人讲义命题1.1.8,这一双侧约定与R的官方Fisher检验文档一致。
| 列1 | 列2 | 合计 | |
|---|---|---|---|
| 行1 | 4 | 0 | 4 |
| 行2 | 0 | 4 | 4 |
| 合计 | 4 | 4 | 8 |
在假设数据中,\(n=8,m=K=4,A_{\mathrm{obs}}=4\)。五个可能计数 \(0,1,2,3,4\) 的条件概率依次为 \(1,16,36,16,1\) 除以七十。正关联的单侧p值为 \(1/70\approx0.01428571\);概率排序双侧p值包含两端,为 \(2/70=1/35\approx0.02857143\)。按双侧5%水平拒绝独立性,1%水平则不拒绝。
未校正的Pearson统计量使用每格拟合期望二,得到 \(Q_{\mathrm{ind}}=4(2^2/2)=8\)。自由度一卡方近似尾概率约 \(0.00467773\),在1%水平会给出不同判断。两种判断的差异说明,小期望频数下的卡方近似可能失准。观测零计数也不等于原假设中概率为零;本例四格模型概率仍均为 \(1/4\)。
18.5.6 稀疏、零概率与分组的适用边界
定理18.5-1固定了类别数和正概率,因此各期望频数随 \(n\) 一起增长。“每格期望至少五”是常见近似诊断规则,并非达到该数就获得精确水平的定理。总数较大,但多数格的期望仍很小,或者类别数也随样本量增长时,不能只凭总数调用固定维证明。已知概率的小问题可枚举多项分布,列联表可考虑适合抽样模型的条件精确分布;数值模拟需说明校准机制与模拟误差。
结构零表示模型规定某类概率为零。如果观察到该类正计数,其事件在原假设下概率为零;如果计数也为零,应在已知有效支持上重新定义类别数,而不是计算零除零。未知边界参数、结构零列联表或秩不足模型,需要重新分析模型维数和极限分布。不能因为某个观测计数恰好为零,就删除该单元格或机械改变自由度。
连续分布分组检验只比较各区间概率,区间内部的不同分布可能给出完全相同的计数模型。分组最好预先指定;若按数据选择分组,选择步骤也应进入校准。不拒绝拟合优度原假设不证明模型正确;拒绝独立性提供分类关联的证据,不单独证明干预的因果效应。
18.5.7 习题与解答
习题18.5-1 相同期望不等于相同分布
令 \(n=1,k=3,p_i^0=1/3\)。计算Pearson统计量,比较其分布与 \(\chi^2_2\)。若把卡方参考尾概率当作精确p值,它在阈值 \(0.4\) 下是否有效?
唯一被观察的类别计数为一,其余为零,所以统计量恒为 \((1-1/3)^2/(1/3)+2(0-1/3)^2/(1/3)=2\)。卡方二自由度的均值也是二,但分布连续且并非常数。对应参考尾概率恒为 \(e^{-1}\approx0.367879\)。因此这个“p值”不超过零点四的概率为一,超过零点四,违反有效p值所需的超均匀性。固定维卡方极限需要样本量趋于无穷;只核对均值与自由度不能获得有限样本保证。
习题18.5-2 一致估计不能自动减一个自由度
在 \(p(\theta)=(\theta,\theta,1-2\theta)\) 模型的真值 \(\theta_0=1/3\) 处,改用 \(\widetilde\theta=\min\{5/12,\max(1/4,N_1/n)\}\)。证明它在该点一致,但代入后的Pearson统计量极限为 \(4\chi^2_1\),而非 \(\chi^2_1\)。
第一类经验比例依概率收敛到三分之一,处于截断区间内部,因此截断以趋于一的概率不生效,估计在该点一致。此时三类比例残差依次为零、\(\widehat p_2-\widehat p_1\)、\(\widehat p_1-\widehat p_2\),两个非零残差的分母都趋于三分之一。故统计量为 \(6n(\widehat p_2-\widehat p_1)^2+o_P(1)\)。单次观察的第二类指示减第一类指示均值为零、方差为三分之二,CLT给出 \(\sqrt n(\widehat p_2-\widehat p_1)\Rightarrow N(0,2/3)\),平方乘六即得到四倍一卡方极限。估计量没有消去有效拟合的正交切空间方向,参数数量为一仍不足以选用标准一卡方参考分布。
习题18.5-3 改变小表后重新求精确尾概率
把例18.5-3的表改为 \(\begin{pmatrix}3&1\\1&3\end{pmatrix}\)。计算正关联的单侧p值、概率排序双侧p值和Pearson统计量。说明相同边际为何允许复用超几何分布。
行列边际仍各为四,总数八,条件分布不含其他表内信息,因而概率仍是 \(1,16,36,16,1\) 除以七十。现在观测左上角为三,右尾取三与四,得到 \(17/70\approx0.242857\)。双侧概率排序包括零、一、三、四,得到 \(34/70=17/35\approx0.485714\)。拟合期望仍每格为二,Pearson统计量为四个单位平方偏差除以二,即二。本例的对称性使双侧值恰为两倍右尾,一般不对称边际下仍须按规定的概率排序求和。
习题18.5-4 观测零计数是否能删除
对假设 \(2\times3\) 表 \(\begin{pmatrix}18&12&0\\6&6&18\end{pmatrix}\),在iid交叉分类模型下计算独立性拟合期望、Pearson统计量与渐近自由度。能否因右上角计数为零而删除这一格?
两行总数都为三十,列总数为二十四、十八、十八,总数六十。拟合期望两行都是十二、九、九。第一行的平方贡献为三、一、九,第二行相同,故统计量为二十六,自由度为 \((2-1)(3-1)=2\)。右上角在拟合独立模型下期望为九,概率并非零,不能把这次观测零值当作结构零删除。删除会改变样本空间和原假设,原有投影证明与自由度计算也不再对应所观察的完整表。
关键术语
- 多项抽样:固定总数的iid分类观察汇总为计数向量,类别计数之间具有总数约束。
- 期望频数:模型类别概率乘样本量;拟合后使用估计概率计算。
- Pearson残差:计数偏差除以期望频数的平方根,其平方和构成Pearson统计量。
- 模型切空间:参数微小变化在标准化概率中产生的一阶方向空间。
- 有效分类拟合:满足式(18.5-6)的估计展开,使残差按信息内积投影到切空间的正交补。
- 自由度:卡方极限中保留下来的独立正态坐标数量,须由模型与拟合条件确定。
- 列联表独立性:两分类变量的联合概率分解为边际概率乘积。
- Fisher条件检验:条件于边际,使用不含未知边际概率的精确表分布进行检验。
- 概率排序双侧p值:累加条件概率不超过观测表概率的所有表,包括并列项。
- 结构零:模型规定不可能出现的类别或单元格,须区别于一次样本的零计数。
参考文献
- University of Maryland,Stat 701,Chi-Square Multinomial Goodness of Fit Test,2002讲义。多项CLT、协方差约束与 \(k-1\) 自由度的完整证明。
- NIST/SEMATECH,e-Handbook:Chi-Square Goodness-of-Fit Test,§1.3.5.15。明确使用分类计数似然估计参数,区别于原始未分组观测的拟合。
- Mathias Drton、Bernd Sturmfels、Seth Sullivant,Lectures on Algebraic Statistics,2008作者版本,§1.1,命题1.1.6与1.1.8。正概率独立模型的渐近卡方检验及条件超几何分布。
- R stats包官方文档,Fisher's Exact Test for Count Data。条件优势比模型、概率排序双侧检验以及精确计算与条件模拟的区别。