专栏首页用户2133719的专栏CS229 课程笔记之九:EM 算法与聚类

CS229 课程笔记之九:EM 算法与聚类

1 K-means聚类

聚类问题是一种「无监督学习」,给定训练集

\{x^{(1)}, \ldots, x^{(m)}\}

,我们希望将其聚合成几个特定的类。k-means 聚类算法的流程如下:

  1. 随机初始化「聚类中心」
\mu_1,\mu_2,\dots, \mu_k \in \mathbb{R}^n
  1. 重复以下步骤直至收敛: 对于每个
i

(训练集大小),令

c^{(i)} := \arg \min_j ||x^{(i)}-\mu_j||^2

对于每个

j

(聚类数量),令

\mu^{(j)} := \frac {\sum _{i=1}^m 1\{c^{(i)}=j\}x^{(i)}} {\sum _{i=1}^m 1\{c^{(i)}=j\}}

该算法的思想为:先将每个训练样本

x^{(i)}

分配到距离其最近的中心

\mu_j

,再将每个聚类中心移动到第一步中分配到该中心的样本的均值。

下图可视化了 k-means 算法的运行流程:

为了证明 k-means 算法能否保证收敛,我们定义「失真函数」(distortion function)为:

J(c,\mu) = \sum_{i=1}^m ||x^{(i)}-\mu_{c^{(i)}}||^2

可以发现 k-means 本质上就是对失真函数进行坐标上升法优化:其内层循环首先保持

\mu

不变关于

c

最小化

J

,然后保持

c

不变关于

\mu

最小化

J

。因此,

J

一定会持续下降,最终达到收敛。一般

c

\mu

也会收敛,但理论上存在同时出现多种聚类组合的可能性,使得失真函数的值一样。

失真函数是一个非凸函数,这意味着坐标上升并不能保证其收敛至全局最优,存在收敛到局部最优的可能性。一般情况下这不会发生,可以通过多次运行 k-means 算法,选择最优解来解决这个问题。

2 混合高斯分布

混合高斯分布可以用于软聚类问题,即输出一个样本属于各个类的概率。我们可以将数据的类别看作一个隐含随机变量

z

,并给出如下假设:

z

服从多项式分布

z^{(i)} \sim \text{Multinomial}(\phi)
  1. 给定不同的
z

x

服从不同的高斯分布

x^{(i)} | z^{(i)} = j \sim \mathcal{N}(\mu_i, \Sigma_j)

使用极大似然法求解该优化问题,可以得到如下的似然函数:

\begin{aligned} \ell(\phi,\mu,\Sigma) &= \sum_{i=1}^m \log p(x^{(i)};\phi,\mu,\Sigma) \\ &= \sum_{i=1}^m \log \sum_{z^{(i)} = 1}^k p(x^{(i)}| z^{(i)};\mu, \Sigma)p(z^{(i)};\phi) \end{aligned}

该问题无法求出闭合解。如果

z^{(i)}

已知,即我们知道每个样本来自于哪个高斯分布,那么极大似然估计的求解是容易的,似然函数如下:

\ell(\phi,\mu,\Sigma) = \sum_{i=1}^m \left(\log p(x^{(i)}| z^{(i)};\mu, \Sigma)+ \log p(z^{(i)};\phi)\right)

求解的结果是:

\begin{aligned} \phi_j &= \frac 1 m \sum_{i=1}^m 1\{z^{(i)} = j\} \\ \mu_j &= \frac{\sum_{i=1}^m 1\{z^{(i)}=j\}x^{(i)}}{\sum_{i=1}^m 1\{z^{(i)}=j\}} \\ \Sigma_j &= \frac{\sum_{i=1}^m 1\{z^{(i)}=j\}(x^{(i)}-\mu_{j})(x^{(i)}-\mu_{j})^T}{\sum_{i=1}^m 1\{z^{(i)}=j\}} \end{aligned}

该结果与之前的高斯判别分析的结论类似(GDA 的协方差矩阵必须相同)。

