R/pred_prob_binom.R

In michelecianfriglia/SampleSizeWass: Wasserstein-based Sample Size Calculations

#' Predictive Probability: Binomial Model
#'
#' Computes the predictive probability for the [PPC_binom()] function.
#'
#' @usage pred_prob_binom(
#'   n,
#'   alpha_1,
#'   beta_1,
#'   alpha_2,
#'   beta_2,
#'   alpha_D,
#'   beta_D,
#'   xi,
#'   q,
#'   M
#' )
#'
#' @param n The sample size.
#' @param alpha_1,beta_1 The parameters of the first beta prior.
#' @param alpha_2,beta_2 The parameters of the second beta prior.
#' @param alpha_D,beta_D The parameters of the design beta prior.
#' @param xi A constant used to compute the predictive probability.
#' @param q The order of the Wasserstein distance.
#' @param M The number of Monte Carlo replications.
#'
#' @return A numeric value.
#'
#' @note This is an internal function.
#'
#' @author Michele Cianfriglia \email{michele.cianfriglia@@uniroma1.it}
#'
#' @keywords internal
#'
#' @importFrom extraDistr rbbinom
#' @importFrom stats integrate qbeta

pred_prob_binom <- function(n, alpha_1, beta_1, alpha_2, beta_2, alpha_D, beta_D, xi, q, M) {

  integrand <- function(x, n, alpha_1, beta_1, alpha_2, beta_2, s_n, q) {
    abs(qbeta(p = x, shape1 = alpha_1 + s_n, shape2 = beta_1 + n - s_n) - qbeta(p = x, shape1 = alpha_2 + s_n, shape2 = beta_2 + n - s_n))^q
  }

  q_wass_dist <- function() {
    s_n <- rbbinom(n = 1, size = n, alpha = alpha_D, beta = beta_D)
    val <- integrate(f = integrand, lower = 0, upper = 1, n = n, alpha_1 = alpha_1, beta_1 = beta_1, alpha_2 = alpha_2, beta_2 = beta_2, s_n = s_n, q = q)$value
    return(val)
  }

  p_n <- mean(replicate(n = M, expr = q_wass_dist()) > xi)
  return(p_n)

}