前情回顾:
医学大数据|R|竞争风险模型:基础、R操作与结果解读-CSDN博客
代码复习,但是大家可见得知道图画的比较丑。
library("survival")
library("cmprsk")
library("mgus2")
data(mgus2)
#预处理
mgus2<-as.data.frame(mgus2)
data<-as.data.frame(mgus2)mgus2$etime <- with(mgus2, ifelse(pstat==0, futime, ptime)) #当pstat==0时,etime=futime,否则etime=ptime
#实际上这个地方,etime当发生竞争事件的时候,比如发生死亡,那么etime等于0-死亡时间
#当没有发生竞争事件的时候,etime等于0-跌落时间mgus2$event <- with(mgus2, ifelse(pstat==0, 2*death, 1)) #当pstat==0时,event=2*death ,否则event=1
# 0 1 2 labels分别对应("censor", "pcm", "death")
event <- factor(event, 0:2)
单因素分析(cuminc)
cumic1<- cuminc(etime,event)
plot(cumic1,xlab = 'Month', ylab = '单因素的Fine-Gray检验',lwd=2,lty=1,
col = c('red','blue'))
legend(0,800,
c("PCM","death"))
print(cumic1)
cumic2<- cuminc(etime,event,mgus2$sex)
print(cumic2)
plot(cumic2,xlab = 'Month', ylab = '单因素的Fine-Gray检验',lwd=2,lty=1,
col = c('red','blue','black','forestgreen'))
#多因素分析
dt<-na.omit(mgus2)
dt<-as.data.frame(dt)
cov <- data.frame(age = dt$age,
sex = ifelse(dt$sex=='M',1,0), ## 设置哑变量
hgb = dt$hgb)
##构建多因素的竞争风险模型。此处需要指定failcode=1, cencode=0,
#分别代表结局事件赋值1与截尾赋值0,其他赋值默认为竞争风险事件2。fit<- cuminc(dt$etime,dt$event,dt$sex)
summary(dt)
crr<-crr(dt$etime, dt$event, cov, failcode=1, cencode=0)
print(crr)
所以我们接下来进一步优化。
在上述的基础上,我们知道我们完成了结果呈现,我们可以从得出单因素结果开始,
把结果的数值导出
event <- data.frame(event_time = cumic1[[1]][[1]], event_c = cumic1[[1]][[2]])
death <- data.frame(death_time = cumic1[[2]][[1]], death_C = cumic1[[2]][[2]])
怎么选择x和y的画布的最小值和最大值呢,可以用summay或者table中的描述性统计的结果
table(event$event_time)
table(event$event_time)
结果如下:
0 2 4 5 6 8 9 10 11 12 13 14 16 17 21 22 23 29 30 33 34 35 36 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 38 39 42 44 45 51 52 56 57 60 61 62 63 67 69 70 73 74 76 79 80 81 83 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 84 86 90 91 93 97 98 101 102 109 111 114 116 118 121 123 124 128 135 138 142 147 150 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 152 153 158 161 165 166 168 179 180 188 190 198 201 228 238 259 275 312 340 373 424 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1
你会发现最大值为424,所以我们x最大值选择450k爱是
发现最大结果为0.16和0.8,因此我们选择最大的结果为0.8387
开始画原始的图
p1<-ggplot(event, aes(event_time, event_c)) +
geom_line(color = "#ff4e00", linewidth = 0.7) +
scale_color_manual(values = "#440154FF") +
labs(x = "time", y = "单因素的Fine-Gray检验") +
scale_y_continuous(limits = c(0,1),expand = c(0,0))+
scale_x_continuous(limits = c(0,450),expand = c(0,0))+
theme_bw() +
theme(legend.position = "top")
print(p1)
summary(death$death_C)
p2<-ggplot(death, aes(death_time, death_C)) +
geom_line(color = "#440154FF", linewidth = 0.7) +
scale_color_manual(values = "#440154FF") +
labs(x = "time", y = "单因素的Fine-Gray检验") +
scale_y_continuous(limits = c(0,1),expand = c(0,0))+
scale_x_continuous(limits = c(0,450),expand = c(0,0))+
theme_bw() +
theme(panel.background = element_rect(fill = NA),#去除背景这一步很关键,否则后续合并图像会导致背景覆盖不显示
legend.position = "top")
print(p2)
首先我们要把这两张图画出来
print(p1)
print(p2)
1.假如你想画在两侧
cowplot::plot_grid(p1,p2,ncol = 2)
结果如下
2》我是想画在同一个画面里
##分别获取基于ggplot绘制的两张图象
g1 <- ggplotGrob(p1)
g2 <- ggplotGrob(p2)
pos <- c(subset(g1$layout, name == "panel", select = t:r))
library(gtable) # Arrange 'Grobs' in Tables
library(grid) # The Grid Graphics Package
g <- gtable_add_grob(g1, g2$grobs[[which(g2$layout$name == "panel")]],
pos$t, pos$l, pos$b, pos$l)
plot(g)
结果如下:
最后再进行美化就好了
整体代码复习如下:
rm(cumic)
library("survival")
library("cmprsk")
library("mgus2")
data(mgus2)
#预处理
mgus2<-as.data.frame(mgus2)
data<-as.data.frame(mgus2)
mgus2$etime <- with(mgus2, ifelse(pstat==0, futime, ptime)) #当pstat==0时,etime=futime,否则etime=ptime
#实际上这个地方,etime当发生竞争事件的时候,比如发生死亡,那么etime等于0-死亡时间
#当没有发生竞争事件的时候,etime等于0-跌落时间
mgus2$event <- with(mgus2, ifelse(pstat==0, 2*death, 1)) #当pstat==0时,event=2*death ,否则event=1
# 0 1 2 labels分别对应("censor", "pcm", "death")
event <- factor(event, 0:2)
单因素分析(cuminc)
cumic1<- cuminc(etime,event)
plot(cumic1,xlab = 'Month', ylab = '单因素的Fine-Gray检验',lwd=2,lty=1,
col = c('red','blue'))
event <- data.frame(event_time = cumic1[[1]][[1]], event_c = cumic1[[1]][[2]])
death <- data.frame(death_time = cumic1[[2]][[1]], death_C = cumic1[[2]][[2]])
table(event$event_time)
summary(event$event_c)
p1<-ggplot(event, aes(event_time, event_c)) +
geom_line(color = "#ff4e00", linewidth = 0.7) +
scale_color_manual(values = "#440154FF") +
labs(x = "time", y = "单因素的Fine-Gray检验") +
scale_y_continuous(limits = c(0,1),expand = c(0,0))+
scale_x_continuous(limits = c(0,450),expand = c(0,0))+
theme_bw() +
theme(legend.position = "top")
print(p1)
summary(death$death_C)
p2<-ggplot(death, aes(death_time, death_C)) +
geom_line(color = "#440154FF", linewidth = 0.7) +
scale_color_manual(values = "#440154FF") +
labs(x = "time", y = "单因素的Fine-Gray检验") +
scale_y_continuous(limits = c(0,1),expand = c(0,0))+
scale_x_continuous(limits = c(0,450),expand = c(0,0))+
theme_bw() +
theme(panel.background = element_rect(fill = NA),#去除背景这一步很关键,否则后续合并图像会导致背景覆盖不显示
legend.position = "top")
print(p2)
cowplot::plot_grid(p1,p2,ncol = 2)
##分别获取基于ggplot绘制的两张图象
g1 <- ggplotGrob(p1)
g2 <- ggplotGrob(p2)
pos <- c(subset(g1$layout, name == "panel", select = t:r))
library(gtable) # Arrange 'Grobs' in Tables
library(grid) # The Grid Graphics Package
g <- gtable_add_grob(g1, g2$grobs[[which(g2$layout$name == "panel")]],
pos$t, pos$l, pos$b, pos$l)
plot(g)
责任编辑:医学大数据刘刘老师:头部医疗大数据公司医学科学部研究员
邮箱:897282268@qq.com
久菜盒子工作室 我们是:985硕博/美国全奖doctor/计算机7年产品负责人/医学大数据公司医学研究员/SCI一区2篇/Nature子刊一篇/中文二区核心一篇/都是我们
主要领域:医学大数据分析/经管数据分析/金融模型/统计数理基础/统计学/卫生经济学/流行与统计学/
擅长软件:R/python/stata/spss/matlab/mySQL
团队理念:从零开始,让每一个人都得到优质的科研教育