On the wikipedia article for Expectation-Maximization it states

Given the statistical model which generates a set $\mathbf{X}$ of observed data, a set of unobserved latent data or missing values $\mathbf{Z}$, and a vector of unknown parameters $\boldsymbol\theta$, along with a likelihood function $L(\boldsymbol\theta; \mathbf{X}, \mathbf{Z}) = p(\mathbf{X}, \mathbf{Z}\mid\boldsymbol\theta)$, the maximum likelihood estimate (MLE) of the unknown parameters is determined by maximizing the marginal likelihood of the observed data $$L(\boldsymbol\theta; \mathbf{X}) = p(\mathbf{X}\mid\boldsymbol\theta) = \int p(\mathbf{X},\mathbf{Z} \mid \boldsymbol\theta) \, d\mathbf{Z} = \int p(\mathbf{X} \mid \mathbf{Z}, \boldsymbol\theta) p(\mathbf{Z} \mid \boldsymbol\theta) \, d\mathbf{Z} $$ However, this quantity is often intractable since $\mathbf{Z}$ is unobserved and the distribution of $\mathbf{Z}$ is unknown before attaining $\boldsymbol\theta$.

I don't understand this last sentence. Surely $\mathbf{Z}$ being unobserved is why we integrate over $\mathbf{Z}$, and we set values of $\boldsymbol{\theta}$ during the optimisation procedure which then allows us to compute $p(\mathbf{Z} \mid \boldsymbol{\theta})$?

I initially thought it was intractable because there are too many settings of $\mathbf{Z}$ to consider, but in the expectation step of the EM algorithm we calculate \begin{align} Q(\boldsymbol{\theta} \mid \boldsymbol{\theta}^{(t)}) &= \mathbb{E}_{\mathbf{Z} \mid \mathbf{X}, \boldsymbol{\theta}^{(t)}}\left[\log L(\boldsymbol\theta; \mathbf{X}, \mathbf{Z})\right]\\ &= \int_Z p(\mathbf{Z} \mid \mathbf{X}, \boldsymbol{\theta}^{(t)})\log L(\boldsymbol\theta; \mathbf{X}, \mathbf{Z})d\mathbf{Z} \end{align} which doesn't look any easier than the first integral. If we can approximate this expectation using some method like monte-carlo, why not just approximate the first integral instead?

  • $\begingroup$ Can't we just iteratively sample given the current estimates and maximize to get a new estimate? Sorry for the new post, I can't post comments yet $\endgroup$
    – Chester
    Jul 23, 2023 at 0:04

2 Answers 2


So I've discussed this with a colleague.

Consider the marginalisation $$p(\mathbf{X} \mid \boldsymbol{\theta}) = \int p(\mathbf{X} \mid \mathbf{Z}, \boldsymbol{\theta})p(\mathbf{Z} \mid \boldsymbol{\theta}) d\mathbf{Z}.$$

This can be rewritten as the expectation $$\mathbb{E}_{\mathbf{Z} \mid \boldsymbol{\theta}}\left[ p(\mathbf{X} \mid \mathbf{Z}, \boldsymbol{\theta}) \right].$$

If this can be calculated exactly, we're fine. However, if the expectation is intractable, we need to compute a numerical approximation. Hence to maximise w.r.t. $\boldsymbol{\theta}$ we will need to use some iterative procedure like gradient-ascent.

Suppose we have the current estimate $\boldsymbol{\hat{\theta}^{(t)}}$ for $\boldsymbol{\theta}$. If we now try and approximate the expectation by sampling from $\mathbf{Z}$, we get $$\frac{1}{n}\sum_{i=1}^n p(\mathbf{X} \mid \mathbf{Z}_i, \boldsymbol{\theta})$$ but the $\mathbf{Z}_i$ were sampled from the distribution $\mathbf{Z} \mid \boldsymbol{\hat{\theta}^{(t)}}$, so this actually approximates the expectation $$\mathbb{E}_{\mathbf{Z} \mid \boldsymbol{\hat{\theta}^{(t)}}}\left[ p(\mathbf{X} \mid \mathbf{Z}, \boldsymbol{\theta}) \right] \approx \frac{1}{n}\sum_{i=1}^n p(\mathbf{X} \mid \mathbf{Z}_i, \boldsymbol{\theta})$$ which is not the same as the original expectation we wanted. Hence if we maximised this new expectation w.r.t. $\boldsymbol{\theta}$ we won't be doing Maximum Likelihood Estimation.

