一般在建立好Cox模型之后,需要對模型進行診斷。診斷內(nèi)容包括模型的前提條件,諸如Cox模型的PH假定(比例風險假定),共線性假定等。本篇我們通過合實際例子講解Cox模型診斷過程,實現(xiàn)軟件R語言。
1.1 COX模型的診斷內(nèi)容 Cox模型的診斷一般包括三方面的內(nèi)容: 對上述三方面的診斷,常見的方法為殘差法。 Schoenfeld 殘差用于檢驗比例風險假定; Deviance 殘差用于影響點(異常值)識別; Martingale殘差用于非線性檢驗;
1.2 R中用于評估Cox模型的包 我們將會用到以下兩個包: install.packages(c('survival','survminer')) 加載包
library('survival') library('survminer')
1.3 建立Cox模型 我們利用survial包中自帶的肺癌數(shù)據(jù)“data(lung)”建立cox模型。 library('survival') res.cox <- coxph(surv(time,="" status)="" ~="" age="" +="" sex="" +wt.loss,="" data="">模型中有三個變量; res.cox#顯示模型結(jié)果 Call: coxph(formula = Surv(time, status) ~ age + sex + wt.loss,data = lung) coefexp(coef) se(coef) z p age 0.02009 1.02029 0.00966 2.08 0.0377 sex -0.52103 0.59391 0.17435 -2.99 0.0028 wt.loss 0.00076 1.00076 0.00619 0.12 0.9024 Likelihood ratio test=14.7 on 3 df, p=0.00212 n= 214, number of events= 152 (14 observationsdeleted due to missingness)
1.4 模型診斷——PH假定 PH假定可通過假設檢驗和殘差圖檢驗。正常情況下,Schoenfeld殘差應該與時間無關,如果殘差與時間有相關趨勢,則違反PH假設的證據(jù)。殘差圖上,橫軸代表時間,如果殘差均勻的分布則,表示殘差與時間相互獨立。 R語言survival包中的函數(shù)cox.zph()可以實現(xiàn)這一個檢驗過程。 test.ph <-> test.ph rhochisq p age -0.0483 0.3780.538 sex 0.1265 2.3490.125 wt.loss 0.0126 0.0240.877 GLOBAL NA 2.8460.416
從上面的結(jié)果可以看出,三個變量的P值都大于0.05,說明每個變量均滿足PH檢驗,而模型的整體檢驗P值0.416,模型整體滿足PH檢驗。 在R語言 survminer中ggcoxzph( )函數(shù)可以畫出Schoenfeld殘差圖。 ggcoxzph(test.ph)
上圖中實線是擬合的樣條平滑曲線,虛線表示擬合曲線上下2個單位的標準差。如果曲線偏離2個單位的標準差則表示不滿足比例風險假定。從上圖中可見,各協(xié)變量滿足PH風險假設。 另一種檢查比例風險假定的圖形方法是繪制log(-log(S(t)))與t或log(t)是非平行,這個方法只能用于協(xié)變量是分類變量的情形。 如果違反比例風險假設可以通過以下方式解決: 模型中添加協(xié)變量與時間的交互相應; 分層分析;
至于如何實現(xiàn),我們后期再做介紹。
1.5 模型診斷——模型影響點(異常值)識別 我們可以通過繪制Deviance殘差圖或者dfbeta值實現(xiàn)上述診斷。在R語言 survminer中ggcoxdiagnostics()函數(shù)可以畫出Deviance殘差圖。 ggcoxdiagnostics(res.cox,type = 'deviance', linear.predictions = FALSE,ggtheme = theme_bw())
殘差值均勻的分布在0上下,表明滿足上述假定。
ggcoxdiagnostics(res.cox,type = 'dfbeta', linear.predictions = FALSE,ggtheme = theme_bw())
影響點的可能來源于數(shù)據(jù)錄入錯誤,樣本中的極值點、協(xié)變量不均衡,數(shù)據(jù)不足等。對本例,上圖顯示,將dfbeta 值大小與回歸系數(shù)比較表明,即使某些dfbeta值非常大,但它們不足以對模型系數(shù)的估計值產(chǎn)生影響。
1.6 模型診斷——非線性診斷 一般情況下,我們假設協(xié)變量與-log(s(t))之間是線性關系。通過繪制Martingale殘差圖可以實現(xiàn)模型協(xié)變量的非線性診斷。非線性診斷一般是針對模型中的連續(xù)型變量。 在R語言 survminer中ggcoxfunctional()函數(shù)可以畫出Martingale殘差圖。 ggcoxfunctional(Surv(time, status) ~ age + log(age) + sqrt(age),data = lung)
圖中顯示年齡局部有非線性趨勢,但整體表現(xiàn)出線性趨勢。
本公眾號精彩歷史文章: 04:如何在R軟件中求一致性指數(shù)( Harrell'concordance index:C-index)? 05:Nomogram 繪制原理及R&SAS實現(xiàn). 06 : Lasso方法簡要介紹及其在回歸分析中的應用 07 : 最優(yōu)模型選擇中的交叉驗證(Cross validation)方法 08 : 用R語言進行分位數(shù)回歸(Quantile Regression) 09 : 樣本數(shù)據(jù)中異常值(Outliers)檢測方法及SPSS & R實現(xiàn) 10 : 原始數(shù)據(jù)中幾類缺失值(Missing Data)的SPSS及R處理方法 11 : [Survival analysis] Kaplan-Meier法之SPSS實現(xiàn) 12 : [Survival analysis] COX比例風險回歸模型在SPSS中的實現(xiàn) 13 : 用R繪制地圖:以疾病流行趨勢為例 14 : 數(shù)據(jù)挖掘方法:聚類分析簡要介紹 及SPSS&R實現(xiàn) 15 : 醫(yī)學研究中的Logistic回歸分析及R實現(xiàn) 16 : 常用的非參數(shù)檢驗(Nonparametric Tests)總結(jié) 17 : 高中生都能看懂的最小二乘法原理 18 : R語言中可實現(xiàn)的常用統(tǒng)計假設檢驗總結(jié)(側(cè)重時間序列) 19 : 如何根據(jù)樣本例數(shù)、均數(shù)、標準差進行T-Test和ANOVA 20 : 統(tǒng)計學中自由度的理解和應用 21 : ROC和AUC介紹以及如何計算AUC 22 : 支持向量機SVM介紹及R實現(xiàn) 23 : SPSS如何做主成分分析? 24 : Bootstrap再抽樣方法簡介 25 : 定量測量結(jié)果的一致性評價及 Bland-Altman 法的應用 26 : 使用R繪制熱圖及網(wǎng)絡圖 27 : 幾種常用的雙坐標軸圖形繪制 28 : 遺失的藝術—諾謨圖(Nomogram) 29 : Nomogram 繪制原理及R&SAS實現(xiàn)(二) 30 : WOE:信用評分卡模型中的變量離散化方法 31 : 結(jié)構(gòu)方程模型(SEM)簡介及教程下載 32 : 重復測量的多因素方差分析SPSS實現(xiàn)操作過程
|