An efficient EM algorithm for both element-wise and structural missingness in matrix-variate normal mixture models

Hanzhang Lu; Jeffrey L. Andrews; Ryan P. Browne

Source-linked AI summary

An efficient EM algorithm for both element-wise and structural missingness in matrix-variate normal mixture models

Hanzhang Lu, Jeffrey L. Andrews, Ryan P. Browne

arXiv:2609.00616v1Sep 1, 2026stat.MEstat.COstat.ML

TL;DR

Standard EM can require vectorizing matrices and dense Gaussian conditioning, making estimation computationally infeasible. The paper introduces a matrix-normal partial EM method that directly updates missing-entry conditional quantities, with a separable specialized version; simulations report substantially improved computational efficiency while maintaining nearly identical performance to exact EM.

Problem
Standard EM may require vectorizing matrices and performing dense Gaussian conditioning, which can remain computationally infeasible.
Method
The paper introduces a partial EM algorithm for the matrix normal distribution that directly updates the conditional mean and covariance of missing entries.
Results
Simulation studies report substantially improved computational efficiency over standard exact EM while maintaining nearly identical results.
Takeaways & Limitations
A specialized version further enhances efficiency by utilizing separable structure.
Takeaways & Limitations
The paper notes that an underlying assumption may be too restrictive in practical applications.

Abstract

from arXiv · show

Matrix-variate data with missing entries arise frequently in applications where observations are naturally organized as two-dimensional arrays. Although the matrix normal distribution provides a parsimonious model through its Kronecker covariance structure, standard EM estimation can be computationally expensive because arbitrary missingness patterns typically destroy this separability in the E-step. In this paper, we propose an efficient partial EM algorithm for matrix-variate normal data with missing entries. The proposed method updates the conditional mean and covariance of the missing component through coordinate-wise approximations, avoiding repeated inversion of pattern-specific covariance matrices and avoiding construction of the full vectorized covariance matrix. We further develop a specialized update for submatrix missingness, where the missing-block precision retains a Kronecker product structure, and the covariance update can be carried out independently in the row and column directions. Simulation studies show that the proposed methods substantially reduce computation time compared with exact EM while preserving nearly identical observed-data likelihood across a range of dimensions and missing proportions. A real-data application to hyperspectral image patches demonstrates that the proposed imputation strategy can be embedded within a matrix-variate mixture model for simultaneous imputation and clustering.

1 Introduction

Matrix-variate normal models use Kronecker covariance structure to represent dependencies across rows and columns with fewer covariance parameters. Missing entries make exact EM computationally difficult, motivating a partial EM method that updates conditional moments efficiently and extends to mixtures for simultaneous imputation and clustering.

Matrix-variate normal modeling: Kronecker covariance structure captures dependencies across rows and columns while reducing the number of free covariance parameters.The parameter count decreases from pq(pq + 1)/2 to p(p+1)/2+q(q +1)/2−1.
Existing EM framework: Standard EM computes missing-entry conditional distributions and then updates mean and covariance factors using conditional sufficient statistics.For complete data, covariance-factor estimation is typically performed using a flip-flop or block-coordinate procedure.
Computational challenge: Arbitrary missingness can destroy Kronecker separability in the observed covariance, requiring dense Gaussian conditioning and pattern-specific matrix inversions.These operations can become computationally intensive in time and memory, making standard EM computationally infeasible in some settings.
Proposed method: The proposed matrix-normal partial EM updates missing-entry conditional means and covariances in matrix form using Kronecker precision factors.It avoids pattern-specific matrix inversions and construction of the full pq × pq covariance matrix.
Submatrix missingness: A specialized submatrix-missingness update preserves an exact Kronecker representation of the missing-block precision matrix.The paper also notes that limited coordinate sweeps provide a computationally efficient partial E-step.
Mixture extension: The procedure extends to mixtures of matrix-variate normal distributions, enabling simultaneous imputation and clustering.The paper reports simulation results for both single-component and mixture settings and applies the methods to hyperspectral image data.

2 Preliminaries

The preliminaries formulate EM imputation through conditional moments of missing entries and explain why arbitrary missingness creates computational difficulty for matrix-normal models. They then introduce partial EM updates based on coordinate descent and precision matrices, alongside the covariance-factor identifiability constraints.