3 EM 算法初步

实际上,我们并不知道

z^{(i)}

的值。我们可以通过 EM 算法进行迭代估计出

z^{i}

从而得到参数,其基本思想如下:

重复以下步骤直至收敛:

  1. 「E-step」对于每个
i,j

,令:

w_j^{(i)} := p(z^{(i)} = j|x^{(i)};\phi,\mu,\Sigma)
  1. 「M-step」更新参数:
\begin{aligned} \phi_j &= \frac 1 m \sum_{i=1}^m w_j^{(i)} \\ \mu_j &= \frac{\sum_{i=1}^m w_j^{(i)}x^{(i)}}{\sum_{i=1}^m w_j^{(i)}} \\ \Sigma_j &= \frac{\sum_{i=1}^m w_j^{(i)}(x^{(i)}-\mu_{j})(x^{(i)}-\mu_{j})^T}{\sum_{i=1}^m w_j^{(i)}} \end{aligned}

关于 EM 算法的几点解释:

「E-step」 中,给定

x^{(i)}

,我们使用当前的参数值来计算

z^{(i)}

的后验概率,即

p(z^{(i)} = j|x^{(i)};\phi,\mu,\Sigma) = \frac {p(x^{(i)}|z^{(i)} = j;\mu,\Sigma)p(z^{(i)}=j;\phi)} {\sum_{l=1}^k p(x^{(i)}|z^{(i)}=l;\mu,\Sigma)p(z^{(i)}=l;\phi)}

该概率代表我们对

z^{(i)}

值的软猜测(即以概率值代替具体的值)。

「M-step」 中,参数的更新公式与之前已知

z^{(i)}

的公式相比,只是把指示函数替换为了概率。

与 K-means 算法相比,EM 算法输出的是样本属于各个类的概率,这是一种软聚类。与 K-means 相似,EM 算法容易陷入局部最优,因此多次尝试不同的初始参数可能是一个好主意。下两节将给出 EM 算法的一般形式,并证明其收敛性。

4 Jensen 不等式

本节将介绍 Jensen 不等式,其在 EM 算法中有着重要的作用。

4.1 函数的凹凸性

在介绍 Jensen 不等式之前,先简单回顾一下函数的凹凸性。对于一个实数域的函数

f

,其为凸函数的条件为

f''(x) \ge 0

。如果输入为向量形式,则该条件可推广为其 Hessian 矩阵半正定(

H \ge 0

)。

如果

f''(x) > 0

H > 0

),则函数「严格」凸。凹函数的判定条件与凸函数完全相反。

4.2 定理

f

是一个凸函数,

X

是一个随机变量,则:

\text{E}[f(X)] \ge f(\text{E}X)

如果

f

严格凸,那么当且仅当

X = \text{E}[X]

时等号成立(即

X

为常量)。可以通过下图对该不等式有一个直观的理解:

f

为凹函数时,不等式方向对调,仍然成立。

5 EM 算法

5.1 算法的导出

假定我们有一个包含 m 个独立样本的训练集,我们希望去拟合一个概率模型

p(x,z)

,其对数似然函数为:

\begin{aligned} \ell(\theta) &= \sum_{i=1}^m \log p(x;\theta) \\ &= \sum_{i=1}^m \log \sum_z p(x,z;\theta) \end{aligned}

这里假定

z

「离散」变量(连续变量需要使用积分)。

直接最大化

\ell(\theta)

是难以求解的。EM 算法的思想是先构建一个

\ell

的下界(E-step),然后去优化这个下界(M-step),达到间接最大化

\ell

的目的。

对于每个

i

,设

Q_i

z

的某种概率分布(

\sum_z Q_i(z) = 1, Q_i(z) \ge 0

)。根据期望的定义以及 Jensen 不等式,我们有:

\begin{align*} \sum_i \log p(x^{(i)};\theta) &= \sum_i \log \sum_{z^{(i)}} p(x^{(i)},z^{(i)};\theta) \tag{1} \\ &= \sum_i \log \sum_{z^{(i)}} Q_i(z^{(i)}) \frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})} \tag{2} \\ &\ge \sum_i \sum_{z^{(i)}} Q_i(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})} \tag{3} \\ \end{align*}
\frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})}

