PaperlayerView exact revision on arXiv

Source-linked AI summary

Austerity in MCMC Land: Cutting the Metropolis-Hastings Budget

Anoop Korattikara, Yutian Chen, Max Welling

arXiv:1304.5299v4cs.LGstat.ML

TL;DR

Large datasets make exact Metropolis-Hastings decisions expensive because they may require evaluating all N likelihood contributions for one binary outcome. The paper replaces this with a sequential hypothesis test that often uses fewer datapoints, introducing controllable bias that can be offset by drawing more samples and reducing variance.

  • Problem

    Exact Metropolis-Hastings testing can require examining all data points for one accept/reject decision, making posterior MCMC costly on very large datasets.

  • Method

    The paper uses a sequential hypothesis test on subsampled log-likelihood differences, increasing the data subset until a prescribed confidence level supports acceptance or rejection.

  • Results

    The approximate MH test often makes confident decisions with n < N datapoints, saving computation while introducing a controllable bias that enables more samples per unit time.

  • Takeaways & Limitations

    The method establishes a bias-variance trade-off for MH sampling in which computational savings can reduce variance enough to offset the introduced bias.

  • Takeaways & Limitations

    The algorithm can behave erratically when the central limit theorem fails, such as with very sparse datasets or extreme outliers.

Abstract

from arXiv · show

Can we make Bayesian posterior MCMC sampling more efficient when faced with very large datasets? We argue that computing the likelihood for N datapoints in the Metropolis-Hastings (MH) test to reach a single binary decision is computationally inefficient. We introduce an approximate MH rule based on a sequential hypothesis test that allows us to accept or reject samples with high confidence using only a fraction of the data required for the exact MH rule. While this method introduces an asymptotic bias, we show that this bias can be controlled and is more than offset by a decrease in variance due to our ability to draw more samples per unit of time.

1 Introduction

Standard Metropolis-Hastings MCMC examines all data for each accept/reject decision, making Bayesian computation costly on large datasets. The paper argues for reconsidering this paradigm under finite computational budgets and introduces approximate decisions based on sequential testing.

  • 1 Introduction: MCMC proposals are accepted or rejected using the Metropolis-Hastings formula after examining the likelihood of all data items.This makes each binary decision potentially expensive when datasets are large.
  • 1 Introduction: Standard MCMC methods may be inefficient for Big Data because they spend substantial computation producing only one accept/reject bit.The paper argues that prior efficiency efforts have focused mainly on parallelization rather than questioning the algorithm itself.
  • 1 Introduction: Finite-time MCMC error reflects both bias and variance, while retained post-burn-in samples reduce variance at O(1/T).The paper frames efficiency as minimizing error within a fixed computational time T.
  • 1 Introduction: Sequential hypothesis tests can decide whether to accept or reject many proposals using a small data fraction with high confidence.This cuts the computational budget of the MH step and permits more sampling within the same time.
  • 1 Introduction: The proposed approximate MH step aims to retain a controllable bias knob while avoiding O(N) computations per iteration.The paper contrasts this with related biased methods that omit MH steps and can produce occasional extremely low-probability samples.

2 The Metropolis-Hastings algorithm

Metropolis-Hastings constructs a Markov chain with the desired stationary distribution by accepting or rejecting proposals. Although the resulting estimator is unbiased after burn-in, expensive acceptance calculations and slow mixing can make its variance prohibitively high.

  • 2 The Metropolis-Hastings algorithm: Metropolis-Hastings transitions from θ_t to a proposal θ′ by computing an acceptance probability and comparing it with a Uniform[0,1] draw.The proposal is accepted when u < P_a; otherwise the chain remains at θ_t.
  • 2 The Metropolis-Hastings algorithm: The chain’s stationary distribution is S_0(θ), so post-burn-in samples provide an unbiased estimator of expectations under S_0.The estimator averages f(θ_t) over collected chain samples.
  • 2 The Metropolis-Hastings algorithm: Estimator variance depends on the variance of f under S_0, the integrated autocorrelation time τ, and the number of samples T.Slow mixing increases τ and reduces the effective number of independent samples.
  • 2 The Metropolis-Hastings algorithm: Expensive evaluation of S_0(θ) can limit T, making the variance of the unbiased estimator prohibitively high.This is especially problematic when chains mix slowly and acceptance probabilities must be evaluated at every step.

