EM and Kullback-Leibler divergence

Question

Let $f$ be a density on $\mathbb{R}^{p}$. Let $f_{\theta} = \sum_{i=1}^{d} \alpha_{i}\mathcal{N}_{p}(\cdot \, ; \, \theta_{i})$ be a mixture of $d$ Gaussian distributions on $\mathbb{R}^{p}$. For each $i$, $\theta_{i}$ is a vector of parameters (mean and covariance) which characterize the $i$-th component of the mixture. I would like to minimize the Kullback-Leibler divergence $K(f||f_{\theta})$. It amounts to finding $\theta^{\ast}$ such that :

$$ \theta \in \mathop{\mathrm{argmax}} \limits_{\theta} \int \log \big( f_{\theta}(x) \big) f(x) \, dx = \int \log \Big( \sum_{i=1}^{d} \alpha_{i} \mathcal{N}_{p}( x \, ; \, \theta_{i} ) \Big) f(x) \, dx $$

How can the EM algorithm be used to find $\theta^{\ast}$ ?

The optimization problem may be rewritten :

$$ \theta^{\ast} \in \mathop{\mathrm{argmax}} \limits_{\theta} \, \mathbb{E}_{f}\left[ \log f_{\theta}(X) \right] $$

If I understand correctly, we have a $n$-sample $(Y_{1},\ldots,Y_{n})$ from $f$ and we know, from Monte Carlo Integration that

$$ \frac{1}{n} \log \big( f_{\theta}(Y_i) \big) $$

is an approximation of $\mathbb{E}_{f}\left[ f_{\theta}(X) \right]$. What I do not really understand is why $\theta^{\ast}$ can be obtained as follows :

$$ \theta^{\ast} = \mathop{\mathrm{argmax}} \limits_{\theta} \frac{1}{n} \sum_{i=1}^{n} \log \big( f_{\theta}(Y_i) \big). $$

Answer 1

This problem is solved in our 2007 paper "Convergence of adaptive mixtures of importance sampling schemes" that considers the optimisation of $\mathbb{E}_{f}\left[ \log f_{\theta}(X) \right]$ using an integrated EM algorithm, that is an EM algorithm with no data. In short, you can run an EM algorithm with an extra integral over the (virtual) data.

If we consider the target $$\arg\max_{\alpha,\theta}\int \log\left( \sum_{d=1}^D\alpha_d q_d(x; \theta_d) \right) \pi(x)\,\text{d}x $$ the EM algorithm relies on introducing the latent variable $z$ such that the joint density $f$ of $x$ and $z$ satisfies $$ f(z)=\alpha_z \quad \text{and} \quad f(x|z)=q_z(x;\theta_z) \, . $$ The expectation corresponding to the E step of the EM algorithm is the expected complete log-likelihood, namely, at iteration $t$ of the algorithm, $$ \mathbb{E}_{\pi}^X \left[ \mathbb{E}_{(\alpha^t,\theta^t)}^Z \left\{ \log\left( \alpha_Z q_Z(X; \theta_Z) \right) | X \right\} \right]\,, $$ where the inner expectation is computed under the conditional distribution of $Z$ in the mixture model given the current value $(\alpha^t,\theta^t)$ of the parameters, i.e. $$ f(z|x)=\alpha^t_z q_z(x;\theta_z^t) \bigg/ \sum_{d=1}^D \alpha^t_d q_d(x;\theta_d^t)\,, $$ while the outer expectation is under the distribution $X \sim \pi$.

The corresponding M step amounts to setting the new parameters $(\alpha^{t+1},\theta^{t+1})$ equal to $$ (\alpha^{t+1},\theta^{t+1})=\arg\max_{(\alpha,\theta)}\,%\left( \mathbb{E}_\pi^X\left[\mathbb{E}_{(\alpha^t,\theta^t)}^Z\left\{\log(\alpha_Z q_Z(X; \theta_Z))|X\right\}\right]\,,%\right)\,. $$ It is straightforward to check that the convexity argument used for the EM algorithm also applies in this setup and, hence, that the Kullback-Leibler distance $(\mathfrak{K}(\pi,q_{(\alpha^t,\theta^t)}))_{t\geq 1}$ is a non-decreasing sequence. Setting $$ \rho_d(X;\alpha,\theta) = {\alpha_d q_d(X;\theta_d)}\bigg/{\sum_{\ell=1}^D \alpha_\ell q_\ell(X;\theta_\ell)} \, , $$ the maximisation program reduces to \begin{align*} & \alpha^{t+1}=\arg\max_{\alpha}\mathbb{E}_\pi^X\left[\sum_{d=1}^D \rho_d(X;\alpha^t,\theta^t)\log(\alpha_d)\right]\,, \\ & \theta^{t+1}=\arg\max_{\theta}\mathbb{E}_\pi^X\left[\sum_{d=1}^D \rho_d(X;\alpha^t,\theta^t)\log(q_d(X;\theta_d))\right]\,. \end{align*} where the first maximisation to be carried out under the constraint that $$\sum_{d=1}^D \alpha^{t+1}_d = 1$$ Hence, $$\begin{align*} & \alpha_d^{t+1}=\mathbb{E}_\pi^X\left[\rho_d(X;\alpha^t,\theta^t) \right]\,, \\ & \theta_d^{t+1}=\arg\max_{\theta_d}\mathbb{E}_\pi^X\left[ \rho_d(X;\alpha^t,\theta^t)\log(q_d(X;\theta_d)) \right]\,. \end{align*}$$