山西醫(yī)科大學(xué)衛(wèi)生統(tǒng)計教研室(030001) 王 慧 高 雪 虞明星 王 彤

【提 要】 目的 本文以logistic回歸為例，介紹一類基于Wald檢驗的樣本量和功效計算方法。方法 推導(dǎo)一般情形下的計算方法，給出兩類簡單情形下的樣本量計算公式：(1)模型僅含有一個暴露因素；(2)模型含有一個暴露因素和一個混雜因素，并給出實例分析。在暴露因素為一個二分類變量時，通過模擬研究比較了教科書上常用計算公式與基于Wald檢驗樣本量的區(qū)別。對于病例對照研究，給定功效時，可通過最小化樣本量來估計最優(yōu)的病例與對照的比例。結(jié)果 理論推導(dǎo)和模擬研究均顯示，教材中的公式可能高估也可能低估樣本量，僅當(dāng)H0成立時，三種方法得到的樣本量估計相同。結(jié)論 研究設(shè)計中計算樣本量的統(tǒng)計量與假設(shè)檢驗采用的統(tǒng)計量建議保持一致。

Wald檢驗是回歸系數(shù)的假設(shè)檢驗中最常用的一種方法，考慮到研究設(shè)計中計算樣本量的統(tǒng)計量應(yīng)與假設(shè)檢驗采用的統(tǒng)計量一致，因此，本文將介紹Demidenko[1]提出的一類基于Wald檢驗的樣本量和功效的計算方法。

對于單個參數(shù)θ的假設(shè)檢驗H0:θ=θ0,H1:θ≠θ0，Wald檢驗通常構(gòu)造如下：

但對于假設(shè)檢驗：H0:β=β0,H1:β≠β0，給定I類錯誤水平α和功效P時，中文教材[2-4]中常用的樣本量計算公式為

(1)

而國外文獻(xiàn)[5-6]中常用的樣本量計算公式為

(2)

(3)

令(2)式中V0=V，則

(4)

采用最大似然法估計參數(shù)及其方差，(4)式對應(yīng)的檢驗統(tǒng)計量為

(5)

本文將以logistic回歸為例，給出一般情形下基于Wald檢驗的樣本量和功效的計算方法。給出兩類簡單情形下的樣本量的計算公式：(1)模型僅含有一個暴露因素；(2)模型含有一個暴露因素和一個混雜因素，并在暴露因素為二分類變量時，通過模擬研究比較樣本量公式(1)、(2)和(4)的差異。

當(dāng)疾病或結(jié)局發(fā)生率較低的時候，OR近似RR，這時可以采用logistic回歸分析隊列數(shù)據(jù)。多項研究[7-8]指出隊列研究和病例對照研究采用logistic回歸建模時，參數(shù)估計和假設(shè)檢驗一樣，只是截距項的含義不同。在隊列研究中，截距表示參考組(未暴露且協(xié)變量取0或參考水平)中結(jié)局的發(fā)生率；而在病例對照研究中，截距項表示參考組中病例與對照的比例，因此本文將探討病例對照研究中最優(yōu)的病例與對照比例問題。除特別說明以外，本文中將不區(qū)分研究類型。

一般情形

假定我們感興趣的結(jié)局為二分類變量(如是否患有某種疾病)，記n個個體的數(shù)據(jù)為y=(y1,…,yn)T,X=(xij)n×m，我們通常采用以下logistic回歸對其建模

(6)

其中，α0為截距項，β=(β1,…,βm)T，為m維列向量，Xi.=(xi1,…,xim)為m維行向量。假定n個個體{Xi.,i=1,…,n}獨立同分布。本文中粗體表示向量或矩陣。

Wald檢驗中最重要的一步就是計算Fisher信息矩陣。對于一維參數(shù)分布族f(x,θ)(θ∈Θ)，稱

為這個分布族的Fisher信息量。在(6)式定義的logistic回歸中，由Fisher信息的可加性[9]，n個獨立同分布的樣本含有的參數(shù)(α0,β)的Fisher信息量I*(α0,β)為(m+1)×(m+1)維方陣，且可以寫成如下分塊結(jié)構(gòu)：

(7)

(8)

