基于拟贝叶斯方法的二项分布总体率的区间估计

安金兵 李 勇

基于二项分布B(n，p)的总体率p的区间估计在生物医学领域有比较广泛的应用，其中很常见的就是Wald区间估计，不过Wald区间估计在覆盖率及稳定性等方面存在比较明显的缺陷〔1，2〕。Agresti等人〔3，4〕给出了与Wald区间估计形式上比较类似的估计方法;而以精确概率为基础的区间估计〔5，6〕以及贝叶斯区间估计〔7〕在实际中也有很多应用。本文基本上沿用贝叶斯统计的思想，通过一种基于Kolmogrov分布距离的拟贝叶斯方法给出了总体率p的一种新的显式解。

原理和方法

传统的贝叶斯方法基于似然函数而构造，但是似然函数容易受到异常值和模型偏差的影响，因此一个比较常用的方法就是用一个对异常取值和模型偏差不敏感的函数代替似然函数构造拟贝叶斯后验分布，进而进行统计分析而获得稳健的性质。假设总体X的分布函数为F(t;p)，p∈Θ为未知参数，通过随机抽样得到的容量为n的样本所对应的经验分布函数为G(t)，则我们可以定义G(t)与F(t;p)之间的Kolmogrov分布距离d(p)以及如下的一个关于参数p函数K(p):

由于Kolmogrov距离d(p)在参数p为真实值时会随样本量的增大而趋于0，因而函数K(p)在参数p为真实值时取较大的数值，反之则取较小的数值，这种性质与似然函数非常类似。那么，如果参数p的先验分布密度为 π(p)，则类似于 Lazar〔8〕、Greco〔9〕等人介绍的方法，利用函数K(p)代替似然函数，就可以将贝叶斯后验分布修正为如下的拟贝叶斯形式:

进而利用上述拟贝叶斯公式就可以对参数p进行点估计或者区间估计。

假设总体X服从二项分布B(n，p)，参数p∈Θ=(0，1)，通过随机抽样得到总体X的一个观测x。区别于以往的方法，本文并不是从B(n，p)分布的样本观测x的分布出发考虑问题，而是将总体率p看作是两点分布B(1，p)的参数，同时将二项分布B(n，p)的观测x看作是一个服从两点分布B(1，p)的取值为0，1的样本序列:y=(y1，…，yn)的总和。

然后，注意到以p为参数的两点分布B(1，p)的分布函数F(t;p)与样本序列y对应的经验分布函数G(t)均为(0，1)上的阶梯函数，并且除去界值0，1之外分别只有一个跃度p和^p=x/n，则经过计算得到:对于任意给定的p∈(0，1)，F(t;p)与 G(t)之间的Kolmogrov分布距离的数值恰好为p－^p。这样，只要给定总体率p的先验分布π(p)，则对于给定的置信水平α，就可以按照上述后验分布密度π(p|x)给出参数p的区间估计［L(x)，U(x)］。

本文选择参数p的先验分布为均匀分布U(0，1)，即π(p)=1，则对于二项分布B(n，p)的一个实际观测x，其样本频率为 ^p=x/n，利用上述基于Kolmogrov分布距离的拟贝叶斯公式，参数p定义在参数空间Θ=(0，1)上的后验分布密度π(p|x)就是一个以x/n为最高点的单峰分布。若给定置信水平α，令:β=α/2，则置信下限L(x)可以通过下式确定:

模拟研究

本文利用随机模拟的方法对总体率p的Wald区间估计和本文基于拟贝叶斯方法所得的区间估计进行比较，本文选择区间估计的覆盖率和平均置信长度来评价区间估计的优劣。假设从二项分布B(n，p)中获得一个观测x，如果按照某种方法给出p的一个置信水平为 α 的区间估计［L(x)，U(x)］，则［L(x)，U(x)］对应的覆盖率P及平均长度L分别为:在下面的部分，我们给定置信水平α=0.05，样本量则分别选择25、50和100三种情形，再从均匀分布U(0，1)中随机获取100个数作为总体率p的取值。对于每一个p，本文分别计算了Wald区间估计和拟贝叶斯区间估计的覆盖率和平均区间长度并进行比较。模拟结果见图1。分别用虚线及实线表示上述两种区间估计的覆盖率和平均置信长度。其中，上面的三个图给出覆盖率的相关结果，而下面的三个图则给出平均置信长度的相关结果。为了清楚起见，在覆盖率的几个子图中均加了一条取值为0.95的水平线以便观察每种区间估计的覆盖率与置信度的差别。

图1 不同样本量下Wald区间估计和拟贝叶斯区间估计的覆盖率和平均置信长度

由图1可以看出，Wald区间估计的平均区间宽度较小，而拟贝叶斯区间估计的平均区间宽度则稍大一点;Wald区间估计在很多情形下根本未达到0.95的覆盖率，而拟贝叶斯区间估计的覆盖率则几乎总在0.95水平线之上，高于Wald区间估计的覆盖率。值得注意的是，Wald区间估计的稳定性比较差，在参数空间Θ=(0，1)上有比较大的振荡，并且比较容易受到总体率p的极端取值的影响，而拟贝叶斯区间估计的稳定性则比较好。所以，本文介绍的基于拟贝叶斯方法所得的区间估计对极端取值的稳健性较好。

讨 论

首先，本文将样本频率 ^p=x/n与来自B(1，p)的样本y对应的经验分布函数G(t)联系起来，再基于G(t)及某个参数p对应的分布函数F(t;p)的Kolmogrov分布距离d(p)构造了一个关于参数p的函数K(p)，并代替似然函数构造拟贝叶斯后验分布，进而给出参数p的区间估计的一个显式解。函数K(p)以分布函数为基础进行定义，因而不容易受到模型偏差以及异常观测的影响，而且随机模拟的结果表明:当参数p在参数空间Θ=(0，1)上取值时，尤其当p接近0，1时，拟贝叶斯区间估计的覆盖率比较高。事实上，生物医学中很多有意义的数值例如治愈率、死亡率等常常在0，1附近取值，因此本方法存在实际的价值。

其次，根据拟贝叶斯后验分布π(p|x)的表达，函数K(p)与似然函数有着类似的功能，因而拟贝叶斯公式不但能代替传统的贝叶斯公式完成类似的功能，而且在一定程度上满足稳健性。此外，根据拟贝叶斯公式的构造，函数K(p)还可以根据实际情况选择其他形式。

最后，由于二项分布的总体率p的后验分布是一个显式表达，因此该方法很容易扩展到多个总体率的函数的估计问题，例如两个总体率的差值、比值等的区间估计以及诊断试验的ROC曲线估计等，这也是我们将来需要探讨的问题。措施。包括对患儿进行积极的隔离治疗，加强重点人群和重点地区的疫情主动搜索，开展行之有效的健康教育宣传，加强医疗机构和托幼机构相关人员的业务培训，加强托幼机构的晨午检制度和消毒制度，通过村医加大对散居儿童的监测和管理等，取得了较好的防控效果，使手足口病发病在4月份出现高峰后开始下降。

手足口病的人群分布与其他一些地区报道的情况基本相似〔4，5〕，男性发病高于女性，可能与男孩喜好活动，相互接触密切有关。3岁以下儿童发病率高，主要原因可能为该年龄组儿童抵抗力低，家庭、社区中隐性感染者为本病的重要传染源，散居儿童居家治疗隔离不规范。

聚集性病例中83.6%发生在托幼机构，提示托幼机构为本病的重点防控单位。通过对疫情的处理，认为只要采取有效的措施即可控制手足口病的发病和流行。处理过程中除对患儿进行严格的隔离治疗外，还对幼儿园进行彻底的消毒(活动室和教室开窗通风，被褥置于阳光下暴晒，餐具一用一消，玩具、教具规范消毒)，并严格落实晨午检制度，发现疑似病例及时隔离，严禁病例带病上课，同时利用专家现场讲座、张贴宣传画、发放手足口病明白纸进行健康宣教，提高幼儿及家长对手足口病的认识，掌握预防知识，取得了较好的防控效果。

综上述结果分析，手足口病虽然传播途径多、传染性强、隐性感染者多又无疫苗、药物等特异性防控手段，但只要严格病例隔离制度、重点地区消毒制度、托幼机构晨午检制度、重点人群和重点地区的疫情主动搜索制度，同时做好手足口病的防治知识宣传，提高群众的防病意识和能力就能取得良好的防控效果。

1．吕冬梅，刑秀暖．46例手足口病发病情况分析．河北医科大学学报，2008，29(5):755-756．

2．周永运，武桂珍．加强疾病监测与动态进行人类肠道病毒71型生物风险评估．中华预防医学杂志，2009，43:328-330．

3．卫生部．手足口病预防控制指南(2009 版)．2009．6．4．http://www．moh．gov．cn．

4．郑惠贞，李剑森，康敏，等．广东省手足口病确诊病例发病及诊疗情况调查．华南预防医学，2008，34(5):10-13．

5．项娜，于海柱，崔兰梅．北京市房山区233例手足口病流行病学分析．中国学校卫生，2008，29(9):817-818．