# 建模基础教学：EM-期望最大化算法

47 篇文章 75 订阅

## python算法实现

#! /usr/bin/env python
#! -*- coding=utf-8 -*-

#模拟两个正态分布的均值估计

from numpy import *
import numpy as np
import random
import copy

SIGMA = 6
EPS = 0.0001
#生成方差相同,均值不同的样本
def generate_data():
Miu1 = 20
Miu2 = 40
N = 1000
X = mat(zeros((N,1)))
for i in range(N):
temp = random.uniform(0,1)
if(temp > 0.5):
X[i] = temp*SIGMA + Miu1
else:
X[i] = temp*SIGMA + Miu2
return X

#EM算法
def my_EM(X):
k = 2
N = len(X)
Miu = np.random.rand(k,1)
Posterior = mat(zeros((N,2)))
dominator = 0
numerator = 0
#先求后验概率
for iter in range(1000):
for i in range(N):
dominator = 0
for j in range(k):
dominator = dominator + np.exp(-1.0/(2.0*SIGMA**2) * (X[i] - Miu[j])**2)
#print dominator,-1/(2*SIGMA**2) * (X[i] - Miu[j])**2,2*SIGMA**2,(X[i] - Miu[j])**2
#return
for j in range(k):
numerator = np.exp(-1.0/(2.0*SIGMA**2) * (X[i] - Miu[j])**2)
Posterior[i,j] = numerator/dominator
oldMiu = copy.deepcopy(Miu)
#最大化
for j in range(k):
numerator = 0
dominator = 0
for i in range(N):
numerator = numerator + Posterior[i,j] * X[i]
dominator = dominator + Posterior[i,j]
Miu[j] = numerator/dominator
print (abs(Miu - oldMiu)).sum()
#print '\n'
if (abs(Miu - oldMiu)).sum() < EPS:
print Miu,iter
break

if __name__ == '__main__':
X = generate_data()
my_EM(X)

## 获取建模资料

• 0
点赞
• 0
评论
• 4
收藏
• 打赏
• 扫一扫，分享海报

04-16 367

09-18 82
04-24 2498
02-22
07-03 3730
11-16 66
12-03 2192
03-15
01-08
06-18 450
08-29 4163
07-08 2456
11-07 9901

mathor_cup

¥2 ¥4 ¥6 ¥10 ¥20

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、C币套餐、付费专栏及课程。