胡純嚴(yán) ,胡良平 ,2*
(1.軍事科學(xué)院研究生院,北京 100850;2.世界中醫(yī)藥學(xué)會(huì)聯(lián)合會(huì)臨床科研統(tǒng)計(jì)學(xué)專業(yè)委員會(huì),北京 100029*通信作者:胡良平,E-mail:lphu927@163.com)
對(duì)高維表資料進(jìn)行差異性分析的基本思路是將高維表降為二維表,降維的重要舉措就是按一個(gè)因素的全部水平或多個(gè)因素的全部水平組合對(duì)資料進(jìn)行分層,從而使每層中的資料都是一個(gè)二維表資料。一種特殊的高維表就是分層后的二維表為2×2表(即含一個(gè)二值的原因變量和一個(gè)二值的結(jié)果變量),簡記為“g×2×2表”。此時(shí),就可依據(jù)“隊(duì)列研究設(shè)計(jì)(含橫斷面研究設(shè)計(jì))”或“病例對(duì)照研究設(shè)計(jì)”,推導(dǎo)出兩套對(duì)此種高維表資料進(jìn)行差異性分析的方法。此法涉及到若干個(gè)具體的分析內(nèi)容,概括地表述為“定性資料的Meta分析[1]”,因篇幅所限,本文只討論g×2×2表資料的齊性檢驗(yàn)問題。
高維表(g×2×2表)資料的表達(dá)模式見表1。
表1 高維表(g×2×2表)的第h層2×2表資料的表達(dá)模式
在高維表資料中,至少有2個(gè)因素(或自變量),1個(gè)定性的結(jié)果變量。除了采用回歸分析可以同時(shí)考察多個(gè)因素對(duì)定性結(jié)果變量的影響之外,差異性分析的思路是將其他因素當(dāng)作分層變量,只研究剩余的一個(gè)因素對(duì)二值定性結(jié)果變量的影響,這被稱為將高維表降維后使其成為二維表。顯然,在分層變量(它可以是1個(gè)因素,也可以是多個(gè)因素的水平組合)的每個(gè)水平下,都有一張二維表。假定分層變量有g(shù)(g≥2)個(gè)水平,則有g(shù)張2×2表(注:本文不考慮g張R×C表)。只有在它們的內(nèi)部構(gòu)成基本一致時(shí),才能將它們按某種規(guī)則合并或壓縮成一張2×2表。于是,人們將“內(nèi)部構(gòu)成基本一致”的g張2×2表資料稱為滿足“齊性的高維表資料”。
如何度量高維表資料是否滿足齊性?若資料取自隊(duì)列研究設(shè)計(jì),則分層后的各張2×2表的相對(duì)危險(xiǎn)度(RR)或危險(xiǎn)率差(RD)之間接近相等,就可以認(rèn)定此種高維表資料滿足齊性;若資料取自病例對(duì)照研究設(shè)計(jì),則分層后的各張2×2表的優(yōu)勢(shì)比(OR,簡稱為優(yōu)比)之間接近相等,就可以認(rèn)定此種高維表資料滿足齊性。由于對(duì)同一個(gè)資料而言,通常OR值與RR值是比較接近的(但必須注意:選擇2×2表資料的第1列還是第2列來計(jì)算RR的數(shù)值是至關(guān)重要的,其中有一個(gè)與OR值接近,而另一個(gè)與OR值的倒數(shù)接近),故在SAS統(tǒng)計(jì)軟件[2]中,僅給出基于OR為效應(yīng)指標(biāo)的“齊性檢驗(yàn)”的計(jì)算公式和SAS軟件實(shí)現(xiàn)方法。然而,在文獻(xiàn)[3-7]中,先給出一個(gè)統(tǒng)一的用于齊性檢驗(yàn)的“Q檢驗(yàn)統(tǒng)計(jì)量”;然后再按照3種效應(yīng)指標(biāo)(分別為“相對(duì)危險(xiǎn)度RR或lnRR”“危險(xiǎn)率差RD”和“優(yōu)勢(shì)比OR或lnOR”)分別展開“Q檢驗(yàn)統(tǒng)計(jì)量”。也就是說,實(shí)際存在3個(gè)不同的用于齊性檢驗(yàn)的“Q檢驗(yàn)統(tǒng)計(jì)量”。
在SAS/STAT的FREQ 過程[2]中,高維表資料(特指g×2×2表)優(yōu)勢(shì)比齊性檢驗(yàn)的方法有以下5種 ,其 中 ,“Breslow-Day檢 驗(yàn) ”“Breslow-Day-Tarone檢驗(yàn)”和“Q檢驗(yàn)”在本質(zhì)上都是χ2檢驗(yàn);“I2度量統(tǒng)計(jì)量及其不確定性限值”在本質(zhì)上是Z檢驗(yàn)(也可被視為自由度為1的χ2檢驗(yàn));而“Zelen's精確檢驗(yàn)”在本質(zhì)上類似于分析二維表資料的Fisher's精確檢驗(yàn),即基于超幾何分布而推導(dǎo)出的計(jì)算方法。
2.1.1 Breslow-Day檢驗(yàn)統(tǒng)計(jì)量與Breslow-Day-Tarone檢驗(yàn)統(tǒng)計(jì)量
Breslow-Day檢驗(yàn)統(tǒng)計(jì)量QBD公式如下,見式(1):
在式(1)中,QBD為服從自由度為q-1的χ2分布的檢驗(yàn)統(tǒng)計(jì)量。其中,“ORMH”“E(.)”和“Var(.)”分別代表“基于Mantel-Haenszel方法計(jì)算的優(yōu)勢(shì)比”“第h層2×2表中第(1,1)網(wǎng)格上的期望頻數(shù)”和“第h層2×2表中第(1,1)網(wǎng)格上的方差”,其計(jì)算公式分別見下式:
應(yīng)用式(1)所需要滿足的前提條件:各層2×2表中理論頻數(shù)>5的格子數(shù)應(yīng)大于80.0%。若資料不滿足前述的前提條件,SAS軟件建議改用Tarone's校正的Breslow-Day檢驗(yàn)統(tǒng)計(jì)量,簡稱為Breslow-Day-Tarone檢驗(yàn)統(tǒng)計(jì)量QBDT,見式(5):
在上式中,QBDT為服從自由度為g-1的χ2分布的檢驗(yàn)統(tǒng)計(jì)量。
2.1.2 Q檢驗(yàn)統(tǒng)計(jì)量
Q檢驗(yàn)統(tǒng)計(jì)量Q的公式如下,見式(6):
在式(6)中,Q漸近地服從自由度為g-1的χ2分布。分別見下式:
在上述各式的計(jì)算過程中,若某一層的表格中有0頻數(shù),就將該層的全部頻數(shù)分別加上0.5。
2.1.3 I2度量統(tǒng)計(jì)量及其不確定性限值
Higgins和Thompson于2002年提出了I2統(tǒng)計(jì)量[2],它被用于度量在分層2×2表中層間非齊性的一種測(cè)度。I2以百分?jǐn)?shù)的形式來表達(dá),它可以被解釋為層間變異占總變異的比例。其計(jì)算公式見式(11):
在式(11)中,“g”代表表的層數(shù),Q的計(jì)算公式見前文式(6)。
Higgins和Thompson于2002年提出,通過構(gòu)建H統(tǒng)計(jì)量的置信區(qū)間,可以轉(zhuǎn)換成構(gòu)建I2的不確定限值(類似于“置信區(qū)間”)[2]。首先,定義H統(tǒng)計(jì)量如下式:
其次,需要給出log(H)[即ln(H)]的標(biāo)準(zhǔn)誤的計(jì)算公式,見式(13)和式(14):
最后,構(gòu)建I2的不確定性限值(類似于“置信區(qū)間”)。
基于式(16),通過轉(zhuǎn)換式(15)間接獲得I2的不確定性限值:
當(dāng)I2為0時(shí),SAS/STAT中FREQ過程將置信限的下限設(shè)置為0,并以顯著性水平α取代α/2來決定置信限的上限(注:在本質(zhì)上,此時(shí)所求的就是單側(cè)置信區(qū)間)。
2.1.4 Zelen's精確檢驗(yàn)
Zelen's精確檢驗(yàn)是Breslow-Day近似檢驗(yàn)的精確版本。用于Zelen's精確檢驗(yàn)的“參考集”包括所有可能的g×2×2表,它們的行、列和層合計(jì)頻數(shù)以及各個(gè)2×2表中第(1,1)網(wǎng)格上的合計(jì)頻數(shù)與觀測(cè)到的高維表是相同的。在固定邊緣合計(jì)的條件下,這個(gè)檢驗(yàn)統(tǒng)計(jì)量是觀測(cè)到的g×2×2表出現(xiàn)的概率(它可被稱為“理論概率”),它是超幾何分布概率的乘積。
Zelen's精確檢驗(yàn)的“P值”是前述提及的“理論概率”中小于等于觀測(cè)到的“表概率”的“理論概率之和”。這個(gè)檢驗(yàn)與二維表時(shí)基于Fisher's精確檢驗(yàn)[8-10]是相同的。
2.2.1 高維表資料相對(duì)危險(xiǎn)度RR齊性檢驗(yàn)的具體算法
基于“g×2×2表資料”估計(jì)共同相對(duì)危險(xiǎn)度的前提條件是高維表資料應(yīng)滿足齊性。針對(duì)“相對(duì)危險(xiǎn)度”的齊性檢驗(yàn)(也稱為異質(zhì)性檢驗(yàn))的具體方法如下[3-5,11]:
在式(20)中,RRMH為共同相對(duì)危險(xiǎn)度,其計(jì)算公式見式(21);Q為服從自由度為df=g-1的χ2分布的檢驗(yàn)統(tǒng)計(jì)量。
基于Mantel-Haenszel估計(jì)量(簡稱MH估計(jì)量)估計(jì)高維表資料共同相對(duì)危險(xiǎn)度[2],見式(21):
2.2.2 高維表資料危險(xiǎn)率差RD齊性檢驗(yàn)的具體算法
基于“g×2×2表資料”估計(jì)共同危險(xiǎn)率差的前提條件是高維表資料應(yīng)滿足齊性。針對(duì)“危險(xiǎn)率差”的齊性檢驗(yàn)(也稱為異質(zhì)性檢驗(yàn))的具體方法如下[3-5,11]:
在式(25)中,共同危險(xiǎn)率差RDMH是基于MH法計(jì)算的,見下式:
在式(26)中,權(quán)重w′h由下式定義:
2.3.1 問題與數(shù)據(jù)
【例1】文獻(xiàn)[3]提供了如下資料,試分析5項(xiàng)研究的OR值之間是否滿足齊性。見表2。
表2 吸煙與肝細(xì)胞癌關(guān)系的5項(xiàng)病例對(duì)照研究結(jié)果
【例2】文獻(xiàn)[3]提供了如下資料,試分析6項(xiàng)研究的OR值之間是否滿足齊性。見表3。
表3 鼻咽癌與EB病毒感染關(guān)系的6項(xiàng)病例對(duì)照研究結(jié)果
2.3.2 多項(xiàng)研究的OR值之間齊性檢驗(yàn)的SAS實(shí)現(xiàn)
【例3】沿用例1中的“問題與數(shù)據(jù)”,試檢驗(yàn)5項(xiàng)研究的OR值之間是否滿足齊性。
【分析與解答】設(shè)所需要的SAS程序如下:
【程序說明】“tables語句”中的選項(xiàng)“cmh”要求對(duì)高維表資料進(jìn)行CMHχ2檢驗(yàn);其后括號(hào)內(nèi)的兩個(gè)選項(xiàng)的含義分別為:“I2”代表要求計(jì)算“I2測(cè)度統(tǒng)計(jì)量及其不確定性限值”的數(shù)值;“QOR”代表要求對(duì)各層OR是否滿足齊性進(jìn)行Q檢驗(yàn)(注意:系統(tǒng)默認(rèn)進(jìn)行“Breslow-Day檢驗(yàn)”)?!癳xact語句”中的選項(xiàng)“eqor”,要求對(duì)各層OR是否滿足齊性進(jìn)行Zelen's精確檢驗(yàn)。
【SAS輸出結(jié)果及解釋】
以上是“優(yōu)比齊性的Breslow-Day檢驗(yàn)”的輸出結(jié)果,QBD=0.6277,df=4,P=0.9599,說明5項(xiàng)研究的OR值之間滿足齊性。
以上是“優(yōu)比齊性的Q檢驗(yàn)”的輸出結(jié)果,Q=0.6273,df=4,P=0.9600,說明5項(xiàng)研究的OR值之間滿足齊性。
以上輸出的是“優(yōu)比齊性檢驗(yàn)I2測(cè)度統(tǒng)計(jì)量及其不確性限值”的計(jì)算結(jié)果,I2=0.00%;其不確定性限值的下限與上限都是0.00%。說明5項(xiàng)研究的OR值之間滿足齊性。
【說明】本例實(shí)際上沒有進(jìn)行Zelen's精確檢驗(yàn),因?yàn)镾AS的“日志窗口”中顯示了一條警告信息,告知用戶未進(jìn)行Zelen's精確檢驗(yàn)的理由(就本例而言,顯示“頻數(shù)太大”)。
【例4】沿用例2中的“問題與數(shù)據(jù)”,試檢驗(yàn)6項(xiàng)研究的OR值之間是否滿足齊性。
【分析與解答】設(shè)所需要的SAS程序如下:
【程序說明】在SAS程序的“tables語句”的CMH選項(xiàng)后的括號(hào)內(nèi),加上選項(xiàng)“BDT”,要求輸出Tarone's校正的Breslow-Day檢驗(yàn)統(tǒng)計(jì)量QBDT的數(shù)值及其P值。
【SAS輸出結(jié)果及解釋】
以上為兩種檢驗(yàn)的輸出結(jié)果:其一,“優(yōu)比齊性的Tarone's校正的Breslow-Day檢驗(yàn)”的結(jié)果,QBDT=21.9424,df=5,P=0.0005,說明6項(xiàng)研究的OR值之間不滿足齊性;其二,“優(yōu)比齊性的Zelen's精確檢驗(yàn)”的結(jié)果中,第1行上的“P<0.0001”是指觀測(cè)到的高維表發(fā)生的概率,而第2行上的“P=0.0004”是指各種符合特定條件的分層表的概率之和。此值與前面的校正結(jié)果“P=0.0005”是比較接近的。
以上為“優(yōu)比齊性的Q檢驗(yàn)”的輸出結(jié)果,Q=19.2304,df=5,P=0.0017,說明6項(xiàng)研究的OR值之間不滿足齊性。
以上輸出的是“優(yōu)比齊性檢驗(yàn)I2測(cè)度統(tǒng)計(jì)量及其不確性限值”的計(jì)算結(jié)果,I2=74.00%>50.00%;其不確定性限值的下限與上限分別是40.69%與88.60%,說明6項(xiàng)研究的OR值之間不滿足齊性。
【說明】因篇幅所限,“高維表資料相對(duì)危險(xiǎn)度或危險(xiǎn)率差齊性檢驗(yàn)的SAS實(shí)現(xiàn)”從略。
本文介紹的統(tǒng)計(jì)分析方法主要適用于g×2×2表,而不適用于g×R×C表(R與C中至少有一個(gè)大于2);當(dāng)分層后的2×2表資料均來自病例對(duì)照研究設(shè)計(jì)時(shí),才適合以“優(yōu)勢(shì)比OR”為效應(yīng)指標(biāo);SAS中的FREQ過程只交待“優(yōu)勢(shì)比OR值的齊性檢驗(yàn)”,而未提及“相對(duì)危險(xiǎn)度RR的齊性檢驗(yàn)”。本文補(bǔ)充介紹了“相對(duì)危險(xiǎn)度RR”和“危險(xiǎn)率差RD”的齊性檢驗(yàn)方法。本文例3的計(jì)算結(jié)果(QBD=0.6277,df=4,P=0.9599)與文獻(xiàn)[3]的結(jié)果(Q=0.60,df=4,P>0.05)非常接近;本文例4的計(jì)算結(jié)果(QBDT=21.9424,df=5,P=0.0005)與文獻(xiàn)[3]的相應(yīng)結(jié)果(Q=12.64,df=5,P<0.05)相差較大,但最終的結(jié)論相同。檢驗(yàn)統(tǒng)計(jì)量Q值上存在的差距,主要原因在于所采取的計(jì)算公式不同所致。
【說明】在對(duì)g×2×2表資料進(jìn)行齊性檢驗(yàn)時(shí),SAS軟件[2]并沒有嚴(yán)格區(qū)分資料所對(duì)應(yīng)的研究設(shè)計(jì)類型是隊(duì)列研究設(shè)計(jì)、橫斷面研究設(shè)計(jì),還是病例對(duì)照研究設(shè)計(jì),只是基于不同的數(shù)學(xué)原理推導(dǎo)出稍有區(qū)別的齊性檢驗(yàn)公式或算法而已,但在輸出結(jié)果中,一律顯示“優(yōu)比齊性檢驗(yàn)”的字樣;而在統(tǒng)計(jì)學(xué)教科書(例如:文獻(xiàn)[3-5])和RevMan軟件[6-7]中,區(qū)分得比較清楚。
本文圍繞g×2×2表資料齊性檢驗(yàn)問題,解釋了進(jìn)行齊性檢驗(yàn)的必要性,介紹了多種齊性檢驗(yàn)的計(jì)算公式,并結(jié)合兩個(gè)實(shí)例,基于SAS軟件給出了優(yōu)勢(shì)比齊性檢驗(yàn)的結(jié)果,對(duì)輸出結(jié)果給出了解釋,并據(jù)此做出了統(tǒng)計(jì)結(jié)論和專業(yè)結(jié)論。