Time varying coefficient in Cox model

Question

I have a model for survival after an injury that is borderline passing the Schoenfeld test for the proportional hazards assumption (cox.zph() in R).

However, suspecting that there would be an increased mortality within the first month I fitted an Aalen additive regression model also hinting at a varying coefficient since the cumulative coefficient drops steeply within the first 4 days, then gradually continues to decline before it stabilizes at a nearly constant level after a year for my variable of interest.

Now - how can I accommodate such a model in R?

• Introduce an interaction term [for my variable of interest] varying with time?
• Create different variables for time 0-4, 4-365 and 365 and then also expand the dataset with three rows per patient?

Answer 1

Just answered my own question with the same problem. Basically you need to do a time-split as you describe and then add an interaction term for the time. As you suggested it sounds like a good idea to perhaps split the time more fine grained in the beginning and then increase the intervals. Using my Greg package you can set the by-time to either a single interval or to a vector, e.g. if we create a dataset with four subjects:

test_data <- data.frame(
  id = 1:4,
  time = c(4, 3.5, 1, 5),
  event = c("censored", "dead", "alive", "dead"),
  age = c(62.2, 55.3, 73.7, 46.3),
  date = as.Date(
    c("2003-01-01", 
      "2010-04-01", 
      "2013-09-20",
      "2002-02-23")))

Looking like this:

| id| time|event    |  age|date       |
|--:|----:|:--------|----:|:----------|
|  1|  4.0|censored | 62.2|2003-01-01 |
|  2|  3.5|dead     | 55.3|2010-04-01 |
|  3|  1.0|alive    | 73.7|2013-09-20 |
|  4|  5.0|dead     | 46.3|2002-02-23 |

We split each subject into several through:

library(Greg)
library(dplyr)
split_data <- 
  test_data %>% 
  select(id, event, time, age, date) %>% 
  timeSplitter(by = c(.1, .5, 2), # The time that we want to split by
               event_var = "event",
               time_var = "time",
               event_start_status = "alive",
               time_related_vars = c("age", "date"))

knitr::kable(split_data)

Gives:

| id|event    |  age|     date| Start_time| Stop_time|
|--:|:--------|----:|--------:|----------:|---------:|
|  1|alive    | 62.2| 2002.999|        0.0|       0.1|
|  1|alive    | 62.3| 2003.099|        0.1|       0.5|
|  1|alive    | 62.7| 2003.499|        0.5|       2.0|
|  1|censored | 64.2| 2004.999|        2.0|       4.0|
|  2|alive    | 55.3| 2010.246|        0.0|       0.1|
|  2|alive    | 55.4| 2010.346|        0.1|       0.5|
|  2|alive    | 55.8| 2010.746|        0.5|       2.0|
|  2|dead     | 57.3| 2012.246|        2.0|       3.5|
|  3|alive    | 73.7| 2013.718|        0.0|       0.1|
|  3|alive    | 73.8| 2013.818|        0.1|       0.5|
|  3|alive    | 74.2| 2014.218|        0.5|       1.0|
|  4|alive    | 46.3| 2002.145|        0.0|       0.1|
|  4|alive    | 46.4| 2002.245|        0.1|       0.5|
|  4|alive    | 46.8| 2002.645|        0.5|       2.0|
|  4|dead     | 48.3| 2004.145|        2.0|       5.0|

As described in my previous answer you now just need to model using Surv(Start_time, End_time, event) together with an additional : interaction term for your variable (note that you should also have the variable without the interaction in the model).