---
title: "Survival Models Part 2 - PSA"
date: "`r Sys.Date()`"
output:
rmarkdown::html_vignette:
toc: true
vignette: >
%\VignetteIndexEntry{Survival Models Part 2 - PSA}
%\VignetteEngine{knitr::rmarkdown}
\usepackage[utf8]{inputenc}
---
```{r setup, include=FALSE}
library(heemod)
library(flexsurv)
```
The aim of this vignette is to show how to perform probabilistic sensitivity
analyses (PSA) on survival data.
The main function to use is `resample_surv()` within `define_psa()`.
# From a user-defined parametric distribution
The `resample_surv()` function uses a random generator with the parameters of the
initially defined distribution. All you need to do is specify the `n` argument to
define the number of draws, and you can control variability in this way (the higher
`n` is, the lower the variability).
Under the hood, an empirical cumulative function is created with these random draws
and a nonlinear model determining Least Squares estimates of the new parameters
is then fitted.
```{r}
surv_dist <- define_surv_dist("gamma", shape = 2, rate = 0.1)
psa <- define_psa(surv_dist ~ resample_surv(n = 500))
```
You can display the distribution and its confidence interval using the `plot()`
function. This is useful to check the expected variability of your model.
```{r}
psa2 <- define_psa(surv_dist ~ resample_surv(n = 50))
plot(surv_dist, psa = psa)
plot(surv_dist, psa = psa2)
```
# From real data
A non-parametric bootstrap (random sampling) of the `data.frame` is performed.
At each iteration of the PSA, a new `data.frame` is created and the model runs
with this new data.
To perform a PSA with real data, `resample_surv()` must not contain any arguments.
```{r}
fit_cov <- flexsurv::flexsurvreg(survival::Surv(rectime, censrec) ~ group,
data = bc,
dist = "exp")|>
define_surv_fit()
psa <- define_psa(fit_cov ~ resample_surv())
plot(fit_cov, times = 1:1000, psa = psa, Nrep = 10)
```
It is possible to carry out a PSA with all objects of class `surv_object`, including
complex objects, created using a sequence of operations.
Taking the previous vignette as an example:
```{r}
fitcov_poor <- set_covariates(fit_cov, group = "Poor")
fitcov_medium <- set_covariates(fit_cov, group = "Medium")
fit_w <- flexsurvreg(
formula = Surv(futime, fustat) ~ 1,
data = ovarian, dist = "weibull"
) |>
define_surv_fit()
fit_cov |>
set_covariates(group = "Good") |>
apply_hr(hr = 2) |>
mix(
fitcov_medium,
weights = c(0.25, 0.75)
) |>
add_hazards(
fit_w
) |>
join(
fitcov_poor,
at = 500
) |>
plot(psa = psa, 1:1000, Nrep = 10)
```