南方醫科大學公共衛生與熱帶醫學學院生物統計學系(510515) 潘建紅 陳 靖 何春拉 陳平雁
病例對照研究中,如果分組和暴露變量均為二分類,則資料的四格表形式見表1。我們通常用比數比(odds ratio,OR)估計暴露的相對危險度。四格表的OR值的估計有兩種方法,即條件最大似然估計(conditional maximum likelihood estimate,CMLE)和非條件最大似然估計(unconditional maximum likelihood estimate,UCMLE)。兩種方法對OR值的點估計究竟孰優孰劣,本文將應用Monte Carlo方法對其進行比較,以期為應用提供理論依據。

表1 病例對照研究資料四格表形式
假定四格表為雙邊固定,即分組和暴露的邊際和固定,則頻數變量a,b,c,d均服從超幾何分布。由于此時的自由度為1,我們只需計算其中某一格子觀測頻數(如 a)的概率分布〔1〕,即

其中 max(0,m1- n0)≤ x≤ min(n1,m1),1≤ m1≤余類推。以P(a)構建似然函數,因為OR=eθ,于是可以得到關于θ的似然函數


對lnL(θ)求參數θ的導數,有

根據定義,a的期望為

對C求參數θ的導數,有

將式(6)帶入式(4)可得


此方程為一元多次方程,對于非線性方程,可以通過Newton-Raphson迭代求解得到c(下標 C表示conditional)值。
Newton-Raphson迭代運算,可以使用SAS中proc iml的call nlpnra函數實現。當然對于整個最大似然估計,也可以使用proc nlp過程實現,在選項tech=中定義迭代方法。其主要語句如下:
proc nlp tech=newrap;/*指定迭代方法*/
max logf;/*定義最大似然估計的對數函數*/
parms or='initial-value';/*給OR值賦迭代初始值*/
logf=log(pdf('hyper',a,N,n1,m1,or));/*超幾何分布對數概率函數*/
run;
proc nlp具體細節可以參考SAS軟件的幫助〔2〕。
假定四格表為單邊固定,即分組的邊際和固定但暴露的邊際和不固定,則病例組和對照組的暴露數A和C分別服從二項分布,且變量A和C相互獨立,A~B(n1,π1),C ~ B(n0,π0)。二者的聯合概率為各自的概率之積〔3〕,即

對于病例組來說,其暴露比數ω1=π1/(1-π1);同理,對照組暴露比數ω0=π0/(1-π0),暴露比數比,進一步整理可得

將式(10)代入(9)有

以聯合概率構建似然函數,又 OR=eθ,m1=a+c,可以得到關于θ的似然函數。
進一步化簡,去掉與θ和π0無關的項,似然函數可以轉化成

其對數似然函數為

分別對θ,π0求一階偏導,并令偏導函數等于0,得似然方程組:


非條件最大似然估計的OR值是由觀測的4個頻數計算而得,故又稱之為觀測OR值〔4〕。
對于上述兩種估計方法關于OR值點估計的準確性比較,我們采用相對誤差進行評價,即×100%,這里表示c或u。采用Monte Carlo方法,用SAS中的ranbin隨機數函數分別產生病例組與對照組的暴露人數。模擬過程考慮4種參數,即樣本量 N(10,20,40,60,100)、OR 值(0.3,0.5,0.7,1,2,4)、對照組暴露概率 π0(0.1,0.3,0.5,0.7)和病例組與對照組比例 n1:n0(1:9,1:4,1:1,4:1,9:1),共600種參數組合,每種組合重復1000次。當a,b,c,d中任一數為0時,每個格子數均加0.5計算u值。本研究采用SAS 9.2統計軟件編程實現。
對于UCMLE方法來說,其相對誤差比較穩定,不隨N以及OR值改變而發生較大幅度波動;且不管n1:n0取何值,UCMLE相對誤差均小于CMLE方法(見圖1和圖2)。對于CMLE方法來說,影響其相對誤差的因素主要為N值的大小,當N不斷增大時,相對誤差逐漸減小且與CMLE相對誤差逐漸接近,這種趨勢在n1:n0取值的5種不同組合下均一致。但是我們知道,當N增大時,迭代運算量也逐漸增大,所以即使在二者相對誤差相差不大時,UCMLE較CMLE方法仍更為可取(見圖1)。

圖1 不同樣本量和分組例數比下兩種方法相對誤差比較

圖2 不同OR值和分組例數比下兩種方法相對誤差比較
UCMLE與CMLE方法的主要區別在于“條件”二字,前者的條件是基于單邊固定,即分組的邊際和固定但暴露的邊際和不固定,故各組的暴露數分別服從二項分布;而后者的條件是基于雙邊固定,即分組和暴露的邊際和均固定,故各組的暴露數服從超幾何分布。就數據的實際情況而言,它們符合單邊固定而不是雙邊固定情形,所以UCMLE方法更符合實際情況。
從兩種方法的模擬結果看,對OR值的點估計,UCMLE方法在我們所設定的任何參數組合情況下均顯示其估計誤差明顯小于CMLE方法。從兩種方法的計算量看,由于CMLE方法使用迭代算法,其計算量遠比UCMLE方法大,且隨著樣本量的增大而更甚。
綜上所述,我們的研究結論認為:對于四格表資料,無論是實際資料的符合情況,還是OR值估計的準確性和計算的便捷性,非條件最大似然估計方法比基于超幾何分布的條件最大似然估計方法更值得提倡。
1.Zelterman D.Models for Discrete Data.New York:Oxford,2006.
2.SAS Institute Inc.SAS/STAT(R)9.2 User's Guide,Second Edition,North Carolina,2009.
3.Newman SC.Biostatistical Methods in Epidemiology.New York:Wiley,2001.
4.Breslow NE,Day NE.Statistical methods in cancer research.Vol.1.The analysis of case-control studies.Lyons:IARC,1980.