可以看做

z^{(i)}

的随机变量,其概率分布为

Q_i

,期望可以通过

\sum_{z^{(i)}} Q_i(z^{(i)}) \frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})}

得到。

f(x) = \log x

是一个凹函数,应用 Jensen 不等式时注意方向对调。

为了执行 EM 算法,我们需要选择合适的

Q_i

以保证在当前的参数设置下取到下界,即目前的

\theta

值能够使得 (3) 式的等号成立。之前的定理表明等号成立的条件是随机变量为「常数」

\frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})} = c

因此只要

Q_i(z^{(i)})

p(x^{(i)},z^{(i)};\theta)

成比例即可:

Q_i(z^{(i)}) \propto p(x^{(i)},z^{(i)};\theta)

实际上,因为

\sum_z Q_i(z^{(i)}) = 1

,将其代入上述公式,可以得到:

\sum_z Q_i(z^{(i)}) = \frac {\sum_z p(x^{(i)},z^{(i)};\theta)} c = 1

因此

c = \sum_z p(x^{(i)},z^{(i)};\theta)

,从而有:

\begin{aligned} Q_i(z^{(i)}) &= \frac {p(x^{(i)},z^{(i)};\theta)} {\sum_z p(x^{(i)},z^{(i)};\theta)} \\ &= \frac {p(x^{(i)},z^{(i)};\theta)} {p(x^{(i)};\theta)} \\ &= p(z^{(i)}|x^{(i)};\theta) \end{aligned}

根据上述推导,我们只需要将

Q_i

设置为给定

x^{(i)}

z^{(i)}

的后验分布即可(以

\theta

为参数)。

综上所述,EM 算法的具体步骤为:

  • 「E-step」:对于每个
i

,令

Q_i(z^{(i)}) := p(z^{(i)}|x^{(i)};\theta)
  • 「M-step」:更新参数
\theta := \arg \max_\theta \sum_i \sum_{z^{(i)}} Q_i(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})}

重复以上两个步骤直至收敛。

5.2 收敛性证明

下面证明该算法的收敛性。假定

\theta^{(t)}

\theta^{(t+1)}

来自于 EM 算法的两次成功的迭代,那么通过证明

\ell(\theta^{(t)}) \le \ell(\theta^{(t+1)})

,就可以表明 EM 算法是单调收敛的(对数似然函数单调递增)。

当参数为

\theta^{(t)}

时,根据算法步骤,我们令

Q_i^{(t)}(z^{(i)}) := p(z^{(i)}|x^{(i)};\theta^{(t)})

,这一选择保证了等号成立,即:

\ell(\theta^{(t)}) = \sum_i \sum_{z^{(i)}} Q_i^{(t)}(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\theta^{(t)})}{Q_i^{(t)}(z^{(i)})}

而参数

\theta^{(t+1)}

是通过最大化上式的右边部分得出的(更新

\theta

),因此有:

\begin{align*} \ell(\theta^{(t+1)}) &\ge \sum_i \sum_{z^{(i)}} Q_i^{(t)}(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\theta^{(t+1)})}{Q_i^{(t)}(z^{(i)})} \tag{4} \\ &\ge \sum_i \sum_{z^{(i)}} Q_i^{(t)}(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\theta^{(t)})}{Q_i^{(t)}(z^{(i)})} \tag{5}\\ &= \ell(\theta^{(t)}) \tag{6} \end{align*}

(4) 式的得出来源于 (3) 式;(5) 式的得出来源于

\theta^{(t+1)}

是上一步的下界最大化的结果;(6) 式的得出之前已经证明(满足等号成立的条件)。

综上所述,EM 算法可以保证似然函数的「单调收敛」。一个可取的收敛条件是

\ell(\theta)

的增长速度小于某个临界值。

5.3 坐标上升法

