Convergence-Guaranteed EM Algorithm

• Convergence-guaranteed EM algorithms are iterative methods with explicit conditions ensuring monotone ascent, stability, and finite-sample accuracy.
• They utilize various techniques such as gradient EM, MAP-EM, and regularized EM, often employing Lyapunov functions and contraction analyses to secure convergence.
• The guarantees depend on factors like component separation, strong initialization, and signal-to-noise ratio, which dictate whether convergence is global or local.

A convergence-guaranteed EM algorithm is an Expectation-Maximization procedure, or an EM-type variant such as gradient EM, MAP-EM, approximate EM, or regularized EM, for which the literature supplies explicit assumptions implying monotone ascent, local asymptotic stability, geometric contraction, or finite-sample statistical accuracy. In the classical theory, EM is an iterative method for incomplete-data likelihood maximization whose limit points are stationary points of the likelihood, and whose local rate is linear at best (Neath, 2012). More recent work replaces this generic stationary-point statement by model-specific guarantees: contraction to the true parameter in mixtures of Gaussians and mixed linear regression, asymptotic stability via Lyapunov functions for MAP-EM, and explicit finite-sample error floors for sample-splitting, regularized, and approximate EM schemes (Romero et al., 2020).

1. Formal convergence notions

The standard EM iteration is defined through the complete-data surrogate

$$Q(\theta\mid\theta^{(t)};y)=E_{U|Y=y;\theta^{(t)}}[\log f(y,U;\theta)],$$

with update

$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$

Under mild regularity, each EM step never decreases the observed-data log-likelihood, and if the EM algorithm converges then it converges to a stationary point of the likelihood (Neath, 2012). In the same line of theory, the local rate is governed by the missing-information fraction, so the convergence is linear at best (Neath, 2012).

A dynamical-systems formulation sharpens this picture for MAP-EM. Writing the iterate as a state vector $\theta_k\in\Theta\subseteq\mathbb R^p$, the update is

$$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$

with

$$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$

This makes EM a proximal-point update,

$$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$

and allows one to use the Lyapunov candidate

$V(\theta)=\ell(\theta)-\ell(\theta^*)$

to prove local asymptotic stability, and under radial unboundedness plus uniqueness of the fixed point, global asymptotic stability (Romero et al., 2020).

A further refinement appears in data-adaptive analyses of sample EM. There, the local convergence rate of sample EM is a random variable $\overline K_n$ derived from the data generating distribution, and $\overline K_n$ concentrates on the population-level optimal convergence rate $\overline\kappa$ (Wu et al., 2016). This framework explains the phenomenon that the finite sample version of the algorithm sometimes converges faster even than the population version (Wu et al., 2016).

2. Sufficient conditions and proof templates

Across the literature, convergence guarantees are attached to explicit sufficient conditions rather than to EM in the abstract. These conditions include separation between components, nondegenerate mixing weights, initialization inside a contraction basin, curvature of the surrogate or base loss, and concentration of the empirical update around the population map. This suggests that “convergence-guaranteed EM” is best understood as a family of conditional results rather than as a model-independent theorem.

Setting	Guarantee	Key conditions
$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$0-component mixture of linear regressions	$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$1	$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$2, moment initialization, basin $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$3
$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$4-component Gaussian mixture, gradient EM	$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$5	$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$6, known weights, initialization radius $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$7
Agnostic mixture of $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$8 parametric functions	$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$9	$\theta_k\in\Theta\subseteq\mathbb R^p$0-strong convexity, $\theta_k\in\Theta\subseteq\mathbb R^p$1-smoothness, non-degeneracy, misspecification bound, separation, sample-splitting
Gaussian latent tree model	convergence to unique interior fixed point	interior initialization, unique non-trivial stationary point, boundary fixed points are higher-order saddles

The proof templates are similarly recurrent. In multi-component Gaussian mixtures, one uses a Gradient Stability condition,