2.1 EM algorithm for imputation: EM replaces unavailable sufficient statistics for missing entries with conditional expectations given observed values.The E-step computes conditional means and covariances, while the M-step updates model parameters from these conditional sufficient statistics.
2.1 EM algorithm for imputation: Different missingness patterns can require many matrix inversions at every iteration, increasing computational complexity.The pattern-specific inversion calculation is identified as a main source of complexity.
2.2 EM algorithm for matrix normal distribution: In matrix-normal models, arbitrary missingness disrupts Kronecker structure after vectorization and restriction to observed entries.The full covariance remains separable, but the observed submatrix is not necessarily representable as a smaller Kronecker product.
2.2 EM algorithm for matrix normal distribution: Matrix-normal covariance estimation projects the conditional covariance onto row and column covariance spaces and updates the two factors separately.The projected terms are used in row- and column-covariance update equations, with flip-flop procedures also available.
2.3 Partial EM algorithm: Partial EM reduces KL divergence rather than minimizing it exactly by updating missing coordinates through coordinate descent using the precision matrix.Observed coordinates remain fixed, and the missing-data covariance is updated through analogous coordinate-wise operations.
2.3 Partial EM algorithm: When the precision matrix is known, the partial E-step avoids calculating matrix inverses for unique missingness patterns.This provides the computational mechanism underlying the partial update strategy.
2.4 Identifiability: Kronecker covariance factors are scale-non-identifiable, so the model introduces a variance scalar and determinant constraints.The constraints are |Σ1| = 1 and |Σ2| = 1, yielding Σ = σ2Σ2 ⊗Σ1.

3 Methods for missingness in matrix-variate normals

The method approximates conditional moments coordinate-wise to preserve matrix structure during partial E-steps, avoiding full vectorized precision construction. For structural missingness, exact separability enables independent row- and column-direction covariance updates; simulations report lower computation time with nearly unchanged likelihood.

Arbitrary missingness disrupts Kronecker separability, preventing the E-step from directly exploiting the matrix-normal covariance structure.
The proposed partial E-step preserves Kronecker structure by updating conditional moments coordinate-wise, reducing computational and memory burden.
Conditional mean: A Gauss-Seidel sweep updates missing-entry conditional means without increasing the objective, while repeated sweeps converge to the exact conditional mean.
Conditional covariance: Sequential covariance-coordinate updates preserve positive definiteness and converge, under the stated conditions, to the exact conditional covariance.
The implementation operates only on missing coordinates and avoids constructing the full pq × pq precision or covariance matrix.
Structural missingness: For structural missingness, the missing-block conditional covariance remains separable, allowing row and column factors to be updated independently.
Simulation results: Simulations report substantially reduced computation time versus exact EM while maintaining nearly the same observed-data likelihood across most settings.
Simulation results: The plug-in method is usually fastest but ignores covariance correction, potentially reducing likelihood or parameter accuracy.

4 Missingness in mixtures of matrix normals

The mixture extension combines matrix-variate normal components with latent membership estimation, component-specific imputation, and covariance structures tailored to spatial data. Simulation results indicate that MPEM and Rect-MPEM substantially reduce computation while preserving clustering and estimation accuracy with small likelihood loss.

Finite matrix-normal mixtures use latent component membership and can jointly perform clustering and component-specific imputation.
In the E-step, the method estimates missing-component moments and latent memberships from the exact marginal density of observed entries.Membership probabilities use the observed-entry density rather than a density evaluated at an imputed complete matrix.
For arbitrary missingness, the approach updates component-specific sufficient statistics and mixture parameters, while submatrix missingness retains Kronecker factorization for simpler density evaluation.The mixture updates include mixing proportions, component means, and row- and column-space covariance corrections.
Mixtures of spatial factor analyzers: The spatial factor-analyzer formulation uses flexible spatial decay with factor-based row and column covariance components, while the sigmoid-decay alternative can be restrictive and numerically unstable.The spatial covariance structure is modeled through I-spline decay functions and factor loadings.
Simulation studies: 63.0× speedup was achieved by Rect-MPEM over exact EM for submatrix missingness while preserving comparable accuracy and producing a smaller likelihood gap than nonstructural MPEM.
Simulation studies: MPEM substantially reduces computational cost while maintaining essentially the same clustering and estimation accuracy with only small likelihood loss.Cholesky closely matches EM in likelihood but provides little computational benefit.

