bgm()

Fits a single-network Bayesian graphical model with bgm() and returns an S7 object of class bgms.

Description

bgm() fits a single network for binary, ordinal, continuous, or mixed data, with optional Bayesian edge selection. This page is the technical reference for function arguments and returned objects. For interpretation and workflow, start with Getting Started, then see Graphical Models and The Bayesian Approach. Methodological background is provided in Marsman et al. (2025).

Usage

bgm(
  x,
  variable_type = "ordinal",
  baseline_category,
  iter = 2000,
  warmup = 2000,
  interaction_prior = cauchy_prior(scale = 1),
  threshold_prior = beta_prime_prior(alpha = 0.5, beta = 0.5),
  means_prior = normal_prior(scale = 1),
  precision_scale_prior = gamma_prior(shape = 1, rate = 1),
  delta = NULL,
  edge_selection = TRUE,
  edge_prior = bernoulli_prior(0.5),
  na_action = c("listwise", "impute"),
  update_method = c("nuts", "adaptive-metropolis"),
  target_accept,
  nuts_max_depth = 10,
  learn_mass_matrix = TRUE,
  chains = 4,
  cores = parallel::detectCores(),
  display_progress = c("per-chain", "total", "none"),
  seed = NULL,
  standardize = FALSE,
  verbose = getOption("bgms.verbose", TRUE),
  progress_callback = NULL
)

Arguments

Data

Argument Description
x A data frame or matrix with n rows and p columns. Columns may contain binary, ordinal, or continuous variables (see variable_type). Discrete variables are automatically recoded to non-negative integers (0, 1, ..., m); for regular ordinal variables, unobserved categories are collapsed, while Blume–Capel variables retain all categories. Continuous variables are column-centered internally so that the GGM likelihood is formulated with a zero-mean assumption.
variable_type Character or character vector. Specifies the type of each variable in x. Allowed values: "ordinal", "blume-capel", or "continuous". A single string applies to all variables. A per-variable vector that mixes discrete ("ordinal" / "blume-capel") and "continuous" types fits a mixed MRF. Binary variables are automatically treated as "ordinal". Default: "ordinal".
baseline_category Integer or vector. Baseline category used in Blume–Capel variables. Can be a single integer (applied to all) or a vector of length p. Required if at least one variable is of type "blume-capel".
na_action Character. Specifies missing data handling. Either "listwise" (drop rows with missing values) or "impute" (impute within Gibbs sampler, propagating uncertainty). Default: "listwise". See Missing Data for details.

Prior specification

Priors are specified by passing prior-object constructors (see Prior Constructors). The legacy scalar arguments listed below the table (pairwise_scale, main_alpha, main_beta, inclusion_probability, beta_bernoulli_alpha/beta, …) are deprecated as of bgms 0.2.0; they still work but emit a lifecycle warning and forward to the equivalent prior object.

Argument Description
interaction_prior A bgms_parameter_prior object for pairwise interaction parameters. Allowed: cauchy_prior() (default), normal_prior(), beta_prime_prior(). Default: cauchy_prior(scale = 1).
threshold_prior A bgms_parameter_prior object for threshold (main effect) parameters. Allowed: beta_prime_prior() (default), cauchy_prior(), normal_prior(). Default: beta_prime_prior(alpha = 0.5, beta = 0.5).
means_prior A bgms_parameter_prior object for continuous-variable means (mixed MRF only). Allowed: normal_prior() (default), cauchy_prior(), beta_prime_prior(). Ignored for pure ordinal or pure GGM models. Default: normal_prior(scale = 1).
precision_scale_prior A bgms_scale_prior object for the diagonal elements of the precision matrix (GGM and mixed MRF). Allowed: gamma_prior() (default), exponential_prior(). Ignored for pure ordinal models. Default: gamma_prior(shape = 1, rate = 1).
delta Non-negative numeric, or NULL for the dimension-adaptive default. Determinant-tilt exponent on the continuous-block precision matrix \(K\) (GGM) or \(K_{yy}\) (mixed MRF): multiplies the prior by \(|K|^{\delta}\), softly repelling the chain from the boundary of the positive-definite cone. NULL (default) auto-resolves to \(0.5 \log(p)\), with \(p\) the number of continuous variables. Pass delta = 0 for the untilted prior. Ignored for pure ordinal models.
edge_selection Logical. Whether to perform Bayesian edge selection. If FALSE, the model estimates all edges. Default: TRUE.
edge_prior A bgms_indicator_prior object specifying the prior for edge inclusion. Allowed: bernoulli_prior() (default), beta_bernoulli_prior(), sbm_prior(). Legacy character strings ("Bernoulli", "Beta-Bernoulli", "Stochastic-Block") are accepted but deprecated. Default: bernoulli_prior(0.5).
standardize Logical. If TRUE, the prior scale for each pairwise interaction is adjusted to the range of the response scores, so that shrinkage is comparable across variable pairs with different numbers of categories. The adjusted scale for the interaction between variables \(i\) and \(j\) is scale * m_i * m_j (with \(m_i\) the maximum score after recoding); scale itself then applies to the binary case. See Parameter Priors for the Blume-Capel and mixed-pair rules. Default: FALSE.
NoteDeprecated scalar prior arguments

