Assessing time dependent categorical variable in cox ph survival model

Question

I am trying to calculate the hazard rate of two treatments, which have more risk early in time compared to later. I've provided an example dataset and analysis below to highlight my situation.

I've been using the time-dependent vignette to learn from (https://cran.r-project.org/web/packages/survival/vignettes/timedep.pdf). But the examples are often for a continuous covariate and I have a categorical. The tt() solutions don't seem to work for a category so I've tried the step fn solution via an interaction between treatment and a time strata. But I get NA's in the results and I can't test assumptions or predict() results due to dimension issues. Looking for some help on syntax to make this work. Thanks!

library(tidyverse)
library(survival)
library(survminer)


##make data
set.seed(1234)
n <- 200
time=round(runif(n,1,365),0)
status=rnorm(n,0.5,0.1)%>%round(0)
df.interim <- data_frame(time, status, treatment=rep(c("w","p"), times=n/2))

##add more morality early on in "w" group
df <- df.interim%>%
  rbind(data_frame(time=round(runif(50,1,50),0), 
                   status=1, 
                   treatment="w"))


##look at curves
survfit2(Surv(time, status)~ treatment, data = df) %>% 
  ggsurvfit(linewidth = 1) +
  add_confidence_interval() +
  add_risktable() +
  add_quantile(y_value = 0.6, color = "gray50", linewidth = 0.75) +
  scale_ggsurvfit()

##model
m1 <- coxph(Surv(time, status) ~ treatment , data=df)
m1

             coef exp(coef) se(coef)     z        p
treatmentw 0.7215    2.0575   0.1662 4.341 1.42e-05

plot(cox.zph(m1, global=TRUE))
cox.zph(m1, global=TRUE)

          chisq df       p
treatment  18.8  1 1.4e-05
GLOBAL     18.8  1 1.4e-05



df.split <- survSplit(Surv(time, status) ~ ., data=df, cut=c(50), episode= "tgroup", zero=0)
m2 <- coxph(Surv(time, status) ~ treatment:strata(tgroup) , data=df.split)
m2

                       coef exp(coef) se(coef)      z        p
treatmentp:strata(tgroup)tgroup=1 -1.8873    0.1515   0.4004 -4.713 2.44e-06
treatmentw:strata(tgroup)tgroup=1      NA        NA   0.0000     NA       NA
treatmentp:strata(tgroup)tgroup=2 -0.2594    0.7715   0.2020 -1.284    0.199
treatmentw:strata(tgroup)tgroup=2      NA        NA   0.0000     NA       NA

plot(cox.zph(m2, global=TRUE))
cox.zph(m2, global=TRUE)

Error in dimnames(tbl) <- list(c(termname, "GLOBAL"), c("chisq", "df",  : 
  length of 'dimnames' [1] not equal to array extent

Assessing time dependent categorical variable in cox ph survival model

Answers (1)

Related Questions