Scalable Bayesian Inference for Discrete Markov Random Fields

Estimation of discrete Markov random fields through the pseudo-likelihood, for both Ising models (Ising (1925) ) and ordinal Markov random fields (Marsman et al. (2025) ). Point estimates, either maximum pseudo-likelihood or maximum a posteriori, are obtained through a trust region optimization algorithm, with robust standard errors from the Godambe-Huber-White (GHW, or sandwich) estimator, for both full and constrained network structures. Bayesian inference corrects the pseudo-posterior through coordinate rescaling (Arena and Marsman (2026) ), which matches its covariance to the GHW covariance, and samples from it with the Fisher adaptive Metropolis-adjusted Langevin algorithm (Titsias (2024) ). For comparison, the posterior based on the full likelihood can also be sampled, exactly for small networks or with the double Metropolis-Hastings algorithm (Liang (2010) ). Savage-Dickey Bayes factors test each pairwise interaction, computed either by sampling importance resampling (Skare et al. (2003) ) or from the posterior draws.