R語言生存分析: 時變競爭風險模型分析淋巴瘤患者

原文連接：http://tecdat.cn/?p=22422 

在本文中，咱們描述了靈活的競爭風險迴歸模型。迴歸模型被指定爲轉移機率，也就是競爭性風險設置中的累積發生率。該模型包含Fine和Gray（1999）的模型做爲一個特例。這能夠用來對次分佈危險的比例假設作擬合度測試（Scheike和Zhang 2008）。還能夠爲預測的累積發病率曲線構建置信區間。咱們將這些方法應用於Pintilie（2007）的濾泡細胞淋巴瘤數據，其中競爭風險是疾病復發和沒有復發的死亡。

工做實例：濾泡細胞淋巴瘤研究

咱們考慮Pintilie（2007）的濾泡細胞淋巴瘤數據。該數據集由541名疾病早期的濾泡細胞淋巴瘤（I或II）患者組成，並接受單純放療（化療=0）或放療和化療的聯合治療（化療=1）。疾病復發或無反應和緩解期死亡是兩個競爭風險。患者的年齡（年齡：平均=57，sd=14）和血紅蛋白水平（hgb：平均=138，sd=15）也被記錄。隨訪時間的中位數是5.5年。首先咱們讀取數據，計算死亡緣由指標並對協變量進行編碼。

R> table(cause)
cause
0 1 2
193 272 76
R> stage <- as.numeric(clinstg == 2)
R> chemo <- as.numeric(ch == "Y")
R> times1 <- sort(unique(time\[cause == 1\]))

有272個（無治療反應或復發）因疾病引發的事件，76個競爭性風險事件（無復發的死亡）和193個刪減的個體。事件時間用dftime表示。變量times1給出了緣由爲 "1 "的事件時間。咱們首先估計非參數累積發病率曲線進行比較。

咱們指定事件時間並刪減變量爲cause == 0。迴歸模型只包含一個截距項（+1）。cause變量給出了與不一樣事件相關的緣由。cause= 1指定咱們考慮類型1的事件。計算/基於估計值的時間能夠由參數times = times1給出。

圖1（a）顯示了估計的兩種緣由的累積發生率曲線。在圖1（b）中，咱們構建95%的置信區間（虛線）和95%的置信帶。

risk(Surv(dftime, cause == 0) ~ + 1, 
causeS = 1, n.sim = 5000, cens.code = 0, model = "additive")

圖1

R> fit <- cum(time, cause, group)
R> plot(fit)

子分佈危險法和直接二項式模型法都是基於反機率的刪減加權技術。在應用這種權重時，關鍵是刪減權重的估計不能有誤差，不然累積發病率曲線的估計也可能有誤差。

在這個例子中，咱們發現刪減分佈明顯取決於協變量血紅蛋白、階段和化療，並能夠由Cox的迴歸模型很好地描述。Cox模型的擬合是經過累積殘差來驗證的，進一步的細節見Martinussen和Scheike（2006）。所以，對剔除權重使用簡單的KaplanMeier估計可能會致使嚴重的誤差估計。所以，咱們在調用中加入了cens.model = "cox "的選項，這就使用了Cox模型中競爭風險模型的全部協變量做爲剔除權數。通常來講，反機率刪減權重的迴歸模型能夠用來提升效率（Scheike等人，2008）。

如今咱們來擬合模型

咱們首先擬合一個通常比例模型，容許全部協變量具備時變效應。在下面的調用中，只有模型（6）中的協變量x被定義。模型（6）中的協變量z是由一個const操做符指定的。

summary(outf)
OUTPUT:
Competing risks Model
Test for nonparametric terms
Test for non-significant effects
Supremum-test of significance p-value H_0: B(t)=0
(Intercept) 3.29 0.0150
stage 5.08 0.0000
age 4.12 0.0002
chemo 2.79 0.0558
hgb 1.16 0.8890
Test for time invariant effects
Kolmogorov-Smirnov test p-value H_0:constant effect
(Intercept) 8.6200 0.0100
stage 1.0400 0.0682
age 0.0900 0.0068
chemo 1.7200 0.0004
hgb 0.0127 0.5040
Cramer von Mises test p-value H_0:constant effect
(Intercept) 3.69e+01 0.0170
stage 2.52e+00 0.0010
age 4.26e-03 0.0014
chemo 1.50e+00 0.0900
hgb 2.64e-04 0.4220

基於非參數檢驗的顯著性檢驗顯示，在非參數模型中，階段和年齡是顯著的，化療是較顯著的（p = 0.056），血紅蛋白是不顯著的（p = 0.889）。

圖2

繪製估計的迴歸係數αj (t)及其95%的置信帶，並分別繪製常數效應的觀察檢驗過程和空值下的模擬檢驗過程。

R> plot(outf, score = 1)

圖2顯示了這些效應並不隨時間變化而變化，在早期的時間段內效應至關明顯。95%的指向性置信區間，以及95%的置信區間。

