Skip to content

Commit

Permalink
si, bf, describe_posterior
Browse files Browse the repository at this point in the history
[skip]
  • Loading branch information
mattansb committed Sep 5, 2024
1 parent e49b705 commit 4b3e792
Show file tree
Hide file tree
Showing 8 changed files with 231 additions and 90 deletions.
54 changes: 34 additions & 20 deletions R/bayesfactor_parameters.R
Original file line number Diff line number Diff line change
Expand Up @@ -85,26 +85,22 @@
#' (Note that by default, `brms::brm()` uses flat priors for fixed-effects;
#' See example below.)
#' \cr\cr
#' It is important to provide the correct `prior` for meaningful results.
#' It is important to provide the correct `prior` for meaningful results,
#' to match the `posterior`-type input:
#'
#' - When `posterior` is a numerical vector, `prior` should also be a numerical vector.
#' - When `posterior` is a `data.frame`, `prior` should also be a `data.frame`, with matching column order.
#' - When `posterior` is a `stanreg`, `brmsfit` or other supported Bayesian model:
#' - `prior` can be set to `NULL`, in which case prior samples are drawn internally.
#' - `prior` can also be a model equivalent to `posterior` but with samples from
#' the priors *only*. See [unupdate()].
#' - **Note:** When `posterior` is a `brmsfit_multiple` model, `prior` **must** be provided.
#' - When `posterior` is an output from a `{marginaleffects}` function, `prior` should also be an an output
#' from a `{marginaleffects}` function equivalent to `posterior` but created
#' with a model of priors samples *only*.
#' - When `posterior` is an `emmGrid` / `emm_list` object:
#' - `prior` should also be an `emmGrid` / `emm_list` object equivalent to `posterior` but
#' created with a model of priors samples *only*. See [unupdate()].
#' - `prior` can also be the original (posterior) *model*. If so, the function will try to
#' update the `emmGrid` / `emm_list` to use the [unupdate()]d prior-model.
#' (*This cannot be done for `brmsfit` models.*)
#' - **Note**: When the `emmGrid` has undergone any transformations (`"log"`, `"response"`, etc.),
#' or `regrid`ing, then `prior` must be an `emmGrid` object, as stated above.
#' - **A numeric vector** - `prior` should also be a _numeric vector_, representing the prior-estimate.
#' - **A data frame** - `prior` should also be a _data frame_, representing the prior-estimates, in matching column order.

Check warning on line 92 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_parameters.R,line=92,col=121,[line_length_linter] Lines should not be more than 120 characters. This line is 122 characters.

Check warning on line 92 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_parameters.R,line=92,col=121,[line_length_linter] Lines should not be more than 120 characters. This line is 122 characters.
#' - If `rvar_col` is specified, `prior` should be _the name of an `rvar` column_ that represents the prior-estimates.
#' - **Supported Bayesian model (`stanreg`, `brmsfit`, etc.)**
#' - `prior` should be _a model an equivalent model with MCMC samples from the priors *only*_. See [unupdate()].
#' - If `prior` is set to `NULL`, [unupdate()] is called internally (not supported for `brmsfit_multiple` model).
#' - **Output from a `{marginaleffects}` function** - `prior` should also be _an equivalent output_ from a `{marginaleffects}` function based on a prior-model

Check warning on line 97 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_parameters.R,line=97,col=121,[line_length_linter] Lines should not be more than 120 characters. This line is 158 characters.

Check warning on line 97 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_parameters.R,line=97,col=121,[line_length_linter] Lines should not be more than 120 characters. This line is 158 characters.
#' (See [unupdate()]).
#' - **Output from an `{emmeans}` function**
#' - `prior` should also be _an equivalent output_ from an `{emmeans}` function based on a prior-model (See [unupdate()]).

Check warning on line 100 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_parameters.R,line=100,col=121,[line_length_linter] Lines should not be more than 120 characters. This line is 124 characters.

Check warning on line 100 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_parameters.R,line=100,col=121,[line_length_linter] Lines should not be more than 120 characters. This line is 124 characters.
#' - `prior` can also be _the original (posterior) model_, in which case the function
#' will try to "unupdate" the estimates (not supported if the estimates have undergone
#' any transformations -- `"log"`, `"response"`, etc. -- or any `regrid`ing).
#'
#' @section Interpreting Bayes Factors:
#' A Bayes factor greater than 1 can be interpreted as evidence against the
Expand Down Expand Up @@ -406,13 +402,31 @@ bayesfactor_parameters.comparisons <- bayesfactor_parameters.emmGrid


