朱 凯李 悦
RPT方法在多元游程检验中的应用
朱 凯1李 悦2△
1.江西师范大学科学技术学院(330027)
2.南昌大学公共卫生学院
游程检验是根据样本标志表现排列所形成的游程的多少进行判断的检验方法,在检验样本的随机性及判断数据的规律性等方面有着重要的作用。本文针对多元小样本随机游程检验问题,介绍可以有效解决这一问题的方法——randomized permutation test(简称RPT),并运用Matlab软件编写程序实现该方法。
设二元序列中0和1的个数分别为n1和n2,n1+ n2=n,一般来说n1和n2都大于0,设R表示总游程数,则R的最小值是2,最大值为n,若游程总数R的值过大(即0和1呈周期性变化的趋势),或游程总数R的值过小(即序列二元0和1的出现过于集中),可认为样本数据受到了某些非随机因素的干扰,它不符合随机抽取的原则[1]。
在实际应用中,随机游程对某一样本观测值序列进行检验时,通常有参数法和确切概率两种方法。当n1,n2都较大时,可使用参数法。利用检验统计量(游程)近似服从正态分布的特性,通过矩法估计求出其参数,再根据现有样本的游程数得出p值[2]。但对小样本数据参数法会产生较大的误差。此外,当样本序列中出现了类似于0,1,2的三元或三元以上的序列时,游程数的确切分布更加难以获得,国内目前尚没有涉及三元及以上游程检验方法的文献报道。显然参数法只能针对游程为二元的大样本数据使用。于是基于确切概率的非参数检验方法可用于游程检验,该方法既可以避免小样本数据的误差,又能够用于分析多元游程数据,故此法较传统的参数法有很多优势。其基本过程为:首先计算现有样本的游程总数,记为Robs,称Robs为游程总数的临界值,然后找出该样本所有可能组合及各组合对应的游程,根据各组合对应的游程求出大于(或小于)Robs的频率,即为p值。该方法为全排列方法,即exact permutation test(简称EPT)。当样本量偏大时往往由于排列较为复杂且组合数太多(当n1=10,n2=18时,组合数为C1028=13123110),难以获得游程数的确切分布,此时基于全排列的EPT方法难以实现。故可以在所有的排列中随机抽取其中一部分排列,得出游程数的近似分布,该方法称之为random permutation test(RPT方法)[4]。RPT的实现过程为:利用软件产生的随机数,对所有可能的排列进行抽样,得到游程总数R的近似抽样分布,即对所有可能的组合进行随机抽样后计算抽样游程总数大于临界值Robs的频率,该方法在游程检验中的关键步骤如下:
1.建立假设,确定检验水准,H0:样本序列符合随机抽取的原则,α=0.05(单侧检验);
2.计算现有样本的游程数Robs,称Robs为R的临界值;
3.利用软件对构成现有的元素进行随机排列,得到与现有样本各元素数均相同的随机样本,计算其游程数R,反复进行该步骤得到检验统计量R的经验抽样分布;
例:在某次临床随机试验中,共入选14人。现用字母代表他们的病情程度,A、B、C分别代表重度、中度、轻度。按照入组时间的先后顺序,这14人实验前的病情程度依次为:A A B C C C B A C B B A A B。
试问:能否认为该实验病人入选顺序随机,α=0.05。
解:H0:病人入选顺序是随机的,H1:病人入选顺序非随机
由题中数据可知总例数n =14,nA=5,nB=5,nC=4,共有3个A游程,4个B游程,2个C游程,游程数观察值为Robs=9。用RPT模拟100000次,即对5 个A,5个B,4个C进行随机排序,重复进行100000次,分别求出每次模拟中总游程数,再求总游程数R大于Robs的频率。笔者编写了相应的Matlab程序实现该方法,程序[3]及注释如下:
a =[0 0 1 2 2 2 1 0 2 1 1 0 0 1];%0 1 2分别代表A B C有三种结果
t =youcheng(a);%计算样本游程数对应的概率
n =100000;%模拟次数
m =length(a);
p0 =0;
for i =1:n
b =matric_randperm(a);%将样本进行随机排列,形成一个模拟抽样
y(i)=youcheng(b);%模拟抽样对应的概率
if y(i)>=t
p0 =p0 +1;end
end
pp =p0/ n;
p =min(pp,1-pp)
%以下是主程序中调用youcheng函数的代码
function[t,b]=youcheng(a)
t =1;b(1)=1;
for i =2:length(a)
if a(i-1)~=a(i)
t =t +1;b(t)=1;
else
b(t)=b(t)+1;
end
end
通过运行该程序,得到p =P(RA≥Robs)=0.30503,P值大于0.05,故接受H0,认为病人出现顺序是随机的。若使用EPT方法进行全排列,虽然所有不同的组合只有种,但寻找这些不同的组合过程十分困难,需从314=4782969种不同组合中将它们筛选出来,整个过程涉及多个复杂的程序及大量计算,而且随着样本量的增加,EPT方法运算次数还将呈几何数增长。故EPT方法在游程检验中是难以实现的。该问题游程总数的确切分布见表1。
表1 总游程数R的EPT分布
根据表1我们可以得到该问题的确切概率(单侧检验)p =P(RA≥Robs)=0.304196,RPT方法与之的相对误差仅有0.2742%。据此我们可知RPT方法较EPT方法简化了计算过程且减少了运算次数,未引起过大的误差。RPT方法可以针对多元数据进行游程检验,突破了参数法只能进行二元游程检验的局限。此外该方法可根据实际情况适当调整抽样次数,以达到增加精度或减少计算量的目的[4]。
随机游程检验是检验一个序列中数据出现是否与顺序无关的常用方法。本文将RPT方法应用于游程检验,解决了多元游程检验的问题。但在使用过程中需注意以下三点:
1.在样本量不大时使用参数方法进行游程检验,p值会有一定的误差。RPT方法可有效地提高p值的精度;
2.本文只讨论了二元和三元序列的游程检验,若遇上三元以上序列的游程检验问题,可直接调用本文给出的程序求出值;
3.由于程序运行中要用到随机数,最终的结果必然产生一定的误差,因此建议在计算机允许的情况下加大模拟次数,并反复运行程序,待结果相对稳定后再下结论[5-6]。
RPT作为EPT的一种近似方法,具有使用方便、误差小、执行效率高的优点,是一种解决多元游程检验并有效提高检验效率的好方法。
参考文献
[1]付云廷.游程检验法在制订业务计划中的应用.中国卫生统计1990,7(4):42-43.
[2]王星.非参数统计.北京:清华大学出版社2009,68-71.
[3]丁元林,孔丹莉.多个样本及其两两比较的秩和检验SAS程序.中国卫生统计,2002,19(5):313-314.
[4]荀鹏程,赵杨,柏建岭,等.Permutation Test在假设检验中的应用.数理统计与管理,2006,26(5):616-621.
[5]Cai JW,Shen Y.Permutation tests for comparing marginal survival functions with clustered failure time data.Statist,Med,2000,19:2963-2973.
[6]朱凯,李悦.RPT对秩和检验的改进及Matlab实现.中国卫生统计,2012,20(4):597-597.
(责任编辑:邓 妍)
·综述·
通信作者:△李悦,Email:liyue0803@126.com