朱 凱李 悅

?

RPT方法在多元游程檢驗中的應(yīng)用

朱 凱1李 悅2△

1.江西師范大學(xué)科學(xué)技術(shù)學(xué)院（330027）

2.南昌大學(xué)公共衛(wèi)生學(xué)院

游程檢驗是根據(jù)樣本標(biāo)志表現(xiàn)排列所形成的游程的多少進(jìn)行判斷的檢驗方法，在檢驗樣本的隨機(jī)性及判斷數(shù)據(jù)的規(guī)律性等方面有著重要的作用。本文針對多元小樣本隨機(jī)游程檢驗問題，介紹可以有效解決這一問題的方法——randomized permutation test（簡稱RPT），并運用Matlab軟件編寫程序?qū)崿F(xiàn)該方法。

隨機(jī)游程檢驗簡介

設(shè)二元序列中0和1的個數(shù)分別為n1和n2，n1+ n2=n，一般來說n1和n2都大于0，設(shè)R表示總游程數(shù)，則R的最小值是2，最大值為n，若游程總數(shù)R的值過大（即0和1呈周期性變化的趨勢），或游程總數(shù)R的值過?。葱蛄卸?和1的出現(xiàn)過于集中），可認(rèn)為樣本數(shù)據(jù)受到了某些非隨機(jī)因素的干擾，它不符合隨機(jī)抽取的原則［1］。

在實際應(yīng)用中，隨機(jī)游程對某一樣本觀測值序列進(jìn)行檢驗時，通常有參數(shù)法和確切概率兩種方法。當(dāng)n1，n2都較大時，可使用參數(shù)法。利用檢驗統(tǒng)計量（游程）近似服從正態(tài)分布的特性，通過矩法估計求出其參數(shù)，再根據(jù)現(xiàn)有樣本的游程數(shù)得出p值［2］。但對小樣本數(shù)據(jù)參數(shù)法會產(chǎn)生較大的誤差。此外，當(dāng)樣本序列中出現(xiàn)了類似于0，1，2的三元或三元以上的序列時，游程數(shù)的確切分布更加難以獲得，國內(nèi)目前尚沒有涉及三元及以上游程檢驗方法的文獻(xiàn)報道。顯然參數(shù)法只能針對游程為二元的大樣本數(shù)據(jù)使用。于是基于確切概率的非參數(shù)檢驗方法可用于游程檢驗，該方法既可以避免小樣本數(shù)據(jù)的誤差，又能夠用于分析多元游程數(shù)據(jù)，故此法較傳統(tǒng)的參數(shù)法有很多優(yōu)勢。其基本過程為：首先計算現(xiàn)有樣本的游程總數(shù)，記為Robs，稱Robs為游程總數(shù)的臨界值，然后找出該樣本所有可能組合及各組合對應(yīng)的游程，根據(jù)各組合對應(yīng)的游程求出大于（或小于）Robs的頻率，即為p值。該方法為全排列方法，即exact permutation test（簡稱EPT）。當(dāng)樣本量偏大時往往由于排列較為復(fù)雜且組合數(shù)太多（當(dāng)n1=10，n2=18時，組合數(shù)為C1028=13123110），難以獲得游程數(shù)的確切分布，此時基于全排列的EPT方法難以實現(xiàn)。故可以在所有的排列中隨機(jī)抽取其中一部分排列，得出游程數(shù)的近似分布，該方法稱之為random permutation test（RPT方法）［4］。RPT的實現(xiàn)過程為：利用軟件產(chǎn)生的隨機(jī)數(shù)，對所有可能的排列進(jìn)行抽樣，得到游程總數(shù)R的近似抽樣分布，即對所有可能的組合進(jìn)行隨機(jī)抽樣后計算抽樣游程總數(shù)大于臨界值Robs的頻率，該方法在游程檢驗中的關(guān)鍵步驟如下：

1.建立假設(shè)，確定檢驗水準(zhǔn)，H0：樣本序列符合隨機(jī)抽取的原則，α=0.05（單側(cè)檢驗）；

2.計算現(xiàn)有樣本的游程數(shù)Robs，稱Robs為R的臨界值；

3.利用軟件對構(gòu)成現(xiàn)有的元素進(jìn)行隨機(jī)排列，得到與現(xiàn)有樣本各元素數(shù)均相同的隨機(jī)樣本，計算其游程數(shù)R，反復(fù)進(jìn)行該步驟得到檢驗統(tǒng)計量R的經(jīng)驗抽樣分布；

實例應(yīng)用與Matlab實現(xiàn)

例：在某次臨床隨機(jī)試驗中，共入選14人?，F(xiàn)用字母代表他們的病情程度，A、B、C分別代表重度、中度、輕度。按照入組時間的先后順序，這14人實驗前的病情程度依次為：A A B C C C B A C B B A A B。