The following arguments are deprecated as of bgms 0.2.0. They still work but emit a lifecycle warning. Use the prior-object equivalents in the table above.

Deprecated Replacement
pairwise_scale = s interaction_prior = cauchy_prior(scale = s)
main_alpha = a, main_beta = b threshold_prior = beta_prime_prior(alpha = a, beta = b)
inclusion_probability = p edge_prior = bernoulli_prior(p)
beta_bernoulli_alpha = a, beta_bernoulli_beta = b edge_prior = beta_bernoulli_prior(a, b)
beta_bernoulli_alpha_between, beta_bernoulli_beta_between, dirichlet_alpha, lambda edge_prior = sbm_prior(alpha_between = …, beta_between = …, dirichlet_alpha = …, lambda = …)
interaction_scale = s interaction_prior = cauchy_prior(scale = s) (deprecated since 0.1.6.0)
threshold_alpha = a, threshold_beta = b threshold_prior = beta_prime_prior(alpha = a, beta = b) (deprecated since 0.1.6.0)
burnin = n warmup = n (deprecated since 0.1.6.0)
save No replacement; raw samples are always stored in raw_samples (deprecated since 0.1.6.0)

Sampler settings

Argument Description
iter Integer. Number of post-warmup iterations (per chain). Default: 2000.
warmup Integer. Number of warmup iterations before collecting samples. Default: 2000. Short warmups trigger progressive warnings; see Warmup Schedule for the adaptation stages.
update_method Character. Specifies how the MCMC sampler updates the model parameters: "adaptive-metropolis" (componentwise adaptive Metropolis–Hastings with Robbins–Monro proposal adaptation) or "nuts" (the No-U-Turn Sampler with dynamically chosen trajectory lengths). Default: "nuts".
target_accept Numeric between 0 and 1. Target acceptance rate for the sampler. Defaults are set automatically if not supplied: 0.44 for adaptive Metropolis and 0.80 for NUTS. Under update_method = "nuts" this governs the step-size adaptation; the edge-move Metropolis proposals then use a fixed 0.44 target regardless of this argument.
nuts_max_depth Integer. Maximum tree depth in NUTS. Must be positive. Default: 10.
learn_mass_matrix Logical. If TRUE, adapt a diagonal mass matrix during warmup (NUTS only). If FALSE, use the identity matrix. Default: TRUE.
chains Integer. Number of parallel chains to run. Default: 4.
cores Integer. Number of CPU cores for parallel execution. Default: parallel::detectCores().
display_progress Character. Controls progress reporting during sampling. Options: "per-chain" (separate bar per chain), "total" (single combined bar), or "none" (no progress). Default: "per-chain".
seed Optional integer. Random seed for reproducibility. Must be a single non-negative integer.
verbose Logical. If TRUE, prints informational messages during data processing (e.g., missing data handling, variable recoding). Defaults to getOption("bgms.verbose", TRUE). Set options(bgms.verbose = FALSE) to suppress messages globally.
progress_callback Optional R function with signature function(completed, total), called at regular intervals during sampling with the number of completed and total iterations across all chains. Useful for external front-ends (e.g., JASP) that supply their own progress reporting. Default: NULL (no callback).

Value

An S7 object of class bgms with posterior summaries, posterior mean matrices, and raw MCMC draws. Supports print(), summary(), coef(), predict(), and simulate().

Main components:

  • posterior_summary_main — Data frame with posterior summaries (mean, sd, MCSE, ESS, Rhat) for main-effect parameters. For OMRF models these are category thresholds; for mixed MRF models these are discrete thresholds and continuous means. NULL for GGM models.
  • posterior_summary_quadratic — Data frame with posterior summaries for the residual variance parameters (GGM and mixed MRF). NULL for OMRF models.
  • posterior_summary_pairwise — Data frame with posterior summaries for partial association parameters.
  • posterior_summary_indicator — Data frame with posterior summaries for edge inclusion indicators (if edge_selection = TRUE).
  • posterior_mean_main — Posterior mean of main-effect parameters. NULL for GGM models. For OMRF: a matrix (p x max_categories) of category thresholds. For mixed MRF: a list with $discrete (threshold matrix) and $continuous (q x 1 matrix of means).
  • posterior_mean_pairwise — Symmetric matrix of posterior mean partial associations (zero diagonal). For continuous variables these are unstandardized partial correlations; for discrete variables these are half the log adjacent-category odds ratio. Use extract_precision(), extract_partial_correlations(), or extract_log_odds() to convert to interpretable scales.
  • posterior_mean_residual_variance — Named numeric vector of posterior mean residual variances \(1/\Theta_{ii}\). Present for GGM and mixed MRF models; NULL for OMRF.
  • posterior_mean_indicator — Symmetric matrix of posterior mean inclusion probabilities (if edge selection was enabled).
  • raw_samples — A list of raw MCMC draws per chain with the following sublists:
    • main — List of main effect samples.
    • pairwise — List of pairwise effect samples.
    • indicator — List of indicator samples (if edge selection enabled).
    • allocations — List of cluster allocations (if SBM prior used).
    • nchains — Number of chains.
    • niter — Number of post-warmup iterations per chain.
    • parameter_names — Named lists of parameter labels.
  • arguments — A list of function call arguments and metadata.

