CMMSE – A Hybrid EM and Bayesian‐Stochastic Approximation Method for Variance Components in Incomplete Mixed Models

ABSTRACT Variance‐component estimation in Gaussian linear mixed models (LMMs) is routinely performed by least squares/analysis of variance (LS/ANOVA) and likelihood‐based methods (ML/REML), but these approaches can become numerically fragile and computationally expensive under incomplete and highly unbalanced designs. We introduce Bayesian‐stochastic approximation expectation‐maximization (BSA‐EM), a hybrid procedure tailored to variance‐component inference with arbitrary missing patterns. The method targets maximum‐a‐posteriori (MAP) estimation under an explicit prior on the variance components. Missing responses are handled through a deterministic EM completion: at each iteration, the conditional mean and variance of the missing block given the observed block are computed analytically from the Gaussian covariance partition, avoiding any Monte Carlo E‐step or latent‐variable simulation. Given these conditional moments, inverse‐gamma regularization yields closed‐form MAP updates for each variance component, and a Robbins‐Monro damping recursion stabilizes the iteration and enforces strictly positive updates in sparse regimes where boundary/singularity issues are common for likelihood optimization. Extensive Monte Carlo experiments in crossed two‐way random‐effects designs with progressively removed cells benchmark BSA‐EM against LS, ANOVA), ML, REML, SAEM, and online EM. Across designs, BSA‐EM remains stable under extreme imbalance and is orders of magnitude faster than iterative likelihood and stochastic‐EM baselines. While it is less accurate than REML, it bridges the gap between spectral speed and REML accuracy in sparse layouts, enabling large‐scale parametric bootstrap uncertainty quantification and practical prior‐sensitivity analyses. Additional experiments confirm that these computational advantages persist as the number of random‐effect levels increases.