The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.
This document reproduces the computations in Rufibach, Jordan, and Abt (2016b) (doi). All the references below to sections, figures, or tables are with respect to that paper. The functions used are all implemented in the R
package Rufibach, Jordan, and Abt (2016a) and all the code below is also available in the help file of the function bpp_1interim
.
Further discussions around Bayesian Predictive Power (BPP), mainly about the choice of a prior, are provided in Rufibach, Burger, and Abt (2016) (doi).
We start with defining the parameters of the pivotal trial and the results from the two external studies:
# load the package:
library(bpp)
# get the code for all the computations below:
?bpp_1interim
# specifications of the the pivotal trial
0.5 # proportion of patients randomized to arm A
propA <- (propA * (1 - propA)) ^ (-1)
fac <- c(0.5, 1) * 1600
nevents <- sqrt(fac / nevents[2])
finalSE <- c(0.001, 0.049)
alphas <- qnorm(1 - alphas / 2)
za <- exp(- za * sqrt(fac / nevents))
hrMDD <- log(hrMDD[2])
successmean <-
# efficacy and futility interim boundary
log(hrMDD[1])
effi <- log(1.025)
futi <-
# specifications of the first external study
0.396
hr1 <- 0.837
sd1 <-
# specifications of the second external study
0.287
hr2 <- 0.658
sd2 <-
# grid to compute densities on
seq(-0.65, 0.3, by = 0.01) thetas <-
Then we define the parameters of the prior distributions, see Sections 4.1 and 4.2.
# ----------------------------------
# mean and sd of Normal prior:
# ----------------------------------
0.85
hr0 <- 0.11
sd0 <-
# all computations are done on the log(hazard ratio) scale, in order to
# be able to approximate the distribution of the log(hazard ratio) via a
# Normal distribution:
log(hr0)
priormean <-
# ----------------------------------
# parameters of pessimistic, or flat, prior:
# ----------------------------------
0.866
hr0flat <- log(hr0flat)
priormeanflat <- 0.21
width1 <- 2.48 height1 <-
Let us compare these two prior distributions:
par(las = 1, mar = c(9, 5, 2, 1), mfrow = c(1, 1))
plot(0, 0, type = "n", xlim = c(-0.6, 0.3), ylim = c(-0.1, 5), xlab = "", ylab = "density", main = "")
title(expression("Normal and pessimistic prior density for "*theta), line = 0.7)
basicPlot(leg = FALSE, IntEffBoundary = NA, IntFutBoundary = NA, successmean = NA, priormean = NA)
lines(thetas, dnorm(thetas, mean = priormean, sd = sd0), col = 2, lwd = 2)
lines(thetas, dUniformNormalTails(thetas, mu = priormeanflat, width = width1, height = height1), lwd = 2, col = 3)
In what follows, we compute all quantities that appear in Table 1 of Rufibach, Jordan, and Abt (2016b). In addition, the densities shown in Figures 1 and 2 are generated as well. All these computations are done for both, the Normal and the pessimistic prior.
For the Normal prior:
# Normal prior probabilities to be below 0.7 or above 1:
c(0.7, 1)
lims <- plnorm(lims[1], meanlog = priormean, sdlog = sd0, lower.tail = TRUE, log.p = FALSE)
pnorm1 <- plnorm(lims[2], meanlog = priormean, sdlog = sd0, lower.tail = FALSE, log.p = FALSE) pnorm2 <-
And the same for the pessimistic prior:
# pessimistic prior probabilities to be below 0.7 or above 1:
pUniformNormalTails(x = log(lims[1]), mu = priormeanflat, width = width1, height = height1)
flat1 <- 1 - pUniformNormalTails(x = log(lims[2]), mu = priormeanflat, width = width1,
flat2 <-height = height1)
Compute BPP based on the prior distributions and the specifications of the pivotal trial.
# ----------------------------------
# Normal prior:
# ----------------------------------
bpp(prior = "normal", successmean = successmean, finalSE = finalSE,
bpp0 <-priormean = priormean, priorsigma = sd0)
# ----------------------------------
# pessimistic prior:
# ----------------------------------
1 <- bpp(prior = "flat", successmean = successmean, finalSE = finalSE,
bpp0_priormean = priormeanflat, width = width1, height = height1)
We have two external studies that define our likelihood with which we update the two prior distributions.
# ----------------------------------
# Normal prior:
# ----------------------------------
NormalNormalPosterior(datamean = log(hr1), sigma = sd1, n = 1, nu = priormean,
up1 <-tau = sd0)
bpp(prior = "normal", successmean = successmean, finalSE = finalSE,
bpp1 <-priormean = up1$postmean, priorsigma = up1$postsigma)
# update prior with second external study (result derived from pooled analysis:
# Cox regression on patient level, stratified by study):
NormalNormalPosterior(datamean = log(hr2), sigma = sd2, n = 1, nu = priormean,
up2 <-tau = sd0)
bpp(prior = "normal", successmean = successmean, finalSE = finalSE,
bpp2 <-priormean = up2$postmean, priorsigma = up2$postsigma)
# ----------------------------------
# pessimistic prior
# ----------------------------------
1 <- integrate(FlatNormalPosterior, lower = -Inf, upper = Inf, successmean = successmean,
bpp1_finalSE = finalSE, interimmean = log(hr1), interimSE = sd1,
priormean = priormeanflat, width = width1, height = height1,
subdivisions = 300)$value
1 <- integrate(FlatNormalPosterior, -Inf, Inf, successmean = successmean,
bpp2_finalSE = finalSE, interimmean = log(hr2),
interimSE = sd2, priormean = priormeanflat,
width = width1, height = height1, subdivisions = 300)$value
The next quantities we are interested in is what happens to BPP if we do not stop at an interim analysis. We assess various scenarios, as in Table 1: interim is set up with futility and efficacy boundary, only one of the two, and we also provide updates using conditional power for the two extreme cases, namely that the estimate at the interim was known to be either equal to the futility or the efficacy boundary. Let us start with the Normal prior:
# ----------------------------------
# compute bpp after not stopping at interim, for Normal prior and various
# assumptions on the amount of information we learn at the interim
# ----------------------------------
# assuming both boundaries:
bpp_1interim(prior = "normal", interimSE = sqrt(fac / nevents[1]),
bpp3.tmp <-finalSE = finalSE, successmean = successmean,
IntEffBoundary = effi, IntFutBoundary = futi, IntFix = 1,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
priorsigma = up2$postsigma)
bpp3.tmp$"BPP after not stopping at interim interval"
bpp3 <- bpp3.tmp$"posterior density interval"
post3 <-
# assuming only efficacy boundary:
bpp_1interim(prior = "normal", interimSE = sqrt(fac / nevents[1]),
bpp3_effi_only <-finalSE = finalSE, successmean = successmean,
IntEffBoundary = effi, IntFutBoundary = log(Inf), IntFix = 1,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
priorsigma =
$postsigma)$"BPP after not stopping at interim interval"
up2
# assuming only futility boundary:
bpp_1interim(prior = "normal", interimSE = sqrt(fac / nevents[1]),
bpp3_futi_only <-finalSE = finalSE, successmean = successmean,
IntEffBoundary = log(0), IntFutBoundary = futi, IntFix = 1,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
priorsigma =
$postsigma)$"BPP after not stopping at interim interval"
up2
# assuming interim efficacy boundary:
bpp_1interim(prior = "normal", interimSE = sqrt(fac / nevents[1]),
bpp4.tmp <-finalSE = finalSE, successmean = successmean,
IntEffBoundary = effi, IntFutBoundary = Inf, IntFix = c(effi, futi),
priormean = up2$postmean, propA = 0.5, thetas = thetas,
priorsigma = up2$postsigma)
bpp4.tmp$"BPP after not stopping at interim exact"[2, 1]
bpp4 <- bpp4.tmp$"posterior density exact"[, 1]
post4 <-
# assuming interim futility boundary:
bpp_1interim(prior = "normal", interimSE = sqrt(fac / nevents[1]),
bpp5.tmp <-finalSE = finalSE, successmean = successmean,
IntEffBoundary = effi, IntFutBoundary = Inf, IntFix = futi,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
priorsigma = up2$postsigma)
bpp5.tmp$"BPP after not stopping at interim exact"[2, 1]
bpp5 <- bpp5.tmp$"posterior density exact" # same as post4[, 2] post5 <-
And the same for the pessimistic prior:
# ----------------------------------
# compute bpp after not stopping at interim, for pessimistic prior and various
# assumptions on the amount of information we learn at the interim
# ----------------------------------
# assuming both boundaries:
1 <- bpp_1interim(prior = "flat", interimSE = sqrt(fac / nevents[1]),
bpp3.tmp_finalSE = finalSE, successmean = successmean,
IntEffBoundary = effi, IntFutBoundary = futi, IntFix = 1,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
width = width1, height = height1)
1 <- bpp3.tmp_1$"BPP after not stopping at interim interval"
bpp3_1 <- bpp3.tmp_1$"posterior density interval"
post3_
# assuming only efficacy boundary:
1_effi_only <- bpp_1interim(prior = "flat", interimSE = sqrt(fac / nevents[1]),
bpp3_finalSE = finalSE, successmean = successmean,
IntEffBoundary = effi, IntFutBoundary = log(Inf), IntFix = 1,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
width = width1,
height = height1)$"BPP after not stopping at interim interval"
# assuming only futility boundary:
1_futi_only <- bpp_1interim(prior = "flat", interimSE = sqrt(fac / nevents[1]),
bpp3_finalSE = finalSE, successmean = successmean,
IntEffBoundary = log(0), IntFutBoundary = futi, IntFix = 1,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
width = width1,
height = height1)$"BPP after not stopping at interim interval"
# assuming interim efficacy boundary:
1.tmp <- bpp_1interim(prior = "flat", interimSE = sqrt(fac / nevents[1]),
bpp4_finalSE = finalSE, successmean = successmean,
IntEffBoundary = log(0), IntFutBoundary = effi, IntFix = effi,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
width = width1, height = height1)
1 <- bpp4_1.tmp$"BPP after not stopping at interim exact"[2, 1]
bpp4_1 <- bpp4_1.tmp$"posterior density exact"
post4_
# assuming interim futility boundary:
1 <- integrate(Vectorize(estimate_toIntegrate), lower = -Inf, upper = Inf, prior = "flat",
bpp5_successmean = successmean, finalSE = finalSE, interimmean = futi,
interimSE = sqrt(fac / nevents[1]), priormean = up2$postmean,
width = width1, height = height1, subdivisions = 300)$value
1.tmp <- bpp_1interim(prior = "flat", interimSE = sqrt(fac / nevents[1]),
bpp5_finalSE = finalSE, successmean = successmean,
IntEffBoundary = log(0), IntFutBoundary = effi, IntFix = futi,
priormean = up2$postmean, propA = 0.5, thetas = thetas,
width = width1, height = height1)
1 <- bpp5_1.tmp$"BPP after not stopping at interim exact"[2, 1]
bpp5_1 <- bpp5_1.tmp$"posterior density exact" post5_
The resulting posteriors after the updates with the external studies are depicted below, compare Figure 1 in Rufibach, Jordan, and Abt (2016b).
par(las = 1, mar = c(9, 5, 2, 1), mfrow = c(1, 2), cex = 0.8)
# ----------------------------------
# Normal prior:
# ----------------------------------
plot(0, 0, type = "n", xlim = c(-0.6, 0.3), ylim = c(-0.1, 5), xlab = "", ylab = "density",
main = "")
title(expression("Normal prior density and corresponding posteriors for "*theta), line = 0.7)
basicPlot(leg = FALSE, IntEffBoundary = effi, IntFutBoundary = futi, successmean = successmean,
priormean = priormean)
lines(thetas, dnorm(thetas, mean = priormean, sd = sd0), col = 2, lwd = 2)
lines(thetas, dnorm(thetas, mean = up1$postmean, sd = up1$postsigma), col = 3, lwd = 2)
lines(thetas, dnorm(thetas, mean = up2$postmean, sd = up2$postsigma), col = 4, lwd = 2)
lines(thetas, post3, col = 1, lwd = 2)
legend(-0.64, 5.2, c("prior", "posterior after Sub1", "posterior after Sub1 & Sub2",
"posterior after Sub1 & Sub2 and not stopping at interim"), lty = 1, col = c(2:4, 1),
bty = "n", lwd = 2)
# ----------------------------------
# pessimistic prior:
# ----------------------------------
# first we have to compute the posteriors after the external updates:
rep(NA, length(thetas))
flatpost1 <- flatpost1
flatpost2 <-for (i in 1:length(thetas)){
estimate_posterior(x = thetas[i], prior = "flat", interimmean = log(hr1),
flatpost1[i] <-interimSE = sd1, priormean = priormeanflat, width = width1,
height = height1)
estimate_posterior(x = thetas[i], prior = "flat", interimmean = log(hr2),
flatpost2[i] <-interimSE = sd2, priormean = priormeanflat, width = width1,
height = height1)
}
plot(0, 0, type = "n", xlim = c(-0.6, 0.3), ylim = c(-0.10, 5), xlab = "", ylab = "density", main = "")
title(expression("Flat prior density and corresponding posteriors for "*theta), line = 0.7)
basicPlot(leg = FALSE, IntEffBoundary = effi, IntFutBoundary = futi, successmean = successmean,
priormean = priormeanflat)
lines(thetas, dUniformNormalTails(thetas, mu = priormeanflat, width = width1, height = height1),
lwd = 2, col = 2)
lines(thetas, flatpost1, col = 3, lwd = 2)
lines(thetas, flatpost2, col = 4, lwd = 2)
lines(thetas, post3_1, col = 1, lwd = 2)
legend(-0.64, 5.2, c("prior", "posterior after Sub1", "posterior after Sub1 & Sub2",
"posterior after Sub1 & Sub2 and not stopping at interim"), lty = 1, col = c(2:4, 1),
bty = "n", lwd = 2)
Next, the posteriors after updating with the interim information, as in Figure 2:
par(las = 1, mar = c(9, 5, 2, 1), mfrow = c(1, 2), cex = 0.8)
# ----------------------------------
# Normal prior:
# ----------------------------------
plot(0, 0, type = "n", xlim = c(-0.6, 0.3), ylim = c(-0.1, 8), xlab = "", ylab = "density",
main = "")
title("Posteriors for after not stopping at interim, Normal prior", line = 0.7)
basicPlot(leg = FALSE, IntEffBoundary = effi, IntFutBoundary = futi, successmean = successmean,
priormean = priormean)
lines(thetas, post3, col = 1, lwd = 2)
lines(thetas, post4, col = 2, lwd = 2)
lines(thetas, post5, col = 3, lwd = 2)
c("interval knowledge", expression(hat(theta)*" = efficacy boundary"),
leg2 <-expression(hat(theta)*" = futility boundary"))
legend(-0.62, 8.2, leg2, lty = 1, col = 1:3, lwd = 2, bty = "n",
title = "posterior after not stopping at interim,")
# ----------------------------------
# pessimistic prior:
# ----------------------------------
plot(0, 0, type = "n", xlim = c(-0.6, 0.3), ylim = c(-0.10, 8), xlab = "", ylab = "density",
main = "")
title("Posteriors after not stopping at interim, pessimistic prior", line = 0.7)
basicPlot(leg = FALSE, IntEffBoundary = effi, IntFutBoundary = futi, successmean = successmean,
priormean = priormeanflat)
lines(thetas, post3_1, col = 1, lwd = 2)
lines(thetas, post4_1, col = 2, lwd = 2)
lines(thetas, post5_1, col = 3, lwd = 2)
c("interval knowledge", expression(hat(theta)*" = efficacy boundary"),
leg.flat <-expression(hat(theta)*" = futility boundary"))
legend(-0.62, 8.2, leg.flat, lty = 1, col = 1:3, lwd = 2, bty = "n",
title = "posterior after not stopping at interim,")
Finally, we collect all the results computed above in one table, reproducing Table 1 in Rufibach, Jordan, and Abt (2016b).
matrix(NA, ncol = 2, nrow = 10)
mat <-1] <- c(pnorm1, pnorm2, bpp0, bpp1, bpp2, bpp3, bpp3_futi_only, bpp3_effi_only, bpp4, bpp5)
mat[, 2] <- c(flat1, flat2, bpp0_1, bpp1_1, bpp2_1, bpp3_1, bpp3_1_futi_only, bpp3_1_effi_only,
mat[, 1, bpp5_1)
bpp4_colnames(mat) <- c("Normal prior", "Flat prior")
rownames(mat) <- c(paste("Probability for hazard ratio to be $\\le$ ", lims[1], sep = ""),
paste("Probability for hazard ratio to be $\\ge$ ", lims[2], sep = ""),
"PoS based on prior distribution", "PoS after Sub1", "PoS after Sub1 and Sub2",
"PoS after not stopping at interim, assuming $\\hat \\theta \\in [\\theta_{futi}, \\theta_{effi}]$",
"PoS after not stopping at interim, assuming $\\hat \\theta \\in [-\\infty, \\theta_{futi}]$",
"PoS after not stopping at interim, assuming $\\hat \\theta \\in [\\theta_{effi}, \\infty]$",
"PoS after not stopping at interim, assuming $\\hat \\theta = \\theta_{effi}$",
"PoS after not stopping at interim, assuming $\\hat \\theta = \\theta_{futi}$")
as.data.frame(format(mat, digits = 2))
mat <-library(knitr); kable(mat)
## Warning: Paket 'knitr' wurde unter R Version 4.1.2 erstellt
Normal prior | Flat prior | |
---|---|---|
Probability for hazard ratio to be \(\le\) 0.7 | 0.039 | 0.039 |
Probability for hazard ratio to be \(\ge\) 1 | 0.070 | 0.147 |
PoS based on prior distribution | 0.702 | 0.612 |
PoS after Sub1 | 0.740 | 0.665 |
PoS after Sub1 and Sub2 | 0.783 | 0.727 |
PoS after not stopping at interim, assuming \(\hat \theta \in [\theta_{futi}, \theta_{effi}]\) | 0.705 | 0.617 |
PoS after not stopping at interim, assuming \(\hat \theta \in [-\infty, \theta_{futi}]\) | 0.822 | 0.782 |
PoS after not stopping at interim, assuming \(\hat \theta \in [\theta_{effi}, \infty]\) | 0.653 | 0.547 |
PoS after not stopping at interim, assuming \(\hat \theta = \theta_{effi}\) | 0.997 | 0.997 |
PoS after not stopping at interim, assuming \(\hat \theta = \theta_{futi}\) | 0.024 | 0.016 |
The theory developed in Rufibach, Jordan, and Abt (2016b) can straightforwardly be extended to accomodate more than one interim analysis. Below, we illustrate the BPP update assuming a trial that did not stop after first a futility and second an efficacy interim analysis. Note that the function bpp_2interim
in Rufibach, Jordan, and Abt (2016a) so far only allows specification of a Normal prior.
# ------------------------------------------------------------------------------------------
# Illustrate the update after two passed interims using the Gallium clinical trial
# ------------------------------------------------------------------------------------------
# mean and sd of Normal prior:
0.9288563
hr0 <- log(hr0)
priormean <- sqrt(4 / 12)
priorsigma <-
# specifications for pivotal study:
0.5
propA <- (propA * (1 - propA)) ^ (-1)
fac <- c(111, 248, 370)
nevents <- sqrt(fac / nevents[1:2])
interimSE <- sqrt(fac / nevents[3])
finalSE <- c(3.9285726330559, 2.5028231888636, 1.9936294555664)
za <- 2 * (1 - pnorm(za))
alphas <- exp(- za * sqrt(fac / nevents))
hrMDD <- log(hrMDD[3])
successmean <-
# 2-d vector of efficacy and futility interim boundaries:
log(c(0, hrMDD[2])) # first interim is for futility only
effi <- log(c(1, Inf)) # second interim is for efficacy only
futi <-
bpp_2interim(prior = "normal", interimSE = interimSE, finalSE = finalSE,
bpp2 <-successmean = successmean, IntEffBoundary = effi, IntFutBoundary = futi,
priormean = priormean, thetas = thetas, priorsigma = priorsigma)
$"initial BPP" bpp2
## [1] 0.41
$"BPP after not stopping at interim interval" bpp2
## [1] 0.3222385
R version and packages used to generate this report:
R version: R version 4.1.1 (2021-08-10)
Base packages: stats / graphics / grDevices / utils / datasets / methods / base
Other packages: knitr / bpp / mvtnorm
This document was generated on 2022-01-13 at 15:00:34.
Rufibach, K., H. U. Burger, and M. Abt. 2016. “Bayesian Predictive Power: Choice of Prior and Some Recommendations for Its Use as Probability of Success in Drug Development.” Pharm. Stat. 15: 438–46.
Rufibach, K., P. Jordan, and M. Abt. 2016a. Bpp: Computations Around Bayesian Predictive Power.
———. 2016b. “Sequentially updating the likelihood of success of a Phase 3 pivotal time-to-event trial based on interim analyses or external information.” J Biopharm Stat 26 (2): 191–201.
These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.