PRML学习总结（9）——Mixture Models and EM

9.1 K-means Clustering

⾸先，我们考虑寻找多维空间中数据点的分组或聚类的问题。假设我们有⼀个数据集 { x 1 , … , x N } \left\{\mathbf{x}_{1}, \ldots, \mathbf{x}_{N}\right\} {x1​,…,xN​}，它由$D$维欧⼏⾥得空间中的随机变量$\mathbf{x}$的$N$次观测组成。我们的⽬标是将数据集划分为$K$个类别。直观上讲，我们会认为由⼀组数据点构成的⼀个聚类中，聚类内部点之间的距离应该⼩于数据点与聚类外部的点之间的距离。引⼊⼀组$D$维向量${\mathbf{\mu }}_k$，其中$k=1,...,K$，且${\mathbf{\mu }}_k$是与第$k$个聚类关联的⼀个代表。正如我们将看到的那样，我们可以认为${\mathbf{\mu }}_k$表⽰了聚类的中⼼。我们的⽬标是找到数据点分别属于的聚类，以及⼀组向量 { μ k } \{\boldsymbol{\mu}_{k}\} {μk​}，使得每个数据点和与它最近的向量k之间的距离的平⽅和最⼩。目标函数如下$$J=\sum _{n=1}^N\sum _{k=1}^K{r}_{nk}{\Vert {\mathbf{x}}_n-{\mathbf{\mu }}_k\Vert }^2$$其中$r_{nk}$代表第$n$个数据是否属于第$k$类，为one-hot向量。我们的目标就是求得 { r n k } \{r_{nk}\} {rnk​}和 { μ k } \{\boldsymbol{\mu}_{k}\} {μk​}，一般采用迭代的方式：（E-step）首先随机初始化 { μ k } \{\boldsymbol{\mu}_{k}\} {μk​}，然后关于 { r n k } \{r_{nk}\} {rnk​}优化$J$；(M-step)之后固定 { r n k } \{r_{nk}\} {rnk​}，在关于 { μ k } \{\boldsymbol{\mu}_{k}\} {μk​}优化$J$。如此迭代知道收敛为止。
在E-step时，$$r_{nk}=\{\begin{array}{ll}1& \text{ if }k=\arg {\min }_j{\Vert {\mathbf{x}}_n-{\mathbf{\mu }}_j\Vert }^2\\ 0& \text{ otherwise }\end{array}$$在M-step时，可以利用求导$$2\sum _{n=1}^N{r}_{nk}\left({\mathbf{x}}_n-{\mathbf{\mu }}_k\right)=0$$$${\mathbf{\mu }}_k=\frac{{\sum }_n{r}_{nk}{\mathbf{x}}_n}{{\sum }_n{r}_{nk}}$$需要说明的是K-means有可能会收敛到局部最优值。当初始化不恰当时，会导致较差的结果。

9.2 Mixtures of Gaussians

概率图模型如下
$$p(\mathbf{z})=\prod _{k=1}^K{\pi }_k^{z_k}$$$$p(\mathbf{x}|\mathbf{z})=\prod _{k=1}^K\mathcal{N}{\left(\mathbf{x}|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)}^{z_k}$$$$p(\mathbf{x})=\sum _{\mathbf{z}}p(\mathbf{z})p(\mathbf{x}|\mathbf{z})=\sum _{k=1}^K{\pi }_k\mathcal{N}\left(\mathbf{x}|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)$$在利用EM算法学习以上模型时，还有一个重要的$$\begin{array}{rl}\gamma \left(z_k\right)\equiv p\left(z_k=1|\mathbf{x}\right)& =\frac{p\left(z_k=1\right)p\left(\mathbf{x}|z_k=1\right)}{{\sum }_{j=1}^{K}p\left(z_j=1\right)p\left(\mathbf{x}|z_j=1\right)}\\ =& \frac{{\pi }_k\mathcal{N}\left(\mathbf{x}|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)}{{\sum }_{j=1}^K{\pi }_j\mathcal{N}\left(\mathbf{x}|{\mathbf{\mu }}_j,{\mathbf{\Sigma }}_j\right)}\end{array}$$这个为一个后验分布。我们可以利用GMM生成数据

9.2.1 Maximum likelihood

现在我们的目标是给定观测数据 { x 1 , … , x N } \left\{\mathbf{x}_{1}, \ldots, \mathbf{x}_{N}\right\} {x1​,…,xN​}，利用GMM对其进行建模。概率图模型如下
似然函数为$$\ln p(\mathbf{X}|\mathbf{\pi },\mathbf{\mu },\mathbf{\Sigma })=\sum _{n=1}^N\ln \left\{\sum _{k=1}^K{\pi }_k\mathcal{N}\left({\mathbf{x}}_n|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)\right\}$$在我们讨论如何最⼤化这个函数之前，有必要强调⼀下由于奇异性的存在造成的应⽤于⾼斯混合模型的最⼤似然框架中的⼀个⼤问题。为了简化起见，我们考虑⼀个⾼斯混合模型，它的分量的协⽅差矩阵为${\mathbf{\Sigma }}_k={\sigma }_k^2\mathbf{I}$，其中$\mathbf{I}$是⼀个单位矩阵，虽然结论对于⼀般的协⽅差矩阵仍然成⽴。假设混合模型的第$j$个分量的均值${\mathbf{\mu }}_j$与某个数据点完全相同，即对于某个$n$值，${\mathbf{\mu }}_j={\mathbf{x}}_n$。这样，这个数据点会为似然函数贡献⼀项，形式为$$\mathcal{N}\left({\mathbf{x}}_n|{\mathbf{x}}_n,{\sigma }_j^2\mathbf{I}\right)=\frac{1}{(2\pi {)}^{1/2}}\frac{1}{{\sigma }_j}$$当${\sigma }_j\to 0$时，似然就会无穷大，而最大化似然时，这就是一个奇异的问题。在GMM中，最大似然时，会很容易造成这种过拟合问题。而对于单高斯模型，则不会，假如退化为一个点，那么对于其他点来说，其似然贡献为0，多个0相差必然导致退化的似然不足以保证为最大。我们后⾯会看到，如果我们使⽤贝叶斯⽅法，那么这种困难之处就不会出现。但是现阶段，我们只需注意，将最⼤似然⽅法应⽤到⾼斯混合模型中时必须避免这种病态解，并且寻找表现较好的似然函数的局部极⼤值。我们可以使⽤合适的启发式⽅法来避免这种奇异性，例如如果检测到⾼斯分量收缩到⼀个点，那么就将它的均值重新设定为⼀个随机选择的值，并且重新将它的⽅差设置为某个较⼤的值，然后继续最优化。
GMM可以直接采用梯度下降的方式进行求解，但下面介绍一种更加优雅的方法——EM算法。

9.2.2 EM for Gaussian mixtures

首先对均值${\mathbf{\mu }}_k$进行求导$$0=-\sum _{n=1}^N\frac{{\pi }_k\mathcal{N}\left({\mathbf{x}}_n|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)}{{\sum }_j{\pi }_j\mathcal{N}\left({\mathbf{x}}_n|{\mathbf{\mu }}_j,{\mathbf{\Sigma }}_j\right)}{\mathbf{\Sigma }}_k\left({\mathbf{x}}_n-{\mathbf{\mu }}_k\right)$$从而得到$${\mathbf{\mu }}_k=\frac{1}{N_k}\sum _{n=1}^N\gamma \left(z_{nk}\right){\mathbf{x}}_n$$其中$N_k={\sum }_{n=1}^N\gamma \left(z_{nk}\right)$
同理对方差${\mathbf{\Sigma }}_k$进行求导$${\mathbf{\Sigma }}_k=\frac{1}{N_k}\sum _{n=1}^N\gamma \left(z_{nk}\right)\left({\mathbf{x}}_n-{\mathbf{\mu }}_k\right){\left({\mathbf{x}}_n-{\mathbf{\mu }}_k\right)}^{\mathrm{T}}$$同理对$\pi$求导$${\pi }_k=\frac{N_k}{N}$$这些结果确实给出了⼀个简单的迭代⽅法来寻找问题的最⼤似然解。正如我们将看到的那样，这个迭代过程是EM算法应⽤于⾼斯混合模型的⼀个实例。

9.3 An Alternative View of EM

本节中，我们介绍EM算法的另⼀种观点，其中潜在变量起着重要的作⽤。我们⾸先使⽤⼀种抽象的⽅式讨论这种⽅法，然后我们再次考虑⾼斯混合模型的例⼦，来具体说明这个模型。
EM算法的⽬标是找到具有潜在变量的模型的最⼤似然解。$$\ln p(\mathbf{X}|\mathbf{\theta })=\ln \left\{\sum _{\mathbf{Z}}p(\mathbf{X},\mathbf{Z}|\mathbf{\theta })\right\}$$
EM算法可以用MAP中，其中E-step保持不变，而M-step只需要$\mathcal{Q}\left(\mathbf{\theta },{\mathbf{\theta }}^{\mathrm{o}\mathrm{l}\mathrm{d}}\right)+\ln p(\mathbf{\theta })$。引入先验能够防止GMM中的奇异问题。

9.3.1 Gaussian mixtures revisited

对于全部都可观测的GMM模型来说
可以直接利用最大似然求解$$p(\mathbf{X},\mathbf{Z}|\mathbf{\mu },\mathbf{\Sigma },\mathbf{\pi })=\prod _{n=1}^N\prod _{k=1}^K{\pi }_k^{z_{nk}}\mathcal{N}{\left({\mathbf{x}}_n|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)}^{z_{nk}}$$ ln ⁡ p ( X , Z ∣ μ , Σ , π ) = ∑ n = 1 N ∑ k = 1 K z n k { ln ⁡ π k + ln ⁡ N ( x n ∣ μ k , Σ k ) } \ln p(\mathbf{X}, \mathbf{Z} | \boldsymbol{\mu}, \boldsymbol{\Sigma}, \boldsymbol{\pi})=\sum_{n=1}^{N} \sum_{k=1}^{K} z_{n k}\left\{\ln \pi_{k}+\ln \mathcal{N}\left(\mathbf{x}_{n} | \boldsymbol{\mu}_{k}, \boldsymbol{\Sigma}_{k}\right)\right\} lnp(X,Z∣μ,Σ,π)=n=1∑N​k=1∑K​znk​{lnπk​+lnN(xn​∣μk​,Σk​)}对于均值和方差的优化，就是$K$个独立的高斯对于属于每个类的单高斯拟合。而${\pi }_k=\frac{1}{N}{\sum }_{n=1}^N{z}_{nk}$。而在实际中$\mathbf{Z}$是不可观测的，则不可能如上有解析解。 E Z [ ln ⁡ p ( X , Z ∣ μ , Σ , π ) ] = ∑ n = 1 N ∑ k = 1 K γ ( z n k ) { ln ⁡ π k + ln ⁡ N ( x n ∣ μ k , Σ k ) } \mathbb{E}_{\mathbf{Z}}[\ln p(\mathbf{X}, \mathbf{Z} | \boldsymbol{\mu}, \boldsymbol{\Sigma}, \boldsymbol{\pi})]=\sum_{n=1}^{N} \sum_{k=1}^{K} \gamma\left(z_{n k}\right)\left\{\ln \pi_{k}+\ln \mathcal{N}\left(\mathbf{x}_{n} | \boldsymbol{\mu}_{k}, \boldsymbol{\Sigma}_{k}\right)\right\} EZ​[lnp(X,Z∣μ,Σ,π)]=n=1∑N​k=1∑K​γ(znk​){lnπk​+lnN(xn​∣μk​,Σk​)}然后会得到之前一样的结果。

9.3.2 Relation to K-means

对⽐⾼斯模型的EM算法与K均值算法，可以看到⼆者有很强的相似性。K均值算法对数据点的聚类进⾏了“硬”分配，即每个数据点只属于唯⼀的聚类，⽽EM算法基于后验概率分布，进⾏了⼀个“软”分配。实际上，我们可以将K均值算法看成⾼斯混合模型的EM算法的⼀个特殊的极限情况，如下所述。
考虑⼀个⾼斯混合模型，其中混合分量的协⽅差矩阵为$ϵ\mathbf{I}$，$ϵ$是⼀个被所有分量共享的⽅差参数，$\mathbf{I}$是单位矩阵，从⽽$$p\left(\mathbf{x}|{\mathbf{\mu }}_k,{\mathbf{\Sigma }}_k\right)=\frac{1}{(2\pi ϵ{)}^{1/2}}\exp \left\{-\frac{1}{2ϵ}{\Vert \mathbf{x}-{\mathbf{\mu }}_k\Vert }^2\right\}$$我们现在考虑$K$个这种形式的⾼斯分布组成的混合模型的EM算法，其中我们将$ϵ$看做⼀个固定的常数，⽽不是⼀个需要重新估计的参数。$$\gamma \left(z_{nk}\right)=\frac{{\pi }_k\exp \left\{-{\Vert {\mathbf{x}}_n-{\mathbf{\mu }}_k\Vert }^2/2ϵ\right\}}{{\sum }_j{\pi }_j\exp \left\{-{\Vert {\mathbf{x}}_n-{\mathbf{\mu }}_j\Vert }^2/2ϵ\right\}}$$当$ϵ\to 0$时，会发现${\Vert {\mathbf{x}}_n-{\mathbf{\mu }}_j\Vert }^2$最小的项会以最慢速度变为0，那么$\gamma \left(z_{nk}\right)$会为0，而$\gamma \left(z_{nj}\right)$为1。从而就跟K-means一样变为硬划分了$\gamma \left(z_{nk}\right)\to r_{nk}$。在同样的极限下$${\mathbb{E}}_{\mathbf{Z}}[\ln p(\mathbf{X},\mathbf{Z}|\mathbf{\mu },\mathbf{\Sigma },\mathbf{\pi })]\to -\frac{1}{2}\sum _{n=1}^N\sum _{k=1}^K{r}_{nk}{\Vert {\mathbf{x}}_n-{\mathbf{\mu }}_k\Vert }^2+\text{ const. }$$跟K-means的目标函数一致。

9.3.3 Mixtures of Bernoulli distributions

对于$D$个二值变量$x_i$$$p(\mathbf{x}|\mathbf{\mu })=\prod _{i=1}^D{\mu }_i^{x_i}{\left(1-{\mu }_i\right)}^{\left(1-x_i\right)}$$其均值和方差 E [ x ] = μ cov ⁡ [ x ] = diag ⁡ { μ i ( 1 − μ i ) } \begin{aligned} \mathbb{E}[\mathbf{x}] &=\mu \\ \operatorname{cov}[\mathbf{x}] &=\operatorname{diag}\left\{\mu_{i}\left(1-\mu_{i}\right)\right\} \end{aligned} E[x]cov[x]​=μ=diag{μi​(1−μi​)}​下面考虑有限混合Bernoulli模型$$p(\mathbf{x}|\mathbf{\mu },\mathbf{\pi })=\sum _{k=1}^K{\pi }_{k}p\left(\mathbf{x}|{\mathbf{\mu }}_k\right)$$where μ = { μ 1 , … , μ K } , π = { π 1 , … , π K } \mu=\left\{\mu_{1}, \ldots, \mu_{K}\right\}, \pi=\left\{\pi_{1}, \ldots, \pi_{K}\right\} μ={μ1​,…,μK​},π={π1​,…,πK​}，对于其中某一个分布有$$p\left(\mathbf{x}|{\mathbf{\mu }}_k\right)=\prod _{i=1}^D{\mu }_{ki}^{x_i}{\left(1-{\mu }_{ki}\right)}^{\left(1-x_i\right)}$$该混合模型的均值和方差$$\begin{array}{rl}\mathbb{E}[\mathbf{x}]& =\sum _{k=1}^K{\pi }_k{\mathbf{\mu }}_k\\ \mathrm{cov}[\mathbf{x}]& =\sum _{k=1}^K{\pi }_k\left\{{\mathbf{\Sigma }}_k+{\mathbf{\mu }}_k{\mathbf{\mu }}_k^{\mathrm{T}}\right\}-\mathbb{E}[\mathbf{x}]\mathbb{E}[\mathbf{x}{]}^{\mathrm{T}}\end{array}$$其中 Σ k = diag ⁡ { μ k i ( 1 − μ k i ) } \boldsymbol{\Sigma}_{k}=\operatorname{diag}\left\{\mu_{k i}\left(1-\mu_{k i}\right)\right\} Σk​=diag{μki​(1−μki​)}。由于协⽅差矩阵$\mathrm{cov}[\mathbf{x}]$不再是对⾓矩阵，因此混合分布可以描述变量之间的相关性，这与单⼀的伯努利分布不同。下面我们考虑一个数据集 X = { x 1 , … , x N } \mathbf{X}=\left\{\mathbf{x}_{1}, \ldots, \mathbf{x}_{N}\right\} X={x1​,…,xN​}，该模型的似然为$$\ln p(\mathbf{X}|\mathbf{\mu },\mathbf{\pi })=\sum _{n=1}^N\ln \left\{\sum _{k=1}^K{\pi }_{k}p\left({\mathbf{x}}_n|{\mathbf{\mu }}_k\right)\right\}$$跟GMM一样，无法直接求解，因此考虑使用EM算法。那么需要引入一个隐变量$\mathbf{z}={\left(z_1,\dots ,z_K\right)}^{\mathrm{T}}$，表示每个$\mathbf{x}$属于第$k$个分量，也就是说其为一个one-hot向量，那么$$p(\mathbf{x}|\mathbf{z},\mathbf{\mu })=\prod _{k=1}^{K}p{\left(\mathbf{x}|{\mathbf{\mu }}_k\right)}^{z_k}$$其先验与GMM一样$$p(\mathbf{z}|\mathbf{\pi })=\prod _{k=1}^K{\pi }_k^{z_k}$$为了推导EM算法，我们⾸先写出完整数据的对数似然函数，形式为 ln ⁡ p ( X , Z ∣ μ , π ) = ∑ n = 1 N ∑ k = 1 K z n k { ln ⁡ π k + ∑ i = 1 D [ x n i ln ⁡ μ k i + ( 1 − x n i ) ln ⁡ ( 1 − μ k i ) ] } \begin{array}{c}{\ln p(\mathbf{X}, \mathbf{Z} | \boldsymbol{\mu}, \boldsymbol{\pi})=\sum_{n=1}^{N} \sum_{k=1}^{K} z_{n k}\left\{\ln \pi_{k}\right.} {+\sum_{i=1}^{D}\left[x_{n i} \ln \mu_{k i}+\left(1-x_{n i}\right) \ln \left(1-\mu_{k i}\right)\right] \}}\end{array} lnp(X,Z∣μ,π)=∑n=1N​∑k=1K​znk​{lnπk​+∑i=1D​[xni​lnμki​+(1−xni​)ln(1−μki​)]}​接下来我们取完整数据对数似然函数关于潜在变量后验概率分布的期望，得 E Z [ ln ⁡ p ( X , Z ∣ μ , π ) ] = ∑ n = 1 N ∑ k = 1 K γ ( z n k ) { ln ⁡ π k + ∑ i = 1 D [ x n i ln ⁡ μ k i + ( 1 − x n i ) ln ⁡ ( 1 − μ k i ) ] } \begin{array}{c}{\mathbb{E}_{\mathbf{Z}}[\ln p(\mathbf{X}, \mathbf{Z} | \boldsymbol{\mu}, \boldsymbol{\pi})]=\sum_{n=1}^{N} \sum_{k=1}^{K} \gamma\left(z_{n k}\right)\left\{\ln \pi_{k}\right.} {+\sum_{i=1}^{D}\left[x_{n i} \ln \mu_{k i}+\left(1-x_{n i}\right) \ln \left(1-\mu_{k i}\right)\right] \}}\end{array} EZ​[lnp(X,Z∣μ,π)]=∑n=1N​∑k=1K​γ(znk​){lnπk​+∑i=1D​[xni​lnμki​+(1−xni​)ln(1−μki​)]}​其中$$\begin{array}{rl}\gamma \left(z_{nk}\right)=\mathbb{E}\left[z_{nk}\right]=& \frac{{\sum }_{z_{nk}}{z}_{nk}{\left[{\pi }_{k}p\left({\mathbf{x}}_n|{\mathbf{\mu }}_k\right)\right]}^{z_{nk}}}{{\sum }_{z_{nj}}{\left[{\pi }_{j}p\left({\mathbf{x}}_n|{\mathbf{\mu }}_j\right)\right]}^{z_{nj}}}\\ =& \frac{{\pi }_{k}p\left({\mathbf{x}}_n|{\mathbf{\mu }}_k\right)}{{\sum }_{j=1}^K{\pi }_{j}p\left({\mathbf{x}}_n|{\mathbf{\mu }}_j\right)}\end{array}$$在M步骤中$${\mathbf{\mu }}_k=\frac{1}{N_k}\sum _{n=1}^N\gamma \left(z_{nk}\right){\mathbf{x}}_n$$$${\pi }_k=\frac{{\sum }_{n=1}^N\gamma \left(z_{nk}\right)}{N}$$利用该模型对手写数字建模如下

9.3.4 EM for Bayesian linear regression

之前在贝叶斯线性回归问题中，对于超参数$\alpha$ and $\beta$是通过最大化evidence$p(\mathbf{t}|\alpha ,\beta )$得到的，其中的隐变量被$\mathbf{w}$积分积掉了。那么可以采取EM的框架求解这个问题：E-step，求得$\mathbf{w}$在当前$\alpha$和$\beta$的后验分布；M-step，求完全数据集下log似然并最大化。$$\ln p(\mathbf{t},\mathbf{w}|\alpha ,\beta )=\ln p(\mathbf{t}|\mathbf{w},\beta )+\ln p(\mathbf{w}|\alpha )$$具体分析可以参见书。

9.4 The EM Algorithm in General

EM算法能很好地解决对于有隐变量的概率模型的最大似然。我们的目标为最大化下式$$p(\mathbf{X}|\mathbf{\theta })=\sum _{\mathbf{Z}}p(\mathbf{X},\mathbf{Z}|\mathbf{\theta })$$引入分布$q(\mathbf{Z})$，则$$\ln p(\mathbf{X}|\mathbf{\theta })=\mathcal{L}(q,\mathbf{\theta })+\mathrm{K}\mathrm{L}(q\Vert p)$$其中$$\begin{array}{rl}\mathcal{L}(q,\mathbf{\theta })& =\sum _{\mathbf{Z}}q(\mathbf{Z})\ln \left\{\frac{p(\mathbf{X},\mathbf{Z}|\mathbf{\theta })}{q(\mathbf{Z})}\right\}\\ \mathrm{K}\mathrm{L}(q\Vert p)& =-\sum _{\mathbf{Z}}q(\mathbf{Z})\ln \left\{\frac{p(\mathbf{Z}|\mathbf{X},\mathbf{\theta })}{q(\mathbf{Z})}\right\}\end{array}$$
对于E-step来说，当参数固定时，$q(\mathbf{Z})=p\left(\mathbf{Z}|\mathbf{X},{\mathbf{\theta }}^{\text{ old }}\right)$
对于M-step来说，$$\begin{array}{rl}\mathcal{L}(q,\mathbf{\theta })& =\sum _{\mathbf{Z}}p\left(\mathbf{Z}|\mathbf{X},{\mathbf{\theta }}^{\mathrm{o}\mathrm{l}\mathrm{d}}\right)\ln p(\mathbf{X},\mathbf{Z}|\mathbf{\theta })-\sum _{\mathbf{Z}}p\left(\mathbf{Z}|\mathbf{X},{\mathbf{\theta }}^{\mathrm{o}\mathrm{l}\mathrm{d}}\right)\ln p\left(\mathbf{Z}|\mathbf{X},{\mathbf{\theta }}^{\text{ old }}\right)\\ & =\mathcal{Q}\left(\mathbf{\theta },{\mathbf{\theta }}^{\text{ old }}\right)+\text{ const }\end{array}$$在参数空间看，整个过程如下我们还可以使用EM算法去求最大后验分布$p(\mathbf{\theta }|\mathbf{X})$，MAP引入了$\mathbf{\theta }$的先验。$$\ln p(\mathbf{\theta }|\mathbf{X})=\ln p(\mathbf{\theta },\mathbf{X})-\ln p(\mathbf{X})$$利用之前的分解，有$$\begin{array}{rl}\ln p(\mathbf{\theta }|\mathbf{X})& =\mathcal{L}(q,\mathbf{\theta })+\mathrm{K}\mathrm{L}(q\Vert p)+\ln p(\mathbf{\theta })-\ln p(\mathbf{X})\\ & ⩾\mathcal{L}(q,\mathbf{\theta })+\ln p(\mathbf{\theta })-\ln p(\mathbf{X})\end{array}$$E-step跟原来一样，M-step则需要考虑$\ln p(\mathbf{\theta })$。