Feature request: plot the survival function with competing events data

[anddis]

With competing events data, besides reporting the Cumulative Incidence Functions (P(T≤t, D=k)) for the different competing events (D=1,...,k), it is often desirable to report the survival function (ie, P(T>t)).

Describe the solution you'd like In a survfit2/survfitms object, the survival function is returned as states = (s0). It would be fantastic if ggcuminc() allowed to plot the survival function as well, instead of limiting the plotting to the CIFs. (I have never used tidycmprsk::cuminc() and I wouldn't know if the survival function is "left behind" there, too).

Describe alternatives you've considered Currently, I create another survfit2 object where all the competing events are combined and considered as the event. I then plot the survival function. This seems so unnecessary, though, since the survival function is already there in the orignial survfit2/survfitms object!

Additional context I see that the survival function is not returned by ggsurvfit::tidy_survfit(). So, I'm not surprised ggcuminc() cannot plot it. tidy_survfit() is also an extremely handy function, and if it returned the survival function, as a "byproduct" of my feature request above, it would be so so great!

Thank you for developing such a great package, Daniel!

[ddsjoberg]

Hej @anddis , thanks for the post!

I've never seen competing risks presented like this before. Can you point to a published example?

[anddis]

Hello Daniel, thank you for the reply.

Well, "desirable" was probably too strong, but sometimes the survival function is a useful complement to the CIFs (but, as always, "it depends").

The Survival function would give you the estimated probability of being still at risk for any competing event (ie, being event-free, or in the initial state) as a function of time. It's the complement to 1 of the sum of the CIFs.

See, for example, Equation 1 in "Andersen, P.K., Abildstrom, S.Z. and Rosthøj, S., 2002. Competing risks as a multi-state model. Statistical methods in medical research, 11(2), pp.203-215." (https://journals.sagepub.com/doi/pdf/10.1191/0962280202sm281ra)

One example is Figure 4 from "Putter, H., Fiocco, M. and Geskus, R.B., 2007. Tutorial in biostatistics: competing risks and multi‐state models. Statistics in medicine, 26(11), pp.2389-2430." (https://onlinelibrary.wiley.com/doi/abs/10.1002/sim.2712).

I would like to plot directly the survival function (by getting Pr(T>t) from ggcuminc()) instead of stacking the CIFs and working visually by difference (as in Figure 4, "event-free" area).

For example:

library(ggsurvfit)
library(tidycmprsk)
library(patchwork)

sf1 <- survfit2(Surv(ttdeath, death_cr) ~ 1, data = trial)
g1 <- ggcuminc(sf1, outcome = "death from cancer") + labs(title = "death from cancer") + ylim(c(0,1))
g2 <- ggcuminc(sf1, outcome = "death other causes") + labs(title = "death from other causes") + ylim(c(0,1))

sf2 <- survfit2(Surv(ttdeath, death_cr %in% c("death other causes", "death from cancer")) ~ 1, data = trial)
g3 <- ggsurvfit(sf2) + labs(title = "Still alive (event-free)") + ylim(c(0,1))

g1 + g2 + g3

# survival function from `sf1`
sf1$pstate[, 1]
# survival function from `sf2`
sf2$surv

> identical(sf1$pstate[, 1], sf2$surv)
[1] TRUE

[ddsjoberg]

Ahhh, thank you for the clarifications!

I think this is something that can be incorporated in a way. But before we go further, can you mock-up a figure like the one you have in mind and post it here?

[anddis]

Sorry if my initial request was unclear, Daniel.

Ideally, one could specify outcome = "(s0)" in ggcuminc() and obtain the figure below:

# Mock-up code/figure
sf1 <- survfit2(Surv(ttdeath, death_cr) ~ 1, data = trial)
ggcuminc(sf1, outcome = "(s0)") # (s0) corresponds to the survival curve

or

# Mock-up code/figure
sf1 <- survfit2(Surv(ttdeath, death_cr) ~ 1, data = trial)
ggcuminc(sf1, outcome = c("death from cancer", "death other causes", "(s0)")) 

Similarly, tidy_survfit() would return, together with the CIFs for all competing events, the survival probabilities too (as outcome = (s0), for example).

Thank you for considering this request!

[ddsjoberg]

The primary reason to use ggsurvfit is for the risktables and the integration with ggplot. Otherwise, it's pretty simple to create figures with the broom::tidy() output.

Implementation will of the survival-like function for competing risks would actually be quite complicated and it's unclear what the risktable should look like in this setting.

Because ggsurvfit figures integrate so well with ggplot, you can simply pipe in an additional geom_step() to add the probability of no events.

library(survival)
library(ggsurvfit)
#> Loading required package: ggplot2

sf <- survfit2(Surv(ttdeath, death_cr) ~ 1, tidycmprsk::trial)

sf |> 
  ggcuminc(outcome = c("death from cancer", "death other causes")) +
  add_confidence_interval() +
  add_risktable() +
  geom_step(
    data = 
      sf |> 
      survival::survfit0() |> 
      broom::tidy() |> 
      dplyr::filter(state %in% "(s0)"),
    mapping = aes(x = time, y = estimate, color = state)
  )

^{Created on 2023-01-25 with reprex v2.0.2}

Is this sufficient for your needs?

[anddis]

I see. Yes, this was useful. Thank you, again, for taking the time to consider my request!

(To get a risktable for the survival-like function alone (with number still at risk/in the initial state and number of events (any competing event)) I'll continue to do something like:)

library(ggsurvfit)
library(tidycmprsk)

sf2 <- survfit2(Surv(ttdeath, death_cr %in% c("death other causes", "death from cancer")) ~ 1, data = trial)
ggsurvfit(sf2, col = "red") + 
  ylim(c(0,1)) + 
  add_risktable()

[ddsjoberg]

Great, thanks for all the details and posts.

I know it's likely a less satisfying outcome than you were hoping for. But for the moment, it would be too much of a time commitment to add and it seems we can create the needed figure.

All the best!