1 Introduction

Futility stopping is a common component of Bayesian adaptive clinical trial designs.

At an interim analysis, the observed data are used to estimate the probability that the trial will eventually meet its prespecified success criterion if enrollment continues to the planned sample size.

This probability is often called the:

Predictive probability of success
Predictive power
Posterior predictive probability

If this probability is sufficiently low, the trial may be stopped early for futility.

2 Trial Scenario

Consider a single-arm clinical trial with a planned total enrollment of 200 patients.

At the interim analysis:

\[ N_{\text{current}} = 50 \]

patients have been enrolled, and

\[ R_{\text{current}} = 15 \]

patients have responded.

The observed interim response rate is therefore:

\[ \frac{15}{50}=0.30. \]

The remaining number of patients is:

\[ N_{\text{remaining}}=200-50=150. \]

The final success criterion is defined as:

\[ \text{Final response rate}>35%. \]

Therefore, the trial is successful if:

\[ \frac{R_{\text{final}}}{200}>0.35. \]

Because the number of responders must be an integer, the trial must have at least 71 total responders:

\[ R_{\text{final}}\geq 71. \]

Since 15 responders have already been observed, at least 56 responses are required among the remaining 150 patients.

3 Method 1: Beta-Binomial Posterior Predictive Simulation

3.1 Prior and Posterior Distribution

Assume a noninformative uniform prior for the response probability:

\[ p \sim \operatorname{Beta}(1,1). \]

After observing 15 responders and 35 nonresponders, the posterior distribution is:

\[ p\mid\text{data} \sim \operatorname{Beta}(1+15,1+35). \]

Therefore:

\[ p\mid\text{data} \sim \operatorname{Beta}(16,36). \]

For every posterior draw of (p), we simulate the number of responders among the remaining 150 patients:

\[ R_{\text{future}}\mid p \sim \operatorname{Binomial}(150,p). \]

The final response rate is:

\[ \frac{15+R_{\text{future}}}{200}. \]

The predictive probability of success is estimated by the proportion of simulated trials in which the final response rate exceeds 35%.

3.2 R Code

set.seed(2026)

# Interim trial information
n_current <- 50
responders_current <- 15
nonresponders_current <- n_current - responders_current

# Planned final sample size
n_total <- 200
n_remaining <- n_total - n_current

# Trial success threshold
success_threshold <- 0.35

# Number of posterior predictive simulations
n_mcmc <- 10000

# Beta(1,1) prior
prior_alpha <- 1
prior_beta <- 1

# Posterior parameters
posterior_alpha <- prior_alpha + responders_current
posterior_beta <- prior_beta + nonresponders_current

# Draw posterior samples of the response probability
p_post_beta <- rbeta(
  n = n_mcmc,
  shape1 = posterior_alpha,
  shape2 = posterior_beta
)

# Simulate future responders
future_responders_beta <- rbinom(
  n = n_mcmc,
  size = n_remaining,
  prob = p_post_beta
)

# Calculate final response rates
total_responders_beta <-
  responders_current + future_responders_beta

final_response_rate_beta <-
  total_responders_beta / n_total

# Predictive probability of trial success
predictive_probability_beta <- mean(
  final_response_rate_beta > success_threshold
)

cat(
  "Beta-binomial predictive probability of success:",
  round(100 * predictive_probability_beta, 2),
  "%\n"
)

## Beta-binomial predictive probability of success: 20.07 %

4 Futility Decision

Suppose the prespecified futility boundary is 5%.

The trial stops for futility if:

\[ P(\text{final success}\mid\text{interim data})<0.05. \]

futility_threshold <- 0.05

if (predictive_probability_beta < futility_threshold) {

  cat(
    "Decision: Stop the trial for futility because",
    "the predictive probability is below 5%.\n"
  )

} else {

  cat(
    "Decision: Do not stop for futility because",
    "the predictive probability is at least 5%.\n"
  )
}

## Decision: Do not stop for futility because the predictive probability is at least 5%.

5 Posterior Distribution of the Response Probability

hist(
  p_post_beta,
  breaks = 40,
  probability = TRUE,
  main = "Posterior Distribution of Response Probability",
  xlab = "Response probability p"
)

abline(
  v = success_threshold,
  lty = 2,
  lwd = 2
)

The dashed vertical line represents the final response-rate success threshold of 35%.

However, the predictive probability is not simply:

\[ P(p>0.35\mid\text{data}). \]

It also accounts for future binomial sampling variation among the remaining patients.

6 Method 2: Logistic-Binomial Model Using rstanarm

The same interim data can also be analyzed using a Bayesian intercept-only logistic regression model.

