EM（期望极大化）算法及其推广

文章目录

1. EM算法
2. EM算法收敛
3. EM算法应用
4. EM算法的推广
5. sklearn.mixture.GaussianMixture

概率模型有时既有观测变量（observable variable），又有隐变量或潜在变量（latent variable）
如果概率模型的变量都是观测变量，那么给定数据，可以直接用极大似然估计法，或贝叶斯估计法估计模型参数。
当模型含有隐变量时，不能简单地使用这些估计方法。EM算法就是含有隐变量的概率模型参数的极大似然估计法，或极大后验概率估计法。

EM 算法与初值的选择有关，选择不同的初值可能得到不同的参数估计值

1. EM算法

EM算法是含有隐变量的概率模型极大似然估计或极大后验概率估计的迭代算法
含有隐变量的概率模型的数据表示为 P(Y,Z∣θ)P(Y,Z| \theta)P(Y,Z∣θ) 。YYY 是观测变量的数据，ZZZ 是隐变量的数据，θ\thetaθ 是模型参数。
EM算法通过迭代求解观测数据的对数似然函数 L(θ)=log⁡P(Y∣θ){L}(\theta)=\log {P}(\mathrm{Y} | \theta)L(θ)=logP(Y∣θ) 的极大化，实现极大似然估计。

每次迭代包括两步：

EEE 步，求期望，即求 log⁡P(Z∣Y,θ)\log P ( Z | Y, \theta )logP(Z∣Y,θ) 关于 P(Z∣Y,θ(i))P (Z | Y, \theta ^{(i)})P(Z∣Y,θ(i)) 的期望：

Q(θ,θ(i))=∑Zlog⁡P(Y,Z∣θ)P(Z∣Y,θ(i))Q\left(\theta, \theta^{(i)}\right)=\sum_{Z} \log P(Y, Z | \theta) P\left(Z | Y, \theta^{(i)}\right)Q(θ,θ(i))=Z∑logP(Y,Z∣θ)P(Z∣Y,θ(i))
称为 QQQ 函数，这里 θ(i)\theta^{(i)}θ(i) 是参数的现估计值

MMM步，求极大，即极大化 QQQ 函数得到参数的新估计值：

θ(i+1)=arg⁡max⁡θQ(θ,θ(i))\theta^{(i+1)}=\arg \max _{\theta} Q\left(\theta, \theta^{(i)}\right)θ(i+1)=argθmaxQ(θ,θ(i))

在构建具体的EM算法时，重要的是定义QQQ函数。每次迭代中，EM算法通过极大化QQQ函数来增大对数似然函数L(θ){L}(\theta)L(θ)

2. EM算法收敛

EM算法在每次迭代后均提高观测数据的似然函数值，单调递增的，即

P(Y∣θ(i+1))⩾P(Y∣θ(i))P\left(Y | \theta^{(i+1)}\right) \geqslant P\left(Y | \theta^{(i)}\right)P(Y∣θ(i+1))⩾P(Y∣θ(i))

在一般条件下EM算法是收敛的，但不能保证收敛到全局最优

3. EM算法应用

EM算法应用极其广泛，主要应用于含有隐变量的概率模型的学习
高斯混合模型的参数估计是EM算法的一个重要应用
下一章的隐马尔可夫模型的非监督学习也是EM算法的一个重要应用

4. EM算法的推广

EM算法还可以解释为 FFF 函数的极大-极大算法
EM算法有许多变形，如 GEM 算法
GEM算法的特点是每次迭代增加 FFF 函数值（并不一定是极大化 FFF 函数），从而增加似然函数值

5. sklearn.mixture.GaussianMixture

sklearn.mixture.GaussianMixture

class sklearn.mixture.GaussianMixture(n_components=1, covariance_type='full', tol=0.001,
reg_covar=1e-06, max_iter=100, n_init=1, init_params='kmeans', weights_init=None,
means_init=None, precisions_init=None, random_state=None, warm_start=False, verbose=0, verbose_interval=10)

参数说明：参考了 https://blog.csdn.net/lihou1987/article/details/70833229

n_components: 混合高斯模型个数，默认为1
covariance_type: 协方差类型，包括{‘full’,‘tied’, ‘diag’, ‘spherical’}四种，分别对应完全协方差矩阵（元素都不为零），相同的完全协方差矩阵（HMM会用到），对角协方差矩阵（非对角为零，对角不为零），球面协方差矩阵（非对角为零，对角完全相同，球面特性），默认‘full’ 完全协方差矩阵
tol：EM迭代停止阈值，默认为1e-3.
reg_covar: 协方差对角非负正则化，保证协方差矩阵均为正，默认为0
max_iter: 最大迭代次数，默认100
n_init: 初始化次数，用于产生最佳初始参数，默认为1
init_params: {‘kmeans’, ‘random’}, defaults to ‘kmeans’.初始化参数实现方式，默认用kmeans实现，也可以选择随机产生
weights_init: 各组成模型的先验权重，可以自己设，默认按照7产生
means_init: 初始化均值，同8
precisions_init: 初始化精确度（模型个数，特征个数），默认按照7实现
random_state :随机数发生器
warm_start :若为True，则fit（）调用会以上一次fit（）的结果作为初始化参数，适合相同问题多次fit的情况，能加速收敛，默认为False。
verbose :使能迭代信息显示，默认为0，可以为1或者大于1（显示的信息不同）
verbose_interval :与13挂钩，若使能迭代信息显示，设置多少次迭代后显示信息，默认10次。

sklearn官方实例

#%%
# ---------sklearn GaussianMixture-----
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
from sklearn import mixture
from mpl_toolkits.mplot3d import Axes3D #3维绘图
n_samples = 300
np.random.seed(0)#%% 以(-10,15) 和(0,20)为中心的高斯分布
shifted_gaussian = np.random.randn(n_samples, 2) + np.array([20, 20])
shifted_gaussian2 = np.random.randn(n_samples, 2) + np.array([-10, 15])#%% 拉伸的（0,0）为中心的高斯分布
C = np.array([[0, -0.7], [3.5, 0.7]])
stretched_gaussian = np.dot(np.random.randn(n_samples, 2), C)#%% 数据堆叠
X_train = np.vstack([shifted_gaussian, stretched_gaussian,shifted_gaussian2])#%% 高斯分布 3个高斯分布
clf = mixture.GaussianMixture(n_components=3, covariance_type='full')#%% 拟合
clf.fit(X_train)#%%
x = np.linspace(-20, 30)
y = np.linspace(-20, 40)#%%
X, Y = np.meshgrid(x, y)
XX = np.array([X.ravel(), Y.ravel()]).T
Z = -clf.score_samples(XX)#%%
Z = Z.reshape(X.shape)
CS = plt.contour(X, Y, Z, norm=LogNorm(vmin=1.0, vmax=1000), levels=np.logspace(1, 4, 18))
CB = plt.colorbar(CS, shrink=0.8, extend='both')
plt.scatter(X_train[:, 0], X_train[:, 1], 0.8)
plt.show()# 3D绘图
fig = plt.figure()
ax = Axes3D(fig)
# ax.scatter(X,Y,Z)
ax.plot_surface(X, Y, Z, cmap=plt.cm.autumn)
plt.show()