[R] Plotting an adjusted survival curve

[Brent Caldwell]

Dear Terry and David
Thank you so much for your rapid and helpful replies!
I've changed tack slightly and decided to focus on the effect modifiers which had the greatest effect on the outcome of the trial (HowmanyZonnic_Q and WISDMPDM)
I've tried to follow your example on page 16 of the manuscript you emailed me, but I get an error which I think is related to trying to do maths on two matrices which are of different dimensions I think???
This is my code:
KM.object.coxph5 <- coxph(Surv(whenfailed,failed) ~ strata(Treatment)+HowmanyZonnic_Q+WISDMPDM,survmat)
asdist <- table(HowmanyZonnic_Q,WISDMPDM)
uage <- as.numeric(dimnames(asdist)[[1]])
tdata <- expand.grid(HowmanyZonnic_Q=uage, WISDMPDM=uage) #This appears to be a table of all the possible combinations of HowmanyZonnic_Q and WISDMPDM
par(mfrow=c(1,2))
sfit4a <- survfit(KM.object.coxph5, newdata=tdata, se.fit=FALSE)
sfit4a$surv <- sfit4a$surv %*% c(asdist/sum(asdist)) #Error in sfit4a$surve %*% c(asdist/sum(asdist)) :   non-conformable arguments
#sfit4a$surv is a table of 10 rows and 25 columns
#c(asdist/sum(asdist)) is a vector with a length of 420
I hope I am not sapping your patience with this question
Thank you - you guys are awesome
Best wishes
Brent

-----Original Message-----
From: Terry Therneau [mailto:therneau at mayo.edu] 
Sent: Tuesday, 27 November 2012 4:03 a.m.
To: r-help at r-project.org; Brent Caldwell
Cc: David Winsemius
Subject: Re: Plotting an adjusted survival curve

First a statistical issue: The survfit routine will produce predicted survival curves for any requested combination of the covariates in the original model.  This is not the same thing as an "adjusted" survival curve.  Confusion on this is prevalent, however.  True adjustment requires a population average over the confounding factors and is closely related to the standardized incidence ratio concept found in epidemiology.

To answer your technical question:
    fit <- coxph(Surv(.........
    mysurv <- survfit(fit, newdata= mydata) This will give a set of predicted curves, one for each observation in mydata.  If we assume 2 treatments and 4 ethnicities, this means that there are 8 possible predicted curves.  You can certainly take the curves for "trt=1, white" and "trt=2, white", plot them together on one graph, and call this your "adjusted" survival curves; the mydata data set would have two observations.  This is not a correct label but is certainly common.

Terry Therneau

On 11/26/2012 05:00 AM, r-help-request at r-project.org wrote:
> Dear R-users
> I am trying to make an adjusted Kaplan-Meier curve (using the Survival 
> package) but I am having difficulty with plotting it so that the plot only shows the curves for the adjusted results.
> My data come from a randomised controlled trial, and I would like the 
> adjusted Kaplan-Meier curve to only show two curves for the adjusted 
> survival: one for those on treatment (Treatment==1) and another curve for those on placebo (Treatment==0).
>
> My problem is that when I plot the survfit of my coxph, I think it 
> displays a curve for every single individual factor in my coxph, 
> whereas I would like it to only display the adjusted curves for when Treatment==1 and Treatment==0.  How can I do this?
>
> A simplified example of my code with only one effect-modifier is:
>
> simple.cox.ethnicity<- coxph(Surv(whenfailed,failed) ~ 
> factor(Treatment) + factor(ethnicity)) #I've my data are attached 
> already
> survfit.simple.cox.ethnicity<- survfit(simple.cox.ethnicity,survmat) 
> #survmat is a data.frame that contains Treatment and ethnicity 
> plot(survfit.simple.cox.ethnicity, col=c("red","black"), 
> main="survfit.simple.cox", xlab="survival time", ylab="propotion 
> surviving")
>
> Thank you so much for your help.
> Yours gratefully,
> Brent Caldwell