零之审判的统计学麻将——检验估计1
前言：这个专栏是我个人基于兴趣得到的统计结果，包括1. 简单的数据统计或随机模拟， 2. 专业的统计检验或参数估计研究，3.一些机器学习模型的应用。 如果出现数据统计，基本来自于2010-2021年的近三百万个四凤南的牌谱，感谢nodocchi.moe网站作者提供资源。
这期的主题是：给出两个玩家的战绩，如何检验其中一个玩家的每局的pt期望是否高于另一个玩家？另外，如何给出每局的pt期望的估计？
为了更好的讨论这个主题，我们以上期的结论为例。
问题1“点5200”的战绩为24个1位，19个2位，33个3位，70个4位。问题2的战绩为1个1位，84个2位，91个3位，45个4位。问题3的战绩为13个1位，60个2位，62个3位，74个4位。我们将三个问题看作三个七段的玩家，如何通过统计检验来确定这三个玩家中两两玩家的pt期望的差异？另外，如何准确地估计这三个玩家的pt期望？
统计检验
1.传统的t检验
如果将战绩看作由90、45、0、-135组成的数据，则玩家的pt期望相当于这个数据的均值。从而比较两个玩家的pt期望的差异相当于比较两组独立数据的均值。我们通常采取独立两样本t检验（Independent two-sample t-test)。检验过程如下：
首先，我们确认了三个玩家对应的样本r1、r2、r3，其之间有显著的方差差异，于是我们采用方差不一致的t检验。检验的结果是玩家1和玩家2有显著差异，玩家2和玩家3有显著差异，而玩家1和玩家3之间检验的p值为0.12，并不显著。那么这是不是意味着我上期的结论“问题1的pt期望显著差于问题3”是错误的呢？并不是这样的。因为我采取了和传统的t检验不一样的检验方法。
在介绍新的检验方法之前，先要明白t检验的关键要求：样本来自正态或近似正态总体。毫无疑问，七段麻将玩家的战绩样本是一个多项分布，由90、45、0、-135组成，与正态分布区别极大。当然，如果对战数相当大时，样本的均值近似为一个正态分布，仍然可以采取t检验。可在我们这个例子中，三个玩家的对战数都在200左右，并不算大样本。因此，我推荐采取下面的检验方法。
2. 贝叶斯假设检验（推荐）
贝叶斯假设检验（Bayesian Hypothesis Testing）不仅避免了传统的假设检验，如t检验下的不足：1.对样本分布有要求； 2.检验值的计算公式使用了近似方法。同时，贝叶斯假设检验还能给出定量的显著指标：贝叶斯因子（Bayes factor）。不同于p值的大小没有多大意义（只需要知道其是否小于0.05或小于0.1），贝叶斯因子的大小可以用来表示假设的靠谱程度。
贝叶斯因子简介（摘自维基百科）：
对于样本D，以及两个备选模型M1，M2。检验得到的贝叶斯因子K，等于样本D在两个模型下出现的概率的比值，使用贝叶斯定理，可以将其转化为两个后验概率Pr(M|D)的比值与两个先验概率Pr(M)的比值的乘积。
具体计算：
回到麻将的例子，我们对玩家1和玩家3的战绩进行贝叶斯假设检验。由于战绩样本是一个多项分布，因此我们通常采用的先验分布为参数空间{p|p1+p2+p3+p4=1}的均匀分布Dir(1,1,1,1)或无信息的Jeffreys先验分布Dir(1/2,1/2,1/2,1/2)。这两种情况得到的K的值区别较小，我们下面给出均匀分布的例子。
现在的样本D相当于两组多项分布数据：一组为多项分布参数p下的数据（24,19,33,70），另一组为多项分布参数q下的数据（13,60,62,74）。模型M1为子参数空间{p,q|90*p1+45*p2-135*p4<90*q1+45*q2-135*q4}，相对应的模型M2为子参数空间{p,q|90*p1+45*p2-135*p4>=90*q1+45*q2-135*q4}。因此M1和M2的先验概率都是1/2。为了计算两者的后验概率的比值，我们使用随机模拟的数值估计方法。参数p和q的后验分布分别为Dir(25,20,34,71)和Dir(14,61,63,75)，于是按这两个分布随机取值并计算（90*p1+45*p2-135*p4)-(90*q1+45*q2-135*q4)，记录其是否小于0。则模型M1的后验概率即为10000次模拟下小于0的记录的比例。
具体的代码如下：
计算可知贝叶斯因子约为15，这说明了数据非常支持模型M1（Strength of evidence is Strong）。因此统计结果告诉我们问题1的pt期望显著差于问题3。
参数估计
在统计学中，与假设检验相对应的还有一个非常重要的领域：参数估计。回到我们一开始的问题，一个玩家目前战绩为24个1位，19个2位，33个3位，70个4位，那么其每局的pt期望是多少？
对于一个没有统计学专业训练的人，可能这是一个愚蠢的问题。难道结果不就是(24*90+19*45-70*135)/146，约为-44吗？然而答案并不这么简单。比如说我们让四个相同的AI对战，有可能一个AI连赢了10把，我们能说它的每局pt期望为135吗？显然不可以，因为其每局pt期望应该是0。因此对每局pt期望的估计是一个不平凡的问题，它是一个完整的统计学问题，包括样本空间、参数空间、似然函数以及损失函数。下面的叙述包括了我的博士研究结果，如需引用请私信我。
我们先给出这个问题的统计学叙述：参数空间为{p|p1+p2+p3+p4=1}，样本空间为{x|x1+x2+x3+x4=n}，其中n为给定的对战数，似然函数为多项分布的似然函数，损失函数为（f(p)-f(p*))^2，其中p*为参数p的估计值，f(p)=90*p1+45*p2-135*p4表示真实的每局pt期望。因此我们的目标是寻找一个好的估计p*，用f(p*)来估计每局pt期望。
一个好的估计应该是怎么样的呢。一个重要的定义是可容许性（Admissibility）。可容许的估计方法p*指的是不存在另外一种估计方法q，使得q的损失函数的期望值在参数空间的所有参数上都小于等于p*的损失函数的期望值，且存在严格小于的情况。那么对于我们上述的问题，有哪些估计是可容许的估计呢？
通过Brown在1980年提出的complete-class theorem，我们可以知道在上述问题中，可容许的估计一定是贝叶斯估计。我们给出两个贝叶斯估计的例子：
先验分布为Dir(0,0,0,0)的贝叶斯估计。通过计算我们可以知道该贝叶斯估计就是最大似然估计，也就是最简单的方法：用战绩中1位到4位所占的比例来估计。例如一个玩家目前战绩为24个1位，19个2位，33个3位，70个4位，那么估计值就是(24/146,19/146,33/146,70/146)。从而每局pt期望的估计值就是(24*90+19*45-70*135)/146，约为-44。然而这个贝叶斯估计在损失函数（f(p)-f(p*))^2下是否是可容许的呢？我的博士研究成果证明了确实是可容许的。
先验分布为Dir(1,1,1,1)的贝叶斯估计。通过计算我们可以知道此时参数的贝叶斯估计值为(25/150,20/150,34/150,71/150)，从而每局pt期望的估计值就是(25*90+20*45-71*135)/150，约为-43。由于这是一个恰当的（proper）贝叶斯估计，因此其必然是可容许的。
但以上问题的结论并不一定适用于真实的玩家，比如四个相同的AI，我们知道其战绩一定是服从(1/4,1/4,1/4,1/4)的多项分布，所以无论其战绩是什么，对其每局pt期望的估计应该为0。因此，问题的关键在于参数空间。上面的问题中参数空间为{p|p1+p2+p3+p4=1}，而对于实际的玩家，其一位率到四位率的参数分布应该为{p|p1+p2+p3+p4=1}的子空间，比如{p|p1+p2+p3+p4=1，0.1<=p1,p2,p3,p4<=0.4}，在这样的参数空间下，统计问题的估计方法会有什么变化吗？
我的博士研究成果证明了，在这种情况下，当对战数n充分大时，最大似然估计是不容许的。因此我们应当使用恰当的贝叶斯估计来计算玩家的每局pt期望。具体说明这里就不展开了。
