需要JoVE订阅才能观看此内容。 请登录或开始免费试用

方法文章

基于R语言的竞争风险模型景观验证

2.6K 次观看

⸱

DOI:

10.3791/64018

⸱

2022年9月16日

 ,  ,  ,  ,  ,  ,  ,  ,  ,  , 

* These authors contributed equally

本文内容

摘要

本方案介绍了用于评估竞争风险模型区分能力和校准能力的 R 语言代码,以及对该模型进行内部和外部验证的代码。

摘要

Cox比例风险模型在临床环境中广泛应用于生存分析,但其无法处理多种生存结局。与传统的Cox比例风险模型不同,竞争风险模型考虑了竞争事件的存在,并结合列线图(nomogram)——一种图形化计算工具,为临床医生提供精确的预后预测。本研究报道了一种建立竞争风险列线图的方法,包括评估其区分能力(即一致性指数和曲线下面积)、校准能力(即校准曲线)以及净收益(即决策曲线分析)。此外,还通过重抽样原始数据集进行内部验证,并利用外部数据集对已建立的竞争风险列线图进行外部验证,以展示其外推能力。该竞争风险列线图可作为临床医生在考虑竞争风险因素下预测预后的有力工具。

引言

近年来,随着精准医学的发展,不断有新的预后因素被发现,结合分子与临床病理因素的预后模型在临床实践中正受到越来越多的关注。然而,诸如Cox比例风险模型等非图形化模型,其结果以系数值形式呈现,临床医生难以直观理解1。相比之下,列线图(nomogram)是一种回归模型(包括Cox回归模型、竞争风险模型等)的可视化工具,是一种用于近似图解计算数学函数的二维图示2。它能够量化临床模型中各变量的不同水平,并计算风险评分(RS),从而预测患者的预后。

模型评估是模型构建中的关键环节,通常有两个被广泛接受的评估特征:区分度和校准度。在临床模型中,区分度指的是模型将发生事件的个体与未发生事件的个体区分开的能力,例如死亡患者与存活患者的区分,通常使用一致性指数(C指数)或受试者工作特征曲线下面积(AUC)来表征3,4。校准度是指将模型预测的概率与实际观察到的概率进行比较的过程,校准曲线已被广泛用于展示这一特性。此外,模型验证(包括内部验证和外部验证)是模型构建中的重要步骤,只有经过验证的模型才能进一步推广应用5。

Cox 比例风险模型是一种在医学研究中用于分析预后因素与生存状态之间关联的回归模型。然而,Cox 比例风险模型仅考虑两种结局状态 [Y (0, 1)],而研究对象常常面临多于两种的状态,此时会出现竞争风险 [Y (0, 1, 2)]1。总生存期(OS)定义为从起始日期(例如治疗开始)到因任何原因死亡的时间,是生存分析中最重要的终点指标。然而,OS 无法区分癌症特异性死亡与非癌症特异性死亡(如心血管事件及其他无关原因),从而忽略了竞争风险的影响6。在这些情况下,应优先采用竞争风险模型来预测考虑了竞争风险的生存状态7。Cox 比例风险模型的构建与验证方法已较为成熟,而关于竞争风险模型的验证报道则较少。

在我们之前的研究中,建立了一个特定的竞争风险列线图,该列线图结合了列线图与竞争风险模型,并基于竞争风险模型构建了风险评分估计方法8。本研究旨在介绍已建立的竞争风险列线图的不同评估与验证方法,该列线图可作为临床医生在考虑竞争风险的前提下预测预后的实用工具。

访问受限。请登录或开始试用以查看此内容。

方案

监测、流行病学和最终结果(SEER)数据库是一个开放获取的癌症数据库,仅包含去标识化的患者数据(SEER ID:12296-Nov2018)。因此,本研究免除了浙江大学医学院附属金华医院伦理审查委员会的审批要求。

