学习笔记——EM算法

EM算法是一种迭代算法，用于含有隐变量(hidden variable)的概率模型参数的极大似然估计，或极大后验概率估计。EM算法的每次迭代由两步组成：E步，求期望(expectation)；M步，求极大（Maximization）。

EM算法的引入

---

给一些观察数据，可以使用极大似然估计法，或贝叶斯估计法估计模型参数。但是当模型含有隐变量时，就不能简单地使用这些方法。有些时候，参数的极大似然估计问题没有解析解，只能通过迭代的方法求解，EM算法就是可以用于求解这个问题的一种迭代算法。

• EM算法
• 输入：观测变量数据Y，隐变量数据Z，联合分布(P(Y, X | heta))，条件分布(P(Z|Y, heta))；
• 输出：模型参数( heta)

1. 选择参数的初值( heta^{(0)})，开始迭代；
2. E步：记( heta^{(0)})为第(i)次迭代参数( heta)的估计值，在第(i+1)次迭代的E步，计算$$Q( heta, heta^{(i)}) = E_Z[log P(Y, Z| heta) | Y, heta^{(i)}]$$ $$= sum_Z log P(Y, Z | heta) P (Z|Y, heta^{(i)})$$ 这里，(P(Z|Y, heta^{(i)}))是在给定观测数据Y和当前的参数估计( heta^{(i)})下隐变量数据Z的条件概率分布；
3. M步：求使(Q( heta, heta^{(i)}))极大化的( heta)，确定第(i+1)次迭代的参数的估计值( heta^{(i+1)}) $$ heta^{(i+1)} = arg max_{ heta} Q( heta, heta^{(i)})$$
4. 重复第(2)步和第(3)步，直到收敛。

其中的函数(Q( heta, heta^{(i)}))是EM算法的核心。称为Q函数。

• Q函数：完全数据的对数似然函数(log P(Y, Z | heta))关于在给定观测数据Y和当前参数( heta^{(i)})下对未观测数据Z的条件概率分布(P(Z|Y, heta^{i}))的期望，$$Q( heta, heta^{(i)}) = E_Z[log P(Y, Z| heta) | Y, heta^{(i)}]$$ $$= sum_Z P(Z | Y, heta^{(i)}) log P(Y, Z | heta)$$

• 可以证明，每次迭代使似然函数增大或达到局部极值。

• EM算法可以用于生成模型的非监督学习，生成模型由联合概率分布(P(X, Y))表示，可以认为非监督学习训练数据是联合概率分布产生的数据，X为观测数据，Y为未观测数据。

• 在应用中，初值的选择变得非常重要。常用的办法是选取几个不同的初值进行迭代，然后对得到的各个估计值加以比较，从中选择最好的。

EM算法在高斯混合模型学习中的应用

---

• 高斯混合模型是指具有以下形式的概率分布模型：$$P(y | heta) = sum_{k = 1}^K alpha_k phi(y| heta_k)$$ 其中，(alpha_k)是系数，(alpha_k geq 0)，(sum_{k=1}^K alpha_k = 1)；(phi)是高斯分布密度，( heta_k = (mu_k, sigma_k^2))，$$phi(y| heta_k) = frac{1}{sqrt{2pi} sigma_k} exp (-frac{(y - mu_k)^2}{2 sigma_k^2})$$称为第k个分模型。

使用EM算法进行高斯混合模型参数估计

假设观测数据(y_1,y_2,...,y_N)由高斯混合模型生成，$$P(y | heta) = sum_{k = 1}^K alpha_k phi(y| heta_k)$$其中，( heta = (alpha_1, alpha_2, ..., alpha_K; heta_1, heta_2,..., heta_K))。我们用EM算法估计高斯混合模型的参数( heta)。

1. 明确隐变量，写出完全数据的对数似然函数

   可以认为观测数据是这样产生的：首先依概率(alpha_k)选择第k个分模型，然后根据这个分模型的概率分布生成一个数据。计算出完全数据的对数似然函数为$$log P(y, gamma | heta) = sum_{k = 1}^K n_k log alpha_k + sum_{j = 1}^N gamma_{jk} [log(frac{1}{sqrt{2pi}}) - log sigma_k - frac{1}{2sigma^2} (y_i - mu_k)^2]$$其中，(gamma_{jk})为1代表第j个观测来自第k个分模型，否则为0，(n_k = sum_{j=1}^N gamma_{jk})，(sum_{k=1}^K n_k = N)。
   观测数据(y_j)及未观测数据(gamma_{jk})，完全数据是((y_j,gamma_{j1},gamma_{j2},...,gamma_{jK}, j = 1,2,...,N)

2. EM算法的E步：确定Q函数
   计算出Q函数为：$$Q( heta, heta^{(i)}) = sum_{k=1}^K n_k log alpha_k + sum_{k=1}^N hat{gamma}{jk} [log (frac{1}{sqrt{2pi}}) - log sigma_k - frac{1}{2sigma_k^2} (y_i - mu_k)^2]$$ 其中，$$hat{gamma}{jk} = E(gamma_{jk} | y, heta) = frac{alpha_k phi (y_i | heta_k)}{sum_{k=1}^K alpha_k phi (y_i | heta_k)}$$ $$n_k = sum_{j = 1}^N Egamma_{jk}$$

3. 确定EM算法中的M步
   迭代的M步是求函数(Q( heta, heta^{(i)}))对( heta)的极大值，即求新一轮迭代的模型参数：$$ heta^{(i+1)} = arg max_{ heta} Q( heta, heta^{(i)})$$
   通过求偏导可以得到计算相应参数的表达式。

• 高斯混合模型参数估计的EM算法
• 输入：观测数据(y_1, y_2,...,y_N)，高斯混合模型
• 输出：高斯混合模型参数。

1. 取参数的初始值开始迭代
2. E步：依据当前模型参数，计算分模型k对观测数据(y_i)的响应度$$hat{gamma}{jk} = frac{alpha_k phi (y_i | heta_k)}{sum{k=1}^K alpha_k phi (y_i | heta_k)}, j = 1,2,...,N; k = 1,2,...,K$$
3. M步：计算新一轮迭代的模型参数：$$hat{mu_k} = frac{sum_{j=1}^N hat{gamma}{jk} y_j}{sum{j = 1}^N hat{gamma}{jk}}, k = 1,2,...,K$$ $$hat{sigma}k^2 = frac{sum{j=1}^N hat{gamma}{jk} (y_j - mu_k)^2}{sum_{k=1}^N hat{gamma}{jk}}, k = 1,2,...,K$$ $$hat{alpha}k = frac{sum{j=1}^N hat{gamma}{jk}}{N}, k = 1,2,...,K$$
4. 重复第（2）步和第（3）步，直到收敛。

• EM算法还可以解释为F函数的极大-极大算法，基于这个解释有若干变形与推广，如广义期望极大（GEM）算法。

(注：本文为读书笔记与总结，侧重算法原理，来源为[《统计学习方法》](http://book.douban.com/subject/10590856/)一书第九章)
作者：[rubbninja](http://www.cnblogs.com/rubbninja/) 出处：[http://www.cnblogs.com/rubbninja/](http://www.cnblogs.com/rubbninja/) 关于作者：目前主要研究领域为机器学习与无线定位技术，欢迎讨论与指正！