\documentclass[12pt]{article}
\usepackage{graphicx}
\usepackage{amssymb,amsmath,amsthm}
\usepackage[round,authoryear]{natbib}
\bibliographystyle{abbrvnat}

\setlength{\parindent}{0in}
\setlength{\parskip}{8pt}
\addtolength{\topmargin}{-.4in}
\setlength{\textheight}{8.7in}
\setlength{\textwidth}{6.5in}
\setlength{\oddsidemargin}{0in}

\renewcommand{\exp}{\mathrm{e}}
\renewcommand{\d}{\mathrm{d}}

\title{A Bayesian Model for the Endpoint Event Incidence Rate in an Interim Analysis of Operational Futility}
\author{}
\date{}

% \VignetteEngine{knitr::knitr}
% \VignetteIndexEntry{Bayesian Model for Incidence Rate}

\begin{document}
\maketitle
Let $n_k$ and $T_k$ denote, respectively, the event count and the observed total person-time at risk at the time of the $k$-th futility analysis, pooling over all treatment arms. Additionally, let $T^{\ast}$ denote the estimated total person-time at risk for the primary efficacy analysis. Let the prior distribution of the treatment arm-pooled incidence rate $p$ be $\mathsf{Ga}(\alpha, \beta)$ parametrized such that the prior mean $E\, p = \alpha / \beta$ (the same Bayesian method applies to treatment arm-specific incidence rates). For the treatment arm-pooled incidence rate, we additionally consider a robust prior distribution described in Section~\ref{s:robustPrior}.

Generally, assuming that, conditional on $p$, the times to event follow $\mathsf{Exp}(p)$, the posterior mean of $p$ at the time of the $k$-th analysis equals
\begin{align}
E [p\,|\,\textrm{data}] &= \frac{\alpha + n_k}{\beta + T_k}\nonumber\\
&= \frac{\alpha}{\beta}\, \frac{\beta}{\beta+T_k} + \frac{n_k}{T_k}\, \frac{T_k}{\beta+T_k},\label{posterior mean}
\end{align}
i.e., the posterior mean can be interpreted as a convex combination of the prior mean and the observed incidence rate. For a given $\beta>0$, the weight on the prior mean at the first analysis depends on the accumulated person-time at risk ($T_1$), and the weight will decrease at subsequent analyses because $\beta/(\beta+T_k)$ is a decreasing function of $T_k$, which is a desirable Bayesian property.

In order to identify $\alpha$ and $\beta$, it is desirable that the prior mean equals the pre-trial assumed treatment arm-pooled incidence rate $p^{\ast}$ (e.g., in a trial in which participants are randomized to treatment and placebo in the 2:1 ratio, assuming the incidence rate of 0.055 endpoints per person-year at risk in the placebo group and $TE=60\%$, $p^{\ast}=(1/3) \times 0.055 + (2/3) \times 0.4 \times 0.055 = 0.033$), i.e.,
\begin{equation}
\frac{\alpha}{\beta} = p^{\ast}.\label{priorMean}
\end{equation}
Furthermore, we propose to consider three values of $\beta$ that correspond to the weights $w=\frac{1}{2}$, $\frac{1}{3}$ and $\frac{1}{4}$ on the prior mean at the time when 50\% of the estimated total person-time at risk has been accumulated, i.e., for each value of $w$, $\beta$ is defined as the solution to the equation
\begin{equation*}
\frac{\beta}{\beta + T^{\ast}/2} = w.
\end{equation*}
It follows that
\begin{equation}
\beta = \beta(w, T^{\ast}) = \frac{w T^{\ast}}{2(1-w)}\, ,\label{beta}
\end{equation}
and the estimation of $T^{\ast}$ is described in Section~\ref{s:estTotalPYRs}. For $w=\frac{1}{2}$, $\frac{1}{3}$ and $\frac{1}{4}$, we obtain $\beta=\frac{T^{\ast}}{2}$, $\frac{T^{\ast}}{4}$, and $\frac{T^{\ast}}{6}$, respectively.

At the $k$-th futility analysis and for each of the three values of $\beta$, we will sample the incidence rate from $\mathsf{Ga}(\alpha+n_k,\beta+T_k)$ for generating future data and report the weight $\frac{\beta}{\beta+T_k}$ on the prior mean in the convex combination \eqref{posterior mean}.

\section{A Robust Mixture Prior Distribution for the Endpoint Event Incidence Rate}\label{s:robustPrior}
The robust prior model \citep{Schmidli2014} is implemented since it is designed to maximize the probability of meeting, e.g., an enrollment expansion guideline for large downward deviations from the protocol-assumed incidence rates, while minimizing a false trigger for protocol-assumed incidence rates.

