1

我希望代码在两者的设置中生成生存曲线

  • 时间相关协变量和
  • 时变系数。

目的是展示计费方法如何影响人寿保险单失效。这很复杂

  1. 客户的计费方式（发票或电子转帐）随时间而变化，
  2. 计费方式对失效的影响会随着时间的推移而逐渐消失，并且
  3. 计费方法对失效的影响取决于其他协变量。

在阅读了有关时间相关协变量的小插曲后，我不知道如何从具有时间相关协变量和时变系数的模型生成生存曲线。

library(survival)

Samp <- data.frame(
  id = c(143,151,680,134),
  time = c(17,16,17,18) ,
  censor= rep(1,4) , 
  covariate = seq(5,20,length.out = 4))
# Lookup provides the values of a tdc
Lookup <- data.frame(
  id =c(rep(134,2),680,143,rep(151,3)) ,
  billing.mode = c("INV",rep("EFT",2),rep("INV",2),"EFT","INV") ,
  switch.time = c(0,3,rep(0,3),2,7))

# create the tdc 
Samp.tdc <- tmerge(data1=Samp,data2=Samp,id=id,
                    lapse=event(time,censor))
Samp.tdc <- tmerge(data1=Samp.tdc,data2=Lookup,id=id,
                    billing.mode=tdc(switch.time,billing.mode))
Samp.tdc$inv = as.numeric(Samp.tdc$billing.mode == "INV")

# the call looks something like this
fit <-coxph(Surv(tstart, tstop, lapse) ~ inv + tt(inv) +
  covariate*inv, data = Samp.tdc, 
            tt = function(x, t, ...) x * t)

当我说我想生成生存曲线时，我的意思是在一组固定的时间和协变量值下预测生存。下面就来说说吧LpsData。 

LpsData <- data.frame(
  tstart = rep(c(0,16,17),times=4),
  tstop = rep(16:18,times=4) ,
  lapse = 0 ,
  covariate = rep(c(10,20),each=3,times=2) ,
  inv = rep(c(0,1),each=6) ,
  curve=rep(c('eft','inv'), each=6)
)
4

1 回答 1

1

这是一个相对复杂的问题，我个人发现survival包的容量在这方面是有限的。例如，您必须预先指定时变的函数形式。另一种方法是使用分段指数加法模型（PAMM），它可以通过估计 mgcv::gam并继承它们的所有灵活性（+ 非线性效应的惩罚估计，包括时变效应）。

一般来说，您必须决定要适合的模型类型。让z成为您的时间相关协变量。比潜在的模型可能

  • 线性协变量效应，线性时变，即你的代码中指定的模型（mgcv公式+ z * t +：）
  • 非线性协变量效应，线性时变（公式：+ s(z, by = t) +）
  • 线性协变量效应，非线性时变（公式：+ s(t, by = z) +）
  • 非线性，非线性时变（公式：+ te(t, z) +）

下面是一个使用包中pbc数据的示例，该数据survival也出现在时间相关协变量的生存小插图中（另请参阅https://adibender.github.io/pammtools/articles/tdcovar.html以与 PAMM 进行比较）：

library(survival)
library(ggplot2)
theme_set(theme_bw())
library(pammtools)
library(mgcv)

数据转换

首先，我将数据转换为分段指数数据 (PED) 格式：

pbc <- pbc %>% filter(id <= 312) %>%
  select(id:sex, bili, protime) %>%
  mutate(status = 1L * (status == 2))

## Transform to piece-wise exponential data (PED) format
pbc_ped <- as_ped(
  data = list(pbc, pbcseq),
  formula = Surv(time, status)~. | concurrent(bili, protime, tz_var = "day"),
  id = "id") %>% ungroup()

pbc_ped <- pbc_ped %>%
  mutate(
    log_bili = log(bili),
    log_protime = log(protime))

拟合分段指数加法模型 (PAM)

在这里，我拟合了一个具有 2 个具有线性协变量效应的时间相关协变量的模型，非线性时变（尽管由于惩罚，估计几乎是线性的）

pbc_pam <- gam(ped_status ~ s(tend, k = 10) + s(tend, by = log_bili) +
  s(tend, by = log_protime),
  data = pbc_ped, family = poisson(), offset = offset)

固定协变量的生存预测

对于预测我

  • 在所有唯一观察到的时间点创建一个新数据集（所有未指定的协变量将设置为均值/模值）
  • log_bili在每个时间点添加一个与时间相关的值
  • 添加生存概率预测 + CI 使用add_surv_prob
ndf <- make_newdata(pbc_ped, tend = unique(tend)) %>%
  mutate(log_bili = runif(n(), min(log_bili), max(log_bili))) %>%
  add_surv_prob(pbc_pam)

绘制预测的生存概率

ggplot(ndf, aes(x = tend, y = surv_prob)) +
  geom_surv() +
  geom_ribbon(aes(ymin = surv_lower, ymax = surv_upper), alpha = 0.3) +
  ylim(c(0, 1))

```

由reprex 包（v0.2.1）于 2018 年 12 月 8 日创建

于 2018-12-08T19:43:59.503 回答