Depicting effect of time-varying covariate in proportional hazards regression by comparing survival curves for specific individuals

Question

I have an dataset in which a time-varying covariate has a strong effect on subsequent risk.

I would like to "illustrate" this effect by comparing the estimated survival curve for a typical individual who does not experience the intervening event to one who does at a particular time.

In the documentation for the rms package in R, there is an example that appears to be very similar in spirit.

(Though I refer to R software here [and paste the example code], I do not believe this is an R-specific question. It may be conceptual.)

In the example, the estimated survival curve for an individual who does not experience the intervening event at 5 days is compared to the curve for one who does.

However, though the two curves differ, they do not differ starting at day 5. They differ all the way back to day 0.

Time-dependent covariates are a challenging topic, but I would have expected the two curves to coincide until day 5, at which point the second individual's estimated survival would be reduced to reflect the newly acquired risk factor.

Perhaps my goal is misguided and I have misunderstood the approach.

Any explanations, pointers, references, etc. would be most appreciated.

Cheers,

Rob

Here is the example copied from help(cph). It assumes you have loaded the rms package.

#Fit a time-dependent covariable representing the instantaneous effect
#of an intervening non-fatal event
rm(age)
set.seed(121)
dframe <- data.frame(failure.time=1:10, event=rep(0:1,5),
                     ie.time=c(NA,1.5,2.5,NA,3,4,NA,5,5,5), 
                     age=sample(40:80,10,rep=TRUE))
z <- ie.setup(dframe$failure.time, dframe$event, dframe$ie.time)
S <- z$S
ie.status <- z$ie.status
attach(dframe[z$subs,])    # replicates all variables

f <- cph(S ~ age + ie.status, x=TRUE, y=TRUE)  
#Must use x=TRUE,y=TRUE to get survival curves with time-dep. covariables


#Get estimated survival curve for a 50-year old who has an intervening
#non-fatal event at 5 days
new <- data.frame(S=Surv(c(0,5), c(5,999), c(FALSE,FALSE)), age=rep(50,2),
                  ie.status=c(0,1))
g <- survfit(f, new)
plot(c(0,g$time), c(1,g$surv[,2]), type='s', 
     xlab='Days', ylab='Survival Prob.')
# Not certain about what columns represent in g$surv for survival5
# but appears to be for different ie.status
#or:
#g <- survest(f, new)
#plot(g$time, g$surv, type='s', xlab='Days', ylab='Survival Prob.')


#Compare with estimates when there is no intervening event
new2 <- data.frame(S=Surv(c(0,5), c(5, 999), c(FALSE,FALSE)), age=rep(50,2),
                   ie.status=c(0,0))
g2 <- survfit(f, new2)
lines(c(0,g2$time), c(1,g2$surv[,2]), type='s', lty=2)
#or:
#g2 <- survest(f, new2)
#lines(g2$time, g2$surv, type='s', lty=2)
detach("dframe[z$subs, ]")
options(datadist=NULL)

Answer 1

There are a couple of "gotchas" here.

The primary one, a general problem in handling time-dependent covariates, is how to specify new data for predictions from your model. You need to let the software know which lines of your data frame correspond to which individuals.

That might be surprising at first, because with a single type of event and no more than 1 event per individual, you don't need to specify individual IDs to build a Cox model. As the time-dependent vignette says in that case at the end of Section 2: "The likelihood equations at any time point use only one copy of any subject, the program picks out the correct row of data at each time."

But for predictions from the model on new data, you need to let the software know which time intervals and associated covariate values correspond to each individual. Otherwise it will assume that each row in the data frame corresponds to a different individual. From the survfit.coxph manual page on the id parameter:

Each group of rows in newdata with the same subject id represents the covariate path through time of a single subject, and the result will contain one curve per subject ... If newid is not present, then each individual row of newdata is presumed to represent a distinct subject.

You need to add an "id" column of some type to each of your data frames for predictions, and specify that column in your call to the prediction software.

Now for the software-specific gotchas: software packages differ in their pickiness with time-dependent covariates. With the survfit function in the rms package I got an error from your model: "Error in survfitcoxph.fit(y, X, weights, X2, risk, newrisk, strata, se.fit, : object 'id' not found," because the original model and data didn't include an id. It also gave a warning: "some aspects of individual=TRUE not yet implemented. Try survfit.coxph." But survfit.coxph in turn didn't seem to play well with the your cph object f.

So I repeated your model with coxph instead, and (as I get confused when I use attached data) I generated a full data frame with your original data to get the Cox model:

f <- coxph(S ~ age + ie.status,data=expandedDF)

Although this model also didn't include an id, survfit.coxph didn't complain about that. Provided you let the software know that multiple lines correspond to the same individual:

new <- data.frame(S=Surv(c(0,5), c(5,999), c(FALSE,FALSE)), age=rep(50,2), ie.status=c(0,1), id=c(1,1))
new2 <- data.frame(S=Surv(c(0,5), c(5, 999), c(FALSE,FALSE)), age=rep(50,2), ie.status=c(0,0),id=c(2,2))
plot(survfit(f,new,id=id),col="red")
lines(survfit(f,new2,id=id),col="blue")

you get the result you want.

Another "gotcha" is in the handling of Surv objects. If you try to combine new and new2 into a single data frame with rbind, their original Surv objects get converted to separate columns S.start, S.stop, S.status. So when survfit tries to use the combined data frame, it complains that "Survival type of newdata does not match the fitted model."

As much as I generally favor cph for most work, sometimes I end up building an identical survival model with coxph, because of some things like this that have different behaviors between the rms and survival packages.