韋新星+趙麗棉
【摘要】本文基于R軟件,對Cox模型的變量選擇進行案例分析。結(jié)果表明,R軟件不僅能很好地對Cox模型進行變量選擇,操作簡便,準確度高,而且還能繪制特定變量的生存函數(shù)圖,供進一步對比分析使用。
【關(guān)鍵詞】R Cox模型 變量選擇 生存函數(shù)圖
【中圖分類號】G642.1 【文獻標識碼】A 【文章編號】2095-3089(2017)39-0217-02
一、引言
Cox模型是英國統(tǒng)計學家Cox于1972年提出的[1]。Cox模型一經(jīng)提出,便被廣泛應用于生存數(shù)據(jù)的處理和分析中。運用Cox模型進行變量選擇,人們不僅可以處理內(nèi)含刪失數(shù)據(jù)的生存數(shù)據(jù),而且還可以處理不服從正態(tài)分布的生存數(shù)據(jù),這為生存數(shù)據(jù)的處理和分析提供了有效的途徑。
然而,在建立Cox模型進行變量選擇的過程中,由于Cox模型隸屬半?yún)?shù)模型,模型表達式較為復雜,若按照公式進行處理,則不僅計算量大,而且還容易出錯。為此,本文在前人研究的基礎上,基于R軟件,對Cox模型的變量選擇進行R軟件的實現(xiàn)研究,以解決上述問題帶來的影響。
二、Cox模型概述
Cox模型,又稱Cox比例風險模型,是一種處理多因素生存數(shù)據(jù)的有效方法。其一般形式為:
其中,為協(xié)變量,為回歸系數(shù),為與無關(guān)的基本風險函數(shù)。之所以又將Cox模型稱為Cox比例風險模型,是因為對于不同的協(xié)變量和,該模型都能滿足一個特點,即風險函數(shù)之比與無關(guān),類似于沒有的情況[2]。于是通常用風險函數(shù)來描述被觀測對象在某一時間點存活的條件下,在未來時間里死亡的概率。另一方面,個體生存時間大于時間的概率被命名為生存函數(shù)。
Cox模型中的回歸系數(shù)通常用偏似然函數(shù)來估計[4]。而生存分析中的很多檢驗也都是基于似然函數(shù)或偏似然函數(shù)的漸進性質(zhì)實現(xiàn)的,如似然比檢驗、Wald檢驗和得分檢驗等,若這些檢驗的統(tǒng)計量值大于臨界值,或檢驗統(tǒng)計量的P值小于檢驗水平,則拒絕原假設,認為回歸系數(shù)有顯著影響。
三、R軟件實現(xiàn)
本節(jié)從典型的生存案例出發(fā),基于R軟件,實現(xiàn)Cox模型變量選擇的過程。
案例采用72個乳腺癌病人的數(shù)據(jù),數(shù)據(jù)來源于文獻[3],共包含4個變量:time,index,age和group,分別表示生存時間、刪失情況(1表示死亡,0表示刪失)、年齡和是否用藥(2表示用了抗癌藥,1表示沒用抗癌藥)。其中,年齡為定量數(shù)據(jù),是否用藥為定性數(shù)據(jù)?,F(xiàn)要探究患者手術(shù)時的年齡與是否用藥對乳腺癌患者術(shù)后生存狀況的影響。
經(jīng)分析,這是一個典型的生存分析問題,宜采用Cox模型來研究age和group這2個協(xié)變量對生存時間time的影響。具體操作如下:
1.數(shù)據(jù)處理
由于文獻[3]的數(shù)據(jù)為sav格式,故需將其另存為逗號分隔符文件(csv格式)或文本文件(txt格式),然后再用read.csv()或read.table()方法進行讀取。在此,我們以另存為到桌面的csv格式為例。步驟:在SPSS的數(shù)據(jù)視圖中,通過選擇“文件”→“另存為”,然后在“保存類型”下拉列中選擇csv,以文件名“案例一”另存為到桌面即可。
2. R代碼編寫
接下來,運用R軟件進行Cox模型變量選擇。
(1)先在R中安裝survival包??赏ㄟ^點擊“Packages”→“Load package”,然后在彈出的select one選擇框中選擇survival,然后點擊“OK”。此時R中會出現(xiàn)以下語句:
> local({pkg <- select.list(sort(.packages(all.available = TRUE)),graphics=TRUE)
+ if(nchar(pkg)) library(pkg, character.only=TRUE)})
Loading required package: splines
接著,用R代碼library(survival)進行加載即可。
(2)然后在R中輸入R代碼。需要指出的是,R中的所有代碼都區(qū)分大小寫,而且其標點符號都需要用半角格式書寫。Cox模型變量選擇的R代碼為:
u=read.csv("C:/Users/Administrator/Desktop/案例一.csv",sep=",")
library(survival)
a=coxph(Surv(time,index)~age+group,data=u)
sa=step(a)
summary(sa)
回車后即可看到Cox模型變量選擇的輸出結(jié)果。
此外,運用R軟件還可以進行生存函數(shù)圖的繪制。在R中輸入以下R代碼:
a0=survfit(Surv(time,index)~group,data=u,type="kaplan-meier")
plot(a0,lty=1:2,xlab="time",ylab="survival fuction")
title("comparison")
legend("topright",c("group=1","group=2"),lty=1:2)
回車后即可看到生存函數(shù)的對比圖。
3.結(jié)果輸出及分析
(1)Cox模型變量選擇的部分輸出結(jié)果如下:
Call:
coxph(formula = Surv(time, index) ~ age + group, data = u)
n= 72, number of events= 36endprint
coef exp(coef) se(coef) z Pr(>|z|)
age 0.2076 1.2307 0.0323 6.423 1.34e-10 ***
group -1.5269 0.2172 0.4086 -3.737 0.000186 ***
---
Signif. codes: 0 ‘*** 0.001 ‘** 0.01 ‘* 0.05 ‘. 0.1 ‘ 1
Concordance= 0.883 (se = 0.055 )
Rsquare= 0.6 (max possible= 0.978 )
Likelihood ratio test= 65.91 on 2 df, p=4.885e-15
Wald test= 47.6 on 2 df, p=4.608e-11
Score (logrank) test = 60.62 on 2 df, p=6.872e-14
可以看到,age和group這2個協(xié)變量的系數(shù)分別為0.2076和-1.5269,其中age的系數(shù)為正,表明高齡乳腺癌患者的死亡率高于低齡的,而group的系數(shù)為負,表明用了抗癌藥的乳腺癌患者的死亡率低于沒用藥的。Cox模型表達式為:
另一方面,為0.6,表明模型的擬合效果較好。而似然比檢驗、Wald檢驗和得分檢驗的P值均小于0.001,可認為回歸系數(shù)有顯著影響。
(2)繪制的生存函數(shù)對比圖為:
圖1 group的生存函數(shù)對比圖
由圖1,我們可以比較用了抗癌藥的乳腺癌患者(group=2)和沒用抗癌藥的乳腺癌患者(group=1)的生存函數(shù)的區(qū)別??梢钥吹剑昧丝拱┧幍娜橄侔┗颊叩纳鏍顩r明顯高于沒用藥的,這再次表明用藥對抗癌有明顯療效。
四、結(jié)論
本文基于R軟件,對Cox模型的變量選擇進行研究分析。從案例的編程及處理上看,盡管數(shù)據(jù)結(jié)構(gòu)不盡相同,但在對Cox模型進行變量選擇時,R軟件操作簡便,實用性強,準確度高,能很好地對Cox模型進行變量選擇。此外,R軟件還可生成特定變量的生存函數(shù)圖,不僅直觀形象,還可供進一步對比分析使用。
參考文獻:
[1] COX D R. Regression Models and Life Tables[J]. Journal of Royal Statistical Society, 1972, 34: 187-220.
[2] 吳喜之. 統(tǒng)計學:從數(shù)據(jù)到結(jié)論[M]. 北京: 中國統(tǒng)計出版社, 2013: 200-207.
[3] 李靜萍, 謝邦昌. 多元統(tǒng)計分析方法與應用[M]. 中國人民大學出版社, 2008: 209-217.
[4] 田俊. Cox回歸模型及其參數(shù)估計[J]. 福建醫(yī)學院學報. 1987, 21(2): 115-117.
基金項目:河池學院碩士專業(yè)學位建設基金課題(2016YT004);廣西大學生創(chuàng)新創(chuàng)業(yè)訓練計劃項目(201610605054)。
作者簡介:韋新星(1990-),女,廣西柳州人,河池學院數(shù)學與統(tǒng)計學院助教,碩士,研究方向:應用統(tǒng)計。endprint