進(jìn)而可以求得參數(shù)(α0,β)的方差協(xié)方差矩陣為V=(vij)(m+1)×(m+1)=I*-1(α0,β)=I-1(α0,β)/n。

(9)

根據(jù)功效的定義，即H1成立時拒絕H0的概率。因此，給定I類錯誤水平α和樣本量n時，功效為|Z|>Z1-α/2的概率。因此，Wald檢驗的功效[10]為

(10)

其中，Φ為標(biāo)準(zhǔn)正態(tài)分布的累積分布函數(shù)。給定I類錯誤水平α和功效P時，所需的樣本量為

(11)

無混雜因素

1.暴露因素為二分類

模型僅含有一個二分類的暴露因素時即m=1，logistic回歸(6)可簡化為

(12)

記A=eα0,B=eβ,px=P(xi=1)，對于假設(shè)檢驗：H0:β=0,H1:β≠0，給定I類錯誤水平α和功效P時拒絕H0所需的樣本量為

(13)

給定I類錯誤水平α和樣本量n時，功效為

(14)

為了得到一個粗略的樣本量估計，我們假設(shè)備擇假設(shè)下β接近于0，即B=1，此時，所需要的樣本量為

(15)

需要注意的是，這里rA類似y的方差，出現(xiàn)在樣本量公式的分母中，而線性回歸模型y=α0+βx+ε中，方差var(y|x)=σ2出現(xiàn)在樣本量公式的分子上。在線性回歸中，樣本量隨著因變量方差的變大而增加。但是在logistic回歸中，樣本量隨著因變量的方差的變大而減少。這是因為，在線性回歸中因變量的方差與回歸系數(shù)是相互獨立的，而在logistic回歸中，因變量的方差是回歸系數(shù)的函數(shù)[11]。

計算樣本量時，需要指定的參數(shù)有：

①I類錯誤水平α和功效P，備擇假設(shè)下暴露因素的效應(yīng)B；

②指定暴露因素和結(jié)局的分布，即人群的暴露率(或比例)px和未暴露人群中結(jié)局的發(fā)生率(或比例)py。

實例1 假定要設(shè)計一項隊列研究分析吸煙與肺癌的關(guān)系，人群中吸煙的比例為px=0.2，不吸煙人群中肺癌的發(fā)生率為py=0.001，吸煙對肺癌的OR為2，我們要計算給定I類錯誤水平α=0.05和功效P=0.8時的樣本量。根據(jù)公式(13)計算得到需要的樣本量為61405。

(1) 三個樣本量公式比較

三個樣本量計算公式(1)、(2)、(4)的差異取決于V0-V，而

因此，公式(1)、(2)高估或者低估樣本量取決于分子(B-1)(1-A2B)的符號。當(dāng)B>1且A2B<1時，即暴露為危險因素，其OR小于未暴露人群中未發(fā)生結(jié)局與發(fā)生結(jié)局的比例的平方，公式(1)、(2)高估樣本量。其他情況下公式(1)、(2)低估樣本量。當(dāng)B=1，即暴露與結(jié)局無關(guān)時，三種方法得到的樣本量估計相同。

圖1展示了三個樣本量計算公式的結(jié)果，其中px=0.01。左側(cè)圖中A=eα0=1/4，因此未暴露組的結(jié)局發(fā)生率為0.02，備擇假設(shè)下OR為2到3，B/16<1，因此，公式(1)、(2)高估樣本量；而右圖中B×0.092>1，低估樣本量。為了檢驗樣本量計算公式(13)，在該公式算得的樣本量下，模擬10000次計算功效，結(jié)果顯示功效接近名義水平80%(見圖1)。

圖1 不同備擇假設(shè)下的樣本量，Ⅰ類錯誤水平α=0.05,功效p=0.8，px=0.01

(2)最優(yōu)病例對照研究設(shè)計

公式(13)用于病例對照研究時，因為參數(shù)A表示未暴露人群中病例與對照的比例，因此，可以通過最小化(13)來計算給定I類錯誤水平α和功效P時未暴露人群中最優(yōu)的病例與對照的比例。令

將f(A)對A求導(dǎo)后令之為0，求解得到

(16)

記n0為總樣本中對照的例數(shù)，n為總樣本例數(shù)，根據(jù)全概率公式

可以求得總樣本中病例對照的比例r和例數(shù)。

