Distribution function of the difference of two milestone survival probabilities

Computes $\Pr(S_t - S_c \le q)$ (or the upper tail $\Pr(S_t - S_c > q)$) for the milestone survival probabilities $S_k = \exp(-\lambda_k t_0)$, where the hazards have independent Gamma distributions, $\lambda_t \sim \textrm{Gamma}(a_t, b_t)$ and $\lambda_c \sim \textrm{Gamma}(a_c, b_c)$ (shape and rate). This is the posterior probability used for a time-to-event endpoint with a concurrent, external or hypothetical control.

Usage

psurvdiff(q, milestone, shape_t, rate_t, shape_c, rate_c, lower.tail = TRUE)

Arguments

q: Numeric vector of quantiles.
milestone: Positive milestone time $t_0$.
shape_t, rate_t: Positive shape and rate of the treatment hazard.
shape_c, rate_c: Positive shape and rate of the control hazard.
lower.tail: Logical; if TRUE (default), probabilities are $\Pr(S_t - S_c \le q)$, otherwise $\Pr(S_t - S_c > q)$.

Value

A numeric vector of probabilities. All arguments except lower.tail are recycled to a common length; each must have length 1 or the common length.

Details

The upper tail is computed as $$\Pr(S_t - S_c > q) = \Pr(S_c \le -q) + \int g_c(\lambda) G_t\left(-\log(q + e^{-\lambda t_0}) / t_0\right) \, d\lambda,$$ where $g_c$ is the density of $\lambda_c$, $G_t$ is the distribution function of $\lambda_t$, the integral runs over the values of $\lambda$ with $0 < q + e^{-\lambda t_0} < 1$, and the first term is zero for $q \ge 0$. The integral is evaluated by adaptive Gauss-Kronrod quadrature (the QUADPACK routine used by integrate), with the range split at quantiles of the control hazard and at the values that correspond to quantiles of the treatment hazard. When the density of $\lambda_c$ is unbounded at 0 ($a_c < 1$), the integral below its median is evaluated in the variable $u = \lambda^{a_c}$, in which the integrand is bounded. The lower tail is computed from the same formula with the roles of the two groups exchanged, so both tails are accurate near 0.

Examples

# Posterior probability that the 3-month survival is higher by more than
# 0.1 with 5 events in 40 months of follow-up versus 9 events in 35 months
# (Gamma(0.1, 1) priors)
psurvdiff(0.1, milestone = 3, shape_t = 5.1, rate_t = 41,
          shape_c = 9.1, rate_c = 36, lower.tail = FALSE)
#> [1] 0.7685102