R语言cmprsk_hr函数进行加权竞争风险模型分析（1）—适用于nhanes复杂加权、逆概率加权、重叠加权等

原创于 2026-02-26 10:17:32 发布 · 527 阅读 ·

本内容遵循CC 4.0 BY-SA版权协议

GEO检测

天桥下的卖艺者

关注

标签

#r语言 #加权竞争风险模型 #nhanes

分类数据科学与机器学习开发

R语言同时被 2 个专栏收录

401 篇文章

订阅专栏

我写的R包

95 篇文章

订阅专栏

竞争风险模型就是指在临床事件中出现和它竞争的结局事件，这是事件会导致原有结局的改变，因此叫做竞争风险模型。比如我们想观察患者肿瘤的复发情况，但是患者在观察期突然车祸死亡，或者因其他疾病死亡，这样我们就观察不到复发情况了，这种情况下不能把缺失数据仅仅当做右删失处理，这样的话会造成数据的估值错误。这是我们应该优先选择竞争风险模型来做数据分析，而不是COX回归。竞争风险模型主要是统计学家想了解真实的事件发生对结局的影响。

在这里插入图片描述
既往我们在多篇文章《手把手教你使用R语言做竞争风险模型并绘制列线图》、《一步到位：手把手教你R语言竞争风险模型建模-列线图-校准曲线-K折验证-外部验证- 决策曲线》、《R语言tidycmprsk包分析竞争风险模型》等文章已经介绍了竞争风险模型的各种应用，

有不少粉丝私信问加权竞争风险模型怎么做，我也花了点时间研究了一下，先把两种竞争风险模型的原理研究了一下，还写了一篇文章《R语言两种方法手搓竞争风险模型（1）》，对比cmprsk包，手搓出来的结果基本一样。
因为能参照的东西不多，加权竞争风险模型函数cmprsk_hr从研究到写出框架前后基本花了一个多月，目前基本把框架写出来了，但是还有很多完善的空间，下面我来演示一下。

先导入数据和R包

library(scitable)
library(ggscitable)
setwd("E:/公众号文章2026年/加权竞争风险模型")
follic1<-read.csv("follic1.csv",sep=',',header=TRUE)

在这里插入图片描述
Trt是咱们研究的变量，单独治疗或者联合治疗，status：结局变量，必须是0，1，2否则就会报错，time是时间，其他是协变量

cmprsk_hr函数做加权竞争风险模型非常简单，就是一句话代码, data就是你的数据，time就是时间，y是结局，x是你的研究变量，weights这里

tab2 <- cmprsk_hr(data=follic1, time = "time",y =  "status",x="Trt",weights = "ipw_ate_stab")

生成结果图下，主要是看out2，下面那个是为后期做列线图的，不用看它
在这里插入图片描述
下图可见生成2个结局的系数，咱们主要看结局1的，这样结果就轻易出来啦

tab2[["out2"]]

在这里插入图片描述
咱们还可以加入协变量，可以看到，协变量对研究结果影响还是挺大的

tab3 <- cmprsk_hr(follic1, time = "time",y =  "status",x="Trt",weights = "ipw_ate_stab",cov = c("age","hgb"))
tab3[["out2"]]

在这里插入图片描述
接下可以绘制加权累积生存曲线图，做法和上面基本一样。因为cif观察的是不同组别生存率的不同，我看了很多资料，这里是不加协变量的，加了也能做，但是算法改变了，生成的图就会有点怪，risk.time=10，这里的意思是比较10年累积生存率有无不同。

out<-ggcmprskcif(data = follic1,x="treatment",y="status",weights = "ipw_ate_stab",time = "time",boot = F,
                 risk.time=10)

在这里插入图片描述
还可以做一些细节的修改，比如该线条颜色

out<-ggcmprskcif(data = follic1,x="treatment",y="status",weights = "ipw_ate_stab",time = "time",boot = F,
                 risk.time=10,col = c("red","blue"))

在这里插入图片描述
修改图例位置

out<-ggcmprskcif(data = follic1,x="treatment",y="status",weights = "ipw_ate_stab",time = "time",boot = F,
                 risk.time=10,col = c("red","blue"),legend.position = c(0.95,0.25))

在这里插入图片描述
如果需要可信区间，可以加上boot=T参数

out<-ggcmprskcif(data = follic1,x="treatment",y="status",weights = "ipw_ate_stab",time = "time",risk.time=10,boot = T)

在这里插入图片描述
如果需要比较两条曲线生存率有无不同，需要定个具体时间，比如我这里定的是10年，第一个是绝对生存率比较，就是对照组减去参考组，第二个是相对生存率比较，就是对照组减去参考组除以参考组

out[["risk"]]

在这里插入图片描述
如果你的绘图能力还行，也可以提取绘图数据，自定义绘图

dt<-out[["cif.data"]]

在这里插入图片描述
Nhanes数据也是能做的，下面我来简单演示一下，导入一个nhanes数据

nhanes_sim<-read.csv("nhanes_sim.csv",sep=',',header=TRUE)
nhanes_sim$status<-ifelse(nhanes_sim$status=="Censored",0,ifelse(nhanes_sim$status=="Event1",1,2))

在这里插入图片描述
数据我就不解释了，咱们按照nhanes的方法来就行了

design <- svydesign(
  ids = ~ psu,          # 聚类变量（初级抽样单元）
  strata = ~ strata,    # 分层变量（可选，但建议提供）
  weights = ~ weight,   # 抽样权重
  data = nhanes_sim,
  nest = TRUE           # 表示 PSU 编号在各层内是独立的（NHANES 需要设为 TRUE）
)

唯一不同的是weights参数改成了svydstr

tab2 <- cmprsk_hr(nhanes_sim, time = "time",y = "status",x="sex",svydstr = design)
tab2[["out2"]]

在这里插入图片描述
绘图

out<-ggcmprskcif(nhanes_sim, time = "time",y = "status",x="sex",svydstr = design,risk.time=10,boot = F)

在这里插入图片描述
本函数目前框架是基本写出来了，后面有时间在慢慢完善了，目前局限是:绘图只能绘制二分类，三分类以上不行，还待开发的有列线图、校准曲线等。

文字苦短，下面视频演示更加详细精彩