Expectation-Maximization Algorithm

Fitting a model to data means adjusting parameters until the observations are as probable as possible under the model. When all the data is observed, this often reduces to setting a derivative to zero and solving. But sometimes part of the data is hidden. A label was never recorded, a grouping was never observed, or some intermediate variable was never observable. The likelihood of the observed data alone involves a sum or integral over all possible values of the hidden variables, and this sum is difficult to optimize directly. The expectation-maximization algorithm, or EM, makes this optimization tractable.

The EM algorithm alternates between two steps. Write $X$ for the observed data, $Z$ for the hidden variables, and $\theta$ for the parameters. The complete data is the pair $(X, Z)$. If $Z$ were known, you could form the complete-data log-likelihood $\log p(X, Z \mid \theta)$ and maximize it directly. Since $Z$ is not known, you replace this log-likelihood with its expected value under the current parameter estimate $\theta^{(t)}$. Specifically, you compute the expected complete-data log-likelihood, averaged over the conditional distribution of $Z$ given the observed data and the current estimate of $\theta$.

\[Q(\theta \mid \theta^{(t)}) = E_{Z \mid X,\, \theta^{(t)}}[\log p(X, Z \mid \theta)]\]

Computing $Q$ is the E step. Finding the $\theta$ that maximizes $Q$ is the M step. That maximizer becomes $\theta^{(t+1)}$, and the cycle repeats.

Each iteration is guaranteed not to lower the observed-data log-likelihood. The log-likelihood $\log p(X \mid \theta)$ can be written as $Q(\theta \mid \theta^{(t)})$ plus a remainder term that depends on the conditional distribution of $Z$ given $X$. The M step raises $Q$, because it picks the $\theta$ that makes $Q$ as large as possible. The remainder term, measured against the fixed reference distribution $p(Z \mid X, \theta^{(t)})$, cannot decrease when $\theta$ moves away from $\theta^{(t)}$. Both effects either increase the observed-data likelihood or leave it unchanged, so the likelihood is nondecreasing along the sequence $\theta^{(0)}, \theta^{(1)}, \theta^{(2)}, \ldots$, and the sequence converges to a stationary point under standard regularity conditions.

Take a concrete case. You flip a coin five times. The coin has an unknown bias $\theta$, the probability of heads on each flip. Three of the outcomes are recorded, and they are H, H, T. The other two outcomes were lost. With all five results in hand, you would count the total number of heads $k$ and set $\hat{\theta} = k/5$. With two results missing, you cannot count directly. Let $Z$ be the number of heads among the two lost flips. The complete data is the full set of five outcomes, with log-likelihood $k \log\theta + (5 - k)\log(1 - \theta)$ where $k = 2 + Z$. In the E step, you compute the expected value of $Z$ given the observed outcomes and the current estimate $\theta^{(t)}$. Each missing flip lands heads independently with probability $\theta^{(t)}$, so $E[Z \mid X, \theta^{(t)}] = 2\theta^{(t)}$. Substituting this expected count into the complete-data log-likelihood and maximizing over $\theta$ gives

\[\theta^{(t+1)} = \frac{2 + 2\theta^{(t)}}{5}.\]

Start with $\theta^{(0)} = 0.5$. The first update gives $\theta^{(1)} = 3/5 = 0.6$. The second gives $\theta^{(2)} = 3.2/5 = 0.64$, and the third gives $\theta^{(3)} = 3.28/5 = 0.656$. At the fixed point, $\theta^\ast = (2 + 2\theta^\ast)/5$, whose solution is $\theta^\ast = 2/3$. The algorithm converges to the same answer you would get by ignoring the missing flips and estimating from the three observed outcomes alone, two heads out of three.

The missing-data setting above is the simplest instance of a much wider pattern. In mixture models, the hidden variable $Z$ indicates which component generated each observation. In factor models, it represents latent variables underlying the observed measurements. In each case, knowing $Z$ would reduce estimation to a standard problem with a closed-form solution. The EM algorithm fills in $Z$ in expectation, solves that standard problem, and repeats. Because no step lowers the observed-data likelihood, the sequence of likelihood values converges whenever it is bounded above. But the likelihood surface may have several peaks, and EM converges to whichever peak its starting point leads to. The guarantee is monotone ascent, not global optimality.