Additional summaries when edge_prior = sbm_prior(...) (see Sekulovski et al., 2025):

  • posterior_summary_pairwise_allocations — Data frame with posterior summaries (mean, sd, MCSE, ESS, Rhat) for the pairwise cluster co-occurrence of the nodes. This serves to indicate whether the estimated posterior allocations, co-clustering matrix and posterior cluster probabilities have converged.
  • posterior_mean_coclustering_matrix — Symmetric matrix of pairwise proportions of occurrence of every variable. This matrix can be plotted to visually inspect the estimated number of clusters and visually inspect nodes that tend to switch clusters.
  • posterior_mean_allocations — A vector with the posterior mean of the cluster allocations of the nodes.
  • posterior_mode_allocations — A vector with the posterior mode of the cluster allocations of the nodes.
  • posterior_num_blocks — A data frame with the estimated posterior inclusion probabilities for all the possible number of clusters.

NUTS diagnostics (tree depth, divergences, energy, E-BFMI) are available in fit$nuts_diag when update_method = "nuts".

Details

For full explanations, see the guide. Model-specific assumptions and parameterization are described in Ordinal MRF, Gaussian Graphical Model, and Mixed MRF. Prior specification and structure learning are discussed in Prior Basics, Edge Selection, and Edge Clustering. To interpret outputs and verify chain quality, see MCMC Output and MCMC Diagnostics.

Examples

Fit an ordinal MRF to the Wenchuan PTSD data (the model below was run in advance; fitting takes a few minutes):

library(bgms)
fit = bgm(Wenchuan, chains = 2, seed = 1234)
summary(fit)
Posterior summaries from Bayesian estimation:

Category thresholds: 
                 mean  mcse    sd   n_eff  Rhat
intrusion (1)   0.457 0.010 0.250 580.592 1.009
intrusion (2)  -2.035 0.024 0.387 259.263 1.006
intrusion (3)  -5.276 0.045 0.661 214.646 1.006
intrusion (4) -10.347 0.066 1.038 248.274 1.003
dreams (1)     -0.817 0.008 0.212 675.187 1.000
dreams (2)     -4.379 0.020 0.420 436.475 1.001
... (use `summary(fit)$main` to see full output)

Pairwise interactions:
                    mean    sd  mcse   n_eff  Rhat
intrusion-dreams   0.369 0.001 0.037 848.100 1.000
intrusion-flash    0.165 0.001 0.034 623.634 1.055
intrusion-upset    0.062 0.055 0.008  46.612 1.058
intrusion-physior  0.007 0.020 0.002 151.748 1.000
intrusion-avoidth  0.000 0.000 0.000  84.121 1.000
intrusion-avoidact 0.000 0.001 0.000 146.722 1.000
... (use `summary(fit)$pairwise` to see full output)

Inclusion probabilities:
                    mean    sd  mcse n0->0 n0->1 n1->0 n1->1    n_eff  Rhat
intrusion-dreams   1.000 0.000           0     0     0  1999               
intrusion-flash    1.000 0.000           0     0     0  1999               
intrusion-upset    0.599 0.490 0.076   781    19    20  1179   41.265 1.116
intrusion-physior  0.106 0.308 0.025  1760    27    27   185  153.394     1
intrusion-avoidth  0.008 0.089 0.002  1969    14    14     2  1578.01  1.03
intrusion-avoidact 0.012 0.109 0.003  1957    18    18     6 1223.503  1.03
... (use `summary(fit)$indicator` to see full output)
Note: NA values are suppressed in the print table; they occur for indicators
that were not updated or whose draws are constant, so ESS/Rhat are undefined.
`summary(fit)$indicator` still contains all computed values.

Use `summary(fit)$<component>` to access full results.
Use `extract_log_odds(fit)` for log odds ratios.
See the `easybgm` package for other summary and plotting tools.

The threshold and pairwise tables report the posterior mean, Monte Carlo standard error (mcse), posterior standard deviation, effective sample size, and R-hat for each parameter. In the pairwise table, an estimate of exactly 0.000 (e.g. intrusion-avoidth) is an edge that edge selection excluded in essentially every draw. The inclusion table (further down the printed output) gives each edge’s posterior inclusion probability in the mean column. See MCMC Output for a walkthrough of every component and MCMC Diagnostics for how to act on n_eff and Rhat.

See also

Getting Started, Methods, Extractor Functions, and bgmCompare().

References

Marsman, M., van den Bergh, D., & Haslbeck, J. M. B. (2025). Bayesian analysis of the ordinal Markov random field. Psychometrika, 90(1), 146–182. https://doi.org/10.1017/psy.2024.4
Sekulovski, N., Arena, G., Haslbeck, J. M. B., Huth, K. B. S., Friel, N., & Marsman, M. (2025). A stochastic block prior for clustering in graphical models. Manuscript Submitted for Publication. https://osf.io/preprints/psyarxiv/29p3m_v3