用诺模图可视化你的模型

花花写于2020.6.2
诺模图出镜率很高,用于多因素cox或logstic模型的可视化展示。

0.输入数据

需要病人临床信息和生存信息表格,如下

rm(list = ls())
load("ph.Rdata")
head(ph)
##                         sample_id gender age T N  M stage event time
## TCGA-MP-A4T4-01A TCGA-MP-A4T4-01A female  68 2 1  0     2     1 2617
## TCGA-05-4250-01A TCGA-05-4250-01A female  79 3 1  0     3     1  121
## TCGA-64-5774-01A TCGA-64-5774-01A   male  60 2 0  0     1     0 2676
## TCGA-97-7937-01A TCGA-97-7937-01A   male  65 2 0 NA     1     0  564
## TCGA-97-A4M6-01A TCGA-97-A4M6-01A female  45 1 0  0     1     0  568
## TCGA-44-A4SU-01A TCGA-44-A4SU-01A female  67 1 0 NA     1     1  409
##                            fp
## TCGA-MP-A4T4-01A -0.005061199
## TCGA-05-4250-01A  0.286351686
## TCGA-64-5774-01A  0.389468347
## TCGA-97-7937-01A  0.098764940
## TCGA-97-A4M6-01A -0.325363900
## TCGA-44-A4SU-01A -0.183674628

1.简约诺模图

使用rms包里的cph函数建模,nomogram函数画图

library(rms)
dd<-datadist(ph)
options(datadist="dd")
mod <- cph(formula =  as.formula(paste("Surv(time, event) ~ ",paste(colnames(ph)[2:7],collapse = "+"))),
           data=ph,x=T,y=T,surv = T)

surv<-Survival(mod) 
surv3<-function(x) surv(1095,x)
surv5<-function(x) surv(1825,x)

x<-nomogram(mod,
            fun = list(surv3,surv5),
            funlabel = c('3-year survival Probability',
                         '5-year survival Probability'))

plot(x)

2.美化版本

regplot画出的图和上面的简约版意义一致,只是可修改的细节多一些。

library(regplot)
library(survival)

mod2 <- coxph(formula =  as.formula(paste("Surv(time, event) ~ ",paste(colnames(ph)[2:7],collapse = "+"))),data=ph)

2.1 简约风

regplot(mod2,
        failtime = c(1095,1825), 
        plots = c("no plot","no plot"),
        points = T,
        prfail = T)

2.2 在图上标出某个病人

regplot(mod2,
        observation=ph[1,], 
        obscol = "#326db1",
        failtime = c(1095,1825), 
        plots = c("no plot","no plot"),
        points = T,
        prfail = T)

2.3. 增加表示分布情况的图形

regplot(mod2,
        observation=ph[1,], 
        failtime = c(1095,1825), 
        plots = c("bars","boxes"),
        points = T,
        prfail = T)

帮助文档中对plots参数的解读:The plots parameter specifies initial plot types. It is length 2. The first item specifies a plot type for non-factor variables as one of: “no plot”, “density”, “boxes”, “spikes”, “ecdf”, “bars”, “boxplot”, “violin” or “bean”. The second item, is for factors and is one of: “no plot”, “boxes”, “bars” or “spikes”.

其他图形也可以试试哦。我觉得简单的就挺好的,已经可以说明问题了。

3.校准曲线

校准曲线的作用是展示模型质量,曲线越贴近0-1对角线越好。

3.1 建模和完成计算

cph的time.inc参数和calibrate的u参数,后面都是写天数,三年就是1095天,5年就是1825天

f3 <- cph(formula =  as.formula(paste("Surv(time, event) ~ ",paste(colnames(ph)[2:7],collapse = "+"))),
           data=ph,x=T,y=T,surv = T, time.inc=1095)
cal3 <- calibrate(f3, cmethod="KM", method="boot", u=1095, m=50, B=1000)
## Using Cox survival estimates at 1095 Days
f5 <- cph(formula =  as.formula(paste("Surv(time, event) ~ ",paste(colnames(ph)[2:7],collapse = "+"))),
           data=ph,x=T,y=T,surv = T,  time.inc=1825)
cal5 <- calibrate(f5, cmethod="KM", method="boot", u=1825, m=50, B=1000)
## Using Cox survival estimates at 1825 Days

3.2 画校准曲线图

plot(cal5)就可以画,但出来的图一言难尽,主要是参考线画的并不是0-1,所以需要下面的代码手动去画咯。

plot(cal3,lwd = 2,lty = 0,errbar.col = c("#2166AC"),
     bty = "l", #只画左边和下边框
     xlim = c(0,1),ylim= c(0,1),
     xlab = "Nomogram-prediced OS (%)",ylab = "Observed OS (%)",
     col = c("#2166AC"),
     cex.lab=1.2,cex.axis=1, cex.main=1.2, cex.sub=0.6)
lines(cal3[,c('mean.predicted',"KM")],
      type = 'b', lwd = 1, col = c("#2166AC"), pch = 16)
mtext("")

plot(cal5,lwd = 2,lty = 0,errbar.col = c("#B2182B"),
     xlim = c(0,1),ylim= c(0,1),col = c("#B2182B"),add = T)
lines(cal5[,c('mean.predicted',"KM")],
      type = 'b', lwd = 1, col = c("#B2182B"), pch = 16)

abline(0,1, lwd = 2, lty = 3, col = c("#224444"))

legend("topleft", #图例的位置
       legend = c("3-year","5-year"), #图例文字
       col =c("#2166AC","#B2182B"), #图例线的颜色,与文字对应
       lwd = 2,#图例中线的粗细
       cex = 1.2,#图例字体大小
       bty = "n")#不显示图例边框

参考:https://blog.csdn.net/weixin_30592887/article/details/112865824 https://www.jianshu.com/p/9085e4e13843 https://www.dxy.cn/bbs/newweb/pc/post/33323120

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

推荐阅读更多精彩内容