Your model is additive so the effect of rx is the same at every age. There is not one survival curve for all ages (unless the beta for age is 0). The curves will shift up and down as you vary age, but they will retain the same relation. A common approach is to use the sample mean of age. Alternatively you can compute the survival curves (or the cumulative hazard curves) for many ages and compute a weighted average of them (weights depending on some reference population perhaps). Finally, if you leave age out of the model you don't even need to specify age to draw the survival curves. It kind of depends on what you are trying to do.

Chris

Ps last time I posted code it all ran together. We'll see what happens this time.

library(survival)

data(kidney)

mod1 <- coxph(Surv(time, status) ~ sex + age, data=kidney)

mod2 <- coxph(Surv(time, status) ~ sex, data=kidney)

meanage <- mean(kidney$age)

plot(survfit(mod1, newdata=data.frame(sex =c(1:2), age=meanage)), xscale=365.25, xlab = "Years", ylab="Survival", col=1:2)

plot(survfit(mod2, newdata=data.frame(sex =c(1:2))), xscale=365.25, xlab = "Years", ylab="Survival", col=1:2)

ttt <- hist(kidney$age)

weights <- ttt$counts/ sum(ttt$counts)

survcurves <- survfit(mod, newdata=expand.grid(sex =c(1:2), age=ttt$mids), se.fit=FALSE)

str(survcurves)

sc2 <- survcurves

sc2$surv <- t(apply(survcurves$surv, 1, function(x) c(mean(x[seq(1,length(x),2)]), mean(x[seq(2,length(x),2)]))))

plot(sc2, xscale=365.25, xlab = "Years", ylab="Survival", col=1:2)

-----Original Message-----

From: Maziar Mohaddes [mailto:

[hidden email]]

Sent: Tuesday, November 27, 2012 12:45 PM

To:

[hidden email]
Subject: [R] Fitting and plotting a coxph with survfit, package(surv)

Hi Dear R-users

I have a database with 18000 observations and 20 variables. I am running cox regression on five variables and trying to use survfit to plot the survival based on a specific variable without success.

Lets say I have the following coxph:

>library(survival)

>fit <- coxph(Surv(futime, fustat) ~ age + rx, data = ovarian) fit

what I am trying to do is plot a survival comparing objects based on rx.

Using this

>plot(survfit(fit, newdata=data.frame(rx =c(1:2), age=c(60)),

xscale=365.25, xlab = "Years", ylab="Survival")) I get the survival for patients at 60, but is there an option to get a survfit for the patients regardless of the value in variable "age"?

Thanks in advance

Maziar Mohaddes

M.D.

Gothenburg, Sweden

[[alternative HTML version deleted]]

**********************************************************

Electronic Mail is not secure, may not be read every day, and should not be used for urgent or sensitive issues

______________________________________________

[hidden email] mailing list

https://stat.ethz.ch/mailman/listinfo/r-helpPLEASE do read the posting guide

http://www.R-project.org/posting-guide.htmland provide commented, minimal, self-contained, reproducible code.