Cox模型變量選擇的R軟件實現(xiàn)

韋新星+趙麗棉

【摘要】本文基于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