$\theta_k\in\Theta\subseteq\mathbb R^p$2

together with the Hessian of $\theta_k\in\Theta\subseteq\mathbb R^p$3, to obtain geometric contraction (Yan et al., 2017). In mixtures of linear regressions, the population M-step is decomposed into a well-conditioned matrix term $\theta_k\in\Theta\subseteq\mathbb R^p$4 and an error term $\theta_k\in\Theta\subseteq\mathbb R^p$5, with $\theta_k\in\Theta\subseteq\mathbb R^p$6 and $\theta_k\in\Theta\subseteq\mathbb R^p$7, yielding the map $\theta_k\in\Theta\subseteq\mathbb R^p$8 (Kwon et al., 2019). In agnostic mixtures, one splits the gradient sum over “good” points in $\theta_k\in\Theta\subseteq\mathbb R^p$9, for which $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$0, and “bad” ones, for which $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$1, so that strong convexity yields contraction on $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$2 and the complement contributes only a small additive term (Ghosh, 7 Apr 2026). In latent Gaussian tree models, global convergence follows from the uniqueness of the non-trivial interior solution of the EM fixed-point equations together with the fact that boundary fixed points are repelling higher-order saddles (Dagan et al., 2022).

3. Gaussian-mixture results

For mixtures of two Gaussians with known covariance, the population theory is unusually sharp. In the symmetric single-mean model, the population EM map has exactly three fixed points, $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$3 and $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$4, and if the initialization is not orthogonal to $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$5, the iterates converge linearly to the nearer of $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$6; the “null-plane” collapses to $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$7 (Xu et al., 2016). In the unconstrained two-mean parameterization $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$8, the center term $$\theta_{k+1}=F(\theta_k)\triangleq \argmax_{\theta\in\Theta}Q(\theta,\theta_k),$$9 contracts linearly and the angle of the difference term $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$0 to the true direction shrinks geometrically, so that population EM from almost every initialization converges to one of the two global maxima or its label-swap (Xu et al., 2016).

An explicit finite-iteration illustration appears in the balanced two-Gaussian model with known equal covariances. There the population EM operator can be written as

$$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$1

and in one dimension with $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$2, ten EM iterations from an infinite initial guess suffice for less than $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$3 error (Daskalakis et al., 2016). In finite samples, the same paper gives $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$4 sample complexity and the optimal rate $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$5 in Mahalanobis distance (Daskalakis et al., 2016).

For general $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$6, the guarantees become local and parameter-sensitive. A gradient-EM analysis for multi-component Gaussian mixtures with known weights and identity covariance shows population contraction under $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$7, where $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$8, and yields a near-optimal local contraction radius $$Q(\theta,\theta')=\log p(\theta\mid y)-d(\theta,\theta'), \qquad d(\theta,\theta')=\mathrm{KL}\bigl(p(x\mid y,\theta')\,\|\,p(x\mid y,\theta)\bigr)\ge 0.$$9 (Yan et al., 2017). An improved analysis later established convergence of both EM and gradient EM from the region

$$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$0

under only $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$1, with sample complexity $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$2 and final statistical error $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$3 up to log-factors (Segol et al., 2021). In full generality, for multi-component Gaussian mixtures with unknown $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$4 and $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$5, a 2025 analysis gives a population contraction factor

$$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$6

in a neighborhood of the truth, and sample-EM rates $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$7 for the means and $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$8 for the covariance, together with minimax lower bounds showing optimality up to $$\theta_{k+1}=\argmin_\theta\{\ell(\theta)+d(\theta,\theta_k)\}, \qquad \ell(\theta)=-\log p(\theta\mid y),$$9 (Bing et al., 10 Sep 2025).

These results also delimit the classical misconception that EM is generically globally convergent for Gaussian mixtures. The two-component case admits global analyses (Xu et al., 2016, Daskalakis et al., 2016), but for $V(\theta)=\ell(\theta)-\ell(\theta^*)$0 the exact-parameterized setting is substantially more difficult, and local rather than global guarantees dominate the literature (Zhou et al., 6 Jun 2025).

