python em算法的实现

2023-07-14 09:00:03 74
'''
数据集：伪造数据集（两个高斯分布混合）
数据集长度：1000
------------------------------
运行结果：
----------------------------
theParameterssetis:
alpha0:0.3,mu0:0.7,sigmod0:-2.0,alpha1:0.5,mu1:0.5,sigmod1:1.0
----------------------------
theParameterspredictis:
alpha0:0.4,mu0:0.6,sigmod0:-1.7,alpha1:0.7,mu1:0.7,sigmod1:0.9
----------------------------
'''

importnumpyasnp
importrandom
importmath
importtime

defloadData(mu0,sigma0,mu1,sigma1,alpha0,alpha1):
'''
初始化数据集
这里通过服从高斯分布的随机函数来伪造数据集
:parammu0:高斯0的均值
:paramsigma0:高斯0的方差
:parammu1:高斯1的均值
:paramsigma1:高斯1的方差
:paramalpha0:高斯0的系数
:paramalpha1:高斯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)

#返回伪造好的数据集
returndataSet

defcalcGauss(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次计算即可将所有计算结果得出，程序上较为简洁
:paramdataSetArr:可观测数据集
:parammu:均值
:paramsigmod:方差
:return:整个可观测数据集的高斯分布密度（向量形式）
'''
#计算过程就是依据式9.25写的，没有别的花样
result=(1/(math.sqrt(2*math.pi)*sigmod**2))*np.exp(-1*(dataSetArr-mu)*(dataSetArr-mu)/(2*sigmod**2))
#返回结果
returnresult


defE_step(dataSetArr,alpha0,mu0,sigmod0,alpha1,mu1,sigmod1):
'''
EM算法中的E步
依据当前模型参数，计算分模型k对观数据y的响应度
:paramdataSetArr:可观测数据y
:paramalpha0:高斯模型0的系数
:parammu0:高斯模型0的均值
:paramsigmod0:高斯模型0的方差
:paramalpha1:高斯模型1的系数
:parammu1:高斯模型1的均值
:paramsigmod1:高斯模型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

#返回两个模型响应度
returngamma0,gamma1

defM_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)

#将更新的值返回
returnmu0_new,mu1_new,sigmod0_new,sigmod1_new,alpha0_new,alpha1_new


defEM_Train(dataSetList,iter=500):
'''
根据EM算法进行参数估计
算法依据“9.3.2高斯混合模型参数估计的EM算法”算法9.2
:paramdataSetList:数据集（可观测数据）
:paramiter:迭代次数
: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
以上就是pythonem算法的实现的详细内容，更多关于pythonem算法的资料请关注毛票票其它相关文章！
python em算法的实现

热门推荐

随机推荐
