R语言与竞争风险模型第一篇

2021-11-25 本文已影响0人小易学统计

转自个人微信公粽号【易学统计】的统计学习笔记：R语言与竞争风险模型
题记：本文是个人的读书笔记，仅用于学习交流使用。

01解决什么问题

研究治疗方案A和肿瘤复发的关系，如果患者在复查的路上出车祸意外死亡，就观察不到肿瘤复发了，也就是说“车祸死亡” 和“复发”存在竞争。这种与复发无关的死亡，失访等，就称为竞争事件。

当存在竞争事件时，传统的做法，通常是将竞争事件的发生视作刪失（censor），来计算结局事件的生存概率，采用cox模型分析，但这会导致结局事件的生存概率估计有偏差，一般会高估累计风险发生率，因此在存在竞争风险时，我们优先使用竞争风险模型，而不是cox等其他风险模型。

竞争风险的单因素分析常用来估计关心终点事件的发生率，多因素分析常用来探索预后影响因素及效应值。

02案例分析

本案例数据来自R自带数据集(survival包的mgus2)，共计1341名单克隆丙种球蛋白病患者，结局事件定义为“发生浆细胞恶性肿瘤”，而某些患者在“发生浆细胞恶性肿瘤”之前因为其他原因死亡，那这些发生其他死亡的患者就无法观察到“发生浆细胞恶性肿瘤”的终点，也就是说“其他死亡”与“发生浆细胞恶性肿瘤”存在竞争风险。故采用竞争风险模型分析。

2.1 加载数据

加载竞争风险模型的包cmprsk，加载数据，并定义结局为因子变量

library(cmprsk)
library(survival)
dt<- mgus2
str(dt)
dt$etime <- ifelse(dt$pstat==0, dt$futime, dt$ptime)
dt$event <- ifelse(dt$pstat==0, 2*dt$death, 1)
dt$event <- factor(dt$event, 0:2)

2.2 数据说明

hgb 血红蛋白；
creat 肌酐；
ptime 直至发展为浆细胞恶性肿瘤（PCM）或最后一次访视的时间（以月为单位）；
pstat 出现PCM(浆细胞恶性肿瘤)：0 =否，1 =是；
futime 直到死亡或最后一次接触的时间，以月为单位；
death 发生死亡：0 =否，1 =是；这里把PCM作为结局事件，death作为竞争事件。

图1.png

2.3 单因素分析

类比两组生存资料log-rank检验，考虑竞争风险事件，同样可以进行单因素分析与多因素分析。下面我们就可以使用cuminc()函数进行单因素的Fine-Gray检验

fit0 <- cuminc(dt$etime,dt$event,dt$sex)
fit0
plot(fit0,xlab = 'Month', ylab = 'CIF',lwd=2,lty=1,
     col = c('red','blue','black','forestgreen'))

图3.png

结果解读：

Test第一行统计量=1.194508, P=0.274422157,表示在控制了竞争风险事件（即第二行计算的统计量和P值）后，男性和女性的累计发展为PCM的风险无统计学差异P=0.274422157。
est表示估计的各时间点男性和女性的累计发展为PCM发生率与累计竞争风险事件发生率（分别用1和2来区分，与第一行第二行一致）。
var表示估计的各时间点男性和女性的累计发展为PCM发生率与累计竞争风险事件发生率的方差（分别用1和2来区分，与第一行第二行一致）。

绘制曲线：

下面我们画出累计PCM发生率与累计竞争风险事件发生率的生存曲线，直观表示上述数字化结果。

图2.png

图形解读：

纵坐标表示累计发生率pcm，横坐标是时间轴。我们从F1对应的红色曲线和M1对应的蓝色曲线可以得出，女性组的发生pcm风险较男性组高，但未达到统计学意义，P=0.27。同理，F2对应的黑色曲线在M2对应的草绿色曲线下方，我们可以得出，女性组的竞争风险事件发生率较男性组低，具有统计学差异，P=0.0064。简单来讲，这个图可以用一句话来概括：在控制了竞争风险事件后，女性组和男性组累计PCM发生率无统计学差异P=0.27。

2.4 竞争风险模型（多因素分析）

下面我们进行考虑竞争风险事件的生存资料的多因素分析方法。在cmprsk包中crr()函数可以很方便的实现多因素分析。

###多因素竞争风险
cov <- data.frame(age = dt$age,
                  sex = ifelse(dt$sex=='M',1,0), ## 设置哑变量
                  hgb = dt$hgb) 

##构建多因素的竞争风险模型。此处需要指定failcode=1, cencode=0, 分别代表结局事件赋值1与截尾赋值0，其他赋值默认为竞争风险事件2。
fit2 <- crr(dt$etime, dt$event, cov, failcode=1, cencode=0)
summary(fit2)

图4.png

结果解读：

在控制了竞争风险事件后，age，即疾病所处阶段是发生pcm的独立影响因素，P=0.0025。其HR及95% CI分别为0.983(0.971,0.994)，说明年龄是个保护因素，年龄越大，发生PCM的风险越低。

2.5 预测

predict_crr <- predict.crr(fit2, data.frame('age'=c(80),'sex'=c(1),'hgb'=c(15.2)))
plot(predict_crr,col='green')

图5.png

2.6 与cox模型对比

图6.png

可以看出，在存在竞争风险时，若不适时使用竞争风险模型，而使用cox等模型，会高估Cumulative event rate。

如果您觉得有用，请点赞，转发哦~
更多统计小知识，请关看 公粽号易学统计

R语言与竞争风险模型第一篇

01解决什么问题

02案例分析

2.1 加载数据

2.2 数据说明

2.3 单因素分析

结果解读：

绘制曲线：

图形解读：

2.4 竞争风险模型（多因素分析）

结果解读：

2.5 预测

2.6 与cox模型对比

猜你喜欢

热点阅读