The EM algorithm therefore takes a different approach.


This is a few years old but I don't believe the current answer really covers it.

What the question boils down is whether $\mathbb{E}_{p(z\mid x;\theta')}[\log p(x, z\,;\, \theta')]$ (in EM lower-bound) is tractable even when $\log p(x ; \theta)=\log(\mathbb{E}_{p(z;\theta)}[ p(x \mid z\,;\,\theta)])$ (standard marginal likelihood) is not. The answer is often Yes. In your question you had left the logarithm out of the marginal likelihood, but that is a key observation; in general $\int \log p(x, z) dz$ may be tractable while $\log (\int p(x, z) dz)=\log p(x)$ is not -- one way to think about this is that sum of logs is much easier than log of sums. The point of EM is to lower-bound this intractable marginal likelihood into a lower-bound that has, possibly, tractable parts.

The EM algorithm is especially popular with models where $p(x, z)$ is in the exponential family, but the marginal $p(x)$ is an intractable non-exponential family distribution. This situation includes household names such as mixture of Gaussians, hidden Markov models, Gaussian state-space models etc. To see the tractability of the lower-bound term in exponential families, consider $\log p(x, z)$ in a general exponential family: \begin{align} \log p(x, z; \theta) = \langle \eta(\theta), t(x, z)\rangle - \log Z(\theta) + constant \end{align}
and therefore: \begin{align} \mathbb{E}_{p(z\mid x;\theta)}[\log p(x, z)] = \langle \eta(\theta), \mathbb{E}_{p(z\mid x;\theta)}[t(x, z)]\rangle - \log Z(\theta) + constant \end{align}
all we need are the expected sufficient statistics under the posterior which is, usually, tractable with conjugate exponential families. In contrast, as mentioned above $\log p(x ; \theta)$ is typically not tractable and nor even exponential family, hence the motivation for the whole EM algorithm.

Addition: the reason (related to the other answer) learning parameters by sampling is tricky:

$\nabla_{\theta}\log \mathbb{E}_{p(z;\theta)}[p(x \mid z; \theta)]=\frac{\nabla_{\theta}\mathbb{E}_{p(z;\theta)}[p(x \mid z; \theta)]}{\mathbb{E}_{p(z;\theta)}[p(x \mid z; \theta)]}$

ignoring the challenge of sampling the denominator, taking gradients of the numerator alone is tricky already as it doesn't only involve gradients $\nabla_{\theta}p(x \mid z;\theta)$, but instead you have to take gradients through the monte-carlo sample: \begin{align} \nabla_{\theta}\mathbb{E}_{p(z;\theta)}[p(x \mid z; \theta)]&=\int p(x \mid z; \theta)\nabla_{\theta}p(z;\theta)dz + \int p(z;\theta)\nabla_{\theta}p(x \mid z; \theta)dz \\ &= \mathbb{E}_{p(z;\theta)}[p(x \mid z; \theta)\nabla_{\theta}\log p(z;\theta)]+\mathbb{E}_{p(z;\theta)}[\nabla_{\theta} p(x \mid z; \theta)] \\ &\approx \frac{1}{M}\sum_{i=1}^M p(x \mid z_i; \theta)\nabla_{\theta}\log p(z_i;\theta) + \frac{1}{M}\sum_{i=1}^M \nabla_{\theta} p(x \mid z_i; \theta) \end{align} Usually these gradients have a very large variance (there are some tricks to alleviate such issues, for example the reparametrization trick). However, I believe if this was easy to get these gradients accurately I would think this approach could work, which would contradict the earlier answer.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.