python实现em聚类算法_EM算法的python实现的方法步骤
前言:前一篇文章大概说了EM算法的整个理解以及一些相关的公式神马的,那些数学公式啥的看完真的是忘完了,那就来用代码记忆记忆吧!接下来将会对python版本的EM算法进行一些分析。
EM的python实现和解析
引入问题(双硬币问题)
假设有两枚硬币A、B,以相同的概率随机选择一个硬币,进行如下的抛硬币实验:共做5次实验,每次实验独立的抛十次,结果如图中a所示,例如某次实验产生了H、T、T、T、H、H、T、H、T、H,H代表正面朝上。
假设试验数据记录员可能是实习生,业务不一定熟悉,造成a和b两种情况
a表示实习生记录了详细的试验数据,我们可以观测到试验数据中每次选择的是A还是B
b表示实习生忘了记录每次试验选择的是A还是B,我们无法观测实验数据中选择的硬币是哪个
问在两种情况下分别如何估计两个硬币正面出现的概率?
以上的针对于b实习生的问题其实和三硬币问题类似,只是这里把三硬币中第一个抛硬币的选择换成了实习生的选择。
对于已知是A硬币还是B硬币抛出的结果的时候,可以直接采用概率的求法来进行求解。对于含有隐变量的情况,也就是不知道到底是A硬币抛出的结果还是B硬币抛出的结果的时候,就需要采用EM算法进行求解了。如下图:
其中的EM算法的第一步就是初始化的过程,然后根据这个参数得出应该产生的结果。
构建观测数据集
针对这个问题,首先采集数据,用1表示H(正面),0表示T(反面):
#硬币投掷结果
observations = numpy.array([[1,0,0,0,1,1,0,1,0,1],
[1,1,1,1,0,1,1,1,0,1],
[1,0,1,1,1,1,1,0,1,1],
[1,0,1,0,0,0,1,1,0,0],
[0,1,1,1,0,1,1,1,0,1]])
第一步:参数的初始化
参数赋初值
第一个迭代的E步
抛硬币是一个二项分布,可以用scipy中的binom来计算。对于第一行数据,正反面各有5次,所以:
#二项分布求解公式
contribution_A = scipy.stats.binom.pmf(num_heads,len_observation,theta_A)
contribution_B = scipy.stats.binom.pmf(num_heads,len_observation,theta_B)
将两个概率正规化,得到数据来自硬币A,B的概率:
weight_A = contribution_A / (contribution_A + contribution_B)
weight_B = contribution_B / (contribution_A + contribution_B)
这个值类似于三硬币模型中的μ,只不过多了一个下标,代表是第几行数据(数据集由5行构成)。同理,可以算出剩下的4行数据的μ。
有了μ,就可以估计数据中AB分别产生正反面的次数了。μ代表数据来自硬币A的概率的估计,将它乘上正面的总数,得到正面来自硬币A的总数,同理有反面,同理有B的正反面。
#更新在当前参数下A,B硬币产生的正反面次数
counts['A']['H'] += weight_A * num_heads
counts['A']['T'] += weight_A * num_tails
counts['B']['H'] += weight_B * num_heads
counts['B']['T'] += weight_B * num_tails
第一个迭代的M步
当前模型参数下,AB分别产生正反面的次数估计出来了,就可以计算新的模型参数了:
new_theta_A = counts['A']['H']/(counts['A']['H'] + counts['A']['T'])
new_theta_B = counts['B']['H']/(counts['B']['H'] + counts['B']['T'])
于是就可以整理一下,给出EM算法单个迭代的代码:
def em_single(priors,observations):
"""
EM算法的单次迭代
Arguments
------------
priors:[theta_A,theta_B]
observation:[m X n matrix]
Returns
---------------
new_priors:[new_theta_A,new_theta_B]
:param priors:
:param observations:
:return:
"""
counts = {'A': {'H': 0, 'T': 0}, 'B': {'H': 0, 'T': 0}}
theta_A = priors[0]
theta_B = priors[1]
#E step
for observation in observations:
len_observation = len(observation)
num_heads = observation.sum()
num_tails = len_observation-num_heads
#二项分布求解公式
contribution_A = scipy.stats.binom.pmf(num_heads,len_observation,theta_A)
contribution_B = scipy.stats.binom.pmf(num_heads,len_observation,theta_B)
weight_A = contribution_A / (contribution_A + contribution_B)
weight_B = contribution_B / (contribution_A + contribution_B)
#更新在当前参数下A,B硬币产生的正反面次数
counts['A']['H'] += weight_A * num_heads
counts['A']['T'] += weight_A * num_tails
counts['B']['H'] += weight_B * num_heads
counts['B']['T'] += weight_B * num_tails
# M step
new_theta_A = counts['A']['H'] / (counts['A']['H'] + counts['A']['T'])
new_theta_B = counts['B']['H'] / (counts['B']['H'] + counts['B']['T'])
return [new_theta_A,new_theta_B]
EM算法主循环
给定循环的两个终止条件:模型参数变化小于阈值;循环达到最大次数,就可以写出EM算法的主循环了
def em(observations,prior,tol = 1e-6,iterations=10000):
"""
EM算法
:param observations :观测数据
:param prior:模型初值
:param tol:迭代结束阈值
:param iterations:最大迭代次数
:return:局部最优的模型参数
"""
iteration = 0;
while iteration < iterations:
new_prior = em_single(prior,observations)
delta_change = numpy.abs(prior[0]-new_prior[0])
if delta_change < tol:
break
else:
prior = new_prior
iteration +=1
return [new_prior,iteration]
调用
给定数据集和初值,就可以调用EM算法了:
print em(observations,[0.6,0.5])
得到
[[0.72225028549925996, 0.55543808993848298], 36]
我们可以改变初值,试验初值对EM算法的影响。
print em(observations,[0.5,0.6])
结果:
[[0.55543727869042425, 0.72225099139214621], 37]
看来EM算法还是很健壮的。如果把初值设为相等会怎样?
print em(observations,[0.3,0.3])
输出:[[0.64000000000000001, 0.64000000000000001], 1]
显然,两个值相加不为1的时候就会破坏这个EM函数。
换一下初值:
print em(observations,[0.99999,0.00001])
输出:[[0.72225606292866507, 0.55543145006184214], 33]
EM算法对于参数的改变还是有一定的健壮性的。
以上是根据前人写的博客进行学习的~可以自己动手实现以下,对于python练习还是有作用的。希望对大家的学习有所帮助,也希望大家多多支持。
python实现em聚类算法_EM算法的python实现的方法步骤相关推荐
- 怎样检查python环境是否安装好_如何搭建pytorch环境的方法步骤
1.conda创建虚拟环境pytorch_gpu conda create -n pytorch_gpu python=3.6 创建虚拟环境还是相对较快的,它会自动为本环境安装一些基本的库,等待时间无 ...
- python环境搭建什么意思_如何搭建Python环境
机器学习,人工智能是最近两年的热门话题,新技术的流行往往驱使开发者们前往一探究竟,那么要掌握机器学习,人工智能方面的技术,俗话说工欲善其事必先利其器,选择语言和IDE是前提,Python作为机器学习比 ...
- 怎么用python画sin函数图像_如何使用python的matplotlib模块画正弦函数图像
python是一个很有趣的语言,可以在命令行窗口运行.python中有很多功能强大的模块,这篇经验告诉你,如何利用python的matplotlib模块,绘制正弦函数y=sin(x)的图像. 工具/原 ...
- em算法python代码_EM 算法求解高斯混合模型python实现
注:本文是对<统计学习方法>EM算法的一个简单总结. 1. 什么是EM算法? 引用书上的话: 概率模型有时既含有观测变量,又含有隐变量或者潜在变量.如果概率模型的变量都是观测变量,可以直接 ...
- 【白话机器学习】算法理论+实战之EM聚类
1. 写在前面 如果想从事数据挖掘或者机器学习的工作,掌握常用的机器学习算法是非常有必要的,常见的机器学习算法: 监督学习算法:逻辑回归,线性回归,决策树,朴素贝叶斯,K近邻,支持向量机,集成算法Ad ...
- 白话机器学习算法理论+实战之EM聚类
1. 写在前面 如果想从事数据挖掘或者机器学习的工作,掌握常用的机器学习算法是非常有必要的,比如我之前写过的一篇十大机器学习算法的小总结,在这简单的先捋一捋, 常见的机器学习算法: 监督学习算法:逻辑 ...
- python 三种聚类算法(K-means,AGNES,DBScan)
python实现鸢尾花三种聚类算法(K-means,AGNES,DBScan) 更新时间:2019年06月27日 14:44:44 作者:weixin_42134141 这篇文章主要介绍了pyth ...
- em算法怎么对应原有分类_EM算法原理
转自:https://www.cnblogs.com/Gabby/p/5344658.html 我讲EM算法的大概流程主要三部分:需要的预备知识.EM算法详解和对EM算法的改进. 一.EM算法的预备知 ...
- Python机器学习---2.聚类算法理论部分
文章目录 1.聚类分析 1.1 无监督学习与聚类算法 1.1.1.旨在理解数据自然结构的聚类 1.1.2 用于数据处理的聚类 1.2 核心概念 1.2.1 聚类分析 1.2.2 簇 1.3 基于原型的 ...
- Python使用K-means聚类算法进行分类案例一则
K-means算法是经典的基于划分的聚类方法,是十大经典数据挖掘算法之一,其基本思想是:以空间中k个点为中心进行聚类,对最靠近它们的对象归类.通过迭代的方法,逐次更新各聚类中心的值,直至得到最好的聚类 ...
最新文章
- Health Check in eShop -- 解析微软微服务架构Demo(五)
- (一).NET SubSonic2.0 的配置
- thingsboard源码结构解析
- BugKuCTF 加密 简单加密
- 在子线程中使用runloop,正确操作NSTimer计时的注意点 三种可选方法
- 4 weekend110的hive入门
- WDK10+VS2015 驱动环境搭建
- JVM内存模型及垃圾回收机制
- HTML网站去色代码
- 【强大的数字设计工具包】Sketch 55.1 for Mac
- 有哪些值得推荐的数据可视化工具?
- 单片机实验六 动态数码管实验
- 基于51单片机的篮球记分牌设计
- 机器学习之Python Sklearn——线性回归
- 【11】Activity的生命周期
- 网络电台管理套件AzuraCast
- 软件保护系统Themida是否需要Internet才能工作?
- 如何使用虚拟机改ip地址?
- Oracle中的within,oracle的 listagg() WITHIN GROUP ()函数使用
- ansible-playbook安装mysql
热门文章
- 硬盘数据恢复:自己在家修复你的硬盘只需要5分钟就够了
- html css 怎么画星形,CSS画各种图形(五角星、吃豆人、太极图等)
- windows XP下openbravo ERP 2.40安装手迹
- 程序员编程艺术第四十一章 四十二章 荷兰国旗 矩阵相乘Strassen算法
- php cookie 注入,LiveZilla 'setCookieValue()'函数PHP对象注入漏洞
- [代码审计]信呼协同办公系统2.2存在文件上传配合云处理函数组合拳RCE
- Indy相关函数用法
- MIP(Mobile instant pages 移动网页加速器)
- [书籍阅读] Spring Persistence with Hibernate
- 使用npm运行react程序报错The 'mode' option has not been set, webpack will fallback to 'production' for th