3. EM算法

3.1. 最大似然

机器学习的一般过程是：

  1. 对数据进行建模。比如 线性回归、逻辑回归

  2. 定义目标函数。比如，最小二乘法、交叉熵、最大似然。

  3. 极大（小）化目标函数，求解参数的值。

最大似然估计是一中普遍用于概率模型中定义目标函数的方法。 假设我们的数据样本集合为 \(\pmb{Y}\) (不是数据特征，是结果变量，比如二分类中的分类结果)，样本数量为N， 模型参数为 \(\pmb{\omega}\) ，似然函数为：

()\[object\ function = L(\mathbf{Y}|\pmb{\omega}) = P(\mathbf{Y}|\pmb{\omega}) = \prod_{i=1}^N P(Y_i|\pmb{\omega}) \tag{1}\]

备注

似然函数的概率意义为，当模型参数为 \(\pmb{\omega}\) 时，所有样本 \(\pmb{Y}\) 同时发生的概率，所以是连乘。

上式中有连乘操作，为了避免计算机浮点数溢出，会增加log把乘法变成加法（这对于极大化求解的参数值时等价的）。

()\[object\ function = ln L(\mathbf{Y}|\pmb{\omega}) = ln P(\mathbf{Y}|\pmb{\omega}) = \sum_{i=1}^N ln\ P(Y_i|\pmb{\omega})\]

我们极大化目标函数，求解出 \(\pmb{\omega}\)

()\[\begin{split}\pmb{\omega} &= \mathop{\arg\max} \ ln\ L(\mathbf{Y}|\pmb{\omega}) \\ &= \mathop{\arg\max}\ ln\ \prod_{i=1}^N l(Y_i|\pmb{\omega}) \\ &= \mathop{\arg\max}\ \sum_{i=1}^N\ ln\ l(Y_i|\pmb{\omega}) \ \ \ \ \text{(2)}\end{split}\]

3.2. EM原理

上文中，我们回顾了最大似然估计，先定义出目标函数，然后极大化目标函数求解，

有些时候我们的模型中，除了参数 \(\pmb{\omega}\) 外，还会包含另外的变量，而这个变量的值是未知的，我们称之为 隐变量 。

从上文中，我们知道目标函数是关于参数 \(\pmb{\omega}\) 的一个函数 \(object\ function = f(\pmb{\omega}|\mathbf{Y})\) , 如果模型中增加了额外的未知参数 \(\pmb{\theta}\) ，目标函数就可以改写成：

()\[object\ function = f(\pmb{\omega},\pmb{\theta}|\mathbf{Y})\]

比如在IRT模型中，除了项目参数 \(\pmb{\omega}\) 外， 学生能力参数 \(\pmb{\theta}\) 也是未知的，这时学生能力参数可以看做是隐变量。

根据边际化理论，我们可以通过对隐变量参数 \(\pmb{\theta}\) 进行边际化，以消除它的影响。

上述 目标函数需要转换成边际似然 。

备注

当不知道一个随机变量的取值时，对一个随机变量进行边际化（边缘化），就相当消去这个随机变量的影响。

()\[object\ function = lnL(\mathbf{Y}|\pmb{\omega}) = ln \sum_{\theta \in \Theta} L(\mathbf{Y},\pmb{\theta}|\pmb{\omega}) \tag{3}\]

备注

回顾前文的最大似然估计， 上式中的函数 \(L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})\) 是似然函数。

备注

在边际化的时候，对于连续值随机变量应该是积分操作，对于离散值随机变量是求和。这里为了简化公式，我们假设 \(\pmb{\theta}\) 是离散随机变量，这样就可以用求和的方式进行边际化。

公式(3)中存在两个问题。第一，其包含对数中求和的操作，这个基本无法处理；第二，包含另外一个未知参数 \(\pmb{\theta}\) ， 也就是前文讲的 缺失数据 。不要慌，我们可以通过一些变换解决这个问题。

定理1 Jensen不等式 :

当 f 是一个下凸函数，并且Y是一个随机变量时，有:

()\[Ef(Y) \ge f(EY)\]

当 f 是一个上凹函数，并且Y是一个随机变量时，有:

()\[Ef(Y) \le f(EY) \tag{4}\]

公式(3)是一个上凹函数（有极大值），结合Jensen不等式可以进行一些变换

()\[\begin{split}lnL(\mathbf{Y}|\pmb{\omega}) &= \ ln \sum_{\theta \in \Theta} L(\mathbf{Y},\pmb{\theta}|\pmb{\omega}) \\ &= ln \sum_{\theta \in \Theta} g(\pmb{\theta}) \frac{L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})}{g(\pmb{\theta})} \\ &\ge \sum_{\theta \in \Theta} g(\pmb{\theta}) ln \frac{L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})}{g(\pmb{\theta})}\end{split}\]

