Online EM for exponential families

Online expectation–maximization replaces a full-data E-step by a running average of posterior sufficient statistics. Cappé and Moulines develop the general method in 2009, On-line expectation–maximization algorithm for latent data models. This note specializes the update to a full exponential family, derives the upstream regret identity, and applies the same statistic order as batch EM for GH.

Read the exponential-family core for natural and expectation coordinates and IG duality for Bregman divergences. Throughout the full-family derivation, assume a regular minimal family, finite posterior statistics, and iterates in an interior region where the expectation-to-natural map exists. Matrix statistic blocks use their independent symmetric coordinates.

Posterior statistics and the sequential update¶

For observed $X$ and hidden $Y$, write

Let $\phi(\eta)=\langle\theta,\eta\rangle-\psi(\theta)$ be the Legendre dual, so that $\theta=\nabla\phi(\eta)$. At observation $x_t$, compute

Starting from $\eta_0=\nabla\psi(\theta_0)$, update

The reciprocal schedule $\rho_t=1/(\tau_0+t)$, with $\tau_0\geq0$, has the exact aggregation identity

Here $\rho_t$ is the step size; $\tau_t=\tau_0+t$ is its reciprocal. The initial statistic has weight $\tau_0$ in pseudo-observation units. Each $\bar t_t$ uses the parameter available when that observation arrives; it is not recomputed under $\theta_T$. Thus this average differs from a batch E-step at the final parameter. With $\tau_0=0$, the first update can lie on the boundary, so the assumed inverse map needs particular care.

For a mini-batch, replace $\bar t_t$ by its within-batch average. A constant $\rho_t$ gives exponential forgetting. Decreasing schedules satisfying $\sum_t\rho_t=\infty$ and $\sum_t\rho_t^2<\infty$ are common in stationary stochastic approximation; for example $(t+t_0)^{-\kappa}$ with $1/2<\kappa\leq1$. These conditions alone do not prove convergence: stability, regularity, and a suitable M-step are also needed. The method’s cost per update uses fewer observations than batch EM, but this does not imply the same convergence rate or likelihood ascent at each observation.

A Bregman regret identity¶

Define the complete-data divergence, with its orientation explicit,

For a sequence $x_1,\ldots,x_T$, compare sequential prediction with a fixed parameter chosen using all observations and this initial penalty:

Suppose this minimum is finite and attained at an interior $\theta_*$, and put $\eta_*=\nabla\psi(\theta_*)$. For the reciprocal schedule above, the upstream decomposition is

This is an identity, not by itself a sublinear regret bound. The posterior KL sum is nonnegative and need not be small. It disappears for fully observed data. The penalty belongs to the comparator in (6); the sequential algorithm need not find its minimizer.

To verify the identity, let $g_t(\theta)=\psi(\theta)-\langle\theta,\bar t_t\rangle$. Parameter-independent carrier terms cancel in likelihood differences. The update gives $\bar t_t=\tau_t\eta_t-\tau_{t-1}\eta_{t-1}$, and Legendre duality gives

Subtract these expressions. The observed-versus-complete likelihood identity adds the posterior KL sum, while $D_\phi(\eta_t\Vert\eta_{t-1})=D_\psi(\theta_{t-1}\Vert\theta_t)$ supplies (7). The calculation is algebraic and requires all displayed expectations and divergences to be finite.

Application to generalized hyperbolic mixtures¶

Use the GH model

The classical tuple is $\vartheta=(\mu,\gamma,\Sigma,p,a,b)$; $\theta$ denotes natural coordinates. In the scalar conditioning note, $W=Y$, $\beta=\gamma$, and $\sigma^2=\Sigma$. For $a,b>0$ and positive-definite $\Sigma$, the posterior is GIG and its power moments are finite. Set $u_t=\mathbb E[Y^{-1}\mid x_t]$, $v_t=\mathbb E[Y\mid x_t]$, and $l_t=\mathbb E[\log Y\mid x_t]$, evaluated at $\vartheta_{t-1}$. In the Batch 1 statistic order,

The first three slots are $(\log Y,Y^{-1},Y)$, a permutation of the source online note’s $(Y^{-1},Y,\log Y)$. The normal parameters are recovered as

The mixing-law update maximizes $L_{\mathrm{GIG}}(p,a,b\mid\eta_{1,t},\eta_{2,t},\eta_{3,t})$. These are the batch M-step formulas applied to running statistics. Their feasibility, nonzero-denominator, and positive-definiteness requirements still apply; boundary mixing laws need their own moment checks. A shrinkage target can regularize the statistics, but changes the update and its objective.

Curved families need a constrained M-step¶

For $f(x,y;u)=h(x,y)\exp\{\langle\theta(u),t(x,y)\rangle-\psi(\theta(u))\}$, the ambient expectation vector need not equal the statistics being averaged. The appropriate map is

At an interior stationary point it satisfies

This is weaker than ambient moment matching. Applying $\nabla\phi$ without the constraint can leave the curved family. This limitation concerns that unconstrained inverse map: it does not exclude online EM for curved families. The Cappé–Moulines construction uses a model-specific maximizing map. Its convergence assumptions must be checked for the particular model; the full-family regret derivation above does not transfer automatically.

GH factor analysis gives a concrete constrained example. The EM update framework separates statistic aggregation, update rules, and parameter recovery. Numerical backends and runnable streaming examples stay in the upstream EM framework.

Source and adaptation¶

Adapted from xshi19/normix, docs/theory/online_em.md, at revision 763bb3608920661a012cf089888d349fbf680aad (2026-09-13 import). Copyright (c) 2020 xshi19. Licensed under MIT. The pinned source records the original version. Notation, mathematical qualifications, and links were adapted for this site. Package interfaces, fitter recipes, and executable cells are omitted; no upstream benchmark or formal-proof verification is claimed.