Summaries from joint posterior draws

Computes either the joint directional statement (the default): for each \(k\), the joint posterior probability that the \(k\) strongest directional claims hold simultaneously, or simultaneous credible intervals (interval = TRUE).

Usage

joint(
  draws,
  interval = FALSE,
  prob = 0.95,
  null.value = 0,
  direction = NULL,
  est.FUN = median
)

Arguments

draws: Numeric matrix or data frame of posterior draws (rows = draws, columns = parameters). At least 1000 rows and 2 columns.
interval: Logical (default FALSE). If FALSE, return the joint directional statement; if TRUE, return simultaneous credible intervals.
prob: Numeric scalar in \((0, 1)\). Joint coverage probability for the simultaneous intervals (default 0.95). Used only when interval = TRUE.
null.value: Numeric scalar or vector giving the null reference value(s) used to define the claimed direction (a scalar is recycled). Used only when interval = FALSE. Defaults to 0.
direction: Character vector of "greater", "less", or "two.sided" (or NULL), the testing mode per parameter. A scalar is recycled; NULL (default) uses "two.sided" (dominant side) for all. Used only when interval = FALSE.
est.FUN: Function returning a point estimate from a numeric vector (default median). Used only when interval = TRUE.

Value

A data.frame. For interval = FALSE, one row per parameter ordered by decreasing pd, with columns median.est, null.value, pd, direction, and joint_cum (the joint posterior probability that the \(k\) strongest claims all hold). For interval = TRUE, one row per parameter with columns lower, est, upper, and prob.

Details

Joint directional statement (interval = FALSE, the default). Each parameter's claimed direction is its dominant posterior side relative to null.value (or the side set by direction). Parameters are ordered by decreasing pd, and the joint_cum column reports, for each \(k\), the joint posterior probability that the \(k\) strongest claims hold simultaneously. It is computed directly from the draws by intersecting the per-draw directional events, so it reflects the posterior dependence among parameters.

Simultaneous credible intervals (interval = TRUE). Equi-tailed intervals calibrated so that all parameters lie within their bounds simultaneously with probability prob, via the quantile-based simultaneous credible band. Coverage is set by examining how extreme each draw is across all parameters jointly: for every draw the minimum (nearer) tail probability across parameters is taken, and the threshold is the \((1 - \texttt{prob})\)-quantile of these minima. An effect can be declared when a parameter's band excludes the null value. A closely related implementation is sim.cred.band in the credsubs package (Schnell et al., 2020).

References

Box, G. E. P., & Tiao, G. C. (1992). Bayesian inference in statistical analysis (1st ed.). Wiley.

Besag, J., Green, P., Higdon, D., & Mengersen, K. (1995). Bayesian computation and stochastic systems. Statistical Science, 10(1), 3–41.

Held, L. (2004). Simultaneous posterior probability statements from Monte Carlo output. Journal of Computational and Graphical Statistics, 13(1), 20–35.

Gelman, A., Tuerlinckx, F. (2000). Type S error rates for classical and Bayesian single and multiple comparison procedures. Computational Statistics 15, 373–390 (2000). https://doi.org/10.1007/s001800000040

Schnell, P. et al. (2020). credsubs: Credible Subsets. R package. https://CRAN.R-project.org/package=credsubs.

Examples

mu <- c(1.5, 0, -2)
Sigma <- matrix(c(1, 0.8, 0.5, 0.8, 1, 0.3, 0.5, 0.3, 1), 3, 3)
draws <- MASS::mvrnorm(2000, mu, Sigma)
colnames(draws) <- c("theta1", "theta2", "theta3")

joint(draws)                    # joint directional statement (default)
#>        median.est null.value     pd direction joint_cum
#> theta3    -1.9596          0 0.9755 two.sided    0.9755
#> theta1     1.5206          0 0.9405 two.sided    0.9160
#> theta2     0.0021          0 0.5005 two.sided    0.4820
joint(draws, interval = TRUE)   # simultaneous credible intervals
#>             lower          est     upper prob
#> theta1 -0.7956087  1.520594092 3.7265003 0.95
#> theta2 -2.4132437  0.002074391 2.2378859 0.95
#> theta3 -4.2133338 -1.959627191 0.3455705 0.95