The model is:

\[ R_{\text{current}} \sim \operatorname{Binomial}(N_{\text{current}},p), \]

with:

\[ \operatorname{logit}(p)=\beta_0. \]

Equivalently:

\[ p=\operatorname{logit}^{-1}(\beta_0). \]

7 Load rstanarm

library(rstanarm)

8 Prepare the Data

dat <- data.frame(
  responders = responders_current,
  nonresponders = nonresponders_current
)

dat

##   responders nonresponders
## 1         15            35

The expression

cbind(responders, nonresponders)

is the grouped binomial outcome.

It tells the model that one row of data contains:

15 successful responses
35 unsuccessful responses

9 Define the Prior

The prior is placed on the logit-scale intercept:

\[ \beta_0 \sim N\left( \operatorname{logit}(0.30), 1^2 \right). \]

Since:

\[ \operatorname{logit}(0.30) = \log\left(\frac{0.30}{0.70}\right), \]

the prior is centered at a response probability of 30%.

prior_probability <- 0.30

prior_logit_mean <- qlogis(prior_probability)

prior_logit_mean

## [1] -0.8472979

The prior standard deviation of 1 is specified on the logit scale, not directly on the probability scale.

10 Fit the rstanarm Model

fit_stan <- stan_glm(
  cbind(responders, nonresponders) ~ 1,
  family = binomial(link = "logit"),
  data = dat,
  prior_intercept = normal(
    location = qlogis(0.30),
    scale = 1,
    autoscale = FALSE
  ),
  chains = 4,
  iter = 5000,
  warmup = 2500,
  seed = 2026,
  refresh = 0
)

11 Model Summary

print(
  fit_stan,
  digits = 3
)

## stan_glm
##  family:       binomial [logit]
##  formula:      cbind(responders, nonresponders) ~ 1
##  observations: 1
##  predictors:   1
## ------
##             Median MAD_SD
## (Intercept) -0.854  0.291
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg

12 Extract Posterior Response Probabilities

The sampled intercept values are initially on the logit scale.

logit_p_samples <- as.matrix(
  fit_stan,
  pars = "(Intercept)"
)[, 1]

The inverse-logit transformation converts each posterior draw to the probability scale:

\[ p = \frac{\exp(\beta_0)} {1+\exp(\beta_0)}. \]

p_post_stan <- plogis(
  logit_p_samples
)

summary(p_post_stan)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  0.1123  0.2591  0.2985  0.3007  0.3411  0.6044

13 Posterior Predictive Simulation Using rstanarm

For each posterior draw of (p), simulate the number of responders among the remaining 150 patients.

set.seed(2026)

future_responders_stan <- rbinom(
  n = length(p_post_stan),
  size = n_remaining,
  prob = p_post_stan
)

total_responders_stan <-
  responders_current + future_responders_stan

final_response_rate_stan <-
  total_responders_stan / n_total

predictive_probability_stan <- mean(
  final_response_rate_stan > success_threshold
)

cat(
  "rstanarm predictive probability of success:",
  round(100 * predictive_probability_stan, 2),
  "%\n"
)

## rstanarm predictive probability of success: 17.03 %

14 Futility Decision Using rstanarm

if (predictive_probability_stan < futility_threshold) {

  cat(
    "Decision: Stop the trial for futility because",
    "the predictive probability is below 5%.\n"
  )

} else {

  cat(
    "Decision: Do not stop for futility because",
    "the predictive probability is at least 5%.\n"
  )
}

## Decision: Do not stop for futility because the predictive probability is at least 5%.

15 Comparison of the Two Methods

comparison <- data.frame(
  Method = c(
    "Beta-binomial",
    "rstanarm logistic model"
  ),
  Prior = c(
    "Beta(1, 1) on p",
    "Normal(logit(0.30), 1) on logit(p)"
  ),
  Posterior_Mean_p = c(
    mean(p_post_beta),
    mean(p_post_stan)
  ),
  Predictive_Probability = c(
    predictive_probability_beta,
    predictive_probability_stan
  )
)

comparison$Posterior_Mean_p <- round(
  comparison$Posterior_Mean_p,
  4
)

comparison$Predictive_Probability <- round(
  comparison$Predictive_Probability,
  4
)

comparison

##                    Method                              Prior Posterior_Mean_p
## 1           Beta-binomial                    Beta(1, 1) on p           0.3082
## 2 rstanarm logistic model Normal(logit(0.30), 1) on logit(p)           0.3007
##   Predictive_Probability
## 1                 0.2007
## 2                 0.1703

16 Visual Comparison of Posterior Distributions

