Binary Probit

Definition

The binary probit model relates a binary outcome $Y_i \in \{0,1\}$ to covariates $x_i$ via $\Pr(Y_i = 1) = \Phi(x_i'\beta)$ , where $\Phi$ is the standard Normal cumulative distribution function (CDF) and $\beta$ is a $k$ -vector of coefficients. The standard parametrization assumes homoskedastic unit-variance latent error, so the coefficient scale is fixed at $\sigma = 1$ for identification.

Key Ideas

Data augmentation (Albert-Chib 1993b): introduce latent $Z_i \sim \mathcal{N}(x_i'\beta, 1)$ with $Y_i = \mathbf{1}[Z_i > 0]$ . Conditional on $Z$ , the model reduces to a standard Normal linear regression.
The latent $Z_i$ representation unifies binary, ordered, and unordered multinomial probit — each case differs only in how the outcome constraint truncates $Z_i$ .
Bayesian Residual Analysis (Albert-Chib 1995): both parametric residuals $r_i = y_i - \Phi(x_i'\beta)$ and latent data residuals $\varepsilon_i = Z_i - x_i^T\beta$ are real-valued with continuous posteriors. The latent residual is $\mathcal{N}(0,1)$ a priori, giving a calibrated reference. Both emerge as free by-products of the Gibbs run.
Alternative links (logit, t-link, complementary log-log) can be embedded via scale mixtures without changing the Gibbs framework.

How It Works

Two-Block Gibbs Sampler (Albert-Chib 1993b, eqs. 5–6)

Under a flat prior on $\beta$ , the joint posterior augmented by $Z = (Z_1,\ldots,Z_n)'$ admits two tractable full conditionals:

Block 1 — $\beta | Z, y$ : $\beta | Z, y \sim \mathcal{N}_k\!\left(\hat\beta_Z,\; (X'X)^{-1}\right), \qquad \hat\beta_Z = (X'X)^{-1}X'Z$

This is exactly the ordinary least squares (OLS) posterior for a Normal linear model with the latent $Z$ as the response.

Block 2 — $Z_i | y_i, \beta$ : $Z_i | y_i, \beta \sim \mathcal{N}(x_i'\beta,\, 1) \text{ truncated to } \begin{cases}(0,\infty) & \text{if } y_i = 1 \\ (-\infty, 0] & \text{if } y_i = 0\end{cases}$

Truncated Normal draws are obtained by inversion of the Normal CDF. The Gibbs cycle alternates between these two blocks. Albert and Chib (1993b) report convergence to the exact posterior at 800 iterations for Finney's (1947) toxicity dataset (39 observations, 3 covariates), verified against numerical integration.

Robust t-Link via Scale Mixture (eqs. 7–11)

To accommodate heavier-tailed link functions, introduce scale mixture weights $\lambda_i$ : $Z_i | \lambda_i, \beta \sim \mathcal{N}(x_i'\beta,\, \lambda_i^{-1}), \qquad \lambda_i \sim \text{Gamma}(\nu/2,\, 2/\nu)$

Marginalizing over $\lambda_i$ gives a $t_\nu$ -link. The four-block Gibbs adds draws of $\beta$ (weighted Normal with $W = \text{diag}(\lambda_i)$ ), $Z_i$ (truncated weighted Normal), $\lambda_i$ (Gamma), and $\nu$ (discrete on $\{4, 8, 16, 32\}$ ). The logistic link is well-approximated by $t(8)$ .

Hierarchical Prior (eqs. 13–14)

Placing $\beta \sim \mathcal{N}(A\beta^0, \sigma^2 I)$ introduces a hyperparameter $\sigma^2$ controlling the degree of shrinkage. The posterior $\beta | Z, \sigma^2$ remains Normal with an augmented design matrix; $\sigma^2 | Z$ is non-standard but can be grid-sampled. This enables Bayesian variable selection: regional effects with small signals are shrunk toward zero.

Ordered Multinomial Extension (eqs. 16–18)

For $J$ ordered categories, $Y_i = j$ iff $\gamma_{j-1} < Z_i \leq \gamma_j$ , with $\gamma_0 = -\infty$ , $\gamma_J = \infty$ , $\gamma_1 = 0$ (identification constraint). Three-block Gibbs:

$\beta | Z$ — same Normal draw as binary case (eq. 5)
$Z_i | \beta, \gamma, y$ — truncated Normal on $(\gamma_{y_i - 1},\, \gamma_{y_i})$
$\gamma_j | Z, \beta, y \sim \text{Uniform}\!\left[\max\{Z_i : Y_i = j\},\; \min\{Z_i : Y_i = j+1\}\right]$

The threshold update (block 3) is exact and requires no tuning parameters.

Unordered Multinomial Probit (eqs. 19–21)

For $J$ unordered categories, $\mathbf{Z}_i = (Z_{i1},\ldots,Z_{iJ})' \sim \mathcal{N}_J(X_i\beta, \Sigma)$ with $Y_i = j$ iff $Z_{ij} > Z_{ik}$ for all $k \neq j$ . The Gibbs draw for $\mathbf{Z}_i$ uses accept-reject: propose from $\mathcal{N}_J(X_i\beta, \Sigma)$ and accept if the observed category is the argmax. Applied to Daganzo's (1977) trivariate probit with 50 observations: posterior mean correlation $\hat\rho = 0.234$ vs. maximum likelihood estimate (MLE) 0.238.

Mixture of Normals Probit Extension (Geweke-Keane 1997)

Geweke and Keane (1999) extend Albert-Chib (1993b) by replacing the fixed $\mathcal{N}(0,1)$ disturbance with a finite $m$ -component Gaussian mixture $p(\varepsilon_t) = \sum_{j=1}^m p_j h_j^{1/2}\phi(h_j^{1/2}(\varepsilon_t - \alpha_j))$ . The Albert-Chib latent utility augmentation $\tilde{y}_t$ is retained; a second latent variable, the component-assignment indicator $L_t \in \{1,\ldots,m\}$ , is added. Conditional on $L_t = j$ , the model reduces to a standard probit with mean $\alpha_j$ and precision $h_j$ , so $\tilde{y}_t$ 's full conditional remains a univariate truncated Normal. The result is a six-block fully conjugate sampler. Three identified variants: scale mixture ( $\alpha_j = 0$ , heavy tails, symmetric), full mixture ( $\alpha_j$ and $h_j$ both free, heavy tails, asymmetric), mean mixture ( $h_j = 1$ , asymmetric). Identification requires one precision fixed ( $h_{j^*}=1$ ) and a labeling restriction on either $\alpha$ or $h$ .

Applied to PSID (Panel Study of Income Dynamics) female labor force participation (T=1,555): Bayes factors of order $10^7$ favor the four-component scale mixture over standard probit. The dominant failure mode of probit in this dataset is thin tails (symmetry without kurtosis), not asymmetry. Predictive distributions for welfare benefit and spouse income effects differ materially between models. See Geweke-Keane (1999) and Mixture of Normals.

Why It Matters

Makes full Bayesian posterior inference for probit models computationally feasible where MLE can fail (near-collinearity, sparse cells, separation problems).
The data augmentation principle extends to any model where the latent variable structure is Gaussian conditional on observables — e.g., Tobit (Chib 1992), selection models, panel probit.
Foundational for modern Bayesian binary regression: every subsequent Gibbs-based probit sampler (hierarchical, random effects, spatial) builds on the same two-block structure.
Bayesian residuals provide a natural, continuous outlier-detection tool unavailable in classical probit.

Open Questions

Efficiency of the accept-reject multivariate truncated Normal sampler degrades for large $J$ (many competitors). Later work (McCulloch-Rossi 1994) proposed more efficient algorithms.
With strong prior-data conflict or near-separation, the truncated Normal draws for $Z_i$ can be extreme; robust priors and reparameterizations help.
Scale mixture representation works for symmetric heavy-tailed links; asymmetric links (complementary log-log) require a different augmentation structure.