3 Approximate MCMC and the Bias-Variance Tradeoff

Approximate MCMC deliberately allows controlled stationary-distribution bias to increase sampling speed and reduce variance within finite computation time. The paper presents this as a bias-variance trade-off governed by a tunable parameter.

  • 3 Approximate MCMC and the Bias-Variance Tradeoff: Allowing a small bias in the stationary distribution can produce Markov chains that are cheaper to simulate.The target estimator then uses samples from S_ϵ rather than exactly from S_0.
  • 3 Approximate MCMC and the Bias-Variance Tradeoff: As ϵ approaches zero, S_ϵ approaches S_0, lowering bias but increasing variance because fewer samples can be collected.The exact-MCMC limit is therefore not necessarily optimal under finite computational time.
  • 3 Approximate MCMC and the Bias-Variance Tradeoff: Moving ϵ away from zero increases the difference between S_ϵ and S_0 while enabling more samples and lower variance.This produces the classical bias-variance trade-off in estimator risk.
  • 3 Approximate MCMC and the Bias-Variance Tradeoff: Risk decomposes as R = B^2 + V, with the optimal ϵ depending on the available computational time.With infinite time, ϵ = 0 is appropriate; with finite time, a biased setting may reduce total risk.
  • 3 Approximate MCMC and the Bias-Variance Tradeoff: The paper extends the bias-knob paradigm to approximate MH tests that can use any proposal distribution for large-data Bayesian posteriors.This provides an alternative to methods such as SGLD and SGFS, which omit MH tests.

4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling

The proposed approximate MH test reformulates acceptance as a statistical decision about the mean log-likelihood difference. It sequentially analyzes data until the decision reaches a prescribed confidence, often using fewer than N observations while trading controlled bias against variance.

  • 4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling: For a dataset of N observations, exact posterior MH testing can require O(N) computation for one accept/reject decision.The paper motivates subsampling because evaluating the posterior for every generated sample becomes costly at very large N.
  • 4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling: The MH accept/reject rule can be reformulated as testing whether the mean log-likelihood difference µ exceeds a threshold µ_0.The proposal is accepted when µ > µ_0 and rejected when µ ≤ µ_0.
  • 4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling: The sequential test computes a sample mean and standard deviation from data drawn without replacement, then evaluates a Student-t statistic when the CLT applies.Its hypotheses are H1: µ > µ_0 versus H2: µ < µ_0.
  • 4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling: The algorithm starts with a mini-batch of size m and adds m datapoints whenever the confidence threshold has not been reached.It terminates when δ < ϵ, deciding according to whether the estimated mean exceeds µ_0.
  • 4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling: Often, confident decisions use n < N datapoints, saving computation while introducing a small stationary-distribution bias that can reduce variance through additional samples.The knob ϵ controls this trade-off: higher ϵ increases bias, whereas ϵ → 0 increases data use and variance.
  • 4 Approximate Metropolis-Hastings Test for Bayesian Posterior Sampling: The method can behave erratically when the central limit theorem does not hold, including very sparse datasets or datasets with extreme outliers.The paper recommends empirically testing the CLT assumption before running the algorithm.

5 Error Analysis and Test Design

