Exact and Approximate MCMC for Doubly-intractable Probabilistic Graphical Models Leveraging the Underlying Independence Model

1. Problem Statement

The paper addresses the challenge of Bayesian inference for doubly-intractable probabilistic graphical models, specifically Pairwise Exponential Family Graphical Models (PEGMs) like the Ising model, Potts model, and Poisson graphical models.

• The Core Issue: The joint probability distribution is defined as $p(x; \theta) = f(x; \theta) / z(\theta)$. While the unnormalized density $f(x; \theta)$ is tractable, the normalizing constant $z(\theta)$ involves a summation (or integration) over all possible configurations of the variables. For high-dimensional models (large $p$), calculating $z(\theta)$ is computationally intractable ($O(2^p)$).
• Inference Barrier: Standard Markov Chain Monte Carlo (MCMC) methods, such as the Metropolis-Hastings (M-H) algorithm, require computing the ratio of normalizing constants $z(\theta)/z(\theta')$ to determine the acceptance probability. Since $z(\theta)$ is unknown, standard M-H cannot be applied directly.
• Limitations of Existing Solutions:
  • Exact Methods (e.g., Exchange Algorithm): Require sampling from the target distribution $p(\cdot; \theta')$ to construct an unbiased estimator of the ratio. This typically relies on perfect samplers (which are often unavailable or slow in high dimensions) or sequential samplers (like Sequential Monte Carlo), which are computationally prohibitive.
  • Approximate Methods: Often rely on noisy estimates that may not preserve the correct stationary distribution or suffer from poor mixing in high dimensions.

2. Methodology

The authors propose a novel framework that leverages the underlying independence model of the intractable graphical model to construct unbiased estimators without requiring perfect sampling.

Key Insight

For a PEGM with parameter matrix $\theta$, if all off-diagonal elements are set to zero (i.e., $\phi = \text{diag}(\theta)$), the model reduces to an independence model.

• The normalizing constant $z(\phi)$ for this independence model is tractable (it is the product of univariate normalizing constants).
• Sampling from the independence model $p(x; \phi)$ is trivial (drawing independent Bernoulli variables).

A. Exact Pseudo-Marginal MCMC (PM-Sampler)

The authors develop an unbiased estimator for the inverse normalizing constant $[z(\theta)]^{-n}$ (where $n$ is the sample size) to use within a Pseudo-Marginal MCMC framework.

1. Taylor Expansion: They express $[z(\theta)]^{-n}$ using a Taylor series expansion around a point related to the independence model. Specifically, they define $\mu = z(\theta)/z(\phi)$ and expand $(1 - \nu\mu)^{-n}$ (where $\nu$ is a tuning parameter).
2. Random Truncation: To estimate the infinite series, they employ a random truncation technique. A random integer $R$ is drawn, and the sum is truncated at $R$.
3. Unbiased Estimator Construction:
   • They construct an unbiased estimator $\hat{\mu}$ for $\mu$ using Importance Sampling:
     $$\hat{\mu} = \frac{1}{N} \sum_{i=1}^N \frac{f(y_i; \theta)}{f(y_i; \phi)}$$
     where $y_i \sim p(\cdot; \phi)$. Since sampling from $\phi$ is easy, this step is computationally efficient.
   • The final estimator $T$ combines the random truncation and the importance sampling estimates to provide an unbiased estimate of $[z(\theta)]^{-n}$.
4. Handling Negative Estimates: Since randomized sum estimators can occasionally be negative, the algorithm runs a chain on the absolute value of the posterior and reweights samples by the sign of the estimator to recover expectations under the true posterior.

B. Noisy Sampler (Approximate MCMC)

To address the computational cost of the exact method in very high dimensions, the authors propose a Noisy Sampler.

• Instead of estimating the M-H ratio unbiasedly, this method approximates the log-acceptance ratio.
• It uses the consistent estimator $\hat{\mu}$ directly in the log-scale acceptance probability.
• Theoretical Guarantee: The paper proves that as the number of Monte Carlo samples ($N$) increases, the transition kernel of the noisy sampler converges to the true kernel, and the chain inherits uniform ergodicity properties from the original chain (under bounded support assumptions).

C. Gradient-Based Proposals

The framework naturally accommodates gradient-based proposals (Langevin and Hamiltonian Monte Carlo) by utilizing the gradient of the log-posterior, which can be estimated using the same independence model sampling strategy.

3. Key Contributions

1. Exact Pseudo-Marginal Approach: A new exact MCMC method for doubly-intractable models that does not require perfect sampling or sequential Monte Carlo. It relies solely on the tractable independence model underlying the complex graph.
2. Scalability: The method decouples the number of auxiliary samples from the number of observed data points, avoiding the bottlenecks of the Exchange Algorithm. It scales well to high dimensions ($p=100+$).
3. Approximate Noisy Sampler: A computationally cheaper alternative that maintains theoretical ergodicity guarantees as $N \to \infty$.
4. Gradient Integration: Demonstrates how to integrate Langevin and Hamiltonian dynamics into this framework for more efficient exploration of the parameter space.

4. Experimental Results

The methods were evaluated on Ising models with varying dimensions ($p = 3$ to $100$) and compared against the Exchange Algorithm (EX).

• Mixing Efficiency:
  • The Pseudo-Marginal (PM) sampler demonstrated significantly better mixing (higher Effective Sample Size per minute) than the Exchange Algorithm, especially in high dimensions.
  • For $p=100$, the PM sampler with Langevin proposals achieved an ESS/minute roughly 3 times higher than the Exchange Algorithm.
• Runtime:
  • The PM sampler is computationally heavier per iteration than the Exchange Algorithm (due to the need for pilot runs to tune $\nu$), but its superior mixing makes it more efficient overall in high dimensions.
  • The Noisy Sampler had runtimes comparable to the Exchange Algorithm but with mixing properties that approached the PM sampler as $N$ increased.
• Parameter Recovery: All methods (PM, Noisy, and Exchange) recovered the true parameters ($\theta_0$) with similar accuracy (Mean Squared Error), provided the chains mixed well.
• Real Data Application: Applied to the MovieLens 32M dataset (modeling user preferences for 50 movies). The proposed methods successfully identified preference structures, showing high agreement in the signs of estimated parameters between the PM and Noisy methods.

5. Significance

This work represents a significant advancement in Bayesian inference for complex graphical models:

• Overcoming the "Perfect Sampler" Bottleneck: By utilizing the underlying independence structure, the authors bypass the need for perfect samplers, which are the primary limitation of the Exchange Algorithm in high dimensions.
• Scalability: The method is the first to offer a viable, exact (or near-exact) MCMC solution for doubly-intractable models in high dimensions ($p \ge 50$) without resorting to approximations that sacrifice theoretical validity.
• Generalizability: While demonstrated on the Ising model, the approach is applicable to any PEGM with an underlying independence structure (e.g., Potts, Poisson graphical models, Boltzmann machines).
• Practical Impact: It enables robust Bayesian inference for large-scale network models in fields like social network analysis, genomics, and recommendation systems, where traditional methods fail due to computational intractability.