Mixture of Normals

Definition

A mixture of normals (finite Gaussian mixture) represents a random variable as a convex combination of $K$ Gaussian components: $f(x) = \sum_{k=1}^K p_k\, \mathcal{N}(x;\, \mu_k,\, \sigma_k^2)$ , where the mixing weights $p_k > 0$ sum to one. Conditional on a latent component indicator $s \in \{1,\ldots,K\}$ , the variable is Gaussian; marginally it can capture skewness, excess kurtosis, multimodality, and heavy tails. In Bayesian time-series econometrics the mixture is used to approximate analytically intractable non-Gaussian distributions or to flexibly model shock distributions without strong parametric assumptions.

Key Ideas

The mixture likelihood is unbounded without regularization: as $\sigma_k^2 \to 0$ for any observation in component $k$ , the likelihood diverges. Proper Bayesian priors (e.g., Gamma priors on precision $h_k = 1/\sigma_k^2$ ) restore boundedness and guarantee a proper posterior; this is the Bayesian counterpart of the Kiefer-Wolfowitz (1956) consistency result for the nonparametric maximum likelihood estimator (NPMLE).
Latent component indicators $s_t \in \{1,\ldots,K\}$ make the model conditionally Gaussian, enabling efficient Gibbs sampling: draw $s_t$ from a discrete distribution, then draw parameters from conjugate Gaussian and Gamma conditionals.
An identification constraint prevents label switching: common choices are $\mu_1 \leq \mu_2 \leq \ldots \leq \mu_K$ , or fixing one mean (e.g., $\mu_2 = 0$ as in Geweke-Keane 2000).
Increasing $K$ increases flexibility but also the risk of spurious components; $K = 3$ and $K = 7$ are the most common choices in applied work.

How It Works

Bayesian Gibbs Sampler for a $K$ -Component Mixture

Let $x_1,\ldots,x_n \overset{\text{iid}}{\sim} \sum_{k=1}^K p_k\, \mathcal{N}(\mu_k,\, h_k^{-1})$ with latent indicators $s_i \sim \text{Multinomial}(p_1,\ldots,p_K)$ .

Full conditionals (standard conjugate forms):

Indicators $s_i \mid x_i, \mu, h, p$ : discrete with $P(s_i = k) \propto p_k\, h_k^{1/2} \exp\{-h_k(x_i - \mu_k)^2/2\}$ .
Weights $p \mid s$ : Dirichlet with updated counts $n_k + \alpha_k$ where $n_k = \#\{i : s_i = k\}$ .
Means $\mu_k \mid s, x, h_k$ : Normal with posterior precision $h_k n_k + h_0$ (subject to ordering constraint via truncated draws).
Precisions $h_k \mid s, x, \mu_k$ : Gamma with updated shape $a + n_k/2$ and rate $b + \sum_{i:s_i=k}(x_i - \mu_k)^2/2$ .

The proper Gamma prior on $h_k$ is essential: without it the conditional for $\sigma_k^2 = h_k^{-1}$ has a pole near zero for any non-empty component ( $n_k \geq 1$ ), causing the Markov Chain Monte Carlo (MCMC) chain to drift toward degenerate components.

7-Component Mixture for $\log \chi^2(1)$ (Kim-Shephard-Chib [KSC] 1998)

In stochastic volatility models, the log-squared return satisfies $\log y_t^2 = h_t + \log \varepsilon_t^2$ , where $\log \varepsilon_t^2 \sim \log \chi^2(1)$ . The KSC approximation replaces this non-Gaussian distribution with a 7-component normal mixture chosen to match moments; the mixture indicators $s_t$ are drawn as a Gibbs block, converting the stochastic volatility (SV) estimation problem into a Gaussian state-space that the simulation smoother can handle directly. See Stochastic Volatility.

Conditional Mixture — MN-GARCH (Haas-Mittnik-Paolella 2004a; Bauwens-Rombouts 2005)

The mixture-of-normals framework extends naturally to time-varying conditional distributions: in the MN-GARCH model each component has its own generalized autoregressive conditional heteroskedasticity (GARCH)(1,1) variance process $h_{k,t} = \omega_k + \alpha_k y_{t-1}^2 + \beta_k h_{k,t-1}$ , so the conditional mixture is

$f(y_t|\mathcal{F}_{t-1}) = \sum_{k=1}^K \pi_k\, \frac{\phi\!\bigl((y_t-\mu_k)/\sqrt{h_{k,t}}\bigr)}{\sqrt{h_{k,t}}}$

The latent state $S_t \in \{1,\ldots,K\}$ is drawn at each $t$ conditional on past data and current GARCH variances; conditional on $S_T$ , the GARCH parameters $(\omega_k,\alpha_k,\beta_k)$ are updated component by component via griddy-Gibbs. The mixing weights follow a Dirichlet conjugate update. Weak stationarity can hold even when individual components are explosive ( $\alpha_k + \beta_k > 1$ ) because the stationarity condition involves all components jointly.

