拓端數據tecdat|R語言中生存分析模型與時間依賴性ROC曲線可視化

原文鏈接：http://tecdat.cn/?p=20650 

 

人們通常使用接收者操作特征曲線（ROC）進行二元結果邏輯回歸。但是，流行病學研究中感興趣的結果通常是事件發生時間。使用隨時間變化的時間依賴性ROC可以更全面地描述這種情況下的預測模型。

時間依賴性ROC定義

令 Mi為用於死亡率預測的基線（時間0）標量標記。 當隨時間推移觀察到結果時，其預測性能取決於評估時間 t。直觀地說，在零時間測量的標記值應該變得不那么相關。因此，ROC測得的預測性能（區分）是時間t的函數 。 

累積病例

累積病例/動態ROC定義了在時間t 處的閾值c處的 靈敏度和特異性，  如下所示。

累積靈敏度將在時間t之前死亡的視為分母（疾病），而將標記值高於 c 的作為真實陽性（疾病陽性）。動態特異性將在時間t仍然活着作為分母（健康），並將標記值小於或等於 c 的那些作為真實陰性（健康中的陰性）。將閾值 c 從最小值更改為最大值會在時間t處顯示整個ROC曲線 。

新發病例

新發病例ROC1在時間t 處以閾值 c定義靈敏度和特異性，  如下所示。

累積靈敏度將在時間t處死亡的人  視為分母（疾病），而將標記值高於 Ç 的人視為真實陽性（疾病陽性）。

數據准備

我們以數據 包中的dataset3survival為例。事件發生的時間就是死亡的時間。Kaplan-Meier圖如下。

1.  
    
2.  
    
3.  
   ## 變成data_frame
4.  
   data <- as_data_frame(data)
5.  
   ## 繪圖
6.  
   plot(survfit(Surv(futime, fustat) ~ 1,
7.  
   data = data)

可視化結果：

在數據集中超過720天沒有發生任何事件。

1.  
    
2.  
   ## 擬合cox模型
3.  
   coxph(formula = Surv(futime, fustat) ~ pspline(age, df = 4) +
4.  
   ##獲得線性預測值
5.  
   predict(coxph1, type = "lp")

累積病例

實現了累積病例

1.  
    
2.  
   ## 定義一個輔助函數，以在不同的時間進行評估
3.  
   ROC_hlp <- function(t) {
4.  
   survivalROC(Stime
5.  
   status
6.  
   marker
7.  
   predict.time = t,
8.  
   method = "NNE",
9.  
   span = 0.25 * nrow(ovarian)^(-0.20))
10.  
    }
11.  
    ## 每180天評估一次
12.  
    ROC_data <- data_frame(t = 180 * c(1,2,3,4,5,6)) %>%
13.  
    mutate(survivalROC = map(t, survivalROC_helper),
14.  
    ## 提取AUC
15.  
    auc = map_dbl(survivalROC, magrittr::extract2, "AUC"),
16.  
    ## 在data_frame中放相關的值
17.  
    df_survivalROC = map(survivalROC, function(obj) {
18.  
     
19.  
    ## 繪圖
20.  
    ggplot(mapping = aes(x = FP, y = TP)) +
21.  
    geom_point() +
22.  
    geom_line() +
23.  
    facet_wrap( ~ t) +

可視化結果：

180天的ROC看起來是最好的。因為到此刻為止幾乎沒有事件。在最后觀察到的事件（t≥720）之后，AUC穩定在0.856。這種表現並沒有衰退，因為高風險分數的人死了。

新發病例

實現新發病例

1.  
    
2.  
   ## 定義一個輔助函數，以在不同的時間進行評估
3.  
    
4.  
   ## 每180天評估一次
5.  
    
6.  
   ## 提取AUC
7.  
   auc = map_dbl(risksetROC, magrittr::extract2, "AUC"),
8.  
   ## 在data_frame中放相關的值
9.  
   df_risksetROC = map(risksetROC, function(obj) {
10.  
    ## 標記欄
11.  
    marker <- c(-Inf, obj[["marker"]], Inf)
12.  
     
13.  
    ## 繪圖
14.  
     
15.  
    ggplot(mapping = aes(x = FP, y = TP)) +
16.  
    geom_point() +
17.  
    geom_line() +
18.  
    geom_label(data = risksetROC_data %>% dplyr::select(t,auc) %>% unique,
19.  
    facet_wrap( ~ t) +

可視化結果：

這種差異在后期更為明顯。最值得注意的是，只有在每個時間點處於風險集中的個體才能提供數據。所以數據點少了。表現的衰退更為明顯，也許是因為在那些存活時間足夠長的人中，時間零點的風險分沒有那么重要。一旦沒有事件，ROC基本上就會趨於平緩。

結論

總之，我們研究了時間依賴的ROC及其R實現。累積病例ROC可能與風險 （累積發生率）預測模型的概念更兼容 。新發病例ROC可用於檢查時間零標記在預測后續事件時的相關性。

參考

1. Heagerty,Patrick J. and Zheng,Yingye,  Survival Model Predictive Accuracy and ROC Curves,Biometrics，61（1）,92-105（2005）. doi：10.1111 / j.0006-341X.2005.030814.x.

---

最受歡迎的見解

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例