4. Mixed linear regression and agnostic mixtures

Mixtures of linear regressions provide one of the clearest examples in which convergence guarantees depend on signal-to-noise ratio, initialization, and whether the analysis is population or finite-sample. For a $V(\theta)=\ell(\theta)-\ell(\theta^*)$1-component mixture of linear regressions with $V(\theta)=\ell(\theta)-\ell(\theta^*)$2 and noise variance $V(\theta)=\ell(\theta)-\ell(\theta^*)$3, let

$V(\theta)=\ell(\theta)-\ell(\theta^*)$4

Under the requirement

$V(\theta)=\ell(\theta)-\ell(\theta^*)$5

equivalently $V(\theta)=\ell(\theta)-\ell(\theta^*)$6 up to poly-log factors, and with initialization satisfying

$V(\theta)=\ell(\theta)-\ell(\theta^*)$7

population EM contracts as

$V(\theta)=\ell(\theta)-\ell(\theta^*)$8

In the sample-splitting setting, after $V(\theta)=\ell(\theta)-\ell(\theta^*)$9 iterations one obtains

$\overline K_n$0

which in original units is $\overline K_n$1, and crucially this bound does not grow with $\overline K_n$2 or $\overline K_n$3; as $\overline K_n$4, it vanishes exactly (Kwon et al., 2019).

For the two-component model, the convergence picture depends sharply on SNR. With $\overline K_n$5, EM attains the standard parametric rate $\overline K_n$6 after $\overline K_n$7 iterations in high SNR, converges to an $\overline K_n$8 neighborhood in medium SNR, and reaches an $\overline K_n$9 neighborhood after $\overline K_n$0 iterations in low SNR (Kwon et al., 2020). The same work shows that these rates match the minimax lower bound $\overline K_n$1 up to constants and mild log-factors (Kwon et al., 2020).

A distinct line of work studies the fully unknown two-component setting with unknown mixing weights and regression parameters. In the noiseless limit, the population recurrence for the sub-optimality angle,

$\overline K_n$2

produces a cycloid trajectory in the plane spanned by $\overline K_n$3 and $\overline K_n$4. The same analysis identifies two phases: linear convergence when the EM estimate is nearly orthogonal to the ground truth, and quadratic convergence when the angle is small (Luo et al., 7 Nov 2025). At the finite-sample level, one-step errors satisfy

$\overline K_n$5

and arbitrary initialization can be handled after a few “easy-EM” steps (Luo et al., 7 Nov 2025).

Beyond generative mixed linear regression, gradient EM has been extended to the agnostic setting. Given samples $\overline K_n$6, one fits $\overline K_n$7 parametric predictors $\overline K_n$8 using a strongly convex and smooth base loss $\overline K_n$9, with soft-assignment weights

$\overline\kappa$0

Under non-degeneracy, misspecification, separation, and suitable initialization, one iteration satisfies

$\overline\kappa$1

and therefore

$\overline\kappa$2

This framework covers ridge-regularized mixed linear regression, mixed logistic regression, mixtures of support-vector machines, and mixtures of generalized linear models (Ghosh, 7 Apr 2026).

5. Approximate, regularized, and structured EM

A convergence guarantee can also be attached to deterministic approximations of EM. In deterministic approximate EM, the exact conditional $\overline\kappa$3 is replaced by a deterministic surrogate $\overline\kappa$4, and the sufficient statistic is updated by

$\overline\kappa$5

Under curved-exponential-family assumptions and an $\overline\kappa$6-consistency condition on the surrogate, the stable approximate EM sequence remains in a compact, $\overline\kappa$7 converges to some value in $\overline\kappa$8, and if $\overline\kappa$9 has empty interior then $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$00 converges to $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$01 (Lartigue et al., 2020). This theory includes Riemann EM, Tempered EM, and Tempered Riemann EM (Lartigue et al., 2020).

