EM 算法-对鸢尾花数据进行聚类
公号：码农充电站pro
主页：https://codeshellme.github.iohtml
以前介绍过K 均值算法，它是一种聚类算法。今天介绍EM 算法，它也是聚类算法，但比K 均值算法更加灵活强大。python
EM 的全称为 Expectation Maximization，中文为指望最大化算法，它是一个不断观察和调整的过程。git
1，和面过程
咱们先来看一下和面的过程。github
一般状况下，若是你事先不知道面与水的比例，和面过程多是下面这样：算法
- 先放入一些面和水。
- 将面团揉拌均匀。
- 观察面团的稀稠程度：若是面团比较稀，则加入少量面；若是面团比较稠，则加入少量水。
- 如此往复第2，3步骤，直到面团的稀稠程度达到预期。
这个和面过程，就是一个EM 过程：shell
- 先加入一些面和水，将面团揉拌均匀，并观察面团的稀稠程度。这是E 过程。
- 不断的加入水和面（调整水和面的比例），直到达到预期面团程度。这是M 过程。
2，再看K 均值算法
在介绍K 均值 聚类算法时，展现过一个给二维坐标点进行聚类的例子。dom
咱们再来看一下这个例子，以下图：函数
上图是一个聚类的过程，共有6 个步骤：优化
- 初始时散点（绿色点）的分布。
- 随机选出两个中心点的位置，
红色x和蓝色x。 - 计算全部散点分别到
红色x和蓝色x的距离，距离红色x近的点标红色，距离蓝色x近的点标蓝色。 - 从新计算两个中心点的位置，两个中心点分别移动到新的位置。
- 从新计算全部散点分别到
红色x和蓝色x的距离，距离红色x近的点标红色，距离蓝色x近的点标蓝色。 - 再次计算两个中心点的位置，两个中心点分别移动到新的位置。中心点的位置几乎再也不变化，聚类结束。
通过以上步骤就完成了一个聚类过程。3d
实际上，K 均值算法也是一个EM 过程：
- 肯定当前各种中心点的位置，并计算各个散点到现有的中心点的距离。这是E 过程。
- 将各个散点归属到各个类中，从新计算各个类的中心点，直到各种的中心点再也不改变。这是M 过程。
3，EM 算法
将二维数据点的聚类过程，扩展为通常性的聚类问题，EM 算法是这样一个模型：对于待分类的数据点，EM 算法让计算机经过一个不断迭代的过程，来构建一个分类模型。
EM 算法分为两个过程：
- E 过程：根据现有的模型，计算各个数据输入到模型中的计算结果，这称为指望值计算过程，即 E 过程。
- M 过程：从新计算模型参数，以最大化指望值，这称为最大化过程，即M 过程。
以二维数据点的聚类过程为例，咱们定义：
- 同一类中各个点到该类中心的平均距离为 d；
- 不一样类之间的平均距离为 D。
那么二维数据点聚类的M 过程，就是寻求最大化的D 和 -d。咱们但愿的聚类结果是，同一类的点距离较近，不一样类之间距离较远。
EM 算法不是单个算法，而是一类算法。只要知足EM 这两个过程的算法均可以被称为EM 算法。常见的EM 算法有GMM 高斯混合模型和HMM 隐马尔科夫模型。
4，最大似然估计
高等数学中有一门课叫作《几率论与数理统计》，其中讲到了参数估计。
统计推断是数理统计的重要组成部分，它是指利用来自整体的样本提供的信息，对整体的某些特征进行估计或推断，从而认识总体。
统计推断分为两大类：参数估计和假设检验。
咱们假设，对于某个数据集，其分布函数的基本形式已知，但其中含有一个或多个未知参数。
参数估计就是讨论如何根据来自整体的样本提供的信息对未知参数作出估计。参数估计包括点估计和区间估计。其中，点估计中有两种方法：矩估计法和最大似然估计法。
最大似然估计是一种经过已知结果，估计未知参数的方法。
5，EM 算法原理
EM 算法使用的是最大似然估计的原理，它经过观察样本，来找出样本的模型参数。
下面经过一个投硬币的例子，来看下EM 算法的计算过程。
这个例子来自《Nature》（天然）期刊的论文《What is the expectation maximization algorithm？》（什么是指望最大化算法？）。
假定有两枚不一样的硬币 A 和 B，它们的重量分布 θA 和 θB 是未知的，则能够经过抛掷硬币，计算正反面各自出现的次数来估计θA 和 θB。
方法是在每一轮中随机抽出一枚硬币抛掷 10 次，一样的过程执行 5 轮，根据这 50 次投币的结果来计算 θA 和 θB 的最大似然估计。
投掷硬币的过程，记录以下：
第1 到5 次分别投掷的硬币是 B，A，A，B，A。H 表明正面，T 表明负面。将上图转化为表格，以下：
| 次数 | 硬币 | 正面数 | 负面数 |
|---|---|---|---|
| 1 | B | 5 | 5 |
| 2 | A | 9 | 1 |
| 3 | A | 8 | 2 |
| 4 | B | 4 | 6 |
| 5 | A | 7 | 3 |
经过这个表格，能够直接计算 θA 和 θB，以下：
显然，若是知道每次投掷的硬币是A 仍是B，那么计算θA 和 θB 是很是简单的。
可是，若是不知道每次投掷的硬币是A 仍是B，该如何计算θA 和 θB 呢？
此时咱们将上面表格中的硬币一列隐藏起来，这时硬币就是隐变量。因此咱们只知道以下数据：
| 次数 | 正面数 | 负面数 |
|---|---|---|
| 1 | 5 | 5 |
| 2 | 9 | 1 |
| 3 | 8 | 2 |
| 4 | 4 | 6 |
| 5 | 7 | 3 |
这时想要计算 θA 和 θB，就要用最大似然估计的原理。
计算过程以下图：
第一步
先为 θA 和 θB 设定一个初始值，好比 θA = 0.6，θB = 0.5。
第二步
咱们知道每一轮投币的正 / 负面的次数：
- 第1轮：5 正 5 负，计算出现这种结果的几率：
- 若是是A 硬币，那么P(H5T5|A) = 0.6^5 * 0.4^5
- 若是是B 硬币，那么P(H5T5|B) = 0.5^5 * 0.5^5
- 将 P(H5T5|A) 和 P(H5T5|B) 归一化处理，可得：
- P(H5T5|A) = 0.45，P(H5T5|B) = 0.55
- 第2轮：9 正 1 负，计算出现这种结果的几率：
- 若是是A 硬币，那么P(H9T1|A) = 0.6^9 * 0.4^1
- 若是是B 硬币，那么P(H9T1|B) = 0.5^9 * 0.5^1
- 将 P(H9T1|A) 和 P(H9T1|B) 归一化处理，可得：
- P(H9T1|A) = 0.8，P(H9T1|B) = 0.2
- 第3轮：8 正 2 负，计算出现这种结果的几率：
- 若是是A 硬币，那么P(H8T2|A) = 0.6^8 * 0.4^2
- 若是是B 硬币，那么P(H8T2|B) = 0.5^8 * 0.5^2
- 将 P(H8T2|A) 和 P(H8T2|B) 归一化处理，可得：
- P(H8T2|A) = 0.73，P(H8T2|B) = 0.27
- 第4轮：4 正 6 负，计算出现这种结果的几率：
- 若是是A 硬币，那么P(H4T6|A) = 0.6^4 * 0.4^6
- 若是是B 硬币，那么P(H4T6|B) = 0.5^4 * 0.5^6
- 将 P(H4T6|A) 和 P(H4T6|B) 归一化处理，可得：
- P(H4T6|A) = 0.35，P(H4T6|B) = 0.65
- 第5轮：7 正 3 负，计算出现这种结果的几率：
- 若是是A 硬币，那么P(H7T3|A) = 0.6^7 * 0.4^3
- 若是是B 硬币，那么P(H7T3|B) = 0.5^7 * 0.5^3
- 将 P(H7T3|A) 和 P(H7T3|B) 归一化处理，可得：
- P(H7T3|A) = 0.65，P(H7T3|B) = 0.35
而后，根据每一轮的 P(HmTn|A) 和 P(HmTn|B)，能够计算出每一轮的正 / 负面次数。
m 为正面次数，n 为负面次数。
对于硬币A，结果以下：
| 轮数 | P(HmTn|A) | m | n | 正面数 | 负面数 |
|---|---|---|---|---|---|
| 1 | 0.45 | 5 | 5 | 0.45*5=2.2 | 0.45*5=2.2 |
| 2 | 0.8 | 9 | 1 | 0.8*9=7.2 | 0.8*1=0.8 |
| 3 | 0.73 | 8 | 2 | 0.73*8=5.9 | 0.73*2=1.5 |
| 4 | 0.35 | 4 | 6 | 0.35*4=1.4 | 0.35*6=2.1 |
| 5 | 0.65 | 7 | 3 | 0.65*7=4.5 | 0.65*3=1.9 |
| 总计 | - | - | - | 21.3 | 8.6 |
对于硬币B，结果以下：
| 轮数 | P(HmTn|B) | m | n | 正面数 | 负面数 |
|---|---|---|---|---|---|
| 1 | 0.55 | 5 | 5 | 0.55*5=2.8 | 0.55*5=2.8 |
| 2 | 0.2 | 9 | 1 | 0.2*9=1.8 | 0.2*1=0.2 |
| 3 | 0.27 | 8 | 2 | 0.27*8=2.1 | 0.27*2=0.5 |
| 4 | 0.65 | 4 | 6 | 0.65*4=2.6 | 0.65*6=3.9 |
| 5 | 0.35 | 7 | 3 | 0.35*7=2.5 | 0.35*3=1.1 |
| 总计 | - | - | - | 11.7 | 8.4 |
第三步
根据上面两个表格，能够得出（第1次迭代的结果） θA 和 θB：
根据这个估计值，再次回到第一步去计算。
如此往复第1、2、三步，通过10次迭代以后，θA 和 θB 的估计值为：
最终，θA 和 θB 将收敛到一个几乎不变的值，此时迭代结束。这样咱们就求解出了θA 和 θB 的最大似然估计值。
咱们将上述过程当中，第一步称为初始化参数，第二步称为观察预期，第三步称为从新估计参数。
第1、二步为E 过程，第三步为M 过程，这就是EM 算法的过程。
若是咱们有一个待聚类的数据集，咱们把潜在的类别当作隐变量，样本当作观察值，这样就能够把聚类问题转化成参数估计问题。这就是EM 聚类的原理。
6，硬聚类与软聚类
与 K 均值算法相比，K 均值算法是经过距离来区分样本之间的差异，且每一个样本在计算的时候只能属于一个分类，咱们称之为硬聚类算法。
而 EM 聚类在求解的过程当中，实际上每一个样本都有必定的几率和每一个聚类相关，这叫作软聚类算法。
7，EM 聚类的缺点
EM 聚类算法存在两个比较明显的问题。
第一个问题是，EM 算法计算复杂，收敛较慢，不太适合大规模数据集和高维数据。
第二个问题是，EM 算法不必定能给出全局最优解：
- 当优化的目标函数是凸函数时，必定能够获得全局最优解。
- 当优化的目标函数不是凸函数时，可能会获得局部最优解，而非全局最优解。
8，GMM 高斯混合模型
上文中介绍过，常见的EM 算法有GMM 高斯混合模型和HMM 隐马尔科夫模型。这里主要介绍GMM 高斯混合模型的实现。
sklearn 库的mixture 模块中的GaussianMixture 类是GMM 算法的实现。
先来看下 GaussianMixture 类的原型：
GaussianMixture( n_components=1, covariance_type='full', tol=0.001, reg_covar=1e-06, max_iter=100, n_init=1, init_params='kmeans', weights_init=None, means_init=None, precisions_init=None, random_state=None, warm_start=False, verbose=0, verbose_interval=10)
这里介绍几个重要的参数：
- n_components：表明高斯混合模型的个数，也就是咱们要聚类的个数，默认值为 1。
- covariance_type：表明协方差类型。一个高斯混合模型的分布是由均值向量和协方差矩阵决定的，因此协方差的类型也表明了不一样的高斯混合模型的特征。协方差类型有 4 种取值：
- full，表明彻底协方差，也就是元素都不为 0；
- tied，表明相同的彻底协方差；
- diag，表明对角协方差，也就是对角不为 0，其他为 0；
- spherical，表明球面协方差，非对角为 0，对角彻底相同，呈现球面的特性。
- max_iter：表明最大迭代次数，默认值为 100。
9，对鸢尾花数据集聚类
在《决策树算法-实战篇-鸢尾花及波士顿房价预测》中咱们介绍过鸢尾花数据集。这里咱们使用GMM 算法对该数据进行聚类。
首先加载数据集：
from sklearn.datasets import load_iris iris = load_iris() # 加载数据集 features = iris.data # 获取特征集 labels = iris.target # 获取目标集
在聚类算法中，只须要特征数据 features，而不须要目标数据labels，但可使用 labels 对聚类的结果作验证。
构造GMM聚类：
from sklearn.mixture import GaussianMixture # 原数据中有 3 个分类，因此这里咱们将 n_components 设置为 3 gmm = GaussianMixture(n_components=3, covariance_type='full')
对数据集进行聚类：
prediction_labels = gmm.fit_predict(features)
查看原始分类：
>>> print(labels) [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2]
查看聚类结果：
>>> print(prediction_labels) [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 2 1 2 1 1 1 1 2 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2]
对比原始分类和聚类结果，聚类结果中只有个别数据分类错误，我用红圈标了出来：
10，评估聚类结果
咱们可使用 Calinski-Harabaz 指标对聚类结果进行评估。
sklearn 库实现了该指标的计算，即 calinski_harabasz_score 方法，该方法会计算出一个分值，分数越高，表明聚类效果越好，也就是相同类中的差别性小，不一样类之间的差别性大。
下面对鸢尾花数据集的聚类结果进行评估，传入特征数据和聚类结果：
>>> from sklearn.metrics import calinski_harabasz_score >>> calinski_harabasz_score(features, prediction_labels) 481.78070899745234
咱们也能够传入特征数据和原始结果：
>>> calinski_harabasz_score(features, labels) 487.33087637489984
能够看到，对于原始结果计算出的分值是487.33，对于预测结果计算出的分值是481.78，相差并很少，说明预测结果仍是不错。
通常状况下，一个须要聚类的数据集并无目标数据，因此只能对预测结果进行评分。咱们须要人工对聚类的含义结果进行分析。
11，总结
本篇文章主要介绍了以下内容：
- EM 算法的过程及原理，介绍了一个投掷硬币的例子。
- 硬聚类与软聚类的区别。
- EM 聚类的缺点：
- 计算复杂度较大。
- 有可能得不到全局最优解。
- 使用GMM 算法对鸢尾花数据进行聚类。
- 使用 Calinski-Harabaz 指标对聚类结果进行评估。
（本节完。）
推荐阅读：
欢迎关注做者公众号，获取更多技术干货。
- 1. KNN分类实战-----对鸢尾花数据进行分类(sklearn)
- 2. 运用scikit-learn库对鸢尾花数据进行分类
- 3. 朴素贝叶斯对鸢尾花数据集进行分类
- 4. 利用KNN对鸢尾花数据进行分类
- 5. EM算法处理鸢尾花数据实战
- 6. [Go]使用Golang对鸢尾花数据集进行k-means聚类
- 7. 【K-Means】鸢尾花的分类聚类
- 8. knn算法--鸢尾花
- 9. 鸢尾花分类
- 10. 实现鸢尾花数据集分类
- 更多相关文章...
- • C# 数据类型 - C#教程
- • SQLite 数据类型 - SQLite教程
- • Flink 数据传输及反压详解
- • TiDB 在摩拜单车在线数据业务的应用和实践
-
每一个你不满意的现在，都有一个你没有努力的曾经。