# ========================== set_user_prior ===================================
#' Set a user-defined prior
#'
#' Constructs a user-defined prior distribution for use as the argument
#' \code{prior} in \code{\link{hef}} or \code{\link{hanova1}}.
#'
#' @param prior An R function returning \strong{the log of} the prior density
#' for of (perhaps a subset of) the hyperparameter vector \eqn{\phi}.
#' @param ... Further arguments giving the names and values of any
#' parameters involved in the function \code{prior}.
#' @param model A character string. Abbreviated name of the model:
#' "beta_binom" for beta-binomial and "gamma_pois" for gamma-Poisson
#' (see \code{\link{hef}}), "anova1" for 1-way ANOVA (see
#' \code{\link{hanova1}}).
#' @param anova_d An integer scalar. Only relevant if \code{model = anova1}.
#' If \code{anova_d = 2} then \code{prior} must return the log-prior
#' density for the standard deviations \eqn{(\sigma_\alpha, \sigma)}
#' and a normal prior with mean \code{mu0} and standard deviation
#' \code{sigma0} is used for \eqn{\mu}. The values of \code{mu0 = 0} and
#' \code{sigma0 = Inf} are set in the call to \code{hanova1}, with default
#' values \code{mu0 = 0} and \code{sigma0 = Inf}.
#' If \code{anova_d = 3} then \code{prior} must return the log-prior
#' density for \eqn{(\mu, \sigma_\alpha, \sigma)}.
#' @details For details of the hyperparameters in \eqn{\phi} see the
#' \strong{Details} section of \code{\link{hef}} for the models
#' \code{beta_binom} and \code{gamma_pois} and of \code{\link{hanova1}}
#' for the model \code{anova1}.
#' @return A list of class \code{"bang_prior"}. Will contain the component
#' \code{prior}, the user-supplied function to evaluate the log of the prior,
#' and any arguments supplied in ....
#' @seealso \code{\link{hef}} for hierarchical exponential family models.
#' @seealso \code{\link{hanova1}} for hierarchical one-way analysis of
#' variance (ANOVA).
#' @examples
#' # User-defined prior, passing parameters
#' # (equivalent to prior = "gamma" with hpars = c(1, 0.01, 1, 0.01))
#' user_prior <- function(x, hpars) {
#' return(dexp(x[1], hpars[1], log = TRUE) + dexp(x[2], hpars[2], log = TRUE))
#' }
#' user_prior_fn <- set_user_prior(user_prior, hpars = c(0.01, 0.01))
#' @export
set_user_prior <- function(prior, ..., model = c("beta_binom", "gamma_pois",
"anova1"), anova_d = 2) {
if (!is.function(prior)) {
stop("prior must be a function")
}
model <- match.arg(model)
# Return a list and add additional arguments from ....
if (model == "anova1") {
temp <- list(prior = prior, ..., anova_d = anova_d)
} else {
temp <- list(prior = prior, ...)
}
return(structure(temp, class = "bang_prior", model = model))
}
# ============================ check_prior ====================================
check_prior <- function(prior, model, hpars, n_groups = NULL) {
# If prior is a character scalar then a default prior is being requested
if (is.character(prior)) {
prior_name <- prior
prior <- list()
# One-way ANOVA
if (model == "anova1") {
if (prior_name %in% c("default", "bda") & n_groups < 3) {
stop("Need >= 3 groups for posterior propriety if bda prior used")
}
prior$prior <- switch(prior_name,
default = anova1_bda_prior,
bda = anova1_bda_prior,
unif = anova1_unif_prior,
cauchy = anova1_cauchy_prior)
prior$anova_d <- 2
if (prior_name == "cauchy") {
prior$hpars <- hpars
if (is.null(hpars)) {
prior$hpars <- anova1_cauchy_hpars()
} else {
if (any(hpars <= 0)) {
stop("The values in 'hpars' must be positive")
}
}
}
}
# Beta-binomial
if (model == "beta_binom") {
prior$prior <- switch(prior_name,
default = beta_bda_prior,
bda = beta_bda_prior,
gamma = beta_gamma_prior)
if (prior_name == "gamma") {
prior$hpars <- hpars
if (is.null(hpars)) {
prior$hpars <- beta_gamma_hpars()
} else {
if (any(hpars <= 0)) {
stop("The values in 'hpars' must be positive")
}
}
}
}
# Poisson-gamma
if (model == "gamma_pois") {
prior$prior <- switch(prior_name,
default = gamma_gamma_prior,
gamma = gamma_gamma_prior)
if (prior_name == "gamma" | prior_name == "default") {
prior$hpars <- hpars
if (is.null(hpars)) {
prior$hpars <- gamma_gamma_hpars()
} else {
if (any(hpars <= 0)) {
stop("The values in 'hpars' must be positive")
}
}
}
}
} else if (inherits(prior, "bang_prior")) {
if (attr(prior, "model") != model) {
stop("model and model set by set_user_prior() don't match")
}
} else {
stop("A user-defined prior must be set using set_user_prior()")
}
#
return(prior)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.