The prior distribution of $p$ is defined as a weighted mixture of two gamma distributions,
\begin{equation*}
(1-w_R)\mathsf{Ga}(\alpha_I,\beta_I)+w_R \mathsf{Ga}(\alpha_V,\beta_V),
\end{equation*}
where we set, e.g., $w_R=0.2$, and $\mathsf{Ga}(\alpha_V,\beta_V)$ and $\mathsf{Ga}(\alpha_I,\beta_I)$ represent the weakly informative and informative component of the mixture prior, respectively. The parameters $\beta_V$ and $\beta_I$ are calculated following \eqref{beta} with, e.g., $w=1/1000$ and $w=1/3$, respectively (and $T^{\ast}$ per Section~\ref{s:estTotalPYRs}). Subsequently, $\alpha_V$ and $\alpha_I$ are calculated following \eqref{priorMean} with $p^{\ast}$ set to the pre-trial assumed treatment arm-pooled incidence rate for both components of the mixture.

The posterior distribution at the time of the $k$-th analysis is derived following the conjugacy principle, as in \eqref{posterior mean}, which results in a mixture of conjugate posteriors with updated weights
\begin{equation*}
(1-\widetilde{w}_{R,k})\mathsf{Ga}(\alpha_I+n_k,\beta_I+T_k)+ \widetilde{w}_{R,k} \mathsf{Ga}(\alpha_V+n_k,\beta_V+T_k),
\end{equation*}
where
\begin{equation*}
\widetilde{w}_{R,k} \propto w_{R,k} f_V/ \left\{ w_{R,k} f_V + ( 1 - w_{R,k}) f_I\right\}
\end{equation*}
with $f_{\cdot}$ equal to
\begin{equation*}
f_{\cdot}=\frac{\Gamma(\alpha_{\cdot}+n_k)/(\beta_{\cdot}+T_k)^{\alpha_{\cdot}+n_k}}{\Gamma(\alpha_{\cdot})/\beta_{\cdot}^{\alpha_{\cdot}}}
\end{equation*}
\citep[see, e.g.,][Section 5.2.3, pages 279--282]{Bernardo2000}.

\section{Estimation of the Total Person-Years at Risk ($T^{\ast}$)}\label{s:estTotalPYRs}
We consider the standard right-censored failure time analysis framework. Denoting the failure and censoring times as $T$ and $C$, respectively, we assume that $T$ is independent of $C$, $T \sim \mathsf{Exp}(p^{\ast})$, and $C \sim \mathsf{Exp}(d^{\ast})$. It follows that $X := \min(T,C) \sim \mathsf{Exp}(p^{\ast}+d^{\ast})$ and
\begin{equation*}
\begin{split}
T^{\ast}&= N \times E[\min(X,\tau)]\\
&= N \times \big\{ E[X\mid X\leq \tau]\, P(X\leq\tau) + \tau\, P(X>\tau)\big\}\\
&= N \times \bigg\{ (p^{\ast}+d^{\ast}) \int_0^{\tau} x\exp^{-(p^{\ast}+d^{\ast})x}\d x + \tau \exp^{-(p^{\ast}+d^{\ast})\tau}\bigg\}\\
&= N \times \frac{1-\exp^{-(p^{\ast}+d^{\ast})\tau}}{p^{\ast}+d^{\ast}}\,.
\end{split}
\end{equation*}

To illustrate, we consider the total target sample size $N=1,500$ with a 2:1 randomization ratio to treatment vs. placebo, the duration of follow-up per participant $\tau=80/52$ years, the pre-trial assumed dropout rate $d^{\ast}=0.1$ dropouts per person-year at risk (PYR), and the pre-trial assumed constant incidence rate of 0.055 endpoints/PYR in the placebo group. Then, in the $TE=60\%$ scenario, the pre-trial assumed treatment arm-pooled endpoint event incidence rate is $p^{\ast}=(1/3) \times 0.055 + (2/3) \times 0.4 \times 0.055 = 0.033$ endpoints/PYR.

These assumptions result in $T^{\ast} = 2086.91$ PYRs. For comparison, if all $N$ participants were followed for $\tau$ years, the total PYRs would be $N\tau = 2307.69$ years.

Subsequently, for $T^{\ast} = 2086.91$ PYRs, if $T_1 = 0.2\, T^{\ast}$, the weights $\frac{\beta}{\beta+T_1}$ on the prior mean at the first futility analysis corresponding to $w=\frac{1}{2}$, $\frac{1}{3}$, and $\frac{1}{4}$ are 0.71, 0.56, 0.45, respectively. If $T_1 = 0.3\, T^{\ast}$, the respective weights on the prior mean are 0.63, 0.45, and 0.36.

\newpage
\bibliography{references}

\end{document}