1. 数据准备和 R 软件包准备

  1. 准备并导入数据。
    > 数据集 <- read.csv("…/Breast cancer Data.xlsx") #导入数据。
    注意:数据已上传至补充文件1。
  2. 安装并加载R软件包。
    > 软件包 <- c("rms","cmprsk","mstate","survival","riskRegression","
    prodlim")
    > req.pcg <- 函数(pcg){
    新 <- pcg[!(pcg %in% installed.packages()[, "Package"])]
    如果 (length(新)) install.packages(新, dependencies = T)
    sapply(pcg, require, ch = T)
     }
    > req.pcg(软件包)

    ​注意:以下操作基于R软件(版本3.6.2)进行,使用rms、cmprsk、mstate、survival、riskRegression和prodlim软件包(http://www.r-projectrg/)。

2. 采用两种不同方法建立竞争风险列线图

  1. 采用直接法建立竞争风险列线图。
    > mod_cph <- cph(Surv(Survivalmonths, status) ~ factor1+ factor2+…,
     x=T, y=T, surv=T, data=Dataset)
    > nom <- nomogram(mod_cph, fun=list(function(x) 1-surv_cph(36, x)…),
     funlabel=c("3-year event1 Prob."…), lp=F)
    #以第36个月为例。
    > mod_crr <- crr(Survivalmonths, fstatus, failcode=1, cov1=cov)
    > score <- log(log((1-real.3y),(1-cif.min36)))/(maxbeta/100)
    > plot(nom)
  2. 采用加权法建立竞争风险列线图。
    > df.w <- crprep("Survivalmonths"," fstatus",
     data=Dataset, trans=c(1,2), cens=0,
     keep=c("factor1"," factor2"…))
    > mod.w <- cph(Surv(Tstart, Tstop, status==1) ~ factor1+factor2+…,
     data=df.w, weight=weight.cens, subset=failcode==1, surv=T)
    ​> nom.w <- nomogram(mod.w…)

3. 竞争风险列线图的区分能力

  1. 区分度C指数
    1. 拟合矩阵 cov 进入竞争风险模型 mod_crr 并获得一个预测矩阵 suv.
      > suv <- predict.crr(mod_crr, cov)
    2. 获取特定月份的累积发病率 suv 并使用该函数计算C指数 rcorr.cens.
      > cif36 <- suv[which(suv[,1]==36),][-1]
      > rcorr <- rcorr.cens(1-cif36, Surv(Dataset$Survivalmonths, Dataset$tumordeath))
      > cindex <- rcorr[1]
  2. 区分的AUC
    1. 使用该函数评估竞争风险模型的预测性能 评分 (riskRegression 包)
      > fgr.w <- FGR(Hist(生存月数, fstatus) ~ 因子1 + 因子2 + …, data=数据集, cause=1)
      > 评分 <- 评分(列表("Fine-Gray" = fgr.w),
    2. 从……中提取AUC "评分".
      > score$AUC

4. 竞争风险模型的校准能力

  1. 具有95%置信区间的竞争风险模型校准曲线
    1. 获取一个数据框,其中包含每个个体在特定失效时间的累积发病率。
      > cif36 <- data.frame(cif36) 以第36个月为例。
      > cif36.36_o 的列名<-c("36m")
    2. 根据估计的累积发病率将队列分为五个亚组,并计算每个亚组的平均预测累积发病率。
      > group36 <- cut(cif36$`36m`,
      quantile(cif36$`36m`, seq(0, 1, 0.2))
      include.lowest = TRUE, labels = 1:5)
      > 平均值36 <- as.vector(by(cif36$`36m`, group36, mean))
    3. 使用该函数计算观察到的累积发病率,即实际累积发病率 姜黄素,然后在特定失效时间点获得带有95%置信区间的观察累积发病率。
      > cum36 <- cuminc(Dataset$Survivalmonths, Dataset$fstatus, group36)
      > obs36 <- 时间点(cum36, Dataset$Survivalmonths)$est[c(1:5),36]
      > obs36var <- timepoints(cum36, Dataset$Survivalmonths)$var[c(1:5), 36]
      > df <- 数据框(mean36,obs36,obs36var)
    4. 使用该函数以预测累积发病率作为 x 轴,观察累积发病率作为 y 轴绘制校准曲线 ggplot.
      > ggplot(df) + geom_point(aes(x=mean36, y=obs36), col="红色")+
       geom_point(aes(x=mean36, y=obs36), col="红色",pch=4)+
       geom_line(col="红色",aes(x=均值36,y=观测值36))+
       geom_errorbar(col="红色",aes(x=均值36,y=观测值36+1.96
      *sqrt(obs36var))
       ymin = obs36 - 1.96 × √(obs36var),ymax = obs36 + 1.96
      *sqrt(obs36var)
      geom_abline(lty=3,lwd=2,col=c(rgb(0,118,192,
      maxColorValue=255)))
  2. 竞争风险模型风险评分的校准曲线
    1. 评估所有变量的每个水平,并获得总RS值。
      > Dataset$factor1[Dataset$factor1==1] <- factor1.scale["Factor1_水平1"]
      >
      … #例如,Dataset$histology[Dataset$histology==1]<-组织学.比例尺["组织学1"]
      > Dataset$rs <- Dataset$factor1 + Dataset$factor2 + Dataset$factor3 + …
      注意:通过将每位患者的各变量得分相加,计算其总RS得分。
    2. 统计不同总风险评分的频次,并计算其观察到的累积发病率。
      > rs.freq <- as.data.frame(table(数据集$rs))
      > obs.36 <- 向量(模式="数值的",长度 = nrow(rs.freq)
      > for (i in 1:nrow(rs.freq)) {
       数据集 <- subset(数据集, 数据集$rs == rs.freq[i,1])
       cif.dataset <- cuminc(数据集$生存月数, 数据集$死亡3)
       cif36.dataset <- 时间点(cif.dataset,36)
       obs.36[i] <- cif36.dataset$est[1]}
    3. 设置x轴的范围,并计算总风险评分的预测累积发病率。
      > RS <- range(nom$total.points)
      > x.36 <- seq(min(RS), max(RS), 0.01)
      > pre.36 <- 1-(1-cif.min36)^exp(x.36*maxbeta/100)
    4. 以风险评分绘制校准曲线。
      > plot(x.36, pre.36, type='l'…)
      > par(new=TRUE)
      > plot(as.vector(rs.freq[,1]), obs.36…)

5. 竞争风险模型的决策曲线分析

  1. 载入 stdca 函数以执行决策曲线分析。
    > source("stdca.R")
  2. 从列线图中提取多项式方程,用于计算生存概率。
    > nomogramEx(nomo = nom)
    > Dataset$predictors <- A * (Dataset$rs ^3) + B * (Dataset$rs ^2) + C * Dataset$rs + D
    #predictors 为通过已建立的列线图计算得出的癌症特异性死亡预测概率
  3. 执行决策曲线分析。
    > stdca(data = Dataset, outcome = "status", ttoutcome = "Survivalmonths", timepoint = 36,
     predictors = "predictors", cmprsk = TRUE, smooth = FALSE, probability = FALSE)

    注:在存在竞争风险的情况下评估结局时,应将 cmprsk 设为 TRUE。

6. 使用自助法进行内部验证

  1. 使用自助法获取平均预测累积发病率。
    1. 重采样原始数据集(数据集)与 替换 用于生成自助数据集(数据集_in)。建立竞争风险模型(mod.in_crr) 与自举数据集一起使用该函数 predict.crr 预测 mod.in_crr 并循环 b 次以生成 suvall.in.
      B=b
      suvall.in <- 列表()
      for(j in 1:B){
       数据集_in <- 数据集[样本(c(1:nrow(数据集)), nrow(数据集),
      replace = TRUE),]
      连接(数据集_输入)
       覆盖范围 <- model.matrix(~因子1 + 因子2 + …)[,-1]
       _crr 中的修饰 <- crr(生存月数, fstatus, failcode=1, cov1=cov.in)
      分离(Dataset. inner)
      suv. in <- predict.crr(mod. in _crr, cov)
      suvall.in[[j]] <- suv.in}
    2. 获取某个月份的平均预测累积发病率。
      cif36all. 内层 <- 向量(模式="数值的",长度=nrow(数据集)
      对于(k 从 1 到 B){
       cif36all. 内部<- cif36all. inner + suvall. inner[[k]][which(suvall. inner[[k]][,1] == 36),][-1]
      }
      cif36.in <- cif36all.in/B
  2. 使用函数通过内部交叉验证计算C指数 rcorr.cens.
    rcorr. 内部 <- rcorr.cens(1-cif36.in, Surv(Dataset$Survivalmonths, Dataset$tumordeath))
    cindex. 内部 <- rcorr. inner[1]
  3. 使用交叉内部验证进行校准。
    ​注意:具有内部验证的竞争风险模型校准曲线的代码与第4节中的代码类似,而 suv 被替换为 suv.in.

7. 竞争风险模型的外部验证

  1. 利用外部数据获取预测的累积发病率。使用外部数据变量矩阵(cov.ex)计算预测的累积发病率。
    suv.ex <- predict.crr(mod_crr,cov.ex)
    cif36.ex <- suv.ex [which(suv.ex $time=="36"),][-1]
  2. 使用外部验证计算C指数。
    rcorr.ex <- rcorr.cens(1-cif36.ex,Surv(Dataset.ex$Survivalmonths,Dataset.ex$tumordeath))
    cindex.ex <- rcorr.ex[1]
  3. 使用外部验证进行校准。
    注意:具有内部验证的竞争风险模型的校准曲线代码与第4节中的代码类似,只是将suv替换为suv.ex。

访问受限。请登录或开始试用以查看此内容。

结果

本研究从SEER数据库中获取了乳腺癌患者的数据,并将其作为示例数据。SEER数据库提供了代表美国约34.6%人口的癌症数据,且已获得访问该数据库的许可(参考编号:12296-Nov2018)。

建立了两个列线图(图1),二者均包含组织学类型、分化程度、T分期和N分期,分别采用直接法和加权法构建。各变量水平的评分及其对应的总评分概率几乎相同,但也可观察到一些细微差异。Zhang等人提出了一种“加权”方法用于构建竞争风险列线图,该方法首先将原始数据转换为加权数据(使用函数crprep),然后利用加权数据构建Cox回归模型(使用函数coxph),最后基于Cox回归模型建立竞争风险列线图(使用函数nomogram)9。相比之下,本研究中的“直接”方法与“加权”方法完全不同。简而言之,本研究使用竞争风险模型(使用函数crr)生成的参数替代了Cox回归模型(使用函数coxph)中的参数,最终利用这些...

访问受限。请登录或开始试用以查看此内容。

讨论

本研究比较了通过两种不同方法建立的竞争风险列线图,并对所建立的列线图进行了评估与验证。具体而言,本研究提供了基于直接法建立列线图的逐步教程,以及计算C指数和绘制校准曲线的方法。

R 软件中的 rms 包被广泛用于构建和评估 Cox 比例风险模型,但不适用于竞争风险模型。对于具有多个结局的模型,Zhang 等人报道了使用 R 软件中的 riskRegression 包对竞争风险模型进行验证的方法,该方法通过计算 AUC 和 Brier 评分来评估模型的区分能力,并绘制校准曲线以验证其校准能力9,10。

然而,riskRegression 软件包存在一些不足之处,例如无法计算C指数。C指数是指在所有可评估的有序对中,一致有序对所占的比例(C指数 = 一致可评估有序对 / 所有可评估有序对)3,...

访问受限。请登录或开始试用以查看此内容。

披露

作者声明不存在任何竞争利益。

致谢

本研究由医学科学基金资助 &浙江省科技计划项目(项目编号:2013KYA212)、浙江省自然科学基金一般项目(项目编号:Y19H160126)以及金华市科学技术局重点项目 & 科技局(项目编号:2016-3-005、2018-3-001d 和 2019-3-013)。

访问受限。请登录或开始试用以查看此内容。

材料

本文使用的材料清单
姓名公司目录编号评论
R 软件无不适用版本 3.6.2 或更高 
计算机系统Microsoft Windows 10 Windows 10 或更高版本

参考文献

  1. Andersen, P. K., Gill, R. D. Cox's regression model for counting processes: A large sample study. The Annals of Statistics. 10 (4), 1100-1120 (1982).
  2. Lubsen, J., Pool, J., vander Does, E. A practical device for the application of a diagnostic or prognostic function. Methods of Information in Medicine. 17 (2), 127-129 (1978).
  3. Harrell, F. E., Lee, K. L., Mark, D. B. Multivariable prognostic models: Issues in developing models, evaluating assumptions and adequacy, and measuring and reducing errors. Statistics In Medicine. 15 (4), 361-387 (1996).
  4. Hung, H., Chiang, C. -T. Estimation methods for time-dependent AUC models with survival data. The Canadian Journal of Statistics / La Revue Canadienne de Statistique. 38 (1), 8-26 (2010).
  5. Moons, K. G. M., et al. Risk prediction models: I. Development, internal validation, and assessing the incremental value of a new (bio)marker. Heart. 98 (9), 683-690 (2012).
  6. Fu, J., et al. Real-world impact of non-breast cancer-specific death on overall survival in resectable breast cancer. Cancer. 123 (13), 2432-2443 (2017).
  7. Fine, J. P., Gray, R. J. A proportional hazards model for the subdistribution of a competing risk. Journal of the American Statistical Association. 94 (446), 496-509 (1999).
  8. Wu, L., et al. Establishing a competing risk regression nomogram model for survival data. Journal of Visualized Experiments. (164), e60684(2020).
  9. Zhang, Z., Geskus, R. B., Kattan, M. W., Zhang, H., Liu, T. Nomogram for survival analysis in the presence of competing risks. Annals of Translational Medicine. 5 (20), 403(2017).
  10. Zhang, Z. H., et al. Overview of model validation for survival regression model with competing risks using melanoma study data. Annals Of Translational Medicine. 6 (16), 325(2018).
  11. Newson, R. Confidence intervals for rank statistics: Somers' D and extensions. Stata Journal. 6 (3), 309-334 (2006).
  12. Davison, A. C., Hinkley, D. V., Schechtman, E. Efficient bootstrap simulation. Biometrika. 73 (3), 555-566 (1986).
  13. Roecker, E. B. Prediction error and its estimation for subset-selected models. Technometrics. 33 (4), 459-468 (1991).
  14. Steyerberg, E. W., Harrell, F. E. Prediction models need appropriate internal, internal-external, and external validation. Journal of Clinical Epidemiology. 69, 245-247 (2016).
  15. Zhang, Z., Chen, L., Xu, P., Hong, Y. Predictive analytics with ensemble modeling in laparoscopic surgery: A technical note. Laparoscopic, Endoscopic and Robotic Surgery. 5 (1), 25-34 (2022).

访问受限。请登录或开始试用以查看此内容。

重印与许可

探索更多文章

Cox 比例风险模型生存分析竞争风险列线图校准曲线C指数计算Bootstrap 验证外部验证决策曲线分析风险回归软件包