上式中， \(g(\pmb{\theta})\) 是 \(\pmb{\theta}\) 的概率密度（质量）函数，公式第2行，乘以一个 \(g(\pmb{\theta})\) ， 然后再除以一个 \(g(\pmb{\theta})\) 等式不变。 \(\sum_{\theta \in \Theta} g(\pmb{\theta})\) 相当于对后面的 部分 \(\frac{L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})}{g(\pmb{\theta})}\) 求期望，这样结合Jensen不等式得出第三行。

公式(5)提供了似然函数的一个下边界，极大化下边界函数，和极大化似然函数可以得到相似的结果。 \(g(\pmb{\theta})\) 我们定义成 \(\pmb{\theta}\) 的 后验概率密度函数 \(g(\pmb{\theta} | \mathbf{Y},\pmb{\omega})\) , 公式(5)就可以转换成

()\[lnL(\mathbf{Y}|\pmb{\omega}) = \ ln \sum_{\pmb{\theta} \in \Theta} L(\mathbf{Y},\pmb{\theta}|\pmb{\omega}) \Leftrightarrow \sum_{\pmb{\theta} \in \Theta} g(\pmb{\theta} | \mathbf{Y},\pmb{\omega}) ln \frac{L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})}{g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})} \tag{6}\]

通过Jensen不等式的转换， 我们把 对数ln内部的求和（积分）操作转化到外部，简化了计算。

备注

\(\pmb{\theta}\) 是学生能力值，我们看作一个随机变量。既然是随机变量就会有其概率分布。 对于连续值随机变量，其概率分布函数称为概率密度函数， 函数值不是概率值 ，需要进行积分才能得出概率值（不懂的，去恶补概率论基础）； 对于离散值随机变量，其概率分布函数称为概率质量函数， 函数值就是概率值 。 公式中的 \(g(\pmb{\theta} | \mathbf{Y},\pmb{\omega})\) 就是随机变量 \(\pmb{\theta}\) 的概率分布函数。 当然，这里 \(\pmb{\theta}\) 学生能力值应该是连续值，但是积分操作过于复杂，所以我们可以把 \(\pmb{\theta}\) 进行离散化， 把 \(\pmb{\theta}\) 转换成离散值，这样就不用积分操作了。离散化的方法以及 \(g(\pmb{\theta} | \mathbf{Y},\pmb{\omega})\) 如何求下文介绍。

依据对数的性质，公式(6)中的除法可以变成减法

()\[lnL(\mathbf{Y}|\pmb{\omega}) = \underbrace {\sum_{\pmb{\theta} \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega}) ln L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})}_{Q(\omega|\omega^{(t)})} - \underbrace {\sum_{\pmb{\theta} \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega}) ln \ g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})}_{H(\omega|\omega^{(t)})} \tag{8}\]

极大化目标函数是一个迭代的过程，我们令t表示迭代的序号，t=0,1,2,…。公式(8)中包含两部分， \(Q(\pmb{\omega}|\pmb{\omega}^{(t)})\) 和 \(H(\pmb{\omega}|\pmb{\omega}^{(t)})\) ， 其中 \(Q(\pmb{\omega}|\pmb{\omega}^{(t)})\) 表示 完全数据的似然函数 ，也就李航<统计学习>中讲的Q函数。同时，也是在EM算法的M步中要极大化的目标函数。

\(H(\pmb{\omega}|\pmb{\omega}^{(t)})\) 表示潜在能力变量的后验概率密度（其实是关于后验概率密度函数的一个函数）。 根据Jensen不等式，对于任意 \(\pmb{\omega} \in \pmb{\Omega}\) , \(H(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)}) \le H(\pmb{\omega}^{(t)}|\pmb{\omega}^{(t)})\) ，所以在极大化过程中可以忽略H部分 。

在EM算法的M步骤中，每次迭代 \(t \rightarrow t+1\) 都是最大化 \(Q(\pmb{\omega}|\pmb{\omega}^{(t)})\) ， 所以可以确保 \(Q(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)}) \ge Q(\pmb{\omega}^{(t)}|\pmb{\omega}^{(t)})\) 。

这样，随着每一次迭代，一定能确保对数似然函数的值（我们的目标是极大化对数似然函数）是增长的

()\[\begin{split}&lnL(\mathbf{Y}|\pmb{\omega}^{(t+1)}) - lnL(\mathbf{Y}|\pmb{\omega}^{(t)}) \\ &= [Q(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)}) - H(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)})] - [ Q(\pmb{\omega}^{(t)}|\pmb{\omega}^{(t)}) - H(\pmb{\omega}^{(t)}|\pmb{\omega}^{(t)})] \\ &= \underbrace{[Q(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)})- Q(\pmb{\omega}^{(t)}|\pmb{\omega}^{(t)})]}_{\ge0} + \underbrace{[H(\pmb{\omega}^{(t)}|\pmb{\omega}^{(t)}) - H(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)})]}_{\ge0} \\ &\ge 0\end{split}\]