#' @rdname bayesfactor_parameters
#' @inheritParams p_direction
#' @export
bayesfactor_parameters.data.frame <- function(posterior,
prior = NULL,
direction = "two-sided",
null = 0,
rvar_col = NULL,
verbose = TRUE,
...) {
if (length(x_rvar <- .possibly_extract_rvar_col(posterior, rvar_col)) > 0L) {

Check warning on line 414 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_parameters.R,line=414,col=14,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.

Check warning on line 414 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_parameters.R,line=414,col=14,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.
cl <- match.call()
cl[[1]] <- bayestestR::bayesfactor_parameters
cl$posterior <- x_rvar
cl$rvar_col <- NULL
if (length(prior_rvar <- .possibly_extract_rvar_col(posterior, prior)) > 0L) {

Check warning on line 419 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_parameters.R,line=419,col=16,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.

Check warning on line 419 in R/bayesfactor_parameters.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_parameters.R,line=419,col=16,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.
cl$prior <- prior_rvar
}
out <- eval.parent(cl)

obj_name <- insight::safe_deparse_symbol(substitute(posterior))
attr(out, "object_name") <- sprintf('%s[["%s"]]', obj_name, rvar_col)

return(.append_datagrid(out, posterior))
}

# find direction
direction <- .get_direction(direction)

Expand Down Expand Up @@ -473,7 +487,7 @@ bayesfactor_parameters.draws <- function(posterior,
...) {
bayesfactor_parameters(
.posterior_draws_to_df(posterior),
prior = prior,
prior = if (!is.null(prior)) .posterior_draws_to_df(prior),
direction = direction,
null = null,
verbose = verbose,
Expand Down
21 changes: 19 additions & 2 deletions R/bayesfactor_restricted.R
Original file line number Diff line number Diff line change
Expand Up @@ -195,7 +195,21 @@ bayesfactor_restricted.predictions <- bayesfactor_restricted.emmGrid
bayesfactor_restricted.comparisons <- bayesfactor_restricted.emmGrid

#' @export
bayesfactor_restricted.data.frame <- function(posterior, hypothesis, prior = NULL, ...) {
#' @rdname bayesfactor_restricted
#' @inheritParams p_direction
bayesfactor_restricted.data.frame <- function(posterior, hypothesis, prior = NULL, rvar_col = NULL, ...) {
if (length(x_rvar <- .possibly_extract_rvar_col(posterior, rvar_col)) > 0L) {

Check warning on line 201 in R/bayesfactor_restricted.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_restricted.R,line=201,col=14,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.

Check warning on line 201 in R/bayesfactor_restricted.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_restricted.R,line=201,col=14,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.
cl <- match.call()
cl[[1]] <- bayestestR::bayesfactor_restricted
cl$posterior <- x_rvar
cl$rvar_col <- NULL
if (length(prior_rvar <- .possibly_extract_rvar_col(posterior, prior)) > 0L) {

Check warning on line 206 in R/bayesfactor_restricted.R

View workflow job for this annotation

GitHub Actions / lint-changed-files / lint-changed-files

file=R/bayesfactor_restricted.R,line=206,col=16,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.

Check warning on line 206 in R/bayesfactor_restricted.R

View workflow job for this annotation

GitHub Actions / lint / lint

file=R/bayesfactor_restricted.R,line=206,col=16,[implicit_assignment_linter] Avoid implicit assignments in function calls. For example, instead of `if (x <- 1L) { ... }`, write `x <- 1L; if (x) { ... }`.
cl$prior <- prior_rvar
}
return(eval.parent(cl))
}


p_hypothesis <- parse(text = hypothesis)

if (is.null(prior)) {
Expand Down Expand Up @@ -251,7 +265,10 @@ bayesfactor_restricted.data.frame <- function(posterior, hypothesis, prior = NUL

#' @export
bayesfactor_restricted.draws <- function(posterior, hypothesis, prior = NULL, ...) {
bayesfactor_restricted(.posterior_draws_to_df(posterior), hypothesis = hypothesis, prior = prior, ...)
bayesfactor_restricted(.posterior_draws_to_df(posterior),
hypothesis = hypothesis,
prior = if (!is.null(prior)) .posterior_draws_to_df(prior),
...)
}

#' @export
Expand Down
53 changes: 52 additions & 1 deletion R/describe_posterior.R
Original file line number Diff line number Diff line change
Expand Up @@ -568,7 +568,58 @@ describe_posterior.double <- describe_posterior.numeric


#' @export
describe_posterior.data.frame <- describe_posterior.numeric
#' @rdname describe_posterior
#' @inheritParams p_direction
describe_posterior.data.frame <- function(posterior,
centrality = "median",
dispersion = FALSE,
ci = 0.95,
ci_method = "eti",
test = c("p_direction", "rope"),
rope_range = "default",
rope_ci = 0.95,
keep_iterations = FALSE,
bf_prior = NULL,
BF = 1,
rvar_col = NULL,
verbose = TRUE,
...) {
if (length(x_rvar <- .possibly_extract_rvar_col(posterior, rvar_col)) > 0L) {
cl <- match.call()
cl[[1]] <- bayestestR::describe_posterior
cl$posterior <- x_rvar
cl$rvar_col <- NULL
if (length(prior_rvar <- .possibly_extract_rvar_col(posterior, bf_prior)) > 0L) {
cl$bf_prior <- prior_rvar
}
out <- eval.parent(cl)

obj_name <- insight::safe_deparse_symbol(substitute(posterior))
attr(out, "object_name") <- sprintf('%s[["%s"]]', obj_name, rvar_col)

return(.append_datagrid(out, posterior))
}


out <- .describe_posterior(
posterior,
centrality = centrality,
dispersion = dispersion,
ci = ci,
ci_method = ci_method,
test = test,
rope_range = rope_range,
rope_ci = rope_ci,
keep_iterations = keep_iterations,
bf_prior = bf_prior,
BF = BF,
verbose = verbose,
...
)

class(out) <- unique(c("describe_posterior", "see_describe_posterior", class(out)))
out
}


#' @export
Expand Down
39 changes: 32 additions & 7 deletions R/si.R
Original file line number Diff line number Diff line change
Expand Up @@ -221,8 +221,31 @@ si.get_predicted <- function(posterior, prior = NULL, BF = 1, use_iterations = F


#' @rdname si
#' @inheritParams p_direction
#' @export
si.data.frame <- function(posterior, prior = NULL, BF = 1, verbose = TRUE, ...) {
si.data.frame <- function(posterior, prior = NULL, BF = 1, rvar_col = NULL, verbose = TRUE, ...) {
if (length(x_rvar <- .possibly_extract_rvar_col(posterior, rvar_col)) > 0L) {
cl <- match.call()
cl[[1]] <- bayestestR::si
cl$posterior <- x_rvar
cl$rvar_col <- NULL
if (length(prior_rvar <- .possibly_extract_rvar_col(posterior, prior)) > 0L) {
cl$prior <- prior_rvar
}
out <- eval.parent(cl)

obj_name <- insight::safe_deparse_symbol(substitute(posterior))
attr(out, "object_name") <- sprintf('%s[["%s"]]', obj_name, rvar_col)

# This doesn't use .append_datagrid because we get a non-grid output
dgrid <- posterior[,vapply(posterior, function(col) !inherits(col, "rvar"), FUN.VALUE = logical(1)), drop = FALSE]
dgrid$Parameter <- unique(out$Parameter)
out_grid <- datawizard::data_join(dgrid, out, by = "Parameter")
class(out_grid) <- class(out)
return(out)
}


if (is.null(prior)) {
prior <- posterior
insight::format_warning(
Expand Down Expand Up @@ -255,7 +278,9 @@ si.data.frame <- function(posterior, prior = NULL, BF = 1, verbose = TRUE, ...)

#' @export
si.draws <- function(posterior, prior = NULL, BF = 1, verbose = TRUE, ...) {
si(.posterior_draws_to_df(posterior), prior = prior, BF = BF, verbose = verbose, ...)
si(.posterior_draws_to_df(posterior),
prior = if (!is.null(prior)) .posterior_draws_to_df(prior),
BF = BF, verbose = verbose, ...)
}

#' @export
Expand All @@ -276,7 +301,7 @@ si.rvar <- si.draws
)
}

out <- data.frame(
data.frame(
Parameter = colnames(posterior),
CI = BF,
CI_low = sis[, 1],
Expand Down Expand Up @@ -310,12 +335,12 @@ si.rvar <- si.draws

f_prior <- .logspline(prior, ...)
f_posterior <- .logspline(posterior, ...)
d_prior <- logspline::dlogspline(x_axis, f_prior)
d_posterior <- logspline::dlogspline(x_axis, f_posterior)
d_prior <- logspline::dlogspline(x_axis, f_prior, log = TRUE)
d_posterior <- logspline::dlogspline(x_axis, f_posterior, log = TRUE)

relative_d <- d_posterior / d_prior
relative_d <- d_posterior - d_prior

crit <- relative_d >= BF
crit <- relative_d >= log(BF)

cp <- rle(stats::na.omit(crit))
if (length(cp$lengths) > 3 && verbose) {
Expand Down
40 changes: 21 additions & 19 deletions man/bayesfactor_parameters.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

50 changes: 30 additions & 20 deletions man/bayesfactor_restricted.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

Loading

0 comments on commit 4b3e792

Please sign in to comment.