R语言ROC曲线绘制01-survivalROC

作者:白介素2

相关阅读:
R语言生存分析04-Cox比例风险模型诊断
R语言生存分析03-Cox比例风险模型
R语言生存分析-02-ggforest
R语言生存分析-01
ggpubr-专为学术绘图而生(二)
ggstatsplot-专为学术绘图而生(一)
生存曲线
R语言GEO数据挖掘01-数据下载及提取表达矩阵
R语言GEO数据挖掘02-解决GEO数据中的多个探针对应一个基因
R语言GEO数据挖掘03-limma分析差异基因
R语言GEO数据挖掘04-功能富集分析

如果没有时间精力学习代码,推荐了解:零代码数据挖掘课程

SurvivalROC包绘制时间依赖的ROC曲线

  • 含有删失数据的生存数据
  • 使用survivalROC包
  • 包括Kaplan-Meier (KM) or Nearest Neighbor Estimation (NNE) 两种方法

假设我们有删失的生存数据与基线marker值,我们希望看到marker如何预测数据集中的受试者的存活时间。特别是,假设我们有几天的生存时间,我们想看看标记如何预测一年的存活(predict.time=365)。该功能roc.km.calc()返回感兴趣的时间点的唯一标记值、TP(真阳性)、FP(假阳性)、对应于感兴趣时间点(predict.time)和AUC(ROC)曲线下面积的Kaplan-Meier生存估计。

  • 返回值:

    • cut.values:由于计算TP和FP的marker值
    • TP: 根据cutoff 判断的TRUE postive 真阳性
    • FP: 根据cutoff判断的假阳性
    • predict.time: 感兴趣的时间截点:可以为5年,3年等等
    • Survival: kaplan-Meier法的预估生存时间
    • AUC: Area under ROC,在时间截点的曲线下面积
  • 实际代码演示

Sys.setlocale('LC_ALL','C')
library(survivalROC)
data(mayo)
head(mayo)
##   time censor mayoscore5 mayoscore4
## 1   41      1  11.251850  10.629450
## 2  179      1  10.136070  10.185220
## 3  334      1  10.095740   9.422995
## 4  400      1  10.189150   9.567799
## 5  130      1   9.770148   9.039419
## 6  223      1   9.226429   9.033388
nobs <- NROW(mayo)
cutoff <- 365
  • MAYOSCORE 4作为marker, NNE法计算
    • marker value可以为表达值,免疫分数,检验参数等任何可以定义为数值的指标
    • Mayo4.1得到的是列表,其内容是每一个marker的cutoff值都计算出相应的TP,FP
    • TP,FP绘图即得到ROC,ROC曲线下面积即AUC
Mayo4.1= survivalROC(Stime=mayo$time,##生存时间
                     status=mayo$censor,## 终止事件    
                     marker = mayo$mayoscore4, ## marker value    
                     predict.time = cutoff,## 预测时间截点
                     span = 0.25*nobs^(-0.20))##span,NNE法的namda
str(Mayo4.1)## list结构
## List of 6
##  $ cut.values  : num [1:313] -Inf 4.58 4.9 4.93 4.93 ...
##  $ TP          : num [1:313] 1 0.997 0.995 0.993 0.99 ...
##  $ FP          : num [1:313] 1 0.997 0.994 0.99 0.987 ...
##  $ predict.time: num 365
##  $ Survival    : num 0.929
##  $ AUC         : num 0.931
## 绘图
plot(Mayo4.1$FP, Mayo4.1$TP, ## x=FP,y=TP
     type="l",col="red", ##线条设置
     xlim=c(0,1), ylim=c(0,1),   
     xlab=paste( "FP", "\n", "AUC = ",round(Mayo4.1$AUC,3)), ##连接
     ylab="TP",
     main="Mayoscore 4, Method = NNE \n  Year = 1")## \n换行符
abline(0,1,col="gray",lty=2)##线条颜色
image.png
  • MAYOSCORE 4作为marker, KM法计算
## MAYOSCORE 4, METHOD = KM
Mayo4.2= survivalROC(Stime=mayo$time,  
    status=mayo$censor,      
    marker = mayo$mayoscore4,     
    predict.time =  cutoff, method="KM")

plot(Mayo4.2$FP, Mayo4.2$TP, 
    type="l",col="red",xlim=c(0,1), ylim=c(0,1),   
    xlab=paste( "FP", "\n", "AUC = ",round(Mayo4.2$AUC,3)), 
    ylab="TP",
    main="Mayoscore 4, Method = KM \n Year = 1")
abline(0,1,col="gray",lty=2)
image.png
  • 将两个ROC曲线绘制到一起
    • lines函数在原有基础上继续绘图
    • legend函数增加legend
    • 这样的基础绘图方式,代码比较复杂,而且并不是很美观
## NNE法
plot(Mayo4.1$FP, Mayo4.1$TP, ## x=FP,y=TP
     type="l",col="red", ##线条设置
     xlim=c(0,1), ylim=c(0,1),   
     xlab=("FP"), ##连接
     ylab="TP",
     main="Time dependent ROC")## \n换行符
abline(0,1,col="gray",lty=2)##线条颜色

## KM法
lines(Mayo4.2$FP, Mayo4.2$TP, type="l",col="green",xlim=c(0,1), ylim=c(0,1))
legend(0.6,0.2,c(paste("AUC of NNE =",round(Mayo4.1$AUC,3)),
                 paste("AUC of KM =",round(Mayo4.2$AUC,3))),
                 x.intersp=1, y.intersp=0.8,
                 lty= 1 ,lwd= 2,col=c("red","green"),
                 bty = "n",# bty框的类型
                 seg.len=1,cex=0.8)# 

image.png
  • ggsci颜色美化
require(ggsci)
library("scales")
pal_nejm("default")(8)
show_col(pal_nejm("default")(8))

## NNE法
plot(Mayo4.1$FP, Mayo4.1$TP, ## x=FP,y=TP
     type="l",col="#BC3C29FF", ##线条设置
     xlim=c(0,1), ylim=c(0,1),   
     xlab=("FP"), ##连接
     ylab="TP",
     main="Time dependent ROC")## \n换行符
abline(0,1,col="gray",lty=2)##线条颜色

## KM法
lines(Mayo4.2$FP, Mayo4.2$TP, type="l",col="#0072B5FF",xlim=c(0,1), ylim=c(0,1))
legend(0.6,0.2,c(paste("AUC of NNE =",round(Mayo4.1$AUC,3)),
                 paste("AUC of KM =",round(Mayo4.2$AUC,3))),
                 x.intersp=1, y.intersp=0.8,
                 lty= 1 ,lwd= 2,col=c("#BC3C29FF","#0072B5FF"),
                 bty = "n",# bty框的类型
                 seg.len=1,cex=0.8)# 

image.png

很显然这样的配色让图上了一个档次

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 204,684评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 87,143评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 151,214评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,788评论 1 277
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,796评论 5 368
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,665评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,027评论 3 399
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,679评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 41,346评论 1 299
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,664评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,766评论 1 331
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,412评论 4 321
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,015评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,974评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,203评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 45,073评论 2 350
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,501评论 2 343

推荐阅读更多精彩内容