5 Real Application

The real-data study applies the proposed imputation framework to Salinas Valley hyperspectral patches for simultaneous imputation and clustering of three land-cover classes. Under missing-data comparisons, the structural method performs especially well for submatrix missingness, while the partial method remains close to exact EM with lower computation.

Data and preprocessing: Homogeneous 8 × 8 pixel patches preserved local spatial dependencies and assigned every patch entirely to one land-cover class.Class labels were used to construct and evaluate the patches.
Model comparison: The MSFA model was fitted after imputation across factor dimensions r ∈{1, . . . , 4}, spline specifications, covariance constraints, and multiple imputation baselines.Comparisons included EM, MPEM, missRanger, and mean imputation.
Results: Exact EM and feature mean imputation sometimes failed during MSFA refitting, preventing accuracy reporting when components collapsed or covariance updates were non-positive-definite.The table marks component collapse with – and non-positive-definite covariance updates with ∗.
Results: Under MCAR, partial EM was very close to exact EM in RMSE and MAE, faster, and higher in downstream BIC and accuracy, whereas missRanger had larger imputation error and lower BIC.For submatrix missingness, the proposed structural method nearly matched exact EM in imputation accuracy, was over twenty times faster, and achieved the largest BIC.
Results: The three land-cover classes showed distinct spatial correlation trajectories, differing in overall scale and functional form.The trajectories were plotted against normalized distance, and the refitted-model BIC was used alongside them to assess the fitted spatial factor analyzer.

6 Summary

The paper addresses the computational burden of exact EM for matrix-variate normal data with missing entries, whose conditioning can disrupt Kronecker covariance structure. It proposes efficient partial updates and a structural submatrix variant, reporting substantial speedups with nearly unchanged observed-data likelihood and successful integration into mixture modeling.

6 Summary: The framework extends partial EM to mixtures of matrix-variate normal distributions and addresses structural missingness when missing entries form a submatrix.The submatrix setting is treated as a specialized structural case.
6 Summary: Exact EM becomes costly because arbitrary missingness requires conditional moments for each pattern and can disrupt the matrix normal covariance’s Kronecker structure.The standard E-step’s main computational burden is associated with the missing block.
6 Summary: The proposed partial E-step uses coordinate-wise approximations for conditional means and covariances without forming the full covariance or repeatedly inverting pattern-specific matrices.The updates preserve the row-column structure of the matrix normal model.
6 Summary: For submatrix missingness, the missing-block precision retains a separable Kronecker form, enabling independent covariance updates along row and column directions.This structure further reduces computational cost.
6 Summary: Simulation studies found substantially improved computational efficiency over exact EM while maintaining nearly identical observed-data likelihood in most scenarios.The proposed method retains the conditional covariance correction through an efficient approximation, unlike the faster plug-in method.
6 Summary: For submatrix missingness, the block version further improves efficiency over the general version by exploiting the separable structure of the missing entries.The comparison distinguishes the structural block method from the generic partial EM algorithm.
6 Summary: The real-data application integrates the imputation method into a matrix-variate mixture model for simultaneous imputation and clustering while preserving spatial and spectral dependence.The authors present the approach as a scalable alternative to exact EM for random and structured missingness.

Appendix A Covariance updates for MSFA

The appendix supplies parameter updates for the MSFA covariance components, including spatial parameter estimation and factor-loading updates. These updates maintain the model’s covariance constraints and identifiability during estimation.

Covariance updates: MSFA imposes spatial constraints on row covariance matrices while modeling the column covariance separately.The appendix provides corresponding updates for the two component covariance matrices.
Covariance updates: Spatial parameters are estimated by generalized least squares using the current estimated covariance component and spline-based design matrix.The design matrix is defined through vectorized basis and parameter matrices.
Covariance updates: The factor-loading update is applied to estimate the column covariance component.The appendix introduces this update after defining the column covariance quantities.
Covariance updates: The same calculations are repeated with the estimated parameters to maintain identifiability.This step follows the covariance and loading updates.