George-McCulloch (1993) Variable Selection Via Gibbs Sampling

Summary

George and McCulloch (1993) introduce Stochastic Search Variable Selection (SSVS), a three-block Gibbs sampler over regression coefficients $\beta$, error variance $\sigma^2$, and a binary inclusion vector $\gamma \in \{0,1\}^p$. The key innovation is a mixture prior that makes each $\beta_i$ effectively zero when $\gamma_i=0$ (spike component $\mathcal{N}(0,\tau_i^2)$ with small $\tau_i$) and unconstrained when $\gamma_i=1$ (slab component $\mathcal{N}(0,c_i^2\tau_i^2)$ with $c_i \gg 1$). Visiting patterns of $\gamma$ across all post-burn-in Gibbs output — the "all post-burn-in" extraction strategy — identify high-posterior-probability subsets without enumerating the $2^p$ model space. The method is demonstrated on examples up to $p=60$ regressors.

Key Claims

• The mixture prior $\beta_i|\gamma_i \sim (1-\gamma_i)\mathcal{N}(0,\tau_i^2) + \gamma_i\mathcal{N}(0,c_i^2\tau_i^2)$ (eq. 2) with $c_i>1$ shrinks excluded coefficients toward zero without hard zeroing; $P(\gamma_i=1)=p_i$ (eq. 3).
• Joint prior: $\beta|\gamma \sim \mathcal{N}_p(0, D_\gamma R D_\gamma)$ (eq. 4), $D_\gamma = \mathrm{diag}[a_1\tau_1,\ldots,a_p\tau_p]$ where $a_i=1$ if $\gamma_i=0$ and $a_i=c_i$ if $\gamma_i=1$ (eq. 5); $\sigma^2|\gamma \sim \mathrm{IG}(\nu_\gamma/2, \nu_\gamma\lambda_\gamma/2)$ (inverse-gamma; eq. 6).
• Three-block Gibbs: (1) $\beta^j \sim \mathcal{N}_p(A_\gamma(\sigma^{j-1})^{-2}X'X\hat\beta_{LS}, A_\gamma)$ (with $\hat\beta_{LS}$ the least-squares estimate) where $A_\gamma = ((\sigma^{j-1})^{-2}X'X + D_\gamma^{-1}R^{-1}D_\gamma^{-1})^{-1}$ (eq. 13); (2) $\sigma^j \sim \mathrm{IG}((n+\nu_\gamma)/2, (\|Y-X\beta^j\|^2+\nu_\gamma\lambda_\gamma)/2)$ (eq. 14); (3) $\gamma_i^j \sim \mathrm{Bernoulli}(a/(a+b))$ componentwise (eqs. 15–16).
• The $\gamma_i$ update does NOT depend on $Y$ directly: $a = f(\beta^j|\gamma_{(i)}, \gamma_i=1)\cdot p_i$ and $b = f(\beta^j|\gamma_{(i)}, \gamma_i=0)\cdot(1-p_i)$ (eqs. 16c–d) — the likelihood of $Y$ given $\beta$ is already absorbed into the current $\beta^j$ draw.
• Calibration recommendations: $(\hat\sigma_{\hat\beta_i}/\tau_i,\, c_i) \in \{(1,5),(1,10),(10,100),(10,500)\}$; ratio $r_i = \sqrt{(\hat\sigma^2_{\hat\beta_i}/\tau_i^2 + c_i^2)/(\hat\sigma^2_{\hat\beta_i}/\tau_i^2 + 1)}$ measures the Bayes factor signal in each coefficient (eq. 9).
• High-frequency $\gamma$ values in Gibbs output identify high-probability subsets; the Gibbs output provides a J-shaped frequency distribution over $2^p$ models — most models unvisited.
• Small example ($p=5$, $n=60$, true model $Y=X_4+1.2X_5+\varepsilon$): true model ranked first at 25.8% frequency. Large example ($p=60$, $n=120$): method navigates $2^{60}$ model space efficiently.
• Correlation structure: $R=I$ (independence prior) or $R \propto (X'X)^{-1}$ (g-prior, Zellner 1986); prior $f(\gamma) = 2^{-p}$ (indifference) or model-size-weighted version (eq. 7).

Concepts Introduced or Extended

• Spike-and-Slab Prior
• Variable Selection
• Gibbs Sampler

Entities Mentioned

• Edward I. George
• Robert E. McCulloch

Quotes

"We propose using the Gibbs sampler to search those regions of the model space with high posterior probability."

"Subset selection... is an important component of regression analysis since it provides parsimonious models and helps guard against overfitting."

My Take

SSVS is elegant because it converts the discrete model-selection problem into a continuous Bayesian computation via the latent $\gamma$ indicators — a pure data-augmentation Gibbs scheme. The three-block structure is fully conjugate given $\gamma$: both $\beta$ and $\sigma^2$ draw from standard families (Normal and Inverse-Gamma). The $\gamma$ block is Bernoulli and requires only evaluating the ratio of two Normal densities at the current $\beta^j$. The computational cost per iteration is dominated by the $p \times p$ matrix inversion in eq. 13, which scales as $O(p^3)$; this can be expensive for large $p$ but was workable for $p=60$ in 1993. The tuning parameters $(\tau_i, c_i)$ require calibration — the suggested choices in Table 1 provide a practical starting point. The method does not produce a point estimate of $\gamma$ directly; instead the researcher reads off high-frequency patterns from the Gibbs output. This "stochastic search" framing was influential: it predates modern spike-and-slab and horseshoe priors by framing variable selection as posterior inference rather than optimization.