0

I have 5 different survfit() plots of different models that I am trying to combine into one plot in the style of a landmark analysis plot, as seen below.

landmark analysis plot

At the moment they are just plot(survfit(model, newdata = )), how could I combine them so that I have the line of days 0-100 of survfit 1, 100-200 of survfit 2 etc.

ikeagoddess
  • 47
  • 4
  • It's easier to help you if you provide a [reproducible example](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) with sample input that can be used to test and verify possible solutions. – MrFlick Sep 02 '22 at 20:14
  • 3
    If you're after a landmark plot, check the `jskm`-package: https://jinseob2kim.github.io/jskm/ – harre Sep 02 '22 at 21:19

1 Answers1

1

Let's create a model using the built-in lung data from the survival package:

library(survival)
library(tidyverse)

mod1   <- survfit(Surv(time, status) ~ sex, data = lung)

This model actually contains all we need to make the plot. We can convert it to a data frame as follows:

df     <- as.data.frame(unclass(mod1)[c(2:7, 15:16)])
df$sex <- rep(c("Male", "Female"), times = mod1$strata)
head(df)
#>   time n.risk n.event n.censor      surv    std.err     lower     upper  sex
#> 1   11    138       3        0 0.9782609 0.01268978 0.9542301 1.0000000 Male
#> 2   12    135       1        0 0.9710145 0.01470747 0.9434235 0.9994124 Male
#> 3   13    134       2        0 0.9565217 0.01814885 0.9230952 0.9911586 Male
#> 4   15    132       1        0 0.9492754 0.01967768 0.9133612 0.9866017 Male
#> 5   26    131       1        0 0.9420290 0.02111708 0.9038355 0.9818365 Male
#> 6   30    130       1        0 0.9347826 0.02248469 0.8944820 0.9768989 Male

With a bit of data manipulation, we can define 100-day periods and renormalize the curves at the start of each period. Then we can plot using geom_step

df %>%
  filter(time < 300) %>%
  group_by(sex) %>%
  mutate(period = factor(100 * floor(time / 100))) %>%
  group_by(sex, period) %>%
  mutate(surv = surv / first(surv)) %>%
  ggplot(aes(time, surv, color = sex, group = interaction(period, sex))) + 
  geom_step(size = 1) +
  geom_vline(xintercept = c(0, 100, 200, 300), linetype = 2) +
  scale_color_manual(values = c("deepskyblue4", "orange")) +
  theme_minimal(base_size = 16) +
  theme(legend.position = "top")

Created on 2022-09-03 with reprex v2.0.2

Allan Cameron
  • 147,086
  • 7
  • 49
  • 87
  • 2
    It doesn’t appear that you have reset the conditional survivals to unity at each landmark. Admittedly the OP had rather loose expression what I was interpreting as his intent so who knows. I got the idea that he want survival curves conditional on survival to a land mark. It did not appear to me anyway that three different datasets was the right way. – IRTFM Sep 03 '22 at 03:56
  • Oh I see @IRTFM - looking at this again I misunderstood the question. I'll have a look at this again later. Thanks – Allan Cameron Sep 03 '22 at 06:49
  • 1
    @IRTM I think it's now what the OP is after, but as you say, the intention wasn't entirely clear. – Allan Cameron Sep 03 '22 at 12:12
  • 1
    Looks good. I think the OP gets more than she deserved. ;-) Hope she checkmarks it so you get what you deserve. – IRTFM Sep 03 '22 at 15:24
  • @AllanCameron, you've saved me thank you so much for your help, sorry I wasn't very clear I was very tired at the time. You have helped so much! – ikeagoddess Sep 04 '22 at 00:30