EM算法实例
  通过实例可以快速了解EM算法的基本思想,具体推导请点文末链接。图a是让我们预热的,图b是EM算法的实例。
  这是一个抛硬币的例子,H表示正面向上,T表示反面向上,参数θ表示正面朝上的概率。硬币有两个,A和B,硬币是有偏的。本次实验总共做了5组,每组随机选一个硬币,连续抛10次。如果知道每次抛的是哪个硬币,那么计算参数θ就非常简单了,如下图所示:
  如果不知道每次抛的是哪个硬币呢?那么,我们就需要用EM算法,基本步骤为:
  1、给θAθ_AθA​和θBθ_BθB​一个初始值;
  2、(E-step)估计每组实验是硬币A的概率(本组实验是硬币B的概率=1-本组实验是硬币A的概率)。分别计算每组实验中,选择A硬币且正面朝上次数的期望值,选择B硬币且正面朝上次数的期望值;
  3、(M-step)利用第三步求得的期望值重新计算θAθ_AθA​和θBθ_BθB​;
  4、当迭代到一定次数,或者算法收敛到一定精度,结束算法,否则,回到第2步。

  计算过程详解:初始值θA(0)θ_A^{(0)}θA(0)​=0.6,θB(0)θ_B^{(0)}θB(0)​=0.5。
由两个硬币的初始值0.6和0.5,容易得出投掷出5正5反的概率是pA=C105∗(0.65)∗(0.45)p_A=C^5_{10}*(0.6^5)*(0.4^5)pA​=C105​∗(0.65)∗(0.45),pB=C105∗(0.55)∗(0.55)p_B=C_{10}^5*(0.5^5)*(0.5^5)pB​=C105​∗(0.55)∗(0.55), pAp_ApA​/(pAp_ApA​+pBp_BpB​)=0.449, 0.45就是0.449近似而来的,表示第一组实验选择的硬币是A的概率为0.45。然后,0.449 * 5H = 2.2H ,0.449 * 5T = 2.2T ,表示第一组实验选择A硬币且正面朝上次数和反面朝上次数的期望值都是2.2,其他的值依次类推。最后,求出θA(1)θ_A^{(1)}θA(1)​=0.71,θB(1)θ_B^{(1)}θB(1)​=0.58。重复上述过程,不断迭代,直到算法收敛到一定精度为止。
  这篇博客对EM算法的推导非常详细,链接如下:

https://blog.csdn.net/zhihua_oba/article/details/73776553

  Python实现

#coding=utf-8
from numpy import *
from scipy import stats
import time
start = time.perf_counter()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 stepfor observation in observations:len_observation = len(observation)num_heads = observation.sum()num_tails = len_observation-num_heads#二项分布求解公式contribution_A = stats.binom.pmf(num_heads,len_observation,theta_A)contribution_B = 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_headscounts['A']['T'] += weight_A * num_tailscounts['B']['H'] += weight_B * num_headscounts['B']['T'] += weight_B * num_tails# M stepnew_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]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 = abs(prior[0]-new_prior[0])if delta_change < tol:breakelse:prior = new_prioriteration +=1return [new_prior,iteration]#硬币投掷结果
observations = 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]])
print (em(observations,[0.6,0.5]))
end = time.perf_counter()
print('Running time: %f seconds'%(end-start))