如果我们定义:

J(Q,\theta) = \sum_i \sum_{z^{(i)}} Q_i(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\theta)}{Q_i(z^{(i)})}

从之前的推导可以知道

\ell(\theta) \ge J(Q,\theta)

那么 EM 算法可以看做是对

J

的坐标上升法:

  • 在 E-step 中,关于
Q

最大化

J

(使等号成立)

  • 在 M-step 中,关于
\theta

最大化

J

6 混合高斯模型复盘

下面将使用 EM 算法的一般形式来对之前混合高斯模型中的公式进行推导,由于篇幅所限,这里只给出

\phi

\mu_j

的推导过程。

之前我们得出的参数更新公式如下:

\begin{aligned} \phi_j &= \frac 1 m \sum_{i=1}^m w_j^{(i)} \\ \mu_j &= \frac{\sum_{i=1}^m w_j^{(i)}x^{(i)}}{\sum_{i=1}^m w_j^{(i)}} \\ \Sigma_j &= \frac{\sum_{i=1}^m w_j^{(i)}(x^{(i)}-\mu_{j})(x^{(i)}-\mu_{j})^T}{\sum_{i=1}^m w_j^{(i)}} \end{aligned}

根据 E-step 的定义,我们可以得到:

w_j^{(i)} = Q_i(z^{(i)} = j) = P(z^{(i)} = j|x^{(i)};\phi,\mu,\Sigma)

在 M-step 中,我们需要通过上述三个参数去最大化下式:

\begin{aligned} \sum_{i=1}^m &\sum_{z^{(i)}} Q_i(z^{(i)}) \log \frac {p(x^{(i)},z^{(i)};\phi,\mu,\Sigma)}{Q_i(z^{(i)})} \\ &= \sum_{i=1}^m\sum_{j=1}^k Q_i(z^{(i)}=j) \log \frac {p(x^{(i)}|z^{(i)}=j;\mu,\Sigma)p(z^{(i)}=j;\phi)}{Q_i(z^{(i)}=j)} \\ &= \sum_{i=1}^m\sum_{j=1}^k w_j^{(i)} \log \frac {\frac 1 {(2\pi)^{n/2}|\Sigma_j|^{1/2}}\exp \left(-\frac 1 2 (x^{(i)}-\mu_j)^T\Sigma_j^{-1}(x^{(i)}-\mu_j)\right)\cdot\phi_j}{w_j^{(i)}} \end{aligned}

我们首先关于

\mu_j

去进行最大化,求导可得:

\begin{aligned} \nabla_{\mu_j}&\sum_{i=1}^m\sum_{j=1}^k w_j^{(i)} \log \frac {\frac 1 {(2\pi)^{n/2}|\sum_j|^{1/2}}\exp \left(-\frac 1 2 (x^{(i)}-\mu_j)^T\Sigma_j^{-1}(x^{(i)}-\mu_j)\right)\cdot\phi_j}{w_j^{(i)}} \\ &= - \nabla_{\mu_j} \sum_{i=1}^m\sum_{j=1}^k w_j^{(i)}\frac 1 2 (x^{(i)}-\mu_j)^T\Sigma_j^{-1}(x^{(i)}-\mu_j) \\ &= \frac 1 2 \sum_{i=1}^m w_j^{(i)} \nabla_{\mu_j} 2\mu_j^T \Sigma_j^{-1}x^{(i)}-\mu_j^T\Sigma_j^{-1}\mu_j \\ &= \sum_{i=1}^m w_j^{(i)} (\Sigma_j^{-1}x^{(i)}-\Sigma_j^{-1}\mu_j) \end{aligned}

上述推导首先去除了不相关的项,然后去除了求和符号(只求当前的

j

的导数),最后合并同类项得到结果。将上式设为 0 求解

\mu_j

,可得:

\mu_j = \frac{\sum_{i=1}^m w_j^{(i)}x^{(i)}}{\sum_{i=1}^m w_j^{(i)}}

下面关于

\phi_j