Regularized EM alters the objective rather than the E-step. For Gaussian mixtures with cluster covariances $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$02, a penalized likelihood is defined by

$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$03

where

$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$04

The E-step is unchanged, while the M-step yields the shrinkage covariance update

$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$05

Starting from $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$06, the penalized log-likelihood never decreases, all iterates remain positive definite, and every convergent subsequence is a stationary point of $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$07 (Houdouin et al., 2023). The guarantee is aimed at low-sample-support settings and structured covariance targets (Houdouin et al., 2023).

Structured latent-variable models admit still different mechanisms. In Gaussian latent tree models, the population EM fixed-point equations reduce to a system of quadratic polynomials in the edge correlations $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$08. The unique non-trivial stationary point of the population log-likelihood is its global maximum, and in the single latent variable case EM is guaranteed to converge to it; trivial boundary fixed points are higher-order saddles (Dagan et al., 2022). A non-statistical but closely related reinterpretation appears in sparse linear systems: an EM-based fixed-point iteration for $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$09 converges geometrically when the solution is unique, and otherwise converges to a minimal Kullback–Leibler divergence point (Chae et al., 2016).

6. Misconceptions, limitations, and emerging directions

A persistent misconception is that EM convergence is synonymous with monotone improvement of the surrogate $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$10. A contrary viewpoint argues that $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$11 may and should decrease for $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$12 to increase, and that slow or local convergence exists only because of small samples and unfair competition (Lu, 2021). In that framework, the Fair Competition Principle recommends initial means near the Voronoi–centroids of an equal-mass partition of the data, equal or mass-proportional initial mixing weights, and non-singular initial covariances, with an initialization map that can vastly save running times for binary Gaussian mixtures (Lu, 2021). This remains a distinct line of argument rather than a replacement for contraction-based analyses.

Another limitation is the role of bad local regions. In over-parameterized Gaussian mixtures, where a student model with $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$13 components learns a single $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$14, gradient EM converges globally at the sublinear rate $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$15, but there exist bad local regions that can trap gradient EM for an exponential number of steps (Xu et al., 2024). The proof proceeds by showing that the likelihood loss satisfies a recurrence of the form

$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$16

which yields $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$17, while symmetry can force exponentially small gradients in certain highly separated initial configurations (Xu et al., 2024).

More recent results show that overspecification can also improve global behavior when additional structure is imposed. In a regular-simplex overspecified model with nonvanishing discrete Fourier transform of the weights, EM is equivalent to gradient descent on a locally strongly convex negative log-likelihood, satisfies a local Polyak–Łojasiewicz inequality, and achieves

$$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$18

so that $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$19 iterations suffice for $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$20-accuracy (Assylbekov et al., 13 Jun 2025). For general well separated $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$21-component isotropic Gaussian mixtures, an over-parameterized learner with $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$22 and random initialization is proved to converge globally to the ground truth at a polynomial rate with polynomial samples, using Hermite polynomials and tensor decomposition (Zhou et al., 6 Jun 2025). These results do not eliminate the dependence on separation, initialization, or model structure, but they indicate that global convergence beyond the special case $$\theta^{(t+1)}=\arg\max_{\theta\in\Theta}Q(\theta\mid\theta^{(t)};y).$$23 becomes accessible once redundancy, geometry, and identifiability are exploited.

Taken together, the literature shows that a convergence-guaranteed EM algorithm is not a single canonical procedure. It is a technically qualified EM map whose guarantee is tied to a specific regime: stationary-point convergence in classical EM, Lyapunov stability in MAP-EM, local contraction under separation and good initialization in Gaussian mixtures and mixed linear regression, exponential convergence under strong convexity and smoothness in agnostic mixtures, and structured or overspecified global recovery in selected modern settings. This suggests that the decisive question is not whether EM converges, but which EM map, under which geometry, and to which notion of target.