張?zhí)鹛?,劉紅偉 ,劉媛媛 ,李長平 ,2,胡良平
(1.天津醫(yī)科大學公共衛(wèi)生學院,天津 300070;2.世界中醫(yī)藥學會聯合會臨床科研統(tǒng)計學專業(yè)委員會,北京 100029;3.軍事科學院研究生院,北京 100850
本期科研方法專題的第三篇文章介紹了生存資料參數回歸模型有關的基礎知識,包括構建三個常見的生存資料參數回歸模型的基本原理、基于圖示法判斷生存時間服從何種概率分布的方法、基于最大似然估計法求解參數回歸模型中的參數和兩個參數回歸模型擬合優(yōu)度的比較。本文通過一個實例,利用SAS軟件的LIFEREG過程介紹生存資料參數回歸模型的SAS實現方法。
【例1】本文所用數據是對149例糖尿病患者17年追蹤調查數據,包括生存結局、生存時間、隨訪開始時年齡、體重指數等[1]。由于存在刪失數據,此資料為生存資料。Diabetic數據集中生存時間(t,年)、隨訪開始時年齡(age1,歲)、體重指數(BMI)、診斷出糖尿病時的年齡(age0,歲)、收縮壓(SBP,mmHg)和舒張壓(DBP,mmHg)是定量自變量;吸煙狀況(smk,0表示不吸煙,1表示曾吸煙,2表示吸煙)、心電圖讀數(ECG,1表示正常,2表示臨界,3表示異常)、是否有冠心病(CHD,0表示無,1表示有)和結局(status,0表示截尾,1表示死亡)都是定性變量。利用以下SAS數據步程序,創(chuàng)建名為Diabetic的數據集:
【SAS程序說明】因篇幅所限,在CARDS語句后的數據僅列出前5個和后5個觀測。因為收縮壓和舒張壓有一定關聯[2],分析時取加權平均血壓(MBP),即令 MBP=SBP*(1/3)+DBP*(2/3),權重系數分別為1/3和2/3。
對數生存圖、對數-對數生存圖和對數失效比生存圖分別見圖1、圖2和圖3。
圖1 對數生存圖
圖2 對數-對數生存圖
圖3 對數失效比生存圖
圖1中的線圖有些彎曲,圖2在整體上呈線性趨勢,與圖3接近。圖示法表明,基于現有數據,該生存資料同時適合Weibull分布模型和Log-logistic分布模型,本文選擇Weibull分布模型進行擬合(圖示法雖然簡單直觀,但不夠精確,若深入探討該數據是否更加適合Weibull分布模型,可參考其他方法[3])。
因為Weibull分布回歸模型嵌套于廣義Gamma分布回歸模型[4],為選擇最優(yōu)的模型,本文將分別擬合Weibull分布回歸模型和廣義Gamma分布回歸模型,根據似然比檢驗結果來確定最優(yōu)模型。利用以下SAS過程步程序,構建Weibull分布回歸模型和廣義Gamma分布回歸模型。
【SAS程序說明】采用LIFEREG過程分別擬合Weibull分布回歸模型和廣義Gamma分布回歸模型,class語句列出離散型自變量,model語句左邊為生存時間和生存結局變量(括號內為截尾值的標記),右邊為預測變量(即自變量),包括離散和連續(xù)自變量。
Weibull分布回歸模型輸出結果:
廣義Gamma分布回歸模型輸出結果:
【SAS結果說明】Fit Statistics表是基于t為響應變量的最大似然估計得到的統(tǒng)計量,可以用來比較不同協變量的模型。Analysis of Maximum Likelihood Parameter Estimates表給出了參數的估計,包括預測變量的回歸系數和參數分布中參數的估計值。其中“Scale”代表“尺度參數”,“Shape”代表“形狀參數”。
因Weibull分布回歸模型包含2個參數,廣義Gamma分布回歸模型包含3個參數,則似然比擬合優(yōu)度檢驗[5]的自由度為 1,χ2統(tǒng)計量為 8.249(=85.899-77.650)(注意:;顯然,8.249>6.635),P<0.01,即兩個分布擬合效果之間差異有統(tǒng)計學意義,故采用-2logL值最小的分布,即廣義Gamma分布(-2logL值為77.650)。廣義Gamma分布回歸模型輸出結果中,只有age1、MBP、CHD、smk和ECG五個自變量有統(tǒng)計學意義。故僅保留這五個自變量,重新采用廣義Gamma分布回歸模型來進行分析。結果如下:
廣義Gamma分布回歸模型的分析結果表明,隨訪開始時年齡的回歸系數為負值(-0.0408),說明隨訪年齡越大,生存時間越短,死亡風險越大;smk變量0水平與2水平比較時,對應的P值小于0.05,其回歸系數為正值(0.4048),說明不吸煙者生存時間長于吸煙者;ECG變量1水平與3水平比較時,對應的P值小于0.001,其回歸系數為正值(1.0879),說明心電圖正常者生存時間長于心電圖異常者。
在現有的統(tǒng)計軟件中,進行生存資料參數回歸模型建模時尚不能篩選自變量,只能依據計算結果,采用手工方法刪除沒有統(tǒng)計學意義的自變量,這在一定程度上影響了最優(yōu)回歸模型的產生。當自變量中含有兩個以上定量自變量時,根據統(tǒng)計分析的經驗,盡可能產生出一些派生自變量(例如平方項、交叉乘積項等)參與構建回歸模型,可能有助于找到擬合優(yōu)度更好的回歸模型。生存資料參數回歸模型的建模策略與logistic回歸模型建模策略類似,對最后的結果而言,保留在回歸模型中的自變量,既要考慮其應具有統(tǒng)計學意義也要考慮其實際意義,即回歸系數正負號的含義在專業(yè)上是可以得到合理解釋的。
本文通過一個實例,比較詳細地介紹了基于SAS軟件實現參數回歸模型的擬合方法;通過圖示法大致判斷生存時間可能服從何種分布類型;最后,還針對所擬合的兩個參數回歸模型,進行了擬合優(yōu)度檢驗,從而可以確定最優(yōu)回歸模型。需注意的是,當基于圖示法確定生存時間符合某種參數分布且其對應的模型屬于嵌套模型(特指其參數取某些特定值時,原先復雜的模型就退化成一個相對簡單的模型,例如,當威布爾分布模型中的形狀參數γ=1時,它就退化成指數分布模型了)時,則需要采用似然比檢驗確定擬合效果最優(yōu)的模型;否則,需要根據所擬合的兩個模型各自的參數數目和其他擬合統(tǒng)計量的數值,綜合考慮后再選定相對較優(yōu)的模型。