plot(
  density(p_post_beta),
  lwd = 2,
  xlim = c(0, 0.7),
  main = "Comparison of Posterior Distributions",
  xlab = "Response probability p"
)

lines(
  density(p_post_stan),
  lwd = 2,
  lty = 2
)

abline(
  v = success_threshold,
  lty = 3,
  lwd = 2
)

legend(
  "topright",
  legend = c(
    "Beta-binomial",
    "rstanarm",
    "Success threshold"
  ),
  lty = c(1, 2, 3),
  lwd = 2,
  bty = "n"
)

17 Why the Results May Differ

The two methods use the same likelihood but do not use the same prior distribution.

17.1 Beta-binomial method

The first method uses:

\[ p\sim\operatorname{Beta}(1,1). \]

This is a uniform prior directly on the probability scale.

17.2 rstanarm method

The second method uses:

\[ \operatorname{logit}(p) \sim N\left( \operatorname{logit}(0.30), 1^2 \right). \]

This is a logistic-normal prior.

A logistic-normal prior is not equivalent to a Beta prior. Therefore, the two approaches can produce different posterior distributions and different predictive probabilities.

The difference is caused by the prior specification, not by MCMC itself.

18 More Comparable rstanarm Prior

To make the two analyses more comparable, one could use a weakly informative prior centered near the same region as the Beta prior.

However, a normal prior on the logit scale cannot exactly reproduce a ((1,1)) prior on the probability scale.

The following is an example of a relatively weak logit-scale prior:

fit_stan_weak <- stan_glm(
  cbind(responders, nonresponders) ~ 1,
  family = binomial(link = "logit"),
  data = dat,
  prior_intercept = normal(
    location = 0,
    scale = 2.5,
    autoscale = FALSE
  ),
  chains = 4,
  iter = 5000,
  warmup = 2500,
  seed = 2026,
  refresh = 0
)

This prior is centered at:

\[ p=0.50, \]

because:

\[ \operatorname{logit}^{-1}(0)=0.50. \]

It is still not identical to a ((1,1)) prior.

19 Exact Beta-Binomial Predictive Probability

Because the Beta prior is conjugate to the binomial likelihood, the predictive probability can also be calculated without simulation.

Under the posterior distribution:

\[ p\mid\text{data} \sim \operatorname{Beta}(16,36), \]

the future number of responders follows a beta-binomial distribution:

\[ R_{\text{future}} \sim \operatorname{BetaBinomial} \left( 150,16,36 \right). \]

The trial requires at least 56 future responders.

minimum_total_responders <-
  floor(success_threshold * n_total) + 1

minimum_future_responders <-
  minimum_total_responders - responders_current

beta_binomial_pmf <- function(k, n, alpha, beta) {

  choose(n, k) *
    beta(k + alpha, n - k + beta) /
    beta(alpha, beta)
}

exact_predictive_probability <- sum(
  sapply(
    minimum_future_responders:n_remaining,
    beta_binomial_pmf,
    n = n_remaining,
    alpha = posterior_alpha,
    beta = posterior_beta
  )
)

cat(
  "Exact beta-binomial predictive probability:",
  round(100 * exact_predictive_probability, 2),
  "%\n"
)

## Exact beta-binomial predictive probability: 19.76 %

The exact result should be close to the Monte Carlo estimate. Small differences are caused by Monte Carlo sampling error.

20 Monte Carlo Error

The posterior predictive probability is estimated as a sample proportion.

Its approximate Monte Carlo standard error is:

\[ \operatorname{MCSE} = \sqrt{ \frac{ \widehat{P}(1-\widehat{P}) }{ M } }, \]

where (M) is the number of posterior predictive simulations.

mcse_beta <- sqrt(
  predictive_probability_beta *
    (1 - predictive_probability_beta) /
    n_mcmc
)

cat(
  "Approximate Monte Carlo standard error:",
  round(mcse_beta, 5),
  "\n"
)

## Approximate Monte Carlo standard error: 0.00401

Increasing n_mcmc will reduce Monte Carlo error.

For example:

n_mcmc <- 100000

21 Important Interpretation

The predictive probability answers the following question:

Given the currently observed data, the specified prior, and the assumed model, what is the probability that the final trial response rate will exceed 35% if enrollment continues to 200 patients?

It is not the same as:

\[ P(p>0.35\mid\text{current data}). \]

The posterior probability (P(p>0.35)) describes uncertainty about the true response probability.

The predictive probability additionally includes the future sampling variability from the remaining patients.

22 Key Conclusions

The Beta-binomial method provides a simple conjugate analysis for a single-arm binary endpoint.
The posterior distribution under a ((1,1)) prior is: