Gibbs sampling

The Gibbs sampler updates each variable in turn, drawing from its full conditional distribution given the most recent values of the other variables. For a target distribution p(x₁,...,x_d), the algorithm cycles through i = 1,...,d, sampling x_i^(t+1) ~ p(x_i | x₁^(t+1),...,x_i-1^(t+1), x_i+1^(t),...,x_d^(t)). Under mild regularity conditions (irreducibility and aperiodicity), the chain converges to the target distribution as t → ∞.¹ Convergence can be slow when variables are highly correlated, leading to poor mixing. Diagnostics such as the Gelman-Rubin statistic assess convergence by comparing multiple chains.²

Gibbs sampling is a cornerstone of Bayesian computation, enabling inference in hierarchical models, mixture models, and latent Dirichlet allocation (LDA).³ In computational physics, it is used to simulate spin systems like the Ising model, where each spin is updated based on its neighbors.⁴ In genetics, it infers population structure from genotype data, as in the STRUCTURE software.⁵ It also underpins probabilistic graphical models, including Bayesian networks and Markov random fields, for tasks such as image denoising and missing data imputation.

Several variants improve efficiency. Blocked Gibbs sampling updates groups of variables jointly, reducing correlation. Collapsed Gibbs sampling integrates out some variables analytically, often accelerating convergence in topic models.³ The Metropolis-within-Gibbs algorithm substitutes a Metropolis-Hastings step when full conditionals are not tractable. For high-dimensional problems, the Gibbs sampler can be combined with Hamiltonian Monte Carlo or used in a tempering scheme to escape local modes. The algorithm also appears in nonparametric Bayesian methods, such as the Dirichlet process mixture model.

The method's name honors Josiah Willard Gibbs, but he never worked on it; the term was coined by Geman and Geman in 1984, who applied it to Bayesian image restoration.¹ Earlier, the idea appeared in the physics literature as the heat bath algorithm, used for lattice spin systems.⁴ A notable edge case: when variables are independent, Gibbs sampling reduces to direct sampling, but the chain still requires burn-in. The algorithm is a special case of the Metropolis-Hastings algorithm with acceptance probability always 1, which is why it is sometimes called a "component-wise Metropolis-Hastings."⁶ In practice, the order of updates (systematic vs. random) can affect convergence speed, though both yield the same stationary distribution.

Full conditional distribution

The probability distribution of one variable given the current values of all other variables.

Mixing

The rate at which the Markov chain explores the support of the target distribution.

Burn-in

The initial set of samples discarded to allow the chain to reach stationarity.