## How MCMC Solves the Posterior MCMC (Markov Chain Monte Carlo) is the core computational technique that made modern Bayesian statistics possible. Here's how it works: ### The Core Problem: Bayes' Theorem's Denominator Bayes' Theorem gives us the posterior distribution: ``` P(θ | D) = P(D | θ) P(θ) / ∫ P(D | θ) P(θ) dθ ``` The **numerator** — likelihood × prior — is easy to compute for any specific θ. The **denominator** (the marginal likelihood / evidence) is an integral over *every possible value* of θ. For all but the simplest conjugate models (e.g., Beta-Binomial), this integral has no closed-form solution. Numerical integration suffers from the **curse of dimensionality**: a grid with 100 points per dimension needs 100 evaluations for 1 parameter, 10,000 for 2, and 10^20 for 10 parameters. ### The MCMC Trick: Cancel the Denominator MCMC sidesteps the denominator entirely. Every MCMC algorithm works with **ratios** of the posterior at two points: ``` P(θ* | D) / P(θ⁽ᵗ⁾ | D) = [P(D | θ*) P(θ*)] / [P(D | θ⁽ᵗ⁾) P(θ⁽ᵗ⁾)] ``` The denominator ∫P(D|θ)P(θ)dθ appears in both the top and bottom and **divides out**. MCMC only ever evaluates the unnormalized posterior — the shape — never the full probability. ### The 6-Step Workflow **Step 1 — Specify the prior** P(θ). Choose a distribution reflecting what is known about each parameter before seeing the data (informative or vague). **Step 2 — Specify the likelihood** P(D|θ). Choose the probability model describing how the data arises (Binomial, Normal, Poisson, etc.). **Step 3 — Form the unnormalized posterior** P(θ|D) ∝ P(D|θ)P(θ). This product is the function every MCMC algorithm evaluates. It has the same shape as the true posterior but doesn't integrate to 1 — which doesn't matter. **Step 4 — Initialize the chain and propose a move.** Start at some value θ⁽⁰⁾. At each step t, generate a candidate θ* from a **proposal distribution** q(θ* | θ⁽ᵗ⁾) — e.g., current value plus a small random Normal step. **Step 5 — Accept or reject the proposal.** Compute an acceptance probability comparing the unnormalized posterior at θ* and θ⁽ᵗ⁾. Accept θ* with that probability; otherwise stay at θ⁽ᵗ⁾ for another step. **Step 6 — Discard burn-in and summarize.** Drop early iterations before the chain reached the high-probability region. Use the remaining samples to compute posterior means, standard deviations, and credible intervals. ### Why It Works: The Key Intuition The chain "wanders" through parameter space. It **accepts moves that go uphill** (higher posterior density) and **sometimes accepts moves that go downhill** (lower density), with a probability proportional to how much lower. Over enough steps, the chain spends more time in regions where the posterior is high and less time where it's low. A histogram of where the chain visited **approximates the posterior distribution** — from which you can read off means, standard deviations, and credible intervals. ### The Main Algorithms | Algorithm | How it proposes moves | Key feature | |---|---|---| | **Metropolis-Hastings** (1953/1970) | Random walk — add a small Normal step | Simple, works for any model; needs tuning of step size (target ~20-50% acceptance) | | **Gibbs Sampling** (1984) | Update one parameter at a time from its full conditional distribution | Every proposal accepted; needs conditionals to be derivable | | **Hamiltonian Monte Carlo** (HMC) | Uses gradient of log-posterior to simulate a physics trajectory | Follows the gradient — moves far while staying in high-probability regions; much less autocorrelation | | **NUTS** (No-U-Turn Sampler) | Automates HMC's trajectory length via binary tree expansion | Removes manual tuning; used in Stan and PyMC | ### Practical Checks - **Burn-in / warmup**: Early samples are discarded because the chain hasn't converged yet - **Trace plots**: Visual check — the chain should look like "fuzzy caterpillars," not trends - **Gelman-Rubin (R-hat)**: Run multiple chains; R-hat near 1.0 means convergence - **Effective sample size**: Accounts for autocorrelation — independent draws are more efficient In summary: **MCMC never computes the posterior directly.** Instead, it constructs a Markov chain whose stationary distribution *is* the posterior, by evaluating only *ratios* of the unnormalized posterior — making the impossible integral vanish. That breakthrough in 1990 (Gelfand & Smith) is what launched the modern era of applied Bayesian statistics.