EM: soft clustering, one half-step at a time

A Gaussian mixture doesn't assign each point to a cluster — it lets every cluster claim a fraction of every point. EM alternates two moves: the E-step freezes the parameters and updates the beliefs p(j|i) (watch the point colors blend), and the M-step freezes the beliefs and re-fits every parameter from soft counts (watch the bells move). Each pair of moves pushes the log-likelihood ℓ uphill — never down.

E: p(j|i) = πⱼN(xᵢ;μⱼ,σⱼ²) / Σₖ πₖN(xᵢ;μₖ,σₖ²)  ·  M: πⱼ, μⱼ, σⱼ² ← soft-count MLEs
cluster 1 cluster 2 cluster 3 mixture density point color = its p(j|i) blend

Responsibilities p(j|i) — each row sums to 1

Controls

Speed
Data dimension
Clusters K
Starting variances — how wide each bell begins
Click the plot to edit points
Saved examples — tuned starting points

This half-step

Iterations — ℓ only ever climbs

tμσ²π

What EM is actually doing

1 · The model: a mixture of bells

Each cluster j is a Gaussian bell N(x; μⱼ, σⱼ²) with a mixture weight πⱼ — the share of all points it claims. The density the model assigns to any location x is the prior-weighted sum over clusters:

p(x; θ) = Σⱼ πⱼ · N(x; μⱼ, σⱼ²)

We want the parameters θ = (π, μ, σ²) that maximize the log-likelihood of the data:

ℓ(D; θ) = Σᵢ log ( Σⱼ πⱼ N(xᵢ; μⱼ, σⱼ²) )

That sum inside the log is the whole problem: it welds the clusters together so no closed-form solution exists. If we knew which cluster made each point, the sum would collapse and every parameter would be a one-line MLE. EM manufactures that missing knowledge — softly — and alternates.

2 · E-step: parameters frozen, beliefs updated

For every point i and cluster j, compute the posterior probability that cluster j generated point i — its responsibility. It's Bayes' rule: the numerator is cluster j's prior-weighted likelihood (how well its bell explains this location, scaled by its share), and the denominator is the evidence — the same quantity summed over all clusters, i.e. the mixture density at xᵢ itself:

p(j|i) = πⱼ · N(xᵢ; μⱼ, σⱼ²)Σₖ πₖ · N(xᵢ; μₖ, σₖ²)

Note the prior πⱼ appears in every numerator and (inside the sum) in the denominator. Dropping it is only harmless while all priors are equal — true at a symmetric start, false the moment the first M-step runs.

3 · M-step: beliefs frozen, parameters updated

Now re-fit every parameter of every cluster, reading only that cluster's responsibility column. Every formula is a familiar MLE with soft counts substituted for hard ones: instead of belonging to one cluster, point i contributes to cluster j with weight p(j|i). First the effective count — a soft head-count of how many points' worth of responsibility cluster j holds:

n̂ⱼ = Σᵢ p(j|i)   (and Σⱼ n̂ⱼ = n exactly, since each row sums to 1)

Then the three updates, in order:

πⱼ = n̂ⱼn     μⱼ = Σᵢ p(j|i) · xᵢn̂ⱼ     σⱼ² = Σᵢ p(j|i) · ‖xᵢ − μⱼ‖²d · n̂ⱼ

The mean is a responsibility-weighted average — each point votes on where the center sits, with voting power p(j|i). The variance uses squared distance from the new mean, so μⱼ must be computed first. In d dimensions the denominator is d·n̂ⱼ, because each d-dimensional point supplies d scalar deviations — that's why the same code runs the 1D, 2D and 3D demos above, with d the only difference.

4 · Why ℓ climbs, and where it stops

Each E-step makes the beliefs consistent with the parameters; each M-step makes the parameters optimal for the beliefs. Formally, the E-step tightens a lower bound on ℓ (Jensen's inequality) until it touches, and the M-step maximizes that bound — so a full iteration can never decrease ℓ. Watch the iterations table: ℓ only ever climbs, then stalls.

But "climbs" only means uphill from where it started. Try the "Both means right of the data" example: EM converges cleanly — to a worse ℓ than the warm start reaches, with one skinny cluster pinned on the right clump and one wide one draped over everything else. EM finds local optima; the initialization chooses which one. That's why practitioners run it from several random starts (🎲) and keep the best ℓ.

A mixture also separates clusters by scale, not just location. Try the "Needle in a haystack" example: a tight clump embedded in a diffuse spread. The best fit is nested — one wide background bell draped over everything, one narrow bell pinned on the clump — and most starting parameters find it, because the E-step assigns points by density ratio: a dense subclump is explained far better by a skinny bell than by the wide one's flat middle. This is how mixtures model "signal + background" in real datasets.

Two footnotes you'd hit in practice: a cluster that locks onto a single point can drive σⱼ² → 0 and ℓ → ∞ (this demo floors the variance at 0.05 to keep runs honest), and perfectly symmetric data with a perfectly symmetric start can leave every responsibility frozen at ½ forever — symmetry has to break before clustering can happen.