The paper analyzes how sequential-test error propagates from individual acceptance decisions to the stationary distribution, and how test parameters can be designed to control data usage and bias.

  • Error analysis: The single-test parameter ϵ bounds individual test error, while complete-test error requires assumptions that t statistics approximate z statistics and mini-batch means are jointly normal.Under these assumptions, the sequential-test statistic follows a Gaussian random walk.
  • Error analysis: The test error E and average data usage ¯π depend on θ, θ′, and u through the standardized mean µstd(u, θ, θ′).This reduction enables dynamic-programming analysis of the test.
  • Error analysis: At µ = µ0, the sequential-test error is largest, so E(0, m, ϵ) provides an upper bound; error decreases sharply as µ moves away from µ0.Figure 1 compares theoretical and simulated error for 1000 sequential tests in the logistic-regression model.
  • Error analysis: Errors across different u values partly cancel, making the actual acceptance-probability error ∆ usually much smaller than the worst-case test error.For fixed (θ, θ′), ∆ can be computed using one-dimensional quadrature.
  • Error analysis: The stationary-distribution error is bounded linearly by ∆max = supθ,θ′ |∆(θ, θ′)| under the stated contraction condition.The theorem bounds the distance between the exact posterior S0 and the approximate chain’s stationary distribution Sϵ.
  • Test design: A practical design uses m ≈500 and the smallest ϵ that maintains low average data usage, while an optimal test minimizes data subject to an acceptance-probability error tolerance.Without a trial run, controlling worst-case error produces a conservative design because worst-case error is usually much higher than average-case error.

6 Experiments

Experiments on logistic regression, ICA, variable selection, and SGLD show that approximate MH sampling can reduce risk or data usage under finite computation, while exposing design and initialization trade-offs.

  • 6.1 Random Walk - Logistic Regression: On MNIST logistic regression, risk remains variance-dominated over time, and the minimum risk occurs with ϵ > 0 rather than exact MH.The experiment uses 12,214 datapoints reduced to 50 dimensions and evaluates predictive-mean risk against a long Hybrid Monte Carlo run.
  • 6.2 Independent Component Analysis: In ICA, the lowest risk after approximately 6400 seconds occurs with ϵ = 0.1 rather than exact MH, although exact MH is expected to win asymptotically.The study mixes 1.95 million samples from four sources and estimates risk using ten chains against a long exact-MH run.
  • 6.3 Variable selection in Logistic Regression: In MiniBooNE variable selection, the lowest predictive-mean risk is obtained with ϵ > 0.The dataset contains 130,065 datapoints; approximate and exact reversible-jump samplers are compared over roughly 3500 seconds.
  • 6.3 Variable selection in Logistic Regression: Both reversible-jump samplers encounter local minima, with learned models and classification error depending strongly on initialization.Starting with one variable yields about 12 features and approximately 15% error, whereas MAP initialization with all features yields about 45 features and approximately 10% error.
  • 6.4 Stochastic Gradient Langevin Dynamics: Uncorrected SGLD can be driven far from the mode by noisy jumps into low-density valleys, whereas approximate MH rejects most such jumps and produces a stationary distribution nearly indistinguishable from the posterior.With ϵ = 0.5, decisions are always made using the first m = 500 datapoints.
  • 6.5 Optimal Design of Sequential Tests: For sequential-test design, average-case design achieves nearly the target test error with much less data than worst-case design, while optimal design can be twice as fast as fixed m = 600 when error is large.The optimal design uses half the data of the fixed-m heuristic in that regime.

7 Conclusions and Future Work

The approximate MH method reduces per-decision computation by sequentially testing on data subsets, introducing controllable bias that can be offset by lower variance from faster sampling.

  • O(N) likelihood evaluations are replaced by progressively larger data subsets until the acceptance decision reaches a prescribed confidence level.The approach increases the subset size only when the current evidence is insufficient.
  • The resulting bias is more than compensated for by reduced variance because the algorithm can draw more samples per unit time.
  • A fixed decision threshold is used here, while adapting that threshold over time is proposed as future work to balance squared bias and variance.

A Distribution of the test statistic

The sequential test is analyzed through approximate Gaussian and random-walk behavior of statistics computed from growing mini-batches, supporting data-efficient decisions with quantified error and usage.

  • Sequential test: A mini-batch t-statistic is computed first, and additional mini-batches are added until the statistic supports a decision.The procedure is guaranteed to terminate.
  • Distributional assumptions: For sufficiently large mini-batches, the central limit theorem supports a normal approximation, and the t-statistic reduces to a z-statistic at high degrees of freedom.The analysis assumes a multivariate normal joint distribution for the sequence of sample means.
  • Distributional assumptions: Under the assumptions, the z-statistic sequence follows a Gaussian random walk whose parameters use the standardized mean and the observed data proportion.
  • Random-walk behavior: When the observed data proportion approaches 1, the mean of z_j approaches infinity while its variance remains 1, yielding a correct decision with all data observed.
  • Error and data usage: The complete sequential-test error and expected data usage are worst when the standardized mean approaches zero, equivalent to the null boundary.
  • Error and data usage: Theoretical error and average data usage from dynamic programming match simulations, and both decline rapidly as the standardized mean moves away from the threshold.

