胡純嚴 ,胡良平 ,2*
(1.軍事科學院研究生院,北京 100850;2.世界中醫(yī)藥學會聯(lián)合會臨床科研統(tǒng)計學專業(yè)委員會,北京 100029*通信作者:胡良平,E-mail:lphu927@163.com)
統(tǒng)計資料可分為無層級結(jié)構(gòu)和有層級結(jié)構(gòu)兩種類型。通常情況下,研究者收集的統(tǒng)計資料為無層級結(jié)構(gòu)的資料。此時,可運用常規(guī)統(tǒng)計分析方法處理資料。然而,在調(diào)查研究中,受試對象常歸屬于不同層級,例如某省某市某縣某區(qū)。又例如,學生常歸屬于某市某中學某年級某班級。當資料中的受試對象歸屬于多層級結(jié)構(gòu)時,一般來說,不適合選用處理無層級結(jié)構(gòu)資料的統(tǒng)計分析方法。本文將結(jié)合一個具有2個層級結(jié)構(gòu)的常見臨床試驗資料,介紹二水平多重Logistic回歸模型的構(gòu)建方法和SAS實現(xiàn)方法。
在統(tǒng)計資料中,有一類變量被稱為“分層變量”。例如,在人口普查資料中,會涉及省、市、縣、區(qū)、街道或村委會,它們都可以被稱為“分層變量”。從全國角度考量,可以把“省”作為最高級別的分層變量,把“市”作為第二級別的分層變量,以此類推,把“街道或村委會”作為最低級別的分層變量。
在研究人的身高與體重關(guān)系的資料中,若受試對象的來源涉及省、市、縣3個行政級別,在統(tǒng)計學上就稱此資料為具有3個層級結(jié)構(gòu)的資料。也就是說,省下面嵌套著市、市下面又嵌套著縣。同理,在一個多中心臨床試驗研究中,常把“試驗中心”和“受試對象”視為2個不同層級上的變量,因為“試驗中心”下面嵌套著“受試對象”。
在一個涉及省、市、縣3個行政級別的關(guān)于居民健康情況的調(diào)查研究中,受調(diào)查者為1水平單位,縣為2水平單位,市為3水平單位,省為4水平單位;另一種等價的表述為:受調(diào)查者為水平1單位,縣為水平2單位,市為水平3單位,省為水平4單位。
所謂個體水平,實際上就是指多層級資料中的最底層,通常為“個體”。若試驗中需要對每個個體進行多次重復觀測,則時間點就成為最底層了,此時,時間點就叫做“個體水平”,而每個個體就被稱為“組水平”。在沒有重復觀測的試驗研究資料中,“個體”就叫做“個體水平”,而它們的上一級就叫做“組水平”。在多層級結(jié)構(gòu)中,必須給“組水平”附上相應的標記,以示區(qū)別。例如,在一個涉及省、市、縣3個行政級別的關(guān)于居民健康情況的調(diào)查研究中,就有3種“組水平”,即“省級組水平”“市級組水平”和“縣級組水平”。
多水平模型是基于層級結(jié)構(gòu)數(shù)據(jù)形成的一種統(tǒng)計模型[1],對于具有層級結(jié)構(gòu)的二值結(jié)果變量的統(tǒng)計資料,需要采用多水平多重Logistic回歸模型分析,它是一般Logistic回歸模型的延伸,通過在模型中納入隨機效應項來處理層級結(jié)構(gòu)數(shù)據(jù)的組內(nèi)相關(guān)[2]。其模型表達見式(1)。
式(1)中,X是固定效應的解釋變量設(shè)計矩陣,Z是隨機效應的解釋變量設(shè)計矩陣,β是水平1固定回歸系數(shù)向量,U是隨機回歸系數(shù)向量,服從均值為0、協(xié)方差為矩陣G的正態(tài)分布。
下面介紹2水平多重Logistic回歸模型的建模過程及分析步驟。假定在一個臨床試驗中,涉及試驗藥種類(drug)和臨床試驗中心(center)兩個主要的原因變量,還涉及一個代表療效的結(jié)果變量(y)。假設(shè):y=0表示治療成功,y=1表示治療失敗,以pij表示個體y=0發(fā)生的概率。建模過程如下:
第一步,建立空模型,計算組內(nèi)相關(guān)系數(shù)(ICC)的值??漳P椭袃H有一個隨機截距而不包含任何自變量,其模型見式(2)。
式(2)中,j代表2水平單位的編號,C代表2水平單位的數(shù)目;i代表第j個2水平單位組中個體的編號,mj代表第j個2水平單位組中個體的數(shù)目,β0j代表第j個2水平單位組中的截距,其表達式見式(3)。
式(2)和式(3)兩個模型可合并成式(4)。
式(4)中,β0為y=0的總平均logit值,μ0j為組水平(本資料為中心)的平均logit變異值,它表示第j個組的平均logit值與總平均logit值之間的差異,且
多水平多重Logistic回歸模型的組間變異也可用ICC進行評估,因Logistic回歸模型的殘差方差為π2/3,所以ICC的計算公式可由式(5)表示。
第二步,建立包含自變量的隨機截距模型,即在隨機截距的基礎(chǔ)上再考察變量drug的固定效應。構(gòu)建的模型見式(6)。
式(6)中的β0j由式(7)給出。
式(6)和式(7)兩個模型可合并成式(8)的形式。
式(8)中,β0+ β1drugij為固定效應,μ0j為隨機效應,且
第三步,建立包含自變量的隨機截距-隨機斜率模型,即截距項和自變量drug的回歸系數(shù)中均含有隨機效應項。構(gòu)建的模型見式(9)。
式(9)中等號右邊的β0j和β1j分別由式(10)和式(11)給出。
式(9)、式(10)、式(11)的3個模型可合并成式(12)的形式。
式(12)中,β0+β1drugij為 固 定 效 應 ,μ0j+μ1jdrugij為隨機效應,且μ0j與μ1j之間的協(xié)方差可能有統(tǒng)計學意義;若無統(tǒng)計學意義,則將它們之間的協(xié)方差設(shè)定為0。
【例1】某臨床研究中,研究者選擇16所醫(yī)院同時開展臨床試驗,每所醫(yī)院均選取受試者120人,在醫(yī)院內(nèi)隨機等分為兩組,分別接受試驗藥物和對照藥物的治療。治療結(jié)果見表1,試比較兩種藥物的療效[3]。
表1 多中心臨床試驗數(shù)據(jù)Table 1 Data from multicenter clinical trials
為了檢驗各中心優(yōu)勢比之間是否具有齊性,可利用 Breslow-Day檢驗[4]。設(shè)所需要的 SAS程序如下:
【SAS主要輸出結(jié)果及解釋】優(yōu)勢比齊性的Breslow-Day檢驗結(jié)果:χ2=37.994,df=15,P=0.000 9。此結(jié)果表明:16個中心的優(yōu)勢比不滿足齊性要求,說明不應將16個中心的資料簡單合并后構(gòu)建多重Logistic回歸模型。
3.3.1 對試驗中心產(chǎn)生啞變量后構(gòu)建多重Logistic回歸模型
設(shè)所需要的SAS程序如下:
【SAS主要輸出結(jié)果及解釋】對二重Logistic回歸模型中兩個參數(shù)進行的“3型效應分析”的結(jié)果見表2。
表2 模型中兩個參數(shù)的“3型效應分析”的結(jié)果Table 2 Results of "Type 3 effect analysis" of two parameters in the model
由以上輸出結(jié)果可知,16個中心療效之間的差異有統(tǒng)計學意義,兩種藥物療效之間的差異也有統(tǒng)計學意義。因模型中參數(shù)估計部分的輸出結(jié)果很多,以下僅簡單呈現(xiàn)與第16中心比較差異有統(tǒng)計學意義的結(jié)果,見表3。
表3 與第16中心比較差異有統(tǒng)計學意義的結(jié)果Table 3 Statistically significant results compared with the 16th center
由以上輸出結(jié)果可知,僅5個中心與第16中心的療效之間的差異有統(tǒng)計學意義;試驗藥與對照藥療效之間的差異有統(tǒng)計學意義。由于程序中設(shè)定:擬合“y=0(治療成功)”的概率模型,并且,drug=0代表試驗藥,drug=1代表對照藥,drug的回歸系數(shù)為-0.826,即對照藥相對于試驗藥的優(yōu)勢比OR=exp(-0.826)=0.438。也就是說,試驗藥相對于對照藥的優(yōu)勢比為1/OR=1/0.438=2.283倍。
3.3.2 將試驗中心視為分層變量構(gòu)建多重Logistic回歸模型
相當于在每個試驗中心內(nèi)部進行配對設(shè)計,再進行條件Logistic回歸分析。設(shè)所需要的SAS程序如下:
【SAS主要輸出結(jié)果及解釋】基于center分層的一重Logistic回歸模型的分析結(jié)果見表4。
表4 條件最大似然估計的分析結(jié)果Table 4 Analysis results of conditional maximum likelihood estimation
由以上輸出結(jié)果可知,參數(shù)drug的回歸系數(shù)為-0.819,即對照藥相對于試驗藥的優(yōu)勢比OR=exp(-0.819)=0.441。也就是說,試驗藥相對于對照藥的優(yōu)勢比為1/OR=1/0.441=2.268倍。
3.3.3 構(gòu)建隨機截距多水平多重Logistic回歸模型
從16個“不同中心”抽取受試對象,同一個中心的個體之間的差異被稱為1水平上的差異;而不同中心的個體之間的差異被稱為2水平上的差異,需要構(gòu)建隨機截距多水平多重Logistic回歸模型。設(shè)所需要的SAS程序如下:
【SAS主要輸出結(jié)果及解釋】與center對應的截距項的計算結(jié)果:截距項中隨機部分V_u0=0.154(即隨機部分的方差)。模型中固定效應的計算結(jié)果見表5。
表5 模型中固定效應的計算結(jié)果Table 5 Calculation results of fixed effects in the model
由以上輸出結(jié)果可知,b0=0.530,b1=-0.814,它們與0之間的差異均有統(tǒng)計學意義。drug的回歸系數(shù)為-0.814,即對照藥相對于試驗藥的優(yōu)勢比OR=exp(-0.814)=0.443。也就是說,試驗藥相對于對照藥的優(yōu)勢比為1/OR=1/0.443=2.257倍。
為了檢驗截距項中隨機部分V_u0=0.154與0之間的差異是否有統(tǒng)計學意義,設(shè)所需要的SAS程序如下:
/*將以上模型估計的3個參數(shù)值代入下面的程序中*/
【SAS主要輸出結(jié)果及解釋】模型中各參數(shù)的估計結(jié)果見表6。
表6 模型中各參數(shù)的估計結(jié)果Table 6 Estimation results of parameters in the model
由以上輸出結(jié)果可知,對3個參數(shù)進行了重新估計,并且,對它們都進行了假設(shè)檢驗,P值均小于0.05。b0=0.533,b1=-0.819,V_u0=0.144。
3.3.4 構(gòu)建隨機截距-隨機斜率多水平多重Logistic回歸模型
還可以嘗試構(gòu)建隨機截距-隨機斜率多水平多重Logistic回歸模型。設(shè)所需要的SAS程序如下:
【SAS主要輸出結(jié)果及解釋】模型中參數(shù)的協(xié)方差參數(shù)估計結(jié)果見表7。
表7 模型中參數(shù)的協(xié)方差參數(shù)估計結(jié)果Table 7 Covariance parameter estimation results of parameters in the model
由以上輸出結(jié)果可知,截距中的隨機部分V_u0=0.432,斜率中的隨機部分V_u1=0.427,它們之間的協(xié)方差Cov_u01=-0.199。模型中參數(shù)的固定效應的分析結(jié)果見表8。
表8 模型中參數(shù)的固定效應的分析結(jié)果Table 8 Analysis results of fixed effects of parameters in the model
由以上輸出結(jié)果可知,截距和斜率的固定部分的估計值分別為b0=0.532,b1=-0.821,它們與0之間的差異均具有統(tǒng)計學意義。drug的回歸系數(shù)為-0.821,即對照藥相對于試驗藥的優(yōu)勢比OR=exp(-0.821)=0.440。也就是說,試驗藥相對于對照藥的優(yōu)勢比為1/OR=1/0.440=2.273倍。
為了檢驗3個隨機部分的估計值與0之間的差異是否有統(tǒng)計學意義,設(shè)所需要的SAS程序如下:
/*將以上模型估計的5個參數(shù)值代入下面的程序中*/
【SAS主要輸出結(jié)果及解釋】模型中參數(shù)的估計結(jié)果見表9。
表9 模型中參數(shù)的估計結(jié)果Table 9 Analysis results of fixed effects of parameters in the model
由以上輸出結(jié)果可知,后三行為隨機部分,它們的標準誤差無法計算出來,因此,無法檢驗它們與0之間的差異是否有統(tǒng)計學意義。從SAS日志可知,無法計算的原因是最終的海森矩陣不是正定矩陣,因此,估計的協(xié)方差矩陣不是滿秩的[4]。
比較“第3.3.1節(jié)”到“第3.3.4節(jié)”的計算結(jié)果,可得出以下結(jié)論:第一,采用4種方法構(gòu)建多重Logistic回歸模型,所獲得的截距和斜率的固定部分的估計值比較接近;第二,采用多水平Logistic回歸模型分析,可以把不同中心資料優(yōu)勢比之間的不齊性以隨機效應定量地呈現(xiàn)出來,即各中心試驗藥與對照藥的成功率之間存在明顯的跨中心效應。也就是說,采用無層級結(jié)構(gòu)的多重Logistic回歸模型[5-7]分析層級之間存在明顯差異的統(tǒng)計資料是不合適的。
由于本例資料中有一個多值名義變量“center(試驗中心)”,使統(tǒng)計分析者有多種可能的分析策略。除了本文已經(jīng)采用的4種建模方法之外,還有兩種誤用的做法:一是忽視變量“center”的存在;二是將變量“center”視為有序變量或定量變量(賦值為1,2,…,16)。
本文介紹了與多水平模型分析有關(guān)的4個基本概念;介紹了構(gòu)建二水平Logistic回歸模型的三個步驟;針對一個多中心藥物臨床試驗的實例,介紹了用SAS實現(xiàn)分析的全過程,包括如何檢驗各試驗中心優(yōu)勢比的齊性,如何將試驗中心變換為啞變量后構(gòu)建多重Logistic回歸模型,如何將試驗中心視為分層變量構(gòu)建多重Logistic回歸模型,如何構(gòu)建隨機截距多重Logistic回歸模型,以及如何構(gòu)建隨機截距和隨機斜率多重Logistic回歸模型。對4種可行的建模方法所獲得的估計結(jié)果進行了比較,作出了結(jié)論。最后,在討論中還指出了兩種可能被誤用的做法。