How do you set a model’s numbers from data? Pick the numbers that make the data most likely. When part of the data is hidden, the EM (expectation-maximization) algorithm does it in two alternating steps.
1 Fit the numbers that make the data likely
Every model has parameters, numbers set from data, such as a coin’s bias or a network’s weights. θ stands for them. Maximum likelihood sets θ to the value that makes the observed data most probable; when part of the data is hidden, EM climbs to a peak of the same score.
Figure 1. With labels seen, take the peak; with labels hidden, EM loops toward it.
2 Same formula, two questions
Flip a coin 10 times and see 7 heads. Probability fixes θ, here the coin’s chance of heads, and asks which flips are likely. Likelihood fixes the flips you saw and asks which θ explains them best.
The likelihood scores each possible θ by how probable it makes the flips you actually saw.
L(θ) =θh(1 − θ)t
h is the number of heads seen. t is the number of tails. Each head adds a factor θ and each tail a factor 1 − θ, so the product is the chance of that exact sequence of flips. Held as a function of θ with h and t fixed, it is the likelihood. A θ near 0 makes the heads factor tiny and one near 1 makes the tails factor tiny, so the curve peaks in between, at the fraction of flips that were heads. (The count of possible orders is left out; it is a constant and does not move the peak.)
3 Maximum likelihood picks the peak
Choosing the θ at the top of the curve is maximum likelihood estimation; for 7 heads in 10 the peak sits at 0.7 (Figure 2).
The peak of the coin likelihood is simply the share of heads.
θ =hh + t
θ (theta hat) is the best θ, the top of the curve. Raising θ makes the heads more probable and the tails less so, and the two pulls balance exactly when θ equals the share of heads. More data with the same share leaves the peak where it is and only sharpens the curve. This counting answer is why counting works in Naive Bayes.
Figure 2. The 70 of 100 button keeps the peak at 0.7 and narrows the curve.
4 Logs turn products into sums
A likelihood multiplies one chance per independent observation, so it shrinks fast: the peak for 70 heads in 100 is near 10−27, and a long audio file underflows to zero. Logs fix that.
Logs turn a product of many small chances into a sum that never underflows and peaks at the same θ.
log L(θ) =Σi log P(xi | θ)
xi is observation number i. Σi adds over all of them. The log of a product is the sum of the logs, so each observation adds its own term and tiny chances become modest negative numbers a computer can add. Log rises whenever L rises, so the peak stays at the same θ.
Figure 3. Taking the log changes the curve’s shape but not where it peaks.
5 Fitting a Gaussian gives the sample mean and variance
A bell curve (see Common Probability Distributions) has two parameters, its mean μ and its standard deviation σ, whose square σ2 is the variance. Its likelihood is the product of each point’s height on the bell.
Figure 4. Each stick is one point’s height; the product peaks at the average.
For a bell curve, maximum likelihood sets the mean to the average of the points and the variance to their average squared gap from that mean.
μ =1N Σi xi
N is the number of points. The best mean μ is their plain average: each point pulls the center toward itself, harder the farther away it is, and the pulls cancel exactly at the average.
σ2 =1N Σi (xi − μ)2
The maximum-likelihood variance is the average squared gap from that mean. Squaring counts gaps on both sides alike and makes far points count much more, so one outlier widens the bell. It divides by N, not N − 1, so on average it runs a little small, by a factor (N − 1)/N (Bishop 2006, sec. 1.2.4).
This snippet fits a bell to 20 points with the card’s two formulas and again with scipy; look for the same mean and variance.
import numpy as np
from scipy.stats import norm
rng = np.random.default_rng(0)
x = rng.normal(loc=5.0, scale=2.0, size=20)
mu_hat = x.mean()
var_hat = np.mean((x - mu_hat) ** 2)
loc, scale = norm.fit(x)
print(mu_hat, var_hat) # sample mean and variance
print(loc, scale ** 2) # the same two numbers
1Draws N = 20 points from a known Gaussian.
2Computes μ and σ2, eq card in Section 5, lines 1 and 2.
3Checks them against scipy's maximum likelihood fit.
6 Losses you know are likelihoods
If each target is the prediction plus Gaussian noise, maximum likelihood is the same as minimizing squared error. For classes it is the same as minimizing cross-entropy, the loss logistic regression trains on.
Figure 5. Squared error is minus the log height of the bump.
Minimizing squared error, or cross-entropy for classes, is the same as maximizing a likelihood.
−log N(y | f(x), σ2) =
y is the true value. f(x) is the model’s prediction for input x. N(y | f(x), σ2) is the bell curve centered on the prediction with standard deviation σ, read at y. The bold upright N names the normal (Gaussian) density; it is not the count N of Section 5. Minus its log turns “make y likely” into a loss to minimize.
(y − f(x))22σ2 + const
Minus the log of a bell is the squared gap between y and f(x), divided by twice the variance, plus a constant. The constant and the scale do not move the best fit, so maximum likelihood under bell-shaped noise is least squares. A bigger σ shrinks the penalty, so noisier data tolerates larger errors.
−log P(c | x)
c is the true class. For categories the model gives each class a chance, and minus the log of the chance it gives c is the cross-entropy loss. It is near 0 when the model is sure and right, and it grows without bound as that chance falls toward 0.
7 When part of the data is hidden
When some labels are hidden, plain counting stops working. Take two coins, A and B, with unknown biases (an example from Do and Batzoglou 2008). Five times, someone picks one and flips it 10 times. With the labels you would just count: in their example A made sets 2, 3 and 5, 24 heads in 30, so θA = 0.80, and B made 9 heads in 20, so θB = 0.45.
Figure 6. Five sets of flips; the coin behind each set is hidden.
Hide the labels and each step needs the other’s answer: counting needs the labels, and guessing the labels needs the biases.
With a hidden label, a sum lands inside the log, and the counting formula breaks.
log L(θ) =
The score is still the log-likelihood of θ. Only the data changed: each set’s coin is now unseen.
Σi log Σz P(xi, z | θ)
z is the hidden label, here which coin made set i. Because z is unseen, each set’s chance adds up both coins before the log is taken. A log of a sum does not split into one term per coin, so a set’s heads cannot be handed to one coin and counted, and the peak has no closed formula. With known labels the inner sum would shrink to one term and plain counting would return.
8 EM guesses, then alternates two easy steps
EM (Dempster, Laird and Rubin 1977) replaces each hidden label with a soft guess, refits θ as if those guesses were data, and repeats. The E step uses Bayes rule to give each set a share for each coin. The M step counts heads with each set weighted by its share, the expected counts, and refits both biases.
Figure 7. From θA = 0.60 and θB = 0.50 the two biases move apart.
Each EM round credits every set to the two coins in proportion to how well each explains it, then refits each coin’s bias from the flips credited to it.
aiA =wA P(xi | θA)
E step. wA is how often coin A is picked, one half here, as in the paper. aiA is the chance of set i’s flips under coin A, times wA. aiB is the same for B. A set with many heads scores higher under whichever coin has the higher bias right now.
riA =aiAaiA + aiB
E step, by Bayes rule. riA, the responsibility, is coin A’s share of the two scores, the chance that set i came from A. Dividing by the total makes the two shares add to 1. When one score is much larger than the other the share goes near 0 or 1, close to a hard label; when they are close, the set is split.
θA =ΣiriA hiΣiriA (hi + ti)
M step. hi is the number of heads in set i. ti is its number of tails. Coin A’s new bias is the share of heads among the flips credited to it, the counting formula from the peak with each set weighted by riA. Coin B uses 1 − riA. A set that coin A barely claims barely moves its bias.
After 10 rounds the biases settle near 0.80 and 0.52, as the paper reports.
9 What EM promises, and what it does not
No EM round lowers the likelihood; the 1977 paper proves it. EM does not promise the highest peak, only the nearest one, so run it from several starts and keep the best.
10 EM in practice
Libraries run EM inside fit. This snippet is the EM loop for the five coin sets of Figure 6; it stops when one round gains less than tol in log-likelihood, or after max_iter rounds, and prints the round it stopped at and the two biases.
import numpy as np
from scipy.stats import binom
heads = np.array([5, 9, 8, 4, 7])
flips, tol, max_iter = 10, 1e-6, 100
theta = np.array([0.6, 0.5])
old_ll = -np.inf
for it in range(max_iter):
a = 0.5 * binom.pmf(heads[:, None], flips, theta)
ll = np.log(a.sum(axis=1)).sum()
if ll - old_ll < tol:
break
old_ll = ll
r = a / a.sum(axis=1, keepdims=True)
theta = r.T @ heads / (r.sum(axis=0) * flips)
print(it, theta.round(2)) # compare with the paper: 0.80 and 0.52
1Sets the five head counts hi of Figure 6 and the starting θ.
2Computes aiA and aiB, eq card in Section 8, line 1.
3Sums the log-likelihood, eq card in Section 7, line 2, and stops when it gains less than tol.
4Computes the E-step shares riA and riB, eq card in Section 8, line 2.
5Computes the M-step update of θ, eq card in Section 8, line 3.
question
answer
when
Once per fit, in rounds, until the gain falls below tol or max_iter rounds pass.
on what
Every round reads the whole training set.
state
The current θ; the shares r are rebuilt each round.
cost
One round costs one likelihood per point per hidden group (see References).
setting
raise it, lower it
scikit-learn default
tol, the average gain per point in log L(θ)
Smaller runs more rounds and lands closer to the peak; too large stops early. The snippet tests the total gain instead.
1e-3
max_iter
Caps the rounds; too few stops before the peak, with a warning.
100
n_init, number of starts
More starts, best log L kept, lower the risk of a poor local peak at a proportional cost.
Equal starts never separate: both coins get half of every set, so they stay equal forever.
Labels can switch; swapping the names A and B leaves the likelihood unchanged. With K hidden groups (K = 2 coins here) there are K! (K factorial, 2 here) equally good answers.
The likelihood can run to infinity. In a Gaussian mixture a bell can shrink onto one point and grow without bound; keep a floor on its variance.
A bad start settles on a poor local peak; use several starts.
With too little data per hidden group, the fit chases noise.
12 What if the hidden labels are clusters in a cloud of points?
Swap coins for bell curves and flips for points, and EM fits a Gaussian mixture model, what the coin paper calls the soft variant of k-means. For a chain of hidden states it becomes Baum-Welch, which trains hidden Markov models.
13 References
Dempster, A. P., et al. (1977). Maximum likelihood from incomplete data via the EM algorithm. Journal of the Royal Statistical Society B, 39(1), 1 to 38.
Do, C. B., and Batzoglou, S. (2008). What is the expectation maximization algorithm? Nature Biotechnology, 26, 897 to 899. doi:10.1038/nbt1406.
Bishop, C. M. (2006). Pattern Recognition and Machine Learning. Springer. Sec. 1.2.4 (Gaussian maximum likelihood), 9.2 and 9.3 (mixtures and EM). PDF.