完整的推导过程是:

()\[\begin{split}ln L(\mathbf{Y}|\pmb{\omega}) &= ln \sum_{\theta \in \Theta} L(\mathbf{Y},\pmb{\theta}|\pmb{\omega}) \\ &\Leftrightarrow \sum_{\theta \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega}) ln \frac{L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})}{g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})} \ \ \text{根据Jensen不等式}\\ &= \underbrace{\sum_{\theta \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})\ ln L(\mathbf{Y},\pmb{\theta}|\pmb{\omega}) }_Q - \underbrace{ \sum_{\theta \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega}) g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})}_H \\ &\Leftrightarrow \underbrace{\sum_{\theta \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})\ ln L(\mathbf{Y},\pmb{\theta}|\pmb{\omega}) }_Q \ \ \text{H部分可忽略}\end{split}\]

根据前文的推导，我们需要极大化Q函数部分，求得参数 \(\pmb{\omega}\) ,

()\[object\ function = Q = \sum_{\theta \in \Theta} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega})\ ln L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})\]

3.3. 隐变量的后验分布

随机变量 \(\pmb{\theta}\) 概率分布函数是什么呢？ 我们既然已经有了学生的作答记录， 那么，我们可以利用学生的作答记录来计算 \(\pmb{\theta}\) 的后验概率分布。所以公式(6)中使用的是 \(\pmb{\theta}\) 的后验概率。 根据贝叶斯定理，\(\pmb{\theta}\) 的后验概率密度函数可以定义为：

()\[g(\pmb{\theta} | \mathbf{Y},\pmb{\omega} ) = \frac{L( \mathbf{Y}|\pmb{\theta},\pmb{\omega})p(\pmb{\theta},\pmb{\omega})} {\sum_{\pmb{\theta}' \in \Theta} L(\mathbf{Y}|\pmb{\theta}',\pmb{\omega})p(\pmb{\theta}',\pmb{\omega}) }\]

备注

什么是先验概率、后验概率，这里不介绍了，自己恶补概率论基础吧。 根据贝叶斯定理，公式(7)的分母部分应该是全概率公式，如果 \(\pmb{\theta}\) 是连续值，这里是要对 \(\pmb{\theta}\) 进行积分的。 之前我们说过，我们会对 \(\pmb{\theta}\) 进行离散化，这样就可以把积分操作转换成求和操作。

3.3.1. 算法总结

以上是简单的推导过程，EM算法的步骤是

  • 首先，第一轮迭代t=0时，随机初始化待求参数 \(\pmb{\omega}^{(t)}\)

  • E步，求解隐藏变量 \(\pmb{\theta}\) 的后验概率分布 \(g(\pmb{\theta} | \mathbf{Y},\pmb{\omega}^{(t)} )\)

    • \(\pmb{\theta}\) 本身是连续值，这时 \(g(\pmb{\theta} | \mathbf{Y},\pmb{\omega}^{(t)} )\) 就是概率密度函数，计算积分比较复杂。

    • 所以可以把 \(\pmb{\theta}\) 离散化，这样 \(g(\pmb{\theta} | \mathbf{Y},\pmb{\omega}^{(t)} )\) 就是概率质量函数，只需要求出其概率分布，然后利用求和的方式计算全概率。

  • M步，极大化Q函数 \(Q(\pmb{\omega}^{(t+1)}|\pmb{\omega}^{(t)})\) 得到新的 \(\pmb{\omega}^{(t+1)}\)

()\[\pmb{\omega}^{(t+1)} = \mathop{\arg\max}_{\pmb{\omega} \in \pmb{\Omega}} \sum_{\pmb{\theta} \in \pmb{\Theta}} g(\pmb{\theta}|\mathbf{Y},\pmb{\omega}^{(t)}) ln L(\mathbf{Y},\pmb{\theta}|\pmb{\omega})\]
  • 重复E步和M步直到满足收敛条件

    • \(\pmb{\omega}\) 不再变化 \(|\pmb{\omega}^{(t+1)} - \pmb{\omega}^{(t)}|<\epsilon\)

    • 对数似然函数不再变化 \(|lnL(\mathbf{Y}|\pmb{\omega}^{(t+1)}) - lnL(\mathbf{Y}|\pmb{\omega}^{(t)})|<\epsilon\)

3.4. 参考内容

　　[1] IRT Parameter Estimation using the EM Algorithm

　　[2] RoutledgeHandbooks-9781315736013-chapter3

　　[3] Optimizing Information Using the Expectation-Maximization Algorithm in Item Response Theory