机器学习入门：最大期望算法-3

机器学习入门：最大期望算法

1、实验描述

• 本实验通过使用EM算法解决部分观测数据的参数估计的难题，利用scikit-learn提供的EM模块，实现EM算法，创建模型，训练模型，然后利用模型计算均值、方差着数据，并最终用3d图示可视化结果。

• 实验时长：45分钟

• 主要步骤：

  • 设定随机种子数

  • 给定均值和方差

  • 生成模拟样本

  • 建立混合高斯模型

  • 利用样本，训练模型

  • 模型评估

  • 结果可视化

2、实验环境

• 虚拟机数量：1
• 系统版本：CentOS 7.5
• scikit-learn版本： 0.19.2
• numpy版本：1.15.1
• matplotlib版本：2.2.3
• python版本：3.5
• IPython版本：6.5.0

3、相关技能

• Python编程
• Sklearn的使用
• Numpy 的使用
• EM算法
• matplotlib 绘图

4、相关知识点

• 混合高斯模型
• 均值，方差，协方差，隐变量
• 准确率
• 多元高斯分布
• 概率密度

5、实现效果

• 使用高斯混合模型完成对样本空间的观察具体效果如下图：

图 1

6、实验步骤

6.1进入Anaconda创建的虚拟环境“ML”

6.1.1在zkpk的家目录下执行如下命令

[zkpk@master ~]$ cd
[zkpk@master ~]$ source activate ML
(ML) [zkpk@ master ML]$ 

6.1.2此时已经进入虚拟环境。键入如下命令，进入ipython交互是编程环境

(ML) [zkpk@ master ML]$ ipython
Python 3.5.4 |Anaconda, Inc.| (default, Nov  3 2017, 20:01:27)
Type 'copyright', 'credits' or 'license' for more information
IPython 6.2.1 -- An enhanced Interactive Python. Type '?' for help.
In [1]:

6.2在Ipython交互式编程环境中开始进行实验

6.2.1导入实验所需的包

In [1]: import numpy as np
   ...: from scipy.stats import multivariate_normal		# 用于生成多元正态分布矩阵
   ...: from sklearn.mixture import GaussianMixture		# 高斯混合模型
   ...: from mpl_toolkits.mplot3d import Axes3D
   ...: import matplotlib as mpl
   ...: import matplotlib.pyplot as plt
   ...: from sklearn.metrics.pairwise import pairwise_distances_argmin
   ...:

6.2.2生成第一部分样本

6.2.2.1设定随机种子数，便于结果复现

In [2]:     np.random.seed(0)  # 设置随机种子数， 便于产生同样的结果

6.2.2.2给定均值和方差

In [3]:     mu1_fact = (0, 0, 0)  # 给定均值
   ...:     cov1_fact = np.diag((1, 2, 3)) # 给定x,y,z轴的方差
   ...:

6.2.2.3生成第一份样本

In [4]:     # 多元高斯分布，给定均值、方差， 生成400个样本；生成多元正态分布矩阵
   ...:     data1 = np.random.multivariate_normal(mu1_fact, cov1_fact, 400)

6.2.3同理生成第二份样本数据

In [5]:     mu2_fact = (2, 2, 1)
   ...:     cov2_fact = np.diag((4., 1., 2.))/3  # 方差
   ...:     # cov2_fact = np.array(((1, 1, 3), (1, 2, 1), (0, 0, 1)))
   ...:     data2 = np.random.multivariate_normal(mu2_fact, cov2_fact, 100)
   ...:     # 根据均值和方差求取100个样本
   ...:

6.2.4两份样本数据叠加

In [6]:     data = np.vstack((data1, data2)) # 将两个样本叠在一起

6.2.5根据样本类别生成对应的标签值（false，true）

In [7]:     y = np.array([True] * 400 + [False] * 100) # 前400个值是True， 后100个是False

6.2.6利用混合高斯模型拟合数据：

6.2.6.1第一个参数，根据假设，指定分组为2

6.2.6.2第二个参数：描述要使用的协方差类别，这里选择full，每一个都保持不同

6.2.6.3第三个参数：迭代阈值

6.2.6.4第四个参数：迭代次数

In [8]:     g = GaussianMixture(n_components=2, covariance_type='full', tol=1e-6, max_iter=1000)
    ...:    g.fit(data)  # fit数据，训练模型。
    ...:

图 2

6.2.7打印出模拟估计的类别概率

In [14]:    print('类别概率:\t', g.weights_[0])

图 3

6.3打印出模型估计的均值和方差，跟我们定义的对比

In [16]:     print('均值:\n', g.means_, '\n')
    ...:     # 计算出的均值，和我们给的均值也比较符合（0，0,0），（2,2，1）
    ...:     print('方差:\n', g.covariances_, '\n')
    ...:     # 估计出的方差也差不多拟合我们给定的数值（1,2,3），（1.3, 0.3, 0.7)

图 4

6.4准确率的计算

6.4.1得到预测的均值和方差

In [17]:     mu1, mu2 = g.means_
    ...:     sigma1, sigma2 = g.covariances_
    ...:

6.4.2利用均值和方差生成对应的两个正态分布

In [18]:     norm1 = multivariate_normal(mu1, sigma1)
    ...:     norm2 = multivariate_normal(mu2, sigma2)

6.4.3基于这两个分布计算在data上的概率密度函数

In [19]:     tau1 = norm1.pdf(data)
    ...:     tau2 = norm2.pdf(data)

6.4.4根据欧式距离，计算属于某个类别，并利用两个概率密度函数的相对值大小，划分属于类别c1,c2.

In [20]:     order = pairwise_distances_argmin([mu1_fact, mu2_fact], [mu1, mu2], me
    ...: tric='euclidean')
    ...:     print(order)
    ...:     # 根据两个概率密度函数的相对大小划分两个类别c1,c2
    ...:     if order[0] == 0:
    ...:         c1 = tau1 > tau2
    ...:     else:
    ...:         c1 = tau1 < tau2
    ...:
    ...:    
    ...:     c2 = ~c1	# 各值取反

6.4.5计算分类准确率

In [21]:     acc = np.mean(y == c1)
    ...:     # 计算准确率
    ...:     print(u'准确率：%.2f%%' % (100 * acc))

图 5

6.5结果可视化

 In [23]:    fig = plt.figure(figsize=(13, 7), facecolor='w')	# 创建figure画布对象
    ...:     ax = fig.add_subplot(121, projection='3d')  # 添加第一个子图, 并以3d形式展示
    ...:     ax.scatter(data[:, 0], data[:, 1], data[:, 2], c='b', s=30, marker='o'
    ...: , depthshade=True)		# marker='o'表示以圆形显示数据坐标点
    ...:     ax.set_xlabel('X')	# 分别设置三个轴的label
    ...:     ax.set_ylabel('Y')
    ...:     ax.set_zlabel('Z')
    ...:     ax.set_title(u'original DATA', fontsize=18)
    ...:
    ...:     ax = fig.add_subplot(122, projection='3d')  # 第二个子图
    ...:
    ...: # 在第二个子图上分别用红色和绿色画出两个类别的分布。和设定相关的图像阴影，提升显示的效果
    ...:     ax.scatter(data[c1, 0], data[c1, 1], data[c1, 2], c='r', s=30, marker=
    ...: 'o', depthshade=True)
    ...:     ax.scatter(data[c2, 0], data[c2, 1], data[c2, 2], c='g', s=30, marker=
    ...: '^', depthshade=True)		# marker='^'表示以顶角向上的三角形显示数据坐标点
    ...:     ax.set_xlabel('X')
    ...:     ax.set_ylabel('Y')
    ...:     ax.set_zlabel('Z')
    ...:     ax.set_title(u'EM Algorithm Classification', fontsize=18)
    ...:     plt.suptitle(u'EM Algorithm ', fontsize=21)	# 添加总标题
    ...:     plt.subplots_adjust(top=0.90)
    ...:     plt.tight_layout()
    ...:     plt.show()
    ...:
    ...:

图 6

7、参考答案

• 代码清单EM.py

图 7

8、总结

本次实验，利用模拟数据创建EM多元高斯模型，训练模型，利用模型做各种数据预测。并使用matplotlib绘制结果。需要掌握及理解的主要概念有：混合（多元）高斯模型、均值、方差、协方差、隐变量、准确率、概率密度等。并且，本实验涉及大量的编码实战，需要大家多做编程练习，并且学会查官方或其它文档的方式解决类或方法或函数不会使用的问题。