試問：能否認(rèn)為該實驗病人入選順序隨機(jī)，α=0.05。

解：H0：病人入選順序是隨機(jī)的，H1：病人入選順序非隨機(jī)

由題中數(shù)據(jù)可知總例數(shù)n =14，nA=5，nB=5，nC=4，共有3個A游程，4個B游程，2個C游程，游程數(shù)觀察值為Robs=9。用RPT模擬100000次，即對5 個A，5個B，4個C進(jìn)行隨機(jī)排序，重復(fù)進(jìn)行100000次，分別求出每次模擬中總游程數(shù)，再求總游程數(shù)R大于Robs的頻率。筆者編寫了相應(yīng)的Matlab程序?qū)崿F(xiàn)該方法，程序［3］及注釋如下：

a =［0 0 1 2 2 2 1 0 2 1 1 0 0 1］；%0 1 2分別代表A B C有三種結(jié)果

t =youcheng（a）；%計算樣本游程數(shù)對應(yīng)的概率

n =100000；%模擬次數(shù)

m =length（a）；

p0 =0；

for i =1：n

b =matric_randperm（a）；%將樣本進(jìn)行隨機(jī)排列，形成一個模擬抽樣

y（i）=youcheng（b）；%模擬抽樣對應(yīng)的概率

if y（i）＞=t

p0 =p0 +1；end

end

pp =p0/ n；

p =min（pp，1-pp）

%以下是主程序中調(diào)用youcheng函數(shù)的代碼

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，認(rèn)為病人出現(xiàn)順序是隨機(jī)的。若使用EPT方法進(jìn)行全排列，雖然所有不同的組合只有種，但尋找這些不同的組合過程十分困難，需從314=4782969種不同組合中將它們篩選出來，整個過程涉及多個復(fù)雜的程序及大量計算，而且隨著樣本量的增加，EPT方法運算次數(shù)還將呈幾何數(shù)增長。故EPT方法在游程檢驗中是難以實現(xiàn)的。該問題游程總數(shù)的確切分布見表1。

表1　總游程數(shù)R的EPT分布

根據(jù)表1我們可以得到該問題的確切概率（單側(cè)檢驗）p =P（RA≥Robs）=0.304196，RPT方法與之的相對誤差僅有0.2742%。據(jù)此我們可知RPT方法較EPT方法簡化了計算過程且減少了運算次數(shù)，未引起過大的誤差。RPT方法可以針對多元數(shù)據(jù)進(jìn)行游程檢驗，突破了參數(shù)法只能進(jìn)行二元游程檢驗的局限。此外該方法可根據(jù)實際情況適當(dāng)調(diào)整抽樣次數(shù)，以達(dá)到增加精度或減少計算量的目的［4］。

結(jié) 論

隨機(jī)游程檢驗是檢驗一個序列中數(shù)據(jù)出現(xiàn)是否與順序無關(guān)的常用方法。本文將RPT方法應(yīng)用于游程檢驗，解決了多元游程檢驗的問題。但在使用過程中需注意以下三點：

1.在樣本量不大時使用參數(shù)方法進(jìn)行游程檢驗，p值會有一定的誤差。RPT方法可有效地提高p值的精度；

2.本文只討論了二元和三元序列的游程檢驗，若遇上三元以上序列的游程檢驗問題，可直接調(diào)用本文給出的程序求出值；

3.由于程序運行中要用到隨機(jī)數(shù)，最終的結(jié)果必然產(chǎn)生一定的誤差，因此建議在計算機(jī)允許的情況下加大模擬次數(shù)，并反復(fù)運行程序，待結(jié)果相對穩(wěn)定后再下結(jié)論［5-6］。

RPT作為EPT的一種近似方法，具有使用方便、誤差小、執(zhí)行效率高的優(yōu)點，是一種解決多元游程檢驗并有效提高檢驗效率的好方法。

參考文獻(xiàn)

［1］付云廷.游程檢驗法在制訂業(yè)務(wù)計劃中的應(yīng)用.中國衛(wèi)生統(tǒng)計1990，7（4）：42-43.

［2］王星.非參數(shù)統(tǒng)計.北京：清華大學(xué)出版社2009，68-71.

［3］丁元林，孔丹莉.多個樣本及其兩兩比較的秩和檢驗SAS程序.中國衛(wèi)生統(tǒng)計，2002，19（5）：313-314.

［4］荀鵬程，趙楊，柏建嶺，等.Permutation Test在假設(shè)檢驗中的應(yīng)用.數(shù)理統(tǒng)計與管理，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對秩和檢驗的改進(jìn)及Matlab實現(xiàn).中國衛(wèi)生統(tǒng)計，2012，20（4）：597-597.

（責(zé)任編輯：鄧 妍）

·綜述·

通信作者：△李悅，Email：liyue0803@126.com