更新:很抱歉再次更新,但我找到了一些可能的分数多项式解决方案和竞争风险包,我需要一些帮助。
问题
我找不到在 R 中进行时间相关系数分析的简单方法。我希望能够将变量系数转换为时间相关系数(非变量),然后绘制随时间变化的变化:
可能的解决方案
1)拆分数据集
我看过这个例子(Se part 2 of the lab session)但是创建一个单独的数据集似乎很复杂,计算成本很高而且不是很直观......
2) 降阶模型 - coxvc 包
coxvc包提供了一种处理问题的优雅方法——这里有一个手册。问题是作者不再开发包(最新版本是自 2007 年 5 月 23 日以来),经过一些电子邮件对话后,我已经让包工作了,但在我的数据集上运行一次需要 5 小时(140 000条目)并在期末给出极端估计。你可以在这里找到一个稍微更新的包——我大部分只是更新了绘图功能。
这可能只是一个调整的问题,但由于该软件不容易提供置信区间并且该过程非常耗时,我现在正在寻找其他解决方案。
3)timereg包
令人印象深刻的timereg 包也解决了这个问题,但我不确定如何使用它,它并没有给我一个流畅的情节。
4)分数多项式时间(FPT)模型
我发现 Anika Buchholz 的优秀论文“评估时间变化的治疗和预后因素的长期影响”在涵盖不同模型方面做得非常出色。她得出结论,Sauerbrei 等人提出的 FPT似乎最适合时间相关系数:
FPT 非常擅长检测时变效应,而 Reduced Rank 方法导致模型过于复杂,因为它不包括时变效应的选择。
这项研究似乎非常完整,但对我来说有点遥不可及。我也有点想知道,因为她恰好和 Sauerbrei 一起工作。不过这似乎很合理,我想可以使用mfp 包进行分析,但我不确定如何。
5) cmprsk 包
我一直在考虑进行竞争风险分析,但计算非常耗时,所以我改用常规的 cox 回归。crr具有时间相关协变量的选项:
....
cov2 matrix of covariates that will be multiplied
by functions of time; if used, often these
covariates would also appear in cov1 to give
a prop hazards effect plus a time interaction
....
有二次示例,但我不太了解时间实际出现的位置,我不确定如何显示它。我还查看了 test.R 文件,但那里的示例基本相同......
我的示例代码
这是我用来测试不同可能性的示例
library("survival")
library("timereg")
data(sTRACE)
# Basic cox regression
surv <- with(sTRACE, Surv(time/365,status==9))
fit1 <- coxph(surv~age+sex+diabetes+chf+vf, data=sTRACE)
check <- cox.zph(fit1)
print(check)
plot(check, resid=F)
# vf seems to be the most time varying
######################################
# Do the analysis with the code from #
# the example that I've found #
######################################
# Split the dataset according to the splitSurv() from prof. Wesley O. Johnson
# http://anson.ucdavis.edu/~johnson/st222/lab8/splitSurv.ssc
new_split_dataset = splitSuv(sTRACE$time/365, sTRACE$status==9, sTRACE[, grep("(age|sex|diabetes|chf|vf)", names(sTRACE))])
surv2 <- with(new_split_dataset, Surv(start, stop, event))
fit2 <- coxph(surv2~age+sex+diabetes+chf+I(pspline(stop)*vf), data=new_split_dataset)
print(fit2)
######################################
# Do the analysis by just straifying #
######################################
fit3 <- coxph(surv~age+sex+diabetes+chf+strata(vf), data=sTRACE)
print(fit3)
# High computational cost!
# The price for 259 events
sum((sTRACE$status==9)*1)
# ~240 times larger dataset!
NROW(new_split_dataset)/NROW(sTRACE)
########################################
# Do the analysis with the coxvc and #
# the timecox from the timereg library #
########################################
Ft_1 <- cbind(rep(1,nrow(sTRACE)),bs(sTRACE$time/365,df=3))
fit_coxvc1 <- coxvc(surv~vf+sex, Ft_1, rank=2, data=sTRACE)
fit_coxvc2 <- coxvc(surv~vf+sex, Ft_1, rank=1, data=sTRACE)
Ft_3 <- cbind(rep(1,nrow(sTRACE)),bs(sTRACE$time/365,df=5))
fit_coxvc3 <- coxvc(surv~vf+sex, Ft_3, rank=2, data=sTRACE)
layout(matrix(1:3, ncol=1))
my_plotcoxvc <- function(fit, fun="effects"){
plotcoxvc(fit,fun=fun,xlab='time in years', ylim=c(-1,1), legend_x=.010)
abline(0,0, lty=2, col=rgb(.5,.5,.5,.5))
title(paste("B-spline =", NCOL(fit$Ftime)-1, "df and rank =", fit$rank))
}
my_plotcoxvc(fit_coxvc1)
my_plotcoxvc(fit_coxvc2)
my_plotcoxvc(fit_coxvc3)
# Next group
my_plotcoxvc(fit_coxvc1)
fit_timecox1<-timecox(surv~sex + vf, data=sTRACE)
plot(fit_timecox1, xlab="time in years", specific.comps=c(2,3))
代码生成以下图表: coxvc 的不同设置以及coxvc 和 timecox图的比较。我想结果还可以,但我认为我无法解释 timecox 图 - 它似乎很复杂......
我的(当前)问题
- 如何在 R 中进行 FPT 分析?
- 如何在 cmprsk 中使用时间协变量?
- 如何绘制结果(最好使用置信区间)?