# Plot of probability of 1-year survival event vs. variable value? Would be derived from a Cox model

I am hoping to plot the probability of a patient surviving to 1 year without progression (status==0) vs. the value of a continuous variable (radiation.dose). This would be based on a univariate cox proportional hazards model. I've Googled this and can't find anyone who has done this before. What I imagine is similar axes to this logistic regression plot: Logistic probability vs. variable

Here is some sample code for you to work with:

library("survival")
require("survival")

days <- rpois(100, 365)
status <- rbinom(100,1,0.34)

mod<-Surv(df$$days,df$$status)


Thank you for your help!

I'm just thinking out loud here, I'm posting this as an answer instead of a comment for better formatting.

According to the help file on predict.coxph, "The survival probability for a subject is equal to exp(-expected)." So we can predict the survival probability

cox.pred<- exp(-predict(cox.mod, type = "expected"))


• Thanks for thinking through this with me! I think this looks like the right function. I'm far from a statistics expert, but I would have thought that since the predicted probabilities are based on an underlying mathematical model there wouldn't be so much variance in the probabilities. I guess I was expecting what looks like a curve sort of like that logistic regression plot I linked to. Does that make sense? What are your thoughts?
– JJM
Commented Dec 18, 2019 at 19:02
• The way the data set is created means both survivors and non-survivors have a 1 sd overlap of treatment between 60 and 80. Status doesn't depend on days, so why are you doing cox at all? My understanding: There are 100 patients. Each patient receives a treatment dose between 0 and 140. Patients then makes a single followup at ~365 days and are either in remission (status = 1) or relapsed (status = 0). Higher doses are more likely to cause remission. I am assuming relapse patients that present too early are likely to be misdiagnosed as in remission and that is how days should be incorporated? Commented Dec 18, 2019 at 19:59

The answer given by @Ron Jensen is incorrect. By simply calling the predict function without specifying the newdata argument, you are predicting the survival probability given the radiation.dose and the observed value of time. You can easily see this because there is no clear relationship in your plot, which has to be the case since you used a cox-regression model.

If you want to plot the survival probability at t as a function of a continuous variable based on a cox-model, you can instead use the contsurvplot package (https://cran.r-project.org/package=contsurvplot) I created. Using your example, you can do it like this:

library(survival)
library(contsurvplot)

set.seed(42)

days <- rpois(100, 365)
status <- rbinom(100, 1, 0.34)
radiation.dose <- sapply(status,function(x){ifelse(x==0, rnorm(1,80,20), rnorm(1,60,20))})

mod <- cox.mod<-coxph(Surv(days, status) ~ radiation.dose, data=df, x=TRUE)

plot_surv_at_t(time="days",
status="status",
data=df,
model=cox.mod,
t=365)


By specifying t = 365, you're telling the plot_surv_at_t function to use the survival probability at 1 year.

Since there are no other independent variables in this cox model, this is equivalent to using:

vals <- seq(min(df$$radiation.dose), max(df$$radiation.dose), 1)