diff --git a/R/predprobDist.R b/R/predprobDist.R index 270807a6..3576bcfb 100644 --- a/R/predprobDist.R +++ b/R/predprobDist.R @@ -5,71 +5,109 @@ NULL #' Compute the predictive probability that the trial will be #' successful, with a prior distribution on the SOC #' +#' @description `r lifecycle::badge("experimental")` +#' #' Compute the predictive probability of trial success given current data. #' Success means that at the end of the trial the posterior probability -#' Pr(P_E > P_S + delta_0 | data) >= thetaT. Then the +#' `Pr(P_E > P_S + delta_0 | data) >= thetaT`. Then the #' predictive probability for success is: -#' pp = sum over i: Pr(Y=i|x,n)*I{Pr(P_E > P_S + delta|x,Y=i)>=thetaT}, -#' where Y is the number of future responses in the treatment group and x is -#' the current number of responses in the treatment group (out of n). -#' Prior is P_E ~ beta(a, b), default uniform which is a beta(1,1). +#' `pp = sum over i: Pr(Y = i | x , n)*I{Pr(P_E > P_S + delta | x,Y=i) >= thetaT}`, +#' where `Y` is the number of future responses in the treatment group and `x` is +#' the current number of responses in the treatment group out of `n`. +#' Prior is `P_E ~ beta(a, b)` and uniform which is a beta(1,1). #' However, also a beta mixture prior can be specified. Analogously -#' for P_S either a classic beta prior or a beta mixture prior can be +#' for `P_S` either a classic beta prior or a beta mixture prior can be #' specified. #' #' Also data on the SOC might be available. Then the predictive probability is #' more generally defined as -#' pp = sum over i, j: Pr(Y=i|x,n)*Pr(Z=j|xS, nS)*I{Pr(P_E > P_S + -#' delta|x,xS,Y=i,Z=j)>=thetaT} -#' where Z is the future number of responses in the SOC group, and xS is the +#' `pp = sum over i, j: Pr(Y = i | x, n)*Pr(Z = j | xS, nS )*I{Pr(P_E > P_S + delta | x,xS, Y=i, Z=j ) >= thetaT}` +#' where `Z` is the future number of responses in the SOC group, and `xS` is the #' current number of responses in the SOC group. #' +#' In the case where `NmaxControl = 0`, a table with the following contents will be included in the return output : +#' - `i`: `Y = i`, number of future successes in `Nmax-n` subjects. +#' - `density`: `Pr(Y = i|x)` using beta-(mixture)-binomial distribution. +#' - `posterior`: `Pr(P_E > P_S + delta | x, Y = i)` using beta posterior. +#' - `success`: indicator `I(b>thetaT)`. +#' #' A table with the following contents will be -#' included in the \code{tables} attribute of the return value +#' included in the return value #' (in the case that NmaxControl is zero): #' i: Y=i (number of future successes in Nmax-n subjects) #' py: Pr(Y=i|x) using beta-(mixture)-binomial distribution #' b: Pr(P_E > P_S + delta | x, Y=i) #' bgttheta: indicator I(b>thetaT) #' -#' If NmaxControl is not zero, i.e., when data on the control treatment +#' If `NmaxControl` is not zero, i.e., when data on the control treatment #' is available in this trial, then a list with will be included with the #' following elements: #' pyz: matrix with the probabilities Pr(Y=i, Z=j | x, xS) #' b: matrix with Pr(P_E > P_S + delta | x, xS, Y=i, Z=j) #' bgttheta: matrix of indicators I(b>thetaT) #' -#' @param x number of successes (in the treatment group) -#' @param n number of patients (in the treatment group) -#' @param xS number of successes in the SOC group (default: 0) -#' @param nS number of patients in the SOC group (default: 0) -##' @param Nmax maximum number of patients at the end of the trial (in the -##' treatment group) -##' @param NmaxControl maximum number of patients at the end of the trial in the -##' SOC group (default: 0) -##' @param delta margin by which the response rate in the treatment group should -##' be better than in the SOC group (default: 0) -##' @param relativeDelta see \code{\link{postprobDist}} -##' @param parE the beta parameters matrix, with K rows and 2 columns, -##' corresponding to the beta parameters of the K components. default is a -##' uniform prior. -##' @param weights the mixture weights of the beta mixture prior. Default are -##' uniform weights across mixture components. -##' @param parS beta prior parameters in the SOC group (default: uniform) -##' @param weightsS weights for the SOC group (default: uniform) -##' @param thetaT threshold on the probability to be used -##' @return The predictive probability, a numeric value. In addition, a -##' list called \code{tables} is returned as attribute of the returned number. -##' -##' @references Lee, J. J., & Liu, D. D. (2008). A predictive probability -##' design for phase II cancer clinical trials. Clinical Trials, 5(2), -##' 93-106. doi:10.1177/1740774508089279 -##' -##' @example examples/predprobDist.R -##' @export +#' @note +#' +#' # Delta from postprobDist : +#' +#' The desired improvement is denoted as `delta`. There are two options in using `delta`. +#' The absolute case when `relativeDelta = FALSE` and relative as when `relativeDelta = TRUE`. +#' +#' 1. The absolute case is when we define an absolute delta, greater than `P_S`, +#' the response rate of the `SOC` group such that +#' the posterior is `Pr(P_E > P_S + delta | data)`. +#' +#' 2. In the relative case, we suppose that the treatment group's +#' response rate is assumed to be greater than `P_S + (1-P_S) * delta` such that +#' the posterior is `Pr(P_E > P_S + (1 - P_S) * delta | data)`. +#' +#' @typed x : number +#' number of successes in the treatment group at interim +#' @typed n : number +#' number of patients in the treatment group at interim +#' @typed xS : number +#' number of successes in the SOC group at interim +#' @typed nS : number +#' number of patients in the SOC group +#' @typed Nmax : number +#' maximum number of patients in the treatment group at the end of the trial +#' @typed NmaxControl : number +#' maximum number of patients at the end of the trial in the +#' SOC group +#' @typed delta : +#' margin by which the response rate in the treatment group should +#' be better than in the SOC group +#' @typed relativeDelta : +#' see `[postprobDist()]` +#' @typed parE : +#' the beta parameters matrix, with K rows and 2 columns, +#' corresponding to the beta parameters of the K components. default is a +#' uniform prior. +#' @typed weights : +#' the mixture weights of the beta mixture prior. Default are +#' uniform weights across mixture components. +#' @typed parS : +#' beta prior parameters in the SOC group +#' @typed weightsS : +#' weights for the SOC group +#' @typed thetaT : +#' threshold on the probability to be used +#' @return A `list` is returned with names `result` for predictive probability and +#' `table` of numeric values with counts of responses in the remaining patients, +#' probabilities of these counts, corresponding probabilities to be above threshold, +#' and trial success indicators. +#' +#' @references Lee, J. J., & Liu, D. D. (2008). A predictive probability +#' design for phase II cancer clinical trials. Clinical Trials, 5(2), +#' 93-106. doi:10.1177/1740774508089279 +#' +#' @example examples/predprobDist.R +#' @export predprobDist <- function(x, n, - xS = 0, nS = 0, - Nmax, NmaxControl = 0, + xS = 0, + nS = 0, + Nmax, + NmaxControl = 0, delta = 0, relativeDelta = FALSE, parE = c(a = 1, b = 1), @@ -77,59 +115,49 @@ predprobDist <- function(x, n, parS = c(a = 1, b = 1), weightsS, thetaT) { - ## ensure reasonable numbers + # ensure reasonable numbers stopifnot( n <= Nmax, nS <= NmaxControl, x <= n, xS <= nS ) - - ## remaining active patients to be seen: + # remaining active patients to be seen: mE <- Nmax - n - - ## if par is a vector => situation where there is only one component + # if par is a vector => situation where there is only one component if (is.vector(parE)) { - ## check that it has exactly two entries + # check that it has exactly two entries stopifnot(identical(length(parE), 2L)) - - ## and transpose to matrix with one row + # and transpose to matrix with one row parE <- t(parE) } - - ## if prior weights of the beta mixture are not supplied + # if prior weights of the beta mixture are not supplied if (missing(weights)) { weights <- rep(1, nrow(parE)) - ## (don't need to be normalized, this is done in h_getBetamixPost) + # (don't need to be normalized, this is done in h_getBetamixPost) } - - ## if parS is a vector => situation where there is only one component + # if parS is a vector => situation where there is only one component if (is.vector(parS)) { - ## check that it has exactly two entries + # check that it has exactly two entries stopifnot(identical(length(parS), 2L)) - - ## and transpose to matrix with one row + # and transpose to matrix with one row parS <- t(parS) } - - ## if prior weights of the beta mixture are not supplied + # if prior weights of the beta mixture are not supplied if (missing(weightsS)) { weightsS <- rep(1, nrow(parS)) } - - ## now compute updated parameters for beta mixture distribution on the - ## treatment proportion + # now compute updated parameters for beta mixture distribution on the + # treatment proportion activeBetamixPost <- h_getBetamixPost(x = x, n = n, par = parE, weights = weights) - - ## now with the beta binomial mixture: + # now with the beta binomial mixture: py <- with( activeBetamixPost, dbetabinomMix(x = 0:mE, m = mE, par = par, weights = weights) ) - if (NmaxControl == 0) { - ## here is the only difference to predprob.R: - ## how to compute the posterior probabilities + # here is the only difference to predprob.R: + # how to compute the posterior probabilities b <- postprobDist( x = x + c(0:mE), n = Nmax, @@ -140,40 +168,29 @@ predprobDist <- function(x, n, parS = parS, weightsS = weightsS ) - - ret <- structure(sum(py * (b > thetaT)), - tables = - round( - cbind( - i = c(0:mE), - py, - b, - bgttheta = (b > thetaT) - ), - 4 - ) + ret <- list( + result = sum(py * (b > thetaT)), + table = data.frame( + counts = c(0:mE), + # cumul_counts = xS + (0:NmaxControl), + density = py, + posterior = b, + success = (b > thetaT) + ) ) } else { - ## in this case also data on the SOC is available! - - ## determine remaining sample size and probabilities of response - ## counts in future SOC patients: + # in this case also data on the SOC is available! + # determine remaining sample size and probabilities of response + # counts in future SOC patients: mS <- NmaxControl - nS - - controlBetamixPost <- h_getBetamixPost( - x = xS, n = nS, par = parS, - weights = weightsS - ) - + controlBetamixPost <- h_getBetamixPost(x = xS, n = nS, par = parS, weights = weightsS) pz <- with( controlBetamixPost, dbetabinomMix(x = 0:mS, m = mS, par = par, weights = weights) ) - - ## determine resulting posterior probabilities: + # determine resulting posterior probabilities: outcomesY <- x + c(0:mE) - outcomesZ <- xS + c(0:mS) - + outcomesZ <- xS + c(0:mS) # 15 more chances to get success counts pyz <- b <- matrix( nrow = 1 + mE, ncol = 1 + mS, @@ -183,11 +200,10 @@ predprobDist <- function(x, n, 0:mS ) ) - - for (i in seq_along(outcomesY)) { + for (i in seq_along(outcomesY)) { # outside? for (j in seq_along(outcomesZ)) { - ## calculate the posterior probability for this combination - ## of counts + # calculate the posterior probability for this combination + # of counts b[i, j] <- postprobDist( x = outcomesY[i], @@ -201,25 +217,22 @@ predprobDist <- function(x, n, parS = parS, weightsS = weightsS ) - - ## what are the joint probabilities of active and control counts? - ## => because they are independent, just multiply them - pyz[i, j] <- py[i] * pz[j] + # what are the joint probabilities of active and control counts? + # => because they are independent, just multiply them + pyz[i, j] <- py[i] * pz[j] # this is the matrix that Daniel was talking about } } - - ## should we print something? + # should we print something? predprob part ret <- structure(sum(pyz * (b > thetaT)), tables = list( pyz = pyz, b = b, - bgttheta = (b > thetaT) + success = (b > thetaT) ) ) } - - return(ret) + ret } -## todo: predprobDistFail +# todo: predprobDistFail diff --git a/design/design-doc_predprobDist.Rmd b/design/design-doc_predprobDist.Rmd new file mode 100644 index 00000000..9bb21c54 --- /dev/null +++ b/design/design-doc_predprobDist.Rmd @@ -0,0 +1,241 @@ +--- +title: "Design Document" +output: html_document +date: "2023-12-04" +--- + +## Example inputs / Sanity check + +```{r setup, include=FALSE} +x <- 16 +n <- 23 +xS <- 5 +nS <- 10 +Nmax <- 40 +NmaxControl <- 20 +delta <- 0 +thetaT <- 0.6 +parE <- c(1, 1) +parS <- c(1, 1) +weights <- c(2, 1) +weightsS <- c(2, 1) +relativeDelta <- FALSE + +# helper test +mE <- Nmax - n +parE <- t(parE) +weights <- rep(1, nrow(parE)) +parS <- t(parS) +weightsS <- rep(1, nrow(parS)) + +activeBetamixPost <- h_getBetamixPost(x = x, n = n, par = parE, weights = weights) +# now with the beta binomial mixture: +density_y <- with(activeBetamixPost, dbetabinomMix(x = 0:mE, m = mE, par = par, weights = weights)) +``` + +# Sanity Check +```{r} +predprobDist( + x = 16, n = 23, Nmax = 40, delta = 0, thetaT = 0.9, + parE = c(0.6, 0.4), parS = c(7, 11) +) +``` + +# Helper function if NmaxControl == 0, or when there is no control +```{r} +h_predprobdist_single_arm <- function(x, + n, + delta, + relativeDelta, + parE, + weights, + parS, + weightsS, + thetaT, + density, + mE) { + posterior_x <- postprobDist( + x = x + c(0:mE), + n = Nmax, + delta = delta, + relativeDelta = relativeDelta, + parE = parE, + weights = weights, + parS = parS, + weightsS = weightsS + ) + list( + result = sum(density * (posterior_x > thetaT)), + table = data.frame( + counts = c(0:mE), + cumul_counts = x + (0:mE), + density = density_y, + posterior = posterior_x, + success = (posterior_x > thetaT) + ) + ) +} +``` + +# Testing helper function if NmaxControl == 0, or when there is no control +```{r} +# testing +h_get_predproblist_nocontrol( + x = x, + n = Nmax, + delta = delta, + relativeDelta = relativeDelta, + parE = parE, + weights = weights, + parS = parS, + weightsS = weightsS, + thetaT = thetaT, + density = density_y, + mE = 17 +) +``` + +# Helper function 1 if NmaxControl != 0 +```{r} +h_getlistpredprob <- function(NmaxControl, nS, xS, par, weightsS, mS, x, mE, density_y) { + mS <- NmaxControl - nS + controlBetamixPost <- h_getBetamixPost( + x = xS, + n = nS, + par = parS, + weights = weightsS + ) + density_z <- with( + controlBetamixPost, + dbetabinomMix(x = 0:mS, m = mS, par = par, weights = weights) + ) + ## determine resulting posterior probabilities: + outcomesY <- x + c(0:mE) + outcomesZ <- xS + c(0:mS) + density_yz <- posterior_yz <- matrix( + nrow = 1 + mE, + ncol = 1 + mS, + dimnames = + list( + 0:mE, + 0:mS + ) + ) + for (i in seq_along(outcomesY)) { + for (j in seq_along(outcomesZ)) { + posterior_yz[i, j] <- + postprobDist( + x = outcomesY[i], + n = Nmax, + xS = outcomesZ[j], + nS = NmaxControl, + delta = delta, + relativeDelta = relativeDelta, + parE = parE, + weights = weights, + parS = parS, + weightsS = weightsS + ) + density_yz[i, j] <- density_y[i] * density_z[j] + } + } + ret <- list( + result = sum(density_yz * (posterior_yz > thetaT)), + table = data.frame( + counts = c(0:mE), + cumul_counts = x + (0:mE), + density = density_yz, + posterior = posterior_yz, + success = (posterior_yz > thetaT) + ) + ) + ret +} +``` + +# Testing helper function if NmaxControl != 0 +```{r} +h_getlistpredprob( + NmaxControl = 20, + nS = 10, + xS = 5, + weightsS = 1, + par = c(1, 1), + x = 10, + mE = 17, + density_y = density_y +) +``` + +# Main function +```{r cars} +predprobDist <- function(x, n, + xS = 0, nS = 0, + Nmax, + NmaxControl = 0, + delta = 0, + relativeDelta = FALSE, + parE = c(a = 1, b = 1), + weights, + parS = c(a = 1, b = 1), + weightsS, + thetaT) { + stopifnot( + n <= Nmax, + nS <= NmaxControl, + x <= n, + xS <= nS + ) + mE <- Nmax - n + if (is.vector(parE)) { + stopifnot(identical(length(parE), 2L)) + parE <- t(parE) + } + if (missing(weights)) { + weights <- rep(1, nrow(parE)) + } + if (is.vector(parS)) { + stopifnot(identical(length(parS), 2L)) + parS <- t(parS) + } + if (missing(weightsS)) { + weightsS <- rep(1, nrow(parS)) + } + activeBetamixPost <- h_getBetamixPost(x = x, n = n, par = parE, weights = weights) + density_y <- with( + activeBetamixPost, + dbetabinomMix(x = 0:mE, m = mE, par = par, weights = weights) + ) + if (NmaxControl == 0) { + ret <- h_get_predproblist_nocontrol( + x = x, + n = Nmax, + delta = delta, + relativeDelta = relativeDelta, + parE = parE, + weights = weights, + parS = parS, + weightsS = weightsS, + thetaT = thetaT, + density = density_y, + mE = mE + ) + } else { + ret <- h_getlistpredprob( + NmaxControl = nS, + nS = nS, xS = xS, par = parS, weightsS = weightsS, + mS = mS, x = x, mE = mE, density_y = density_y + ) + } + ret +} +``` + +# An example +```{r pressure, echo=FALSE} +predprobDist( + x = 16, n = 23, Nmax = 40, delta = 0.1, thetaT = 0.9, + parE = rbind(c(1, 1), c(50, 10)), + weights = c(2, 1), parS = c(6, 11) +) +``` diff --git a/examples/predprobDist.R b/examples/predprobDist.R index 44b35769..90ce5ca9 100644 --- a/examples/predprobDist.R +++ b/examples/predprobDist.R @@ -1,53 +1,53 @@ -## The original Lee and Liu (Table 1) example: -## Nmax=40, x=16, n=23, beta(0.6,0.4) prior distribution, -## thetaT=0.9. The control response rate is 60%: +# The original Lee and Liu (Table 1) example: +# Nmax=40, x=16, n=23, beta(0.6,0.4) prior distribution, +# thetaT=0.9. The control response rate is 60%: predprob( x = 16, n = 23, Nmax = 40, p = 0.6, thetaT = 0.9, parE = c(0.6, 0.4) ) -## So the predictive probability is 56.6%. 12 responses -## in the remaining 17 patients are required for success. +# So the predictive probability is 56.6%. 12 responses +# in the remaining 17 patients are required for success. -## Now to the modified example, where the control response -## rate p has a beta(7, 11) distribution. -## For example if the 60% threshold stems from a trial -## with 10 patients, of which 6 were observed with response. +# Now to the modified example, where the control response +# rate p has a beta(7, 11) distribution. +# For example if the 60% threshold (control) stems from a trial +# with 10 patients, of which 6 were observed with response. predprobDist( x = 16, n = 23, Nmax = 40, delta = 0, thetaT = 0.9, parE = c(0.6, 0.4), parS = c(7, 11) ) -## Now only 11 responses are needed in the remaining 17 patients -## in order to have a successful trial. -## Also, the predictive probability for a Go decision -## is now 70.8% instead of merely 56.6%. +# Now only7 responses are needed in the remaining 17 patients +# in order to have a successful trial. +# Also, the predictive probability for a Go decision +# is now 97.78% instead of merely 56.6%. -## If the response threshold of 60% stems from an SOC estimate (e.g. 50% = 5/10) -## and a margin of e.g. 10%, then the margin can be specified with -## the delta argument: +# If the response threshold of 60% stems from an SOC estimate (e.g. 50% = 5/10) +# and a margin of e.g. 10%, then the margin can be specified with +# the delta argument: predprobDist( x = 16, n = 23, Nmax = 40, delta = 0.1, thetaT = 0.9, parE = c(1, 1), parS = c(6, 11) ) -## This again changes the result: Only 10 future responses are required, -## the predictive probability of success is now 80%. +# This again changes the result: Only 10 future responses are required, +# the predictive probability of success is now 80%. -## Next we will use a beta mixture prior, with a weight of 1/3 for an -## informative beta(50, 10) prior, on the experimental response rate: +# Next we will use a beta mixture prior, with a weight of 1/3 for an +# informative beta(50, 10) prior, on the experimental response rate: predprobDist( x = 16, n = 23, Nmax = 40, delta = 0.1, thetaT = 0.9, parE = rbind(c(1, 1), c(50, 10)), weights = c(2, 1), parS = c(6, 11) ) -## Now still 10 future responses are required for success, but the predictive -## success probability is higher (87.2%). +# Now still 10 future responses are required for success, but the predictive +# success probability is higher (87.2%). -## We can also have additional controls in our trial. For example -## assume that in the interim analysis we had 10 control patients observed -## with 5 responses out of 20 in total at the end. -## We also had historical control data (say 60 patients of -## which 20 responded), but in order to be robust against changes in the response +# We can also have additional controls in our trial. For example +# assume that in the interim analysis we had 10 control patients observed +# with 5 responses out of 20 in total at the end. +# We also had historical control data (say 60 patients of +# which 20 responded), but in order to be robust against changes in the response # rate over time, we will again use a beta mixture, putting only 1/3 of weight -## on the historical controls: +# on the historical controls: predprobDist( x = 16, n = 23, xS = 5, nS = 10, Nmax = 40, NmaxControl = 20, @@ -57,8 +57,8 @@ predprobDist( parS = rbind(c(1, 1), c(20, 40)), weightsS = c(2, 1) ) -## Now we only have 59.9% predictive probability of success. The bgttheta matrix -## lists the cases where we would have success, for all possible combinations -## of future responses in the control and experimental arms. We see that the -## success threshold for experimental responses depends on the number of -## control responses. +# Now we only have 59.9% predictive probability of success. The bgttheta matrix +# lists the cases where we would have success, for all possible combinations +# of future responses in the control and experimental arms. We see that the +# success threshold for experimental responses depends on the number of +# control responses. diff --git a/inst/WORDLIST b/inst/WORDLIST index 1db7dca6..b5d955a2 100644 --- a/inst/WORDLIST +++ b/inst/WORDLIST @@ -125,6 +125,7 @@ PointMass pointmass postL postprob +postprobDist postU ppL PProbust diff --git a/man/predprobDist.Rd b/man/predprobDist.Rd index 9dffe3f2..6389a81e 100644 --- a/man/predprobDist.Rd +++ b/man/predprobDist.Rd @@ -22,129 +22,152 @@ predprobDist( ) } \arguments{ -\item{x}{number of successes (in the treatment group)} +\item{x}{(\code{number}):\cr number of successes in the treatment group at interm} -\item{n}{number of patients (in the treatment group)} +\item{n}{(\code{number}):\cr number of patients in the treatment group at interim} -\item{xS}{number of successes in the SOC group (default: 0)} +\item{xS}{(\code{number}):\cr number of successes in the SOC group at interim} -\item{nS}{number of patients in the SOC group (default: 0)} +\item{nS}{(\code{number}):\cr number of patients in the SOC group} -\item{Nmax}{maximum number of patients at the end of the trial (in the -treatment group)} +\item{Nmax}{(\code{number}):\cr maximum number of patients in the treatment group at the end of the trial} -\item{NmaxControl}{maximum number of patients at the end of the trial in the -SOC group (default: 0)} +\item{NmaxControl}{(\code{number}):\cr maximum number of patients at the end of the trial in the +SOC group} -\item{delta}{margin by which the response rate in the treatment group should -be better than in the SOC group (default: 0)} +\item{delta}{(``):\cr margin by which the response rate in the treatment group should +be better than in the SOC group} -\item{relativeDelta}{see \code{\link{postprobDist}}} +\item{relativeDelta}{(``):\cr see \verb{[postprobDist()]}} -\item{parE}{the beta parameters matrix, with K rows and 2 columns, +\item{parE}{(``):\cr the beta parameters matrix, with K rows and 2 columns, corresponding to the beta parameters of the K components. default is a uniform prior.} -\item{weights}{the mixture weights of the beta mixture prior. Default are +\item{weights}{(``):\cr the mixture weights of the beta mixture prior. Default are uniform weights across mixture components.} -\item{parS}{beta prior parameters in the SOC group (default: uniform)} +\item{parS}{(``):\cr beta prior parameters in the SOC group} -\item{weightsS}{weights for the SOC group (default: uniform)} +\item{weightsS}{(``):\cr weights for the SOC group} -\item{thetaT}{threshold on the probability to be used} +\item{thetaT}{(``):\cr threshold on the probability to be used} } \value{ -The predictive probability, a numeric value. In addition, a -list called \code{tables} is returned as attribute of the returned number. +A \code{list} is returned with names \code{result} for predictive probability and +\code{table} of numeric values with counts of responses in the remaining patients, +probabilities of these counts, corresponding probabilities to be above threshold, +and trial success indicators. } \description{ +\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}} + Compute the predictive probability of trial success given current data. Success means that at the end of the trial the posterior probability -Pr(P_E > P_S + delta_0 | data) >= thetaT. Then the +\code{Pr(P_E > P_S + delta_0 | data) >= thetaT}. Then the predictive probability for success is: -pp = sum over i: Pr(Y=i|x,n)*I{Pr(P_E > P_S + delta|x,Y=i)>=thetaT}, -where Y is the number of future responses in the treatment group and x is -the current number of responses in the treatment group (out of n). -Prior is P_E ~ beta(a, b), default uniform which is a beta(1,1). +\verb{pp = sum over i: Pr(Y = i | x , n)*I\{Pr(P_E > P_S + delta | x,Y=i) >= thetaT\}}, +where \code{Y} is the number of future responses in the treatment group and \code{x} is +the current number of responses in the treatment group out of \code{n}. +Prior is \code{P_E ~ beta(a, b)} and uniform which is a beta(1,1). However, also a beta mixture prior can be specified. Analogously -for P_S either a classic beta prior or a beta mixture prior can be +for \code{P_S} either a classic beta prior or a beta mixture prior can be specified. -} -\details{ + Also data on the SOC might be available. Then the predictive probability is more generally defined as -pp = sum over i, j: Pr(Y=i|x,n)*Pr(Z=j|xS, nS)*I{Pr(P_E > P_S + -delta|x,xS,Y=i,Z=j)>=thetaT} -where Z is the future number of responses in the SOC group, and xS is the +\verb{pp = sum over i, j: Pr(Y = i | x, n)*Pr(Z = j | xS, nS )*I\{Pr(P_E > P_S + delta | x,xS, Y=i, Z=j ) >= thetaT\}} +where \code{Z} is the future number of responses in the SOC group, and \code{xS} is the current number of responses in the SOC group. +In the case where \code{NmaxControl = 0}, a table with the following contents will be included in the return output : +\itemize{ +\item \code{i}: \code{Y = i}, number of future successes in \code{Nmax-n} subjects. +\item \code{density}: \code{Pr(Y = i|x)} using beta-(mixture)-binomial distribution. +\item \code{posterior}: \code{Pr(P_E > P_S + delta | x, Y = i)} using beta posterior. +\item \code{success}: indicator \code{I(b>thetaT)}. +} + A table with the following contents will be -included in the \code{tables} attribute of the return value +included in the return value (in the case that NmaxControl is zero): i: Y=i (number of future successes in Nmax-n subjects) py: Pr(Y=i|x) using beta-(mixture)-binomial distribution b: Pr(P_E > P_S + delta | x, Y=i) bgttheta: indicator I(b>thetaT) -If NmaxControl is not zero, i.e., when data on the control treatment +If \code{NmaxControl} is not zero, i.e., when data on the control treatment is available in this trial, then a list with will be included with the following elements: pyz: matrix with the probabilities Pr(Y=i, Z=j | x, xS) b: matrix with Pr(P_E > P_S + delta | x, xS, Y=i, Z=j) bgttheta: matrix of indicators I(b>thetaT) } +\note{ +Delta from postprobDist : + +The desired improvement is denoted as \code{delta}. There are two options in using \code{delta}. +The absolute case when \code{relativeDelta = FALSE} and relative as when \code{relativeDelta = TRUE}. +\enumerate{ +\item The absolute case is when we define an absolute delta, greater than \code{P_S}, +the response rate of the \code{SOC} group such that +the posterior is \code{Pr(P_E > P_S + delta | data)}. +\item In the relative case, we suppose that the treatment group's +response rate is assumed to be greater than \code{P_S + (1-P_S) * delta} such that +the posterior is \code{Pr(P_E > P_S + (1 - P_S) * delta | data)}. +} +} \examples{ -## The original Lee and Liu (Table 1) example: -## Nmax=40, x=16, n=23, beta(0.6,0.4) prior distribution, -## thetaT=0.9. The control response rate is 60\%: +# The original Lee and Liu (Table 1) example: +# Nmax=40, x=16, n=23, beta(0.6,0.4) prior distribution, +# thetaT=0.9. The control response rate is 60\%: predprob( x = 16, n = 23, Nmax = 40, p = 0.6, thetaT = 0.9, parE = c(0.6, 0.4) ) -## So the predictive probability is 56.6\%. 12 responses -## in the remaining 17 patients are required for success. +# So the predictive probability is 56.6\%. 12 responses +# in the remaining 17 patients are required for success. -## Now to the modified example, where the control response -## rate p has a beta(7, 11) distribution. -## For example if the 60\% threshold stems from a trial -## with 10 patients, of which 6 were observed with response. +# Now to the modified example, where the control response +# rate p has a beta(7, 11) distribution. +# For example if the 60\% threshold (control) stems from a trial +# with 10 patients, of which 6 were observed with response. predprobDist( x = 16, n = 23, Nmax = 40, delta = 0, thetaT = 0.9, parE = c(0.6, 0.4), parS = c(7, 11) ) -## Now only 11 responses are needed in the remaining 17 patients -## in order to have a successful trial. -## Also, the predictive probability for a Go decision -## is now 70.8\% instead of merely 56.6\%. - -## If the response threshold of 60\% stems from an SOC estimate (e.g. 50\% = 5/10) -## and a margin of e.g. 10\%, then the margin can be specified with -## the delta argument: +# Now only7 responses are needed in the remaining 17 patients +# in order to have a successful trial. +# Also, the predictive probability for a Go decision +# is now 97.78\% instead of merely 56.6\%. + +# If the response threshold of 60\% stems from an SOC estimate (e.g. 50\% = 5/10) +# and a margin of e.g. 10\%, then the margin can be specified with +# the delta argument: predprobDist( x = 16, n = 23, Nmax = 40, delta = 0.1, thetaT = 0.9, parE = c(1, 1), parS = c(6, 11) ) -## This again changes the result: Only 10 future responses are required, -## the predictive probability of success is now 80\%. +# This again changes the result: Only 10 future responses are required, +# the predictive probability of success is now 80\%. -## Next we will use a beta mixture prior, with a weight of 1/3 for an -## informative beta(50, 10) prior, on the experimental response rate: +# Next we will use a beta mixture prior, with a weight of 1/3 for an +# informative beta(50, 10) prior, on the experimental response rate: predprobDist( x = 16, n = 23, Nmax = 40, delta = 0.1, thetaT = 0.9, parE = rbind(c(1, 1), c(50, 10)), weights = c(2, 1), parS = c(6, 11) ) -## Now still 10 future responses are required for success, but the predictive -## success probability is higher (87.2\%). - -## We can also have additional controls in our trial. For example -## assume that in the interim analysis we had 10 control patients observed -## with 5 responses out of 20 in total at the end. -## We also had historical control data (say 60 patients of -## which 20 responded), but in order to be robust against changes in the response +# Now still 10 future responses are required for success, but the predictive +# success probability is higher (87.2\%). + +# We can also have additional controls in our trial. For example +# assume that in the interim analysis we had 10 control patients observed +# with 5 responses out of 20 in total at the end. +# We also had historical control data (say 60 patients of +# which 20 responded), but in order to be robust against changes in the response # rate over time, we will again use a beta mixture, putting only 1/3 of weight -## on the historical controls: +# on the historical controls: predprobDist( x = 16, n = 23, xS = 5, nS = 10, Nmax = 40, NmaxControl = 20, @@ -154,11 +177,11 @@ predprobDist( parS = rbind(c(1, 1), c(20, 40)), weightsS = c(2, 1) ) -## Now we only have 59.9\% predictive probability of success. The bgttheta matrix -## lists the cases where we would have success, for all possible combinations -## of future responses in the control and experimental arms. We see that the -## success threshold for experimental responses depends on the number of -## control responses. +# Now we only have 59.9\% predictive probability of success. The bgttheta matrix +# lists the cases where we would have success, for all possible combinations +# of future responses in the control and experimental arms. We see that the +# success threshold for experimental responses depends on the number of +# control responses. } \references{ Lee, J. J., & Liu, D. D. (2008). A predictive probability diff --git a/tests/test-predprob.R b/tests/testthat/test-predprob.R similarity index 100% rename from tests/test-predprob.R rename to tests/testthat/test-predprob.R