實例2 假定設(shè)計一項病例對照研究分析吸煙與肺癌的關(guān)系，人群中吸煙的比例為px=0.2，吸煙的OR為2，我們要計算給定I類錯誤水平α=0.05和功效P=0.08時的樣本量。

如果按照病例與對照1∶1設(shè)計，那么A=1，根據(jù)公式(13)計算得到需要的總樣本量為

其中，病例和對照各225例。

如果采用最優(yōu)設(shè)計，根據(jù)公式(16)計算得到不吸煙人群中的最優(yōu)肺癌患者和非肺癌患者比例為

這意味著不吸煙人群中需要調(diào)查的非肺癌患者比肺癌患者多73%。代入公式(13)計算得到需要的總樣本量為

其中非肺癌患者的例數(shù)為

肺癌患者的例數(shù)為166?？倶颖玖繙p少33例，約7%，病例組樣本量減少59例，約26%。

2.暴露因素為有序多分類

模型僅含有一個有序多分類暴露因素時，logistic回歸(6)可簡化為

(17)

其中，xi=l1,…,ls(i=1,…,n)。記A=eα0,B=eβ,fi=P(x=li)(i=1,…,s)，對于假設(shè)檢驗：H0:β=0,H1:β≠0，給定I類錯誤水平α和功效P時拒絕H0所需的樣本量為

(18)

給定I類錯誤水平α和樣本量n時，功效為

(19)

當(dāng)s=2,l1=0,l2=1,f1=1-px,f2=px時，得到的結(jié)果與暴露為二分類變量時的結(jié)果一致。當(dāng)s=3,l1=0,l2=1,l3=2,f1=1-f2-f3,f2=px1,f3=px2時，給定I類錯誤水平α和功效P時拒絕H0所需的樣本量為

(20)

計算樣本量時，需要指定的參數(shù)有：

① I類錯誤水平α和功效P，備擇假設(shè)下暴露因素的效應(yīng)(OR值)B；

②指定暴露因素和結(jié)局的分布，即人群的暴露率(或比例)px1和px2，未暴露人群中結(jié)局的發(fā)生率(或比例)py。

實例3 假定設(shè)計一項隊列研究分析每日吸煙量與肺癌的關(guān)系，每日吸煙量分為3類：不吸煙、≤1包、>1包，人群中的比例分別為0.8,0.1,0.1，不吸煙人群中肺癌的發(fā)生率為py=0.001，吸煙的OR為2，根據(jù)公式(19)，給定I類錯誤水平α=0.05和功效P=0.8時的樣本量為15117。

3.二分類暴露因素和二分類混雜

模型僅含有一個二分類暴露因素和一個二分類混雜因素時，即m=2，logistic回歸(6)可簡化為

(21)

暴露與混雜因素的關(guān)系定義為

(22)

記A=eα0,B=eβ,C=ec,D=ed,G=eγ,px=P(xi=1)，pz=P(zi=1)，對于假設(shè)檢驗：H0:β=0,H1:β≠0，給定I類錯誤水平α和功效P時拒絕H0所需的樣本量為

(23)

其中，

給定I類錯誤水平α和樣本量n時，功效為

(24)

計算樣本量時，需要指定的參數(shù)有：

① I類錯誤水平α和功效P，備擇假設(shè)下暴露因素的效應(yīng)B；

②指定暴露因素和混雜因素的分布及暴露與混雜的關(guān)系：人群的暴露率(或比例)px，人群中混雜因素的發(fā)生率(或比例)pz，混雜因素與暴露因素的OR值D；或者指定暴露因素和混雜因素的聯(lián)合分布：π00、π10、π01、π11；

③混雜因素的效應(yīng)值G，非暴露人群中結(jié)局的發(fā)生率(或比例)py。

實例4 假定要設(shè)計一項隊列研究分析吸煙與肺癌的關(guān)系，考慮的混雜因素為性別，人群中吸煙的比例為px=0.2，不吸煙人群中肺癌的發(fā)生率為py=0.001，吸煙的OR為2，男女為1∶1，假定性別與吸煙是獨立的(D=1)，男女患肺癌的概率沒有差別(G=1)，則給定I類錯誤水平α=0.05和功效P=0.8時的樣本量為61405。與實例1中不考慮性別這個混雜因素時計算的樣本量一致。因而，與結(jié)局和暴露都無關(guān)的變量不影響樣本量。

實例5 在實例4中，如果假定男性吸煙的比例是女性的2倍(D=2)，那么，給定I類錯誤水平α=0.05和功效P=0.8時的樣本量為63021。因而，與結(jié)局無關(guān)，而與暴露有關(guān)的變量也會影響樣本量。

實例6 在實例4中，如果假定男性患肺癌的風(fēng)險是女性的2倍(G=2)，那么，給定I類錯誤水平α=0.05和功效P=0.8時的樣本量為41027。因而，與結(jié)局有關(guān)，而與暴露無關(guān)的變量會影響樣本量。

實例7 在實例4中，如果假定男性吸煙的風(fēng)險是女性的2倍(D=2)，男性患肺癌的風(fēng)險是女性的2倍(G=2)，那么，給定I類錯誤水平α=0.05和功效P=0.8時的樣本量為39930。因而，與結(jié)局和暴露都有關(guān)的混雜因素會影響樣本量。

與線性回歸不同，在logistic回歸中增加一個與暴露因素獨立的混雜變量也會改變暴露因素的回歸系數(shù)的方差，從而影響樣本量和功效的計算。Robinson和Jewell[12]也指出線性回歸中一些結(jié)論并不能直接外推到logistic回歸中。圖2描述了給定I類錯誤水平和功效時的樣本量隨著暴露的效應(yīng)和暴露與混雜的效應(yīng)的變化，以及與不含混雜變量時樣本量的比較。

討 論

logistic回歸中，現(xiàn)有的樣本量和功效計算的方法主要基于率差的檢驗[10]、似然比檢驗和得分檢驗，盡管有模擬研究[13]比較這些方法，但并沒有得到一致的結(jié)果。Wald檢驗是logistic回歸中回歸系數(shù)的常用的檢驗方法，但在實際應(yīng)用中存在樣本量計算時采用的統(tǒng)計量和統(tǒng)計分析時采用的統(tǒng)計量不一致的情況。為此，本文介紹了一類基于Wald檢驗的樣本量和功效的計算，并推導(dǎo)了兩類簡單情形下的樣本量和功效公式。目前公式(13)和(23)在PASS 15.0中可以實現(xiàn)。Wald檢驗是回歸系數(shù)的假設(shè)檢驗中最常用的一種檢驗方法，因而這類方法可以推廣到其他模型的回歸系數(shù)的檢驗中，現(xiàn)有文獻(xiàn)中已經(jīng)用于中介分析[14]、縱向數(shù)據(jù)[15]和結(jié)構(gòu)方程模型[16]。

當(dāng)暴露因素為k(k≥3)個分類無序多分類變量時，通常設(shè)定k-1個二分類啞變量引入模型中分析，因此，無序多分類暴露因素的分析可以轉(zhuǎn)換為多個二分類暴露進(jìn)行分析。但是當(dāng)k>3時，信息矩陣的逆沒有通用公式，因而不能得到顯示的公式。當(dāng)暴露因素為連續(xù)變量時，信息矩陣的計算涉及積分，因而也得不到通用公式。此外，本方法計算樣本量時需要指定的參數(shù)py是參考人群中結(jié)局的發(fā)生比例，而不是一般人群中的結(jié)局的發(fā)生率，當(dāng)考慮協(xié)變量時，在實際操作中不易獲得py。

確定最優(yōu)病例對照比例是病例對照研究設(shè)計中的一個問題，因為logistic回歸的截距項表示參考人群中病例與對照的比例，而樣本量計算公式是這個比例的函數(shù)，因此可以求得參考人群中病例與對照的比例。本文中推導(dǎo)了僅含有暴露因素時的最優(yōu)病例對照比例，對于其他情況下，不能得到顯性表達(dá)式，需要通過迭代計算。

在logistic回歸中存在與結(jié)局有關(guān)的協(xié)變量時，無論該變量是否與感興趣的暴露因素有關(guān)，都會影響樣本量?，F(xiàn)有方法僅能得到存在一個混雜變量時的計算公式，實際情況中可能多個混雜因素同時存在，這時如何計算樣本量是需要進(jìn)一步解決的問題。