当前位置: 首页 > 工具软件 > IndexR > 使用案例 >

cox计算C-index及两模型进行C-index的比较,survcomp的安装及使用

谢奕
2023-12-01

计算COX模型的C-index

两个COX模型的比较

PART1 计算COX的C-index

方法1:基于cox直接生成系数以及计算的C-index值

library("survminer")
library("survival")
######### 生成COX模型
cox_model <- coxph()
# 直接查看coxph
summary(cox_model)
#  call之后: Concordance= ***  (se = *** ) 即是C-index
#  查看cindex - concordance
# 可以看出这种方法输出了C指数,也输出了标准误,那么95%可信区间就可以通过加减1.96*se得到。

# 生成cox回归模型的系数以及95%置信区间,查看AIC
library(tableone)
ShowRegTable(data)
ShowRegTable(data, digits = 3)
AIC(data)
summary(data) #可直接查看C-index以及se

tip:通过relevel或者level函数中,设置label的顺序,可以选择参照组,即排名第一个即为参照组。

方法2:基于rms包的rcorr.cens()直接生成C-index,
但是做下来好像C-index偏低?不知道为什么

library(rms)
CstatisticCI <- function(x) {
  se <- x["S.D."]/sqrt(x["n"])
  Low95 <- x["C Index"] - 1.96*se 
  Upper95 <- x["C Index"] + 1.96*se 
  cbind(x["C Index"], Low95, Upper95) 
}
cindex <- rcorr.cens(data$prediction,data$event)
cindex
print(CstatisticCI(cindex))
# 该方法可以获得Dxy,但是这个跟cox直接出来的Cindex相比要少,是因为Dxy的原因吗?
# COX的计算
# 根据survival包,可以计算Dxy
# 把-值去掉, |Dxy|/2+0.5

COX的计算
根据survival包,可以计算Dxy
把-值去掉, |Dxy|/2+0.5

方法3:基于Hmisc的rcorrp.cens语法进行

  library(Hmisc)
  library(survival)
  # rcorrp.cens可以获得Dxy以及Cindex值
  rcorrp.cens(data$model1_prediciton, 
              data$model2_prediction, 
              Surv(data$time, data$event))
  # call的结果中
  # C-index1 = 1- C X1
  # C-index2 = 1- C X2
  # Dxy1 和 Dxy2 分别是model1和model2的Dxy,可用于计算Cindex,与上述一致
  # 但是这里的P值没有学会是啥意思
  

方法4、基于survcomp包进行concordance.index
详见最后,这里先讨论一下安装的问题
R升级之后,该包就被移除了CRAN的库,需要通过Bioconductor进行安装
可运行如下代码

To install this package, start R (version "4.0") and enter:

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

# For older versions of R, please refer to the appropriate Bioconductor release. 

具体可参见如下网址:
https://www.bioconductor.org/packages/release/bioc/html/survcomp.html

各个C-index进行比较,除了在TNM分期的时候,survcomp所计算的C-index较高之外,其余均无明显差别。
rcorr.cens() 与 cox自带的cindex没问题,但se或者S.D是否均可用于计算置信区间仍需要讨论。
----------------------------------╮(╯▽╰)╭---------------------------------------

PART2:模型比较方法:

方法1
似然比检验

  # 方法1
  # model1 vs model2 (直接使用model名称即可)
  anova(model1,model2) 
  # 方法2
  library(rms)
all.X <-data.frame(x.T=data.T, x.N=data.N, x.S=data.S, x.G=data.G, x.V=data.V, x.P=data.P, x.CEA2=data.CEA2, x.CA1992=data.CA1992)
TN.model <- cph(Surv(survival.time,survival.status)~ x.T+x.N, 
                data=all.X, na.action=na.omit )
TNC.model <- cph(Surv(survival.time,survival.status)~ x.T+x.N+x.CEA2, 
                 data=all.X, na.action=na.omit )
TN2TNC <- lrtest(TN.model, TNC.model)
  # 本质上就是对LIKEHOOD进行卡方检验
  # 经过检验发现,这两者的结果是一致的,均可采用。

方法2
基于survcomp包进行concordance.index
survcomp如何进行p值的提取呢?暂时没明白。

library("survival")
library("prodlim")
library("survcomp")
C_index1 <- concordance.index(x=data$model1_prediction, surv.time=data$time, surv.event=data$event,method="noether")
# 往上翻可以直接查看到C_index以及置信区间
C_index1
C_index2 <- concordance.index(x=data$model2_prediction, surv.time=data$time, surv.event=data$event, method="noether")
C_index2
cindex.comp(C_index1, C_index2)

该方法与cox自带的C-index的稍高一些,置信区间更小一些。
但是在做TNM分期的时候,不知道为甚C-index突然就高了。

方法3
compareC包

方法4
nricens包
PredictABEL包

方法5 - 生存数据 cox
survIDINRI包

表1.R中可计算IDI的包

包的名称下载地址分类资料结局生存资料结局
PredictABELCRANreclassification函数不支持
survIDINRICRAN不支持IDI.INF函数
 类似资料: