BayesSurv: Compute the posterior mean and a credible band for the...
In BayesSurvival: Bayesian Survival Analysis for Right Censored Data

Description Usage Arguments Details Value References See Also Examples

This is the main function of this package, computing relevant quantities for a Bayesian survival analysis of (possibly right-censored) time-to-event-data. Starting with a piecewise exponential prior with dependent or independent Gamma heights (details below) on the hazard function, the function computes the posterior mean for the hazard, cumulative hazard and survival function, serving as an estimator for the true functions. In addition, for the cumulative hazard and survival function, the radius for a fixed-width credible band is computed. The interpretation of this credible band as a confidence band is justified in Castillo and Van der Pas (2020).

BayesSurv(
  df,
  time = "time",
  event = "event",
  prior = c("Dependent", "Independent"),
  K = ceiling((dim(df)[1]/log(dim(df)[1]))^(1/2)),
  time.max = max(df[[time]]),
  alpha = 0.05,
  N = 1000,
  alpha.dep = 1,
  alpha0.dep = 1.5,
  beta0.dep = 1,
  alpha.indep = 1.5,
  beta.indep = 1,
  surv.factor = 10,
  surv.epsilon = 1e-10
)

`df`	A dataframe, containing at minimum a column with follow-up times and a column with a status indicator (event observed or censored).
`time`	The name of the column in the dataframe containing the (possibly right-censored) follow-up times, that is, the minimum of the time of the event and the time of censoring. Input the name as character/string.
`event`	The name of the column in the dataframe containing the status indicator, which must be coded as: 0 = censored, 1 = event observed. Input the name as character/string.
`prior`	Select either dependent or independent Gamma heights for the piecewise exponential prior on the hazard. Dependent heights (with the Markov structure described below) is default.
`K`	The number of intervals to be used in the piecewise exponential (histogram) prior. Default is set to K = (n / \log{n})^{1/2}, with n the number of observations, as recommended by Castillo and Van der Pas (2020).
`time.max`	The maximum follow-up time to consider, corresponding to the parameter τ in Castillo and Van der Pas (2020).
`alpha`	The function will compute (1-`alpha`)100% credible bands for the cumulative hazard and survival function.
`N`	The number of samples to draw from the posterior.
`alpha.dep`	For the dependent Gamma prior only. The main parameter α for the dependent Gamma prior, as described below. It is recommended to take `alpha.dep` smaller than `alpha0.dep`.
`alpha0.dep`	For the dependent Gamma prior only. The shape parameter for the Gamma prior on the histogram height for the first interval. It is recommended to take `alpha.dep` smaller than `alpha0.dep`.
`beta0.dep`	For the dependent Gamma prior only. The rate parameter for the Gamma prior on the histogram height for the first interval.
`alpha.indep`	For the independent Gamma prior only. The shape parameter for the Gamma prior on the histogram height for each interval.
`beta.indep`	For the independent Gamma prior only. The rate parameter for the Gamma prior on the histogram height for each interval.
`surv.factor`	The survival function is computed on an equispaced grid consisting of `K x surv.factor` (the number of intervals times this factor).
`surv.epsilon`	The survival function is computed on the interval [0, `time.max` - `surv.epsilon`].

There are two options for the prior: a piecewise exponential (histogram) prior with dependent Gamma heights and a piecewise exponential (histogram) prior with independent Gamma heights. Both priors are described in detail in Castillo and Van der Pas (2020). The dependent prior has a Markov structure, where the height of each interval depends on the height of the previous interval. It implements the autoregressive idea of Arjas and Gasbarra (1994). With λ_k the histogram height on interval k and α a user-selected parameter, the structure is such that, with K the number of intervals:

E[λ_k | λ_{k-1}, …, λ_1] = λ_{k-1}, k = 2, …, K.

Var(λ_k | λ_{k-1}, …, λ_1) = (λ_{k-1})^2/α, k = 2, …, K.

In the independent Gamma prior, the prior draws for the λ_k's are independent of each other and are taken from a Gamma distribution with user-specified shape and rate parameters.

The guideline for the number of intervals K suggested by Castillo and Van der Pas (2020) is

K = (n / \log{n})^{1/(1 + 2γ)},

where n is the number of observations and γ is related to the smoothness of the true hazard function. In the absence of information about the smoothness, a default value of γ = 1/2 is recommended and this is implemented as the default in this package. If this choice leads to many intervals with zero events, it is recommended to decrease the number of intervals.

The samplers used for the dependent and independent Gamma priors are described in the Supplement to Castillo and Van der Pas (2020).

`haz.post.mean`	The posterior mean for the hazard, given as the value on each of the K intervals.
`cumhaz.post.mean`	The posterior mean for the cumulative hazard, given as the value at the end of each of the K intervals. The cumulative hazard can be obtained from this by starting at 0 and linearly interpolating between each of the returned values.
`cumhaz.radius`	The radius for the credible set for the cumulative hazard.
`surv.post.mean`	The posterior mean for the survival, given at each value contained in the also returned `surv.eval.grid`.
`surv.radius`	The radius for the credible set for the survival.
`surv.eval.grid`	The grid on which the posterior mean for the survival has been computed. A finer grid can be obtained by increasing `surv.factor` in the function call.
`time.max`	The maximum follow-up time considered.

Castillo and Van der Pas (2020). Multiscale Bayesian survival analysis. <arXiv:2005.02889>.

Arjas and Gasbarra (1994). Nonparametric Bayesian inference from right censored survival data, using the Gibbs sampler. Statistica Sinica 4(2):505-524.

PlotBayesSurv for a function that takes the result from BayesSurv() and produces plots of the posterior mean of the hazard, the posterior mean and credible band for the cumulative hazard, and the posterior mean and credible band for the survival. To obtain direct samples from the posterior for the hazard, see SamplePosteriorDepGamma and SamplePosteriorIndepGamma.

#Demonstration on a simulated data set
library(simsurv)
library(ggplot2)
hazard.true <- function(t,x, betas, ...){1.2*(5*(t+0.05)^3 - 10*(t+0.05)^2 + 5*(t+0.05) ) + 0.7}
sim.df <- data.frame(id = 1:1000)
df <- simsurv(x = sim.df, maxt = 1, hazard = hazard.true)

bs <- BayesSurv(df, "eventtime", "status")
PlotBayesSurv(bs, object = "survival")
PlotBayesSurv(bs, object = "cumhaz")
PlotBayesSurv(bs, object = "hazard")