Skip to content

Feature Request: Support time dependent covariates with counting process, please #349

@jkylearmstrong

Description

@jkylearmstrong
library('tidyverse')
library('survival')  
  
# example from https://cran.r-project.org/web/packages/survival/vignettes/timedep.pdf

dim(cgd0)
#> [1] 128  20

newcgd <- tmerge(data1=cgd0[, 1:13], data2=cgd0, id=id, tstop=futime)
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime1))
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime2))
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime3))
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime4))
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime5))
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime6))
newcgd <- tmerge(newcgd, cgd0, id=id, infect = event(etime7))
newcgd <- tmerge(newcgd, newcgd, id, enum=cumtdc(tstart))

dim(newcgd)
#> [1] 203  17

newcgd[1:5,c(1, 4:6, 13:17)]
#>   id treat sex age futime tstart tstop infect enum
#> 1  1     1   2  12    414      0   219      1    1
#> 2  1     1   2  12    414    219   373      1    2
#> 3  1     1   2  12    414    373   414      0    3
#> 4  2     0   1  15    439      0     8      1    1
#> 5  2     0   1  15    439      8    26      1    2

summary(newcgd)
#> Call:
#> tmerge(data1 = newcgd, data2 = newcgd, id = id, enum = cumtdc(tstart))
#> 
#>        early late gap within boundary leading trailing tied missid
#> infect     0    0   0     44        0       0        0    0      0
#> infect     0    0   0     16        0       0        1    0      0
#> infect     0    0   0      8        0       0        0    0      0
#> infect     0    0   0      3        0       0        0    0      0
#> infect     0    0   0      2        0       0        0    0      0
#> infect     0    0   0      1        0       0        0    0      0
#> infect     0    0   0      1        0       0        0    0      0
#> enum       0    0   0      0       75     128        0    0      0

newcgd <- newcgd %>%
  mutate(outcome = Surv(tstart, tstop, infect)) 

set.seed(1234)

library('tidymodels')
tidymodels_prefer()
library('censored')

# split the data
ini_split <- group_initial_split(newcgd, group = id)

train <- training(ini_split)
test <- testing(ini_split)

# we can train cox model
coxph_model <- coxph(Surv(tstart, tstop, infect) ~ treat + inherit + steroids,
               data =train, 
               cluster = id)

coxph_model %>%
  broom::tidy()
#> # A tibble: 3 × 6
#>   term     estimate std.error robust.se statistic p.value
#>   <chr>       <dbl>     <dbl>     <dbl>     <dbl>   <dbl>
#> 1 treat     -1.17       0.309     0.361    -3.23  0.00123
#> 2 inherit   -0.0762     0.265     0.355    -0.215 0.830  
#> 3 steroids  -0.902      0.528     0.386    -2.34  0.0195

# we can use glmnet

library('glmnet')
#> Loading required package: Matrix
#> 
#> Attaching package: 'Matrix'
#> The following objects are masked from 'package:tidyr':
#> 
#>     expand, pack, unpack
#> Loaded glmnet 4.1-9

X <- train %>%
  select(treat, inherit, steroids) %>%
  as.matrix()

Y_M <- train %>%
  select(outcome) %>%
  as.matrix()

Y <- Surv(Y_M[,'outcome.start'], Y_M[,'outcome.stop'], Y_M[,'outcome.status'])   

glmnet_model <- cv.glmnet(x=X, y=Y, family = 'cox', nfolds = 5)
#> Warning: cox.fit: algorithm did not converge

plot(glmnet_model)

Image

glmnet_model <- glmnet(x=X, y=Y, family = 'cox', lambda = 10^(-3))
#> Warning: cox.fit: algorithm did not converge

glmnet_model %>%
  broom::tidy()
#> # A tibble: 3 × 5
#>   term      step estimate lambda dev.ratio
#>   <chr>    <dbl>    <dbl>  <dbl>     <dbl>
#> 1 treat        1 -1.18     0.001    0.0377
#> 2 inherit      1  0.0395   0.001    0.0377
#> 3 steroids     1 -0.00850  0.001    0.0377

# now let's try tidymodel's approach:

library('workflowsets')

folds <- group_vfold_cv(train, v = 5, group = id)

rec <- recipe(outcome ~ treat + inherit + steroids, data = train)

cox_spec <- proportional_hazards(
  mode = "censored regression",
  engine = "glmnet",
  penalty = tune(),
  mixture = tune()
)

wflowset <- workflow_set(
  preproc = list(recipe = rec), 
  models = list(cox = cox_spec)
)

grid_ctrl <-
  control_grid(
    save_pred = TRUE,
    parallel_over = "everything",
    save_workflow = TRUE
  )

# does not work, expects Surv(time, status) type of outcome:

grid_results <-
  wflowset %>%
  workflow_map(
    seed = 3110,
    resamples = folds,
    grid = 35,
    control = grid_ctrl,
    eval_time = c(15,30,60,90,120),
    metrics = metric_set(brier_survival, roc_auc_survival, brier_survival_integrated)
  )
#> → A | warning: cox.fit: algorithm did not converge
#> There were issues with some computations   A: x1                                                 → B | error:   For this usage, the allowed censoring type is "right".
#> Warning: All models failed. Run `show_notes(.Last.tune.result)` for more
#> information.

Created on 2025-06-05 with reprex v2.1.1

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions