Introduction
This vignette illustrates how to inspect convergence diagnostics and how to interpret spike-and-slab summaries in bgms models. For some of the model variables spike-and-slab priors introduce binary indicator variables that govern whether the effect is included or not. Their posterior distributions can be summarized with inclusion probabilities and Bayes factors.
Example fit
We use a subset of the Wenchuan dataset:
library(bgms)
?Wenchuan
data = Wenchuan[, 1:5]
fit = bgm(data, seed = 1234)
#> Warning: There were 7 rows with missing observations in the input matrix x.
#> Since na_action = listwise these rows were excluded from the analysis.Note: During fitting, progress bars are shown in interactive sessions. In this vignette, they are suppressed for clarity. Sampling can take a while; the progress bars usually help track progress.
Convergence diagnostics
The quality of the Markov chain can be assessed with common MCMC diagnostics:
summary(fit)$pairwise
#> mean sd mcse n_eff
#> intrusion-dreams 0.630465218 0.001505695 0.0645952597 1840.4640
#> intrusion-flash 0.336939066 0.001279900 0.0615246748 2310.7156
#> intrusion-upset 0.198665557 0.067666296 0.0035365783 366.0819
#> intrusion-physior 0.191577877 0.071227938 0.0054596215 170.2063
#> dreams-flash 0.499277762 0.001398210 0.0596431449 1819.6018
#> dreams-upset 0.230190853 0.055082965 0.0013903146 1569.6703
#> dreams-physior 0.005329141 0.021891564 0.0008067065 736.4147
#> flash-upset 0.007130963 0.025475978 0.0010058723 641.4695
#> flash-physior 0.309372051 0.001159184 0.0533283561 2116.4689
#> upset-physior 0.712758036 0.001639670 0.0602202440 1348.8769
#> Rhat
#> intrusion-dreams 1.000830
#> intrusion-flash 1.001655
#> intrusion-upset 1.011155
#> intrusion-physior 1.009193
#> dreams-flash 1.001358
#> dreams-upset 1.004926
#> dreams-physior 1.003559
#> flash-upset 1.001278
#> flash-physior 1.003235
#> upset-physior 1.002455- R-hat values close to 1 (typically below 1.01) suggest convergence (Vehtari et al., 2021).
- The effective sample size (ESS) reflects the number of independent samples that would provide equivalent precision. Larger ESS values indicate more reliable estimates.
- The Monte Carlo standard error (MCSE) measures the additional variability introduced by using a finite number of MCMC draws. A small MCSE relative to the posterior standard deviation indicates stable estimates, whereas a large MCSE suggests that more samples are needed.
Advanced users can inspect traceplots by extracting raw samples and
using external packages such as coda or
bayesplot. Here is an example using the coda
package to create a traceplot for a pairwise effect parameter.
library(coda)
param_index = 1
chains = lapply(fit$raw_samples$pairwise, function(mat) mat[, param_index])
mcmc_obj = mcmc.list(lapply(chains, mcmc))
traceplot(mcmc_obj, col = c("firebrick", "steelblue", "darkgreen", "goldenrod"),
main = "Traceplot of pairwise[1]") 
Spike-and-slab summaries
The spike-and-slab prior yields posterior inclusion probabilities for edges:
coef(fit)$indicator
#> intrusion dreams flash upset physior
#> intrusion 0.00000 1.0000 1.0000 0.93975 0.0570
#> dreams 1.00000 0.0000 0.9665 1.00000 0.0745
#> flash 1.00000 0.9665 0.0000 0.99800 1.0000
#> upset 0.93975 1.0000 0.9980 0.00000 1.0000
#> physior 0.05700 0.0745 1.0000 1.00000 0.0000- Values near 1.0: strong evidence the edge is present.
- Values near 0.0: strong evidence the edge is absent.
- Values near 0.5: inconclusive (absence of evidence).
Bayes factors
When the prior inclusion probability for an edge is equal to 0.5
(e.g., using a Bernoulli prior with
inclusion_probability = 0.5 or a symmetric Beta prior,
main_alpha = main_beta), we can directly transform
inclusion probabilities into Bayes factors for edge presence vs
absence:
# Example for one edge
p = coef(fit)$indicator[1, 5]
BF_10 = p / (1 - p)
BF_10
#> [1] 0.06044539Here the Bayes factor in favor of inclusion (H1) is small, meaning that there is little evidence for inclusion. Since the Bayes factor is transitive, we can use it to express the evidence in favor of exclusion (H0) as
1 / BF_10
#> [1] 16.54386This Bayes factor shows that there is strong evidence for the absence
of a network relation between the variables intrusion and
physior.
Notes on runtime
- Sampling with spike-and-slab priors can take longer.
- In interactive sessions, progress bars are displayed. In this vignette, they are suppressed for readability.
Next steps
- See Getting Started for a simple one-sample workflow.
- See Model Comparison for group differences.