van Dyk (2000) Fitting Mixed-Effects Models Using Efficient EM-Type Algorithms

Summary

Develops a family of fast, numerically stable Expectation-Maximization (EM)-type algorithms for maximum likelihood (ML), restricted maximum likelihood (REML), and MAP (posterior mode) estimation in Gaussian linear mixed-effects models. Combines three methodological advances — ECME (Expectation-Conditional-Maximization-Either), conditional data augmentation, and PXEM (Parameter-Expanded EM) — to build algorithms that converge uniformly fast across signal-to-noise ratios while retaining monotone convergence guarantees absent from standard implementations (SAS proc mixed, S-Plus lme). Also extends the Hobert-Casella (1996) posterior-propriety result to the general multivariate mixed-effects model.

Key Claims

• Model. $Y_i = X_i\beta + Z_i b_i + e_i$ with $b_i \overset{iid}{\sim} N(0, \sigma^2 D)$ and $e_i \overset{iid}{\sim} N(0, \sigma^2 R_i)$, $i=1,\ldots,m$. Conjugate priors $(\beta|\sigma^2)\sim N_p$, $\sigma^2\sim\chi^{-2}$, $D\sim$ inverse-Wishart.
• Fraction-of-missing-information. The convergence rate of EM is controlled by the spectral radius of $\mathcal{D}_\text{EM}(\alpha)=I-I_\text{obs}\,I_\text{aug}^{-1}(\alpha)$, where $\alpha$ is a working parameter in the data-augmentation scheme $p(Y_\text{aug}|\theta,\alpha)$ (with $p(Y_\text{obs}|\theta)$ unchanged). Choosing $\alpha$ to minimise this eigenvalue is the conditional augmentation principle (Meng-van Dyk 1997).
• ECME algorithms. The M-step for $\sigma^2$ is replaced by a conditional maximisation of the observed-data log-posterior (no data augmentation), yielding a closed-form update that is faster per iteration. Three parameterisations of the random-effects covariance give ECMEo ($\alpha=0$, augment with $b_i$), ECME1 ($\alpha=1$, augment with $L^{-1}b_i$, Cholesky-scaled), and ECME2 (arbitrary lower-triangular $L$). Relative efficiency depends on $\Lambda^* = \sigma^{2*}/\text{mean eigenvalue of }Z_iT^*Z_i^T$.
• Parameter-expanded ECME2 (PX-ECME2). For REML and improper-prior MAP, expand $L$ to an arbitrary invertible $q\times q$ working matrix, fit $L$ inside the EM iteration (Liu-Rubin-Wu 1998 PXEM idea), then discard it. This algorithm is uniformly fast across all values of $\Lambda^*$ and competitive with SAS/S-Plus in speed while maintaining monotone convergence.
• PXEM failure under proper priors. When conjugate proper priors are used, the PXEM objective $Q_{px}(\theta,\alpha|\theta',\alpha_0)$ becomes unbounded for improper-like parameter values; optimising it without nested iterations is infeasible. Conditional augmentation (not parameter expansion) is required for proper Bayesian MAP.
• Posterior propriety (Theorem 1). Under the REML improper prior, with $r=\text{rank}(P_XZ)$, sufficient conditions for a proper posterior are: (a) $\eta<1-q$, (b) $m>q-1-r$, and (c) $n>p-qr-\nu$. Extends Hobert-Casella (1996) from the univariate random-intercept model to the general multivariate mixed-effects model.
• Simulation results. 500 datasets with $m=30$ subjects and $q=3$ random effects. Standard ECME algorithms can be slow for extreme $\Lambda^*$; PX-ECME2 converges quickly across the full range and is comparable to S-Plus lme in wall-clock time with superior numerical stability.

Concepts Introduced or Extended

• Mixed-Effects Regression
• EM Algorithm

Entities Mentioned

• (none with existing wiki pages)

Quotes

"We present efficient algorithms for maximum likelihood (ML), restricted maximum likelihood (REML), and computing posterior modes with conjugate proper and improper priors. These algorithms are not only useful in their own right, but also illustrate how parameter expansion, conditional data augmentation, and the ECME algorithm can be used in conjunction to form efficient algorithms."

"PXEM brings attention to this difficulty with the posterior [under an improper prior] — the mode is not a sufficient summary of the posterior, at least on the original scale."

My Take

A practically important paper that solves a real numerical problem: EM for mixed-effects models was well-known to be slow or unreliable in standard software. The theoretical contribution — grounding the algorithm choices in the fraction-of-missing-information framework — turns what could have been ad hoc tuning into a principled design exercise. The clearest practical take-away is the asymmetry between ML/REML (PXEM works beautifully) and proper Bayesian MAP (PXEM fails; use conditional augmentation). This paper focuses on posterior modes, not full posterior sampling; in a Bayesian time-series context, MCMC (Gibbs sampling) is more commonly used for inference, but EM modes are valuable as starting values and for approximate marginal likelihood calculations.