圖3顯示了相關的檢驗過程，用於決定時變效應是否具備顯著的時變性，或者是否能夠接受H0 : αj (t) = βj。這些圖的摘要在輸出中給出，咱們看到階段和化療顯然是時變的，所以與Fine-Gray模型不一致。Kolmogorov-Smirnov和Cramer von Mises檢驗統計數字對檢驗過程的兩種不一樣總結是一致的，總的結論是三個變量都沒有比例的Cox類型效應。咱們看到血紅蛋白被常數很好地描述，所以咱們考慮用血紅蛋白具備常數效應，其他協變量具備時變效應的模型。

圖3

R> summary(outf1)
OUTPUT:
Competing risks Model
Test for nonparametric terms
Test for non-significant effects
Supremum-test of significance p-value H_0: B(t)=0
(Intercept) 5.46 0
stage 5.18 0 age 4.20 0
chemo 3.89 0
Test for time invariant effects
Kolmogorov-Smirnov test p-value H_0:constant effect
(Intercept) 10.100 0.000
stage 1.190 0.048
age 0.101 0.004
chemo 1.860 0.000
Cramer von Mises test p-value H_0:constant effect
(Intercept) 79.90000 0.000
stage 1.84000 0.006
age 0.00583 0.000
chemo 2.53000 0.000
Parametric terms :
Coef. SE Robust SE z P-val
const(hgb) 0.00195 0.00401 0.00401 0.486 0.627

Competing risks Model
Test for nonparametric terms
Test for non-significant effects
Supremum-test of significance p-value H_0: B(t)=0
(Intercept) 6.32 0
Test for time invariant effects
Kolmogorov-Smirnov test p-value H_0:constant effect
(Intercept) 1.93 0
Cramer von Mises test p-value H_0:constant effect
(Intercept) 14.3 0
Parametric terms :
Coef. SE Robust SE z P-val
const(stage) 0.45200 0.13500 0.13500 3.340 0.000838
const(age) 0.01450 0.00459 0.00459 3.150 0.001610
const(chemo) -0.37600 0.18800 0.18800 -2.000 0.045800
const(hgb) 0.00249 0.00401 0.00401 0.622 0.534000

咱們注意到，血紅蛋白的影響與更合適模型（如上圖所示）的影響幾乎相等。但因爲模型中其餘協變量的不適合，估計值可能有嚴重的誤差，所以可能誤導了數據的重要特徵。最後，咱們將FG模型的預測與半參數模型的預測進行比較，後者對效應的描述更爲詳細。咱們考慮對下面由新數據分配定義的兩種不一樣的病人進行預測。患者類型I：疾病I期（階段=0），40歲，沒有化療治療（化療=0），患者類型II：疾病II期（階段=1），60歲，放療加化療聯合治療（化療=1）。

R> newdata <- data.frame(stage = c(0, 1), age = c(40, 60), chemo = c(0, 1),
+ hgb = c(138, 138))
R> predict(out, newdata)

爲了指定計算預測的數據，咱們能夠指定一個newdata參數。

基於該模型的預測可能不是單調的。咱們繪製了沒有點狀置信區間（se = 0）和沒有置信帶（uniform = 0）的預測。圖4（a）中的預測是基於靈活的模型，而圖4（b）中的預測是基於FG模型的。I型和II型病人的復發累積發生率曲線分別用實線和虛線表示。圖5（a）比較了基於靈活模型和FG模型對I型患者的預測結果。一樣地，圖5（b）比較了對II型病人的預測。兩個預測值周圍的斷線表明了基於靈活模型的置信區。

圖4

R> par(mfrow = c(1, 2))
R> plot(f1, se = 0, uniform = 1, col = 1, lty = 1
R> plot(fg, new = 0, se = 0, uniform = 0, col = 2, lty = 2,

較高的疾病階段、較高的年齡和聯合治療會致使較高的累積發病率，其影響在時間段的早期更爲明顯（圖4（a）和圖2）。另外一方面，化療在時間段的最初增長了累積發病率，隨後下降了發病率（圖4（a）和圖2）。圖5顯示，FG模型不能準確地模擬時變效應。儘管有這些差別，在這種狀況下，整體預測有些相似，特別是當考慮到估計的不肯定性。然而，協變量的時變行爲顯然是重要的。

圖5

4. 討論

本文實現了累積發病率曲線的靈活競爭風險迴歸模型，能夠詳細分析協變量效應如何預測累積發病率，並容許協變量的時間變化效應。能夠檢查較簡單的模型的擬合度，同時能夠產生帶有置信區間和置信帶的預測結果，這對研究人員頗有用。

---

最受歡迎的看法

1.R語言繪製生存曲線估計|生存分析|如何R做生存曲線圖

2.R語言生存分析可視化分析

3.R語言如何在生存分析與Cox迴歸中計算IDI，NRI指標

4.r語言中使用Bioconductor 分析芯片數據

5.R語言生存分析數據分析可視化案例

6.r語言ggplot2偏差棒圖快速指南

7.R 語言繪製功能富集泡泡圖

8.R語言如何找到患者數據中具備差別的指標？（PLS—DA分析）

9.R語言中的生存分析Survival analysis晚期肺癌患者4例