去进行最大化,将与

\phi_j

相关的部分提取出来,可以得到需要最大化的项为:

\sum_{i=1}^m\sum_{j=1}^k w_j^{(i)}\log \phi_j

这里我们还有一个额外的约束条件:

\sum_{j=1}^k\phi_j = 1

因此,我们需要构建拉格朗日算子:

\mathcal{L}(\phi) = \sum_{i=1}^m\sum_{j=1}^k w_j^{(i)}\log \phi_j + \beta (\sum_{j=1}^k \phi -1)

其中

\beta

是拉格朗日乘数,对上式求导可得:

\frac {\partial}{\partial\phi_j} \mathcal{L}(\phi) = \sum_{i=1}^m \frac{w_j^{(i)}}{\phi_j} + \beta

将其设为 0,可得:

\phi_j = \frac {\sum_{i=1}^m w_j^{(i)}} {-\beta}

由于

\sum_j \phi_j = 1

,因此

-\beta = \sum_{i=1}^m \sum_{j=1}^k w_j^{(i)} = \sum_{i=1}^m 1 = m

,所以:

\phi_j = \frac 1 m \sum_{i=1}^m w_j^{(i)}

7 思维导图

本文分享自微信公众号 - 口仆(roito33),作者:口仆

原文出处及转载信息见文内详细说明,如有侵权,请联系 yunjia_community@tencent.com 删除。

原始发表时间:2020-06-02

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • CS229 课程笔记之十五:强化学习与控制

    本章将开始介绍「强化学习」与适应性控制。在监督学习中,对于训练集我们均有明确的标签,算法只需要模仿训练集中的标签来给出预测即可。但对于某些情况,例如序列性的决策...

    口仆
  • CS229 课程笔记:机器学习的实用建议

    对于一个学习算法,有着各种各样的调试手段,不同的调试手段可以解决不同的问题,需要根据实际情况进行选择。学习算法的问题大致可以分为两类:「高偏差/方差」问题以及「...

    口仆
  • Playing Atari with Deep Reinforcement Learning

    本文是对 DQN 原始论文 Playing Atari with Deep Reinforcement Learning 的详细解读。

    口仆
  • 数据百问系列:数据开发需要了解机器学习算法吗?

    对于这个问题,有些群友认为是需要的,也有些群友认为是不需要的,本文根据大家的观点及作者的一些认知,对这个话题进行一个总结。

    木东居士
  • 建议收藏!GitHub标星近10万,用Python实现所有算法合集

    这个项目的算法也是按照字典 A-Z 分类排列的,比如第一个大类就是 Arithmetic Analysis,这个大类里面包括了常见的对分法、高斯消元、交叉法、牛...

    大数据文摘
  • 一文了解强化学习

    虽然是周末,也保持充电,今天来看看强化学习,不过不是要用它来玩游戏,而是觉得它在制造业,库存,电商,广告,推荐,金融,医疗等与我们生活息息相关的领域也有很好的应...

    杨熹
  • 人工智能进行连续决策的关键——强化学习入门指南

    用户1737318
  • 如何得到稳定可靠的强化学习算法?微软两篇顶会论文带来安全的平滑演进

    AI 科技评论按:强化学习最常见的应用是学习如何做出一系列决策,比如,如何一步步攀登上三千英尺高的岩壁。有机会用到强化学习并做出高水准结果的领域包括机器人(以及...

    AI科技评论
  • 学界 | 如何得到稳定可靠的强化学习算法?微软两篇顶会论文带来安全的平滑演进

    AI 科技评论按:强化学习最常见的应用是学习如何做出一系列决策,比如,如何一步步攀登上三千英尺高的岩壁。有机会用到强化学习并做出高水准结果的领域包括机器人(以及...

    AI研习社
  • 强化学习是如何解决问题的?

    什么是强化学习算法呢?要回答这个问题,必须先回答强化学习可以解决什么问题,强化学习如何解决这些问题。

    博文视点Broadview

扫码关注云+社区

领取腾讯云代金券