''' 数据集:伪造数据集(两个高斯分布混合) 数据集长度:1000 ------------------------------ 运行结果: ---------------------------- the Parameters set is: alpha0:0.3, mu0:0.7, sigmod0:-2.0, alpha1:0.5, mu1:0.5, sigmod1:1.0 ---------------------------- the Parameters predict is: alpha0:0.4, mu0:0.6, sigmod0:-1.7, alpha1:0.7, mu1:0.7, sigmod1:0.9 ---------------------------- ''' import numpy as np import random import math import time def loadData(mu0, sigma0, mu1, sigma1, alpha0, alpha1): ''' 初始化数据集 这里通过服从高斯分布的随机函数来伪造数据集 :param mu0: 高斯0的均值 :param sigma0: 高斯0的方差 :param mu1: 高斯1的均值 :param sigma1: 高斯1的方差 :param alpha0: 高斯0的系数 :param alpha1: 高斯1的系数 :return: 混合了两个高斯分布的数据 ''' # 定义数据集长度为1000 length = 1000 # 初始化第一个高斯分布,生成数据,数据长度为length * alpha系数,以此来 # 满足alpha的作用 data0 = np.random.normal(mu0, sigma0, int(length * alpha0)) # 第二个高斯分布的数据 data1 = np.random.normal(mu1, sigma1, int(length * alpha1)) # 初始化总数据集 # 两个高斯分布的数据混合后会放在该数据集中返回 dataSet = [] # 将第一个数据集的内容添加进去 dataSet.extend(data0) # 添加第二个数据集的数据 dataSet.extend(data1) # 对总的数据集进行打乱(其实不打乱也没事,只不过打乱一下直观上让人感觉已经混合了 # 读者可以将下面这句话屏蔽以后看看效果是否有差别) random.shuffle(dataSet) #返回伪造好的数据集 return dataSet def calcGauss(dataSetArr, mu, sigmod): ''' 根据高斯密度函数计算值 依据:“9.3.1 高斯混合模型” 式9.25 注:在公式中y是一个实数,但是在EM算法中(见算法9.2的E步),需要对每个j 都求一次yjk,在本实例中有1000个可观测数据,因此需要计算1000次。考虑到 在E步时进行1000次高斯计算,程序上比较不简洁,因此这里的y是向量,在numpy 的exp中如果exp内部值为向量,则对向量中每个值进行exp,输出仍是向量的形式。 所以使用向量的形式1次计算即可将所有计算结果得出,程序上较为简洁 :param dataSetArr: 可观测数据集 :param mu: 均值 :param sigmod: 方差 :return: 整个可观测数据集的高斯分布密度(向量形式) ''' # 计算过程就是依据式9.25写的,没有别的花样 result = (1 / (math.sqrt(2*math.pi)*sigmod**2)) * np.exp(-1 * (dataSetArr-mu) * (dataSetArr-mu) / (2*sigmod**2)) # 返回结果 return result def E_step(dataSetArr, alpha0, mu0, sigmod0, alpha1, mu1, sigmod1): ''' EM算法中的E步 依据当前模型参数,计算分模型k对观数据y的响应度 :param dataSetArr: 可观测数据y :param alpha0: 高斯模型0的系数 :param mu0: 高斯模型0的均值 :param sigmod0: 高斯模型0的方差 :param alpha1: 高斯模型1的系数 :param mu1: 高斯模型1的均值 :param sigmod1: 高斯模型1的方差 :return: 两个模型各自的响应度 ''' # 计算y0的响应度 # 先计算模型0的响应度的分子 gamma0 = alpha0 * calcGauss(dataSetArr, mu0, sigmod0) # 模型1响应度的分子 gamma1 = alpha1 * calcGauss(dataSetArr, mu1, sigmod1) # 两者相加为E步中的分布 sum = gamma0 + gamma1 # 各自相除,得到两个模型的响应度 gamma0 = gamma0 / sum gamma1 = gamma1 / sum # 返回两个模型响应度 return gamma0, gamma1 def M_step(muo, mu1, gamma0, gamma1, dataSetArr): # 依据算法9.2计算各个值 # 这里没什么花样,对照书本公式看看这里就好了 mu0_new = np.dot(gamma0, dataSetArr) / np.sum(gamma0) mu1_new = np.dot(gamma1, dataSetArr) / np.sum(gamma1) sigmod0_new = math.sqrt(np.dot(gamma0, (dataSetArr - muo)**2) / np.sum(gamma0)) sigmod1_new = math.sqrt(np.dot(gamma1, (dataSetArr - mu1)**2) / np.sum(gamma1)) alpha0_new = np.sum(gamma0) / len(gamma0) alpha1_new = np.sum(gamma1) / len(gamma1) # 将更新的值返回 return mu0_new, mu1_new, sigmod0_new, sigmod1_new, alpha0_new, alpha1_new def EM_Train(dataSetList, iter=500): ''' 根据EM算法进行参数估计 算法依据“9.3.2 高斯混合模型参数估计的EM算法” 算法9.2 :param dataSetList:数据集(可观测数据) :param iter: 迭代次数 :return: 估计的参数 ''' # 将可观测数据y转换为数组形式,主要是为了方便后续运算 dataSetArr = np.array(dataSetList) # 步骤1:对参数取初值,开始迭代 alpha0 = 0.5 mu0 = 0 sigmod0 = 1 alpha1 = 0.5 mu1 = 1 sigmod1 = 1 # 开始迭代 step = 0 while (step < iter): # 每次进入一次迭代后迭代次数加1 step += 1 # 步骤2:E步:依据当前模型参数,计算分模型k对观测数据y的响应度 gamma0, gamma1 = E_step(dataSetArr, alpha0, mu0, sigmod0, alpha1, mu1, sigmod1) # 步骤3:M步 mu0, mu1, sigmod0, sigmod1, alpha0, alpha1 = M_step(mu0, mu1, gamma0, gamma1, dataSetArr) # 迭代结束后将更新后的各参数返回 return alpha0, mu0, sigmod0, alpha1, mu1, sigmod1 if __name__ == '__main__': start = time.time() # 设置两个高斯模型进行混合,这里是初始化两个模型各自的参数 # 见“9.3 EM算法在高斯混合模型学习中的应用” # alpha是“9.3.1 高斯混合模型” 定义9.2中的系数α # mu0是均值μ # sigmod是方差σ # 在设置上两个alpha的和必须为1,其他没有什么具体要求,符合高斯定义就可以 alpha0 = 0.3 # 系数α mu0 = -2 # 均值μ sigmod0 = 0.5 # 方差σ alpha1 = 0.7 # 系数α mu1 = 0.5 # 均值μ sigmod1 = 1 # 方差σ # 初始化数据集 dataSetList = loadData(mu0, sigmod0, mu1, sigmod1, alpha0, alpha1) #打印设置的参数 print('---------------------------') print('the Parameters set is:') print('alpha0:%.1f, mu0:%.1f, sigmod0:%.1f, alpha1:%.1f, mu1:%.1f, sigmod1:%.1f' % ( alpha0, alpha1, mu0, mu1, sigmod0, sigmod1 )) # 开始EM算法,进行参数估计 alpha0, mu0, sigmod0, alpha1, mu1, sigmod1 = EM_Train(dataSetList) # 打印参数预测结果 print('----------------------------') print('the Parameters predict is:') print('alpha0:%.1f, mu0:%.1f, sigmod0:%.1f, alpha1:%.1f, mu1:%.1f, sigmod1:%.1f' % ( alpha0, alpha1, mu0, mu1, sigmod0, sigmod1 )) # 打印时间 print('----------------------------') print('time span:', time.time() - start)
以上就是python em算法的实现的详细内容,更多关于python em算法的资料请关注其它相关文章!
广告合作:本站广告合作请联系QQ:858582 申请时备注:广告合作(否则不回)
免责声明:本站资源来自互联网收集,仅供用于学习和交流,请遵循相关法律法规,本站一切资源不代表本站立场,如有侵权、后门、不妥请联系本站删除!
免责声明:本站资源来自互联网收集,仅供用于学习和交流,请遵循相关法律法规,本站一切资源不代表本站立场,如有侵权、后门、不妥请联系本站删除!
暂无评论...
更新日志
2024年12月25日
2024年12月25日
- 小骆驼-《草原狼2(蓝光CD)》[原抓WAV+CUE]
- 群星《欢迎来到我身边 电影原声专辑》[320K/MP3][105.02MB]
- 群星《欢迎来到我身边 电影原声专辑》[FLAC/分轨][480.9MB]
- 雷婷《梦里蓝天HQⅡ》 2023头版限量编号低速原抓[WAV+CUE][463M]
- 群星《2024好听新歌42》AI调整音效【WAV分轨】
- 王思雨-《思念陪着鸿雁飞》WAV
- 王思雨《喜马拉雅HQ》头版限量编号[WAV+CUE]
- 李健《无时无刻》[WAV+CUE][590M]
- 陈奕迅《酝酿》[WAV分轨][502M]
- 卓依婷《化蝶》2CD[WAV+CUE][1.1G]
- 群星《吉他王(黑胶CD)》[WAV+CUE]
- 齐秦《穿乐(穿越)》[WAV+CUE]
- 发烧珍品《数位CD音响测试-动向效果(九)》【WAV+CUE】
- 邝美云《邝美云精装歌集》[DSF][1.6G]
- 吕方《爱一回伤一回》[WAV+CUE][454M]