EM算法实例(Python)相关推荐

  1. EM算法实例及python实现

    现在有两个硬币A和B,要估计的参数是它们各自翻正面(head)的概率.观察的过程是先随机选A或者B,然后扔10次.以上步骤重复5次.如果知道每次选的是A还是B,那可以直接估计(见下图a).如果不知道选 ...

  2. em算法 实例 正态分布_EM算法解GMM

    看了很多介绍EM算法的文章,但是他们都没有代码,所以在这里写出来. Jensen 不等式 参考期望最大算法 Jensen不等式在优化理论中大量用到,首先来回顾下凸函数和凹函数的定义.假设 是定义域为实 ...

  3. em算法 实例 正态分布_Petuum提出序列生成学习算法通用框架

    近日,来自人工智能创业公司 Petuum 的研究人员发表论文,提出序列生成学习算法的通用框架--广义的熵正则化策略优化框架(Generalized Entropy-Regularized Policy ...

  4. em算法 实例 正态分布_人人都能看懂的EM算法推导

    ↑ 点击蓝字 关注极市平台作者丨August@知乎(已授权)来源丨https://zhuanlan.zhihu.com/p/36331115编辑丨极市平台 极市导读 EM算法到底是什么,公式推导怎么去 ...

  5. EM算法及python简单实现

    目录 1.摘要 2.EM算法简介 3.预备知识 3.1贝叶斯公式 3.2 极大似然估计 3.2.1似然函数 3.2.2问题描述 3.2.3极大似然函数估计值的求解步骤 3.3Jensen不等式(琴生不 ...

  6. em算法 实例 正态分布_【机器学习】EM算法详细推导和讲解

    今天不太想学习,炒个冷饭,讲讲机器学习十大算法里有名的EM算法,文章里面有些个人理解,如有错漏,还请读者不吝赐教. 众所周知,极大似然估计是一种应用很广泛的参数估计方法.例如我手头有一些东北人的身高的 ...

  7. 统计学习 EM算法 Python实现

    1.EM算法是什么 EM算法可以用于有监督学习,也可以用于无监督学习.这个算法是根据观测结果求得对含有隐变量的模型的参数的估计.包含E步骤和M步,E步是求期望,M步是求极大似然估计,极大参数估计是对模 ...

  8. 王小草【机器学习】笔记--EM算法

    原文地址:http://blog.csdn.net/sinat_33761963/article/details/53520898 EM算法的英文全称是Expectation Maximization ...

  9. 计算1至1000间的合数c语言,输出1000以内的素数的算法(实例代码)

    输出1000以内的素数的算法(实例代码) 代码如下所示: 复制代码 代码如下: #include "stdafx.h" #include #include bool IsSushu ...

最新文章

  1. Codeforces 895C - Square Subsets
  2. what???现在的研究生和导师普遍都没有真正理解科研的本质
  3. Linux操作系统用户登录失败次数过多被锁定的解决方法
  4. VM pow 函数 :undefined reference to `pow'
  5. AndroidManifest Intent-Filter Action android:name属性
  6. 对话阿里云:解锁视频云的新技术、新场景
  7. 途牛windows转linux,在 Windows 中通过 VirtualBox 启动物理硬盘上的 Linux 操作系统...
  8. Spring源码学习笔记:Spring设计模式对比和Spring的OOB,BOP,AOP,IOC,DI/DL
  9. bae mysql_获取BAE上的MySQL相关信息
  10. 了解Hadoop数据类型,输入输出格式及用户如何自定义。
  11. 简单mysql主从配置
  12. 基于灰度的模板匹配算法
  13. 认识计算机听课记录20篇,【中学信息技术听课记录】 信息技术听课记录15篇及评析_初中信息技术听课记录_高中信息技术听课记录20篇_东城教研...
  14. canvas绘制出货单
  15. LCD液晶显示屏颜色显示波长研究与总结?
  16. 数据链路层(一、二)——差错控制
  17. 兰州交通大学php,航拍兰州交通大学校园∣让我再看你一遍 从南到北
  18. 牛客网SQL--某东篇
  19. HTML学习笔记——列表标签
  20. 上课笔记--商法(上)

热门文章

  1. HDU - 2027 统计元音
  2. 一师一优课课堂实录下载与合并
  3. hill图matlab代码,Hill密码的加密论文(内含matlab程序代码).doc
  4. Hone C# III
  5. 这种让你肚子疼的分子机制找到了!西湖大学Cell论文揭示病菌入侵人体“新大门”...
  6. BVS智能视频分析—智慧城管解决方案
  7. easyUI -datagrid表格数据不显示
  8. 【物联网/智慧城市】2021年物联网与智慧城市国际学术会议 (IoTSC2021)
  9. 毫米波太赫兹电路板设计--高频电路板布线损耗分析
  10. 2023复试题预测:送你一堆冰墩墩,拿走不谢~文都管联院