Key insight (Bauwens-Rombouts 2005): Applied to S&P 500 returns, $K=2$ captures a dominant near-stationary component ( $\pi \approx 0.79$ , $\hat\alpha+\hat\beta \approx 0.98$ ) and a minority explosive component ( $\pi \approx 0.21$ , $\hat\alpha+\hat\beta \approx 1.14$ ). The standard GARCH(1,1) persistence estimate of $\approx 0.996$ is a weighted average of these two levels — near-integrated GARCH (IGARCH) is a misspecification artifact. See Bauwens-Rombouts (2005) and GARCH and BEKK-GARCH.

Scale-Mixture Representation: Student-t Distribution (Geweke 1993)

Geweke (1993) establishes the foundational application of the scale-mixture device to Bayesian regression with leptokurtic disturbances. If $y_i \sim t(x_i'\beta, \sigma^2; \nu)$ (i.i.d. Student-t with $\nu$ degrees of freedom), this is identical to

$y_i \mid x_i, \omega_i \sim \mathcal{N}(x_i'\beta,\; \sigma^2 \omega_i), \qquad \nu/\omega_i \sim \chi^2(\nu)$

with $n$ latent scale weights $\omega_i > 0$ . This is a scale mixture of normals: each observation is Gaussian conditional on its own variance scale $\sigma^2 \omega_i$ , and the marginal distribution integrating out $\omega_i$ is exactly Student-t. Augmenting the parameter space with $\omega = (\omega_1,\ldots,\omega_n)$ converts a non-conjugate problem into one where every full conditional is conjugate:

$\omega_i \mid \beta, \sigma, \nu$ : Inverse-Gamma — $[(\sigma^{-2}u_i^2 + \nu)/\omega_i] \mid (\beta,\sigma) \sim \chi^2(\nu+1)$ (eq. 15).
$\sigma^2 \mid \beta, \omega$ : Inverse-Gamma — $[\sum_i u_i^2/\omega_i]/\sigma^2 \sim \chi^2(n)$ (eq. 14).
$\beta \mid \sigma, \omega$ : Gaussian — generalized least squares (GLS) form $\hat\beta(\omega) = (X'\Omega^{-1}X)^{-1}X'\Omega^{-1}y$ (eq. 12).
$\nu \mid \beta, \sigma, \omega$ : non-standard kernel (eq. 18), requiring specialized sampling or a fixed grid.

This four-block Gibbs sampler is the earliest published conjugate-augmentation scheme for heavy-tailed regression, and the template that Jacquier-Polson-Rossi (1994) carried into the time-varying SV context. The key difference from later applications is that the $\omega_i$ are cross-sectional (one per observation, i.i.d.) rather than time-varying. See Geweke (1993) and Stochastic Volatility.

Mixture of Normals Probit (Geweke-Keane 1997)

Geweke and Keane (1999) apply the finite Gaussian mixture to the disturbance term of a binary probit model, obtaining a flexible semiparametric-style generalization that retains exact Bayesian finite-sample inference. The disturbance is $p(\varepsilon_t) = \sum_{j=1}^m p_j h_j^{1/2}\phi(h_j^{1/2}(\varepsilon_t - \alpha_j))$ . Three variants differ by which parameters are free: the scale mixture ( $\alpha_j = 0$ ) imposes symmetry but allows heavy tails; the full mixture allows both asymmetry and heavy tails; the mean mixture ( $h_j = 1$ ) allows asymmetry only.

The Gibbs sampler augments the Albert-Chib (1993b) latent utility $\tilde{y}_t$ with a component-assignment indicator $L_t \in \{1,\ldots,m\}$ . Conditional on $L_t = j$ : (1) $\tilde{y}_t$ is a univariate truncated Normal with mean $\alpha_j + \beta'x_t$ and precision $h_j$ ; (2) $L_t$ is drawn from a multinomial with weights $\propto p_j \exp[-0.5 h_j (\tilde{y}_t - \alpha_j - \beta'x_t)^2]$ ; (3) $\beta$ is Normal; (4) $\alpha$ is Normal with ordering constraints (Geweke 1991/1995); (5) $h_j \sim$ inverse-chi-squared; (6) $p \sim$ Dirichlet. All six blocks are conjugate; the structure is the direct discrete-choice analogue of the MN-GARCH component-assignment block in Bauwens-Rombouts (2005).

Marginal likelihood is evaluated via the Gelfand-Dey (1994) estimator on the observed-data parameter space $(\beta, \alpha, h, p)$ , marginalizing over $(\tilde{y}, L)$ analytically.

Panel Study of Income Dynamics (PSID) female labor force participation (T=1,555, 17 covariates): best model is the scale mixture of 4 normals ( $\log M=-548.3$ ); conventional probit $\log M=-566.6$ ; Bayes factor (BF) $\approx 10^7$ against conventional probit. The dominant shock distribution is symmetric and leptokurtic (not skewed), which affects counterfactual welfare-benefit and spouse-income elasticities.

Connection to Geweke (1993). The scale mixture probit ( $\alpha_j = 0$ ) with $m\to\infty$ converges to the Student-t probit — the same representation studied statically in Geweke (1993). The $m$ -component scale mixture is a finite-dimensional parametric approximation to that limiting case. See Geweke-Keane (1999) and Binary Probit.

Application to Earnings Dynamics (Geweke-Keane 2000)

Geweke and Keane (2000) use a 3-component mixture with $K = 3$ , parameterized as $\varepsilon \sim \sum_{k=1}^3 p_k\, \mathcal{N}(a_k,\, h_k^{-1})$ for first-period shocks $e_{i1}$ and transitory increments $g_{it}$ separately. The normalization fixes $a_2 = 0$ and imposes $a_1 \leq a_2 \leq a_3$ . Latent component indicators $s_{it}$ are added as Gibbs blocks alongside the model parameters $(\beta, c, \rho, \phi, q_i)$ . Data augmentation fills unobserved panel periods with latent earnings. The resulting posterior shocks are left-skewed and leptokurtic. Crucially, these shapes substantially alter policy-relevant quantities: a 3% vs. 12% estimated education return, and a $94K vs. $150K racial lifetime earnings gap, depending on whether shocks are Gaussian or mixture.

Competitive EM for Mixture Model Selection (Zhang-Zhang-Yi 2004)

Zhang, Zhang, and Yi (2004) address three simultaneous deficiencies of standard Expectation-Maximization (EM) for finite mixture models: local maxima traps, boundary convergence of degenerate components, and the need to pre-specify $k$ . Their Competitive EM (CEM) algorithm replaces uniform operation-type sampling with competitive probabilities derived from local Kullback divergence: $J(m;\theta) = \int f_m(x;\theta)\log[f_m(x;\theta)/p_m(x;\theta_m)]\,dx$ , where $f_m$ is the empirical distribution weighted by posterior component assignments. A component with high $J$ is poorly explained by its kernel and should split; two components with low post-merge $J$ should merge.

CEM operation probabilities (eqs. 9–10):

$P_\text{split}(m;\theta) \propto J(m;\theta)$ — split probability proportional to local fit divergence
$P_\text{merge}(m,l;\theta) \propto \beta/J(m';\theta')$ — merge probability inversely proportional to merged-component divergence
$P_a = \min(\exp\{[L(\theta(t+1),\mathbf{X}) - L(\theta(t),\mathbf{X})]/\gamma\}, 1)$ with $\gamma=10$ — Metropolis acceptance allows uphill moves

Annihilation: component $m$ dies when $n\alpha_m < N$ (fewer samples than parameters); CEM² sequential updating (Celeux et al. 1999) prevents cascading multi-component death by updating components one at a time so survivors can absorb orphaned samples.

Minimum Message Length (MML)-based criterion (from Figueiredo-Jain 2002): $L(\theta,\mathbf{X}) = \log p(\mathbf{X}|\theta) - \frac{N}{2}\sum_m \log(n\alpha_m/12) - \frac{k}{2}\log(n/12) - k(N+1)/2$ .

Key experimental result — unequal-prior 4-component Gaussian mixture model (GMM) ( $\alpha_4=0.1$ , $n=1000$ , 1000 runs): CEM $P_s=90.2\text{–}96.9\%$ , AEM $P_s=60.3\text{–}71.0\%$ , EM $P_s=19.0\text{–}59.2\%$ (Table 2). AEM fails because its global MML penalty is not sensitive to per-component fit quality; CEM's $J(m;\theta)$ correctly targets the under-represented small-prior component for merging/annihilation. With equal priors (8-GMM), CEM and AEM both achieve $\approx 100\%$ . Application: GMM fitted to 2D pixel coordinates of Chinese characters extracts strokes including overlapping radicals. See Zhang-Zhang-Yi (2004).

Why It Matters

Allows flexible non-Gaussian modeling while retaining conjugate computational structure: conditionally Gaussian given indicators.
Critical for stochastic volatility estimation (7-component KSC approximation reduces SV sampling to a Gaussian state-space problem) and for panel earnings models where distributional assumptions affect policy-relevant statistics.
The regularization via proper priors is a general solution to the label-switching and unbounded-likelihood problems that afflict maximum-likelihood mixture estimation.
Demonstrates that normality assumptions in panel data are not innocuous: Geweke-Keane (2000) shows that distributional shape matters for poverty persistence, lifetime earnings gaps, and returns to education.

Open Questions

Choosing $K$ is unsolved in general: information criteria can over-count components; Bayesian marginal likelihood comparisons (e.g., Laplace approximation as in Bauwens-Rombouts 2005, or reversible-jump MCMC) are the principled alternatives but computationally demanding.
Dirichlet process (nonparametric Bayesian) mixtures sidestep the choice of $K$ but are harder to interpret and less tractable.
Label switching is benign for posterior summaries over the full support of a single ergodic chain but can masquerade as multimodality when multiple chains are analyzed separately.
For conditional mixtures (MN-GARCH), griddy-Gibbs for non-conjugate GARCH parameters scales poorly to $K>3$ or large $T$ .