B Error in One Metropolis-Hastings Step

The approximate MH test can have large conditional errors for individual uniform thresholds, yet marginalizing over the MH random variable can substantially reduce the actual acceptance-probability error.

  • The approximate MH test draws a uniform random variable and then applies the sequential test, with conditional error depending on that draw.
  • When the exact acceptance probability is near 0.5, the average absolute conditional error can be large because typical standardized means lie near zero.
  • Cancellation between errors for cases above and below the threshold allows the approximate acceptance probability to closely approximate the exact acceptance probability.
  • The acceptance-probability error can therefore be estimated more accurately using knowledge of the mean and variance than through worst-case analysis.

C.1 Upper Bound Based on One Step Error

The paper bounds stationary-distribution error through one-step acceptance error and designs sequential tests to reduce average data usage under a practical acceptance-probability tolerance.

  • Upper-bound analysis: A contraction condition for the exact transition kernel and a uniform one-step kernel error bound imply a bound on the distance between stationary distributions.
  • Upper-bound analysis: The approximate chain enters a neighborhood of the exact stationary distribution in finitely many steps and remains there thereafter under the stated conditions.
  • Test design: Bounding acceptance-probability error is more practical than bounding stationary-distribution error directly because the latter depends on the difficult-to-compute contraction parameter.
  • Test design: Given an acceptance-error bound, the design minimizes average data usage over sequential-test parameters, but cannot retune the test at every iteration because error depends on the proposed states.
  • Test design: Worst-case error yields a conservative design that wastes the sequential test’s efficiency, motivating optimization against the expected error over parameter distributions.
  • Test design: The method can use flexible threshold sequences beyond the constant Pocock bound, including designs parameterized by π_1, G_0, and α.

F Application to Gibbs Sampling

The approximate Gibbs sampler applies sequential testing to binary-variable updates, replacing exact conditional-probability evaluation with an error-controlled approximation. Under a Dobrushin condition, the stationary-distribution error is bounded.

  • Sequential testing is applied to Gibbs updates for binary variables, with extension to multi-valued variables also possible.
  • Each update computes P(X_i = 1|x_-i), draws u from Uniform[0, 1], and assigns X_i according to the resulting threshold test.Here x_-i denotes all variables except the ith variable.
  • When the joint distribution factors into multiple terms, the sequential test can speed Gibbs sampling by approximating the conditional probability from those factors.
  • For a Gibbs sampler with Dobrushin coefficient η ∈ [0, 1), the distance between exact and approximate stationary distributions is upper bounded.

F.1 Experiments on Markov Random Fields

Experiments on a densely connected synthetic binary Markov random field compare exact Gibbs sampling with approximate variants across error tolerances. Smaller tolerances initially reduce error faster through more frequent sampling, but exact Gibbs eventually performs best.

  • The experiment uses D = 100 densely connected binary variables with D(D −1)(D −2)/6 three-variable potential functions.Each potential function has 2^3 = 8 values.
  • Approximate methods process mini-batches of 500 potential-function pairs and are compared at ϵ ∈ {0.01, 0.05, 0.1, 0.15, 0.2, 0.25}.
  • Performance is evaluated using average L1 distance over M = 1600 randomly selected subsets of 5 variables because the full joint distribution has 2^100 probabilities.
  • At ϵ = 0.01, the observed maximum conditional-probability error is within 0.01.
  • For small ϵ, approximate methods initially reduce error faster by generating more samples per unit time, but exact Gibbs outperforms them at around 1000 seconds.
Loading 1304.5299v4…