在 r 中使用 predict.coxph 预测基线累积危害

Predicting baseline cumulative hazard using predict.coxph in r

我的目标是预测(从下面的拟合模型预测新观察的累积风险)从拟合模型的时间尺度 0 到开始时间的累积风险值。

我已经使用2次(不等于0的开始时间和结束时间)拟合了cox模型。这样我就可以找到结束时间的累积危险(即从 0 到结束时间的累积危险,我已经使用相同的拟合模型计算过)和开始时间的累积危险(即从 0 到累积危险)之间的差异结束时间,我想在这里计算)这将最终给出每次观察的开始和结束时间之间的 cum haz。

因此,为了获得预期的事件数,我使用了 predict(coxph(), newdata, type= "expected") .

我用过的数据如下:

N <- 10^4 # population
H <- within(data.frame(start_time=runif(N, 0, 50), x1=rnorm(N, 2, 1), x2=rnorm(N, -2, 1)), {
  lp <-   0.05*x1 + 0.2*x2 
  Tm <- qweibull(runif(N,pweibull(start_time,shape = 7.5, scale = 84*exp(-lp/7.5)),1), shape=7.5, scale=84*exp(-lp/7.5))
  Cens1 <- 100
  event_time <- pmin(Tm,Cens1)
  status <- as.numeric(event_time == Tm)})  

预测代码为:

H$X <- rep(1,nrow(H))
D = coxph(Surv(start_time, event_time, status) ~ X, data =  H, x = TRUE )
pred2 <- predict(D, newdata = data.frame(start_time = rep(0,nrow(H)),event_time = H$start_time, status = rep(0,nrow(H)), X = rep(1, nrow(H))), type = "expected")

但是 pred2 只会产生“NA”值。谁能指出我的想法或代码是否有错误

如果需要进一步说明,请告诉我。

有两个问题。首先,您 运行 遇到了一个问题,因为当您指定 ~1 时,这意味着拟合一个没有系数的仅截距模型。所以你所有的预测都会给你一个值?

library(survival)
D <- coxph(Surv(H$start_time, H$event_time, H$status) ~ 1, x = TRUE )
names(D)
 [1] "loglik"            "linear.predictors" "method"           
 [4] "residuals"         "n"                 "nevent"           
 [7] "terms"             "assign"            "concordance"      
[10] "x"                 "y"                 "timefix"          
[13] "formula"           "call"  

table(predict(D))

    0 
10000

我认为这没有多大意义,因此您 运行 犯了所有错误。因此,您需要使用回归中使用的自变量进行预测,例如:

D <- coxph(Surv(start_time,event_time,status) ~ x1+x2, data=H )
pred2 <- predict(D, newdata = data.frame(t_0 = rep(0,nrow(H)),time = H$start_time, event_M = rep(0,nrow(H)), X = rep(1, nrow(H))), type = "expected")

predict(D,newdata=data.frame(x1=runif(10,0,1),x2=runif(10,-1,1)))
        1         2         3         4         5         6         7         8 
0.3033206 0.4213120 0.3952827 0.3879701 0.4798670 0.2170032 0.3385253 0.4141698 
        9        10 
0.3690579 0.4128084 

当你用所有 X=1 拟合一个模型时,这会给你所有的 NA,因为已经有一个截距,这使得这个变量变得多余。您可以查看:

H$X = 1
D <- coxph(Surv(start_time, event_time, status) ~ X,data=H)

Call:
coxph(formula = Surv(start_time, event_time, status) ~ X, data = H)

  coef exp(coef) se(coef)  z  p
X   NA        NA        0 NA NA

它仅在 X 是拟合数据中的实际变量时才有效,因此我使用一个包含 2 个协变量的示例:

H$X = runif(nrow(H))
D <- coxph(Surv(start_time, event_time, status) ~ X + x1,data=H)

您可以通过将 X 固定为 1 并改变 x1 来进行预测:

predict(D,newdata=data.frame(X=1,x1=c(0.1,0.2,0.3)))
         1          2          3 
-0.1132548 -0.1084592 -0.1036637 

或 X 在 2:

predict(D,newdata=data.frame(X=2,x1=c(0.1,0.2,0.3)))
                 1          2          3 
-0.1579480 -0.1531524 -0.1483568

我自己找到了答案,这只是一个小技巧,但我不确定是否总能奏效。 如果我在 predict() 函数之前使用以下行:

D$coefficients["X"] <- 0

但是,我得到了使用不接受开始时间(或一次两个变量)的 nelsonaalen() 函数检查的正确值

如果有任何其他正确的解决方法,请告诉我。