zoukankan      html  css  js  c++  java
  • 混合高斯分布与 EM 算法

    概率论中的 Jensen 不等式

    对于 Jensen 不等式,通常情况下是这样的:对于 (f^{prime prime}(x) geq 0) 也就是对于凸函数而言,这个可以用中值定理来证明,同时这个公式还可以用于证明多项式的均值不等式与调和不等式,这里就不给出证明了,在概率的条件下就是:实际上,这就是均值不等式来的啊, E 表示数学期望

    [mathrm{E}[f(X)] geq f(mathrm{E} X) ]

    EM 算法讲了什么

    对于下面的第一行的式子, 本质上这个式子是函数隐变量的, 将这些隐变量显示的写出来就是下面的第二行的式子:

    [egin{aligned} ell( heta) &=sum_{i=1}^{N} log p(x ; heta) \ &=sum_{i=1}^{N} log sum_{z} p(x, z ; heta) end{aligned} ]

    现在,我们假设对于每一个 $ i $ , (Q_{i}) 表示 $ z $ 的分布,那么我们有 (sum_{z} Q_{i}(z)=1, Q_{i}(z) geq0) 。然后我们使用 Jensen 不等式将上面的式子进行放缩,写成下面的这样形式,

    [egin{aligned} ell( heta) =sum_{i} log pleft(x^{(i)} ; heta ight) &=sum_{i} log sum_{z^{(i)}} pleft(x^{(i)}, z^{(i)} ; heta ight) \ &=sum_{i} log sum_{z^{(i)}} Q_{i}left(z^{(i)} ight) frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)} \ & geq sum_{i} sum_{z^{(i)}} Q_{i}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)} end{aligned} ]

    这里的 (Q_{i}) 应该就是 EM 算法不断的迭代的关键。

    上面的式子表明 (ell( heta)) 有一个下界,从极大似然的角度考虑,我们就是要最大化这个下界,得到 $ heta $ 的取值,但是其中的 (Q_{i}) 是一个隐变量的分布,我们并不知道这个分布是什么样子的。

    上面使用 Jensen 不等式关键的一步是:

    [fleft(mathrm{E}_{z^{(i)} sim Q_{i}}left[frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)} ight] ight) geq mathrm{E}_{z^{(i)} sim Q_{i}}left[fleft(frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)} ight) ight] ]

    这个式子取等式的时候,对于数学期望而言,只有常数可以满足数学期望中的 Jensen 不等式相等。这里不做具体的证明,我们可以考虑从均值不等式来理解这个问题,假设这个常数是 $ c$ :

    [frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)}=c ]

    也就是说:

    [Q_{i}left(z^{(i)} ight) propto pleft(x^{(i)}, z^{(i)} ; heta ight) ]

    我们知道:(sum_{z} Q_{i}left(z^{(i)} ight)=1) ,那么我可以考虑这样的情况:

    [1 = sum{Q_ileft( z^{left( i ight)} ight)} propto sum{pleft( x^{left( i ight)},z^{left( i ight)}; heta ight) } = sum_{z} pleft(x^{(i)}, z^{left(i ight)} ; heta ight) ]

    这样的话,就有下面的公式:

    [egin{aligned} Q_{i}left(z^{(i)} ight) &=frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{sum_{z} pleft(x^{(i)}, z^{left(i ight)} ; heta ight)} \ &=frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{pleft(x^{(i)} ; heta ight)} \ &=pleft(z^{(i)} | x^{(i)} ; heta ight) end{aligned} ]

    这一步我们称为 E 步骤,可以得到等式成立的条件, (Q_{i}left(z^{(i)} ight) :=pleft(z^{(i)} | x^{(i)} ; heta ight)) ,下界成立的条件,我们将这个条件带入 $ell( heta) $ 得到的结果就是下界,也就是说,这个是最大似然函数的条件之一, 我们需要最大似然的就是:那么我们可以用下面的步骤来计算 $ heta$ ,

    [ heta :=arg max _{ heta} sum_{i} sum_{z^{(i)}} Q_{i}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)} ]

    那么这个 $ heta$ 是不是最优的呢?接下来我们来证明 EM 算法中一个很关键的问题,那就是上述的 E 步骤与 M 步骤不断放缩与收敛,最终得到最优的 $ heta$ 的值。

    所以 EM 算法的步骤可以表示成下面这样:

    1. 函数声明 (operatorname{EM}left(p_{Y, C}(y, c ; heta), p_{C | Y}(c | y ; heta), heta^{(0)} ight))
    2. for iteration (t in 1,2, ldots) do
    3. (Q_{i}^{(t)} leftarrow Pleft(z_i | x_i ; heta^{(t-1)} ight) quad( ext { E-step }))
    4. ( heta^{(t)} leftarrow operatorname{argmax}_{ heta} mathbb{E}_{Q_{i}^{(t)}}left[P(y, C ; heta) ight] quadleft(mathrm{M}{-mathrm{Step}} ight))
    5. if ( heta^{(t)} approx heta^{(t-1)}) then
    6. return $ heta^{(t)} $

    EM 算法的收敛问题证明

    这个收敛的思想是这样的:我们的 $ heta $ 是直接通过极大似然函数算出的来的。那么 EM 算法迭代的步骤就是,我们不断地最大化极大似然估计函数,也就是说 我们每次都最大化了 (ellleft( heta^{(t)} ight)) 的下界,只有保证 (ellleft( heta^{(t)} ight) leq ellleft( heta^{(t+1)} ight)) , 这样就会不断地逼近最优解, 说明算法是正确的。在 EM 迭代的过程中,我们不断改变是 (Q_{i}) 也就是分布函数的选择(我们本质是更新 ( heta) , 改变 ( heta) 等价于改变 (Q_i) ),

    为了更好的说明,假设上一步我们已经得到了一个最优解 (ellleft( heta^{(t)} ight)) 以及 $ heta^{(t)} $ 那么在这一步,我们满足下面的情况:

    [Q_{i}^{(t)}left(z^{(i)} ight) :=pleft(z^{(i)} | x^{(i)} ; heta^{(t)} ight) Jensen不等式条件\ 用于计算 ellleft( heta^{(t+1)} ight) ]

    [ellleft( heta^{(t)} ight)=sum_{i} sum_{z^{(i)}} Q_{i}^{(t)}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta^{(t)} ight)}{Q_{i}^{(t)}left(z^{(i)} ight)} ]

    而更新的参数 ( heta^{(t+1)}) 就是来自最大化上面右边的式子, 那么我们可以推出下面的不等式:

    [egin{aligned} ellleft( heta^{(t+1)} ight) & geq sum_{i} sum_{z^{(i)}} Q_{i}^{(t)}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta^{(t+1)} ight)}{Q_{i}^{(t)}left(z^{(i)} ight)} \ & geq sum_{i} sum_{z^{(i)}} Q_{i}^{(t)}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta^{(t)} ight)}{Q_{i}^{(t)}left(z^{(i)} ight)} \ &=ellleft( heta^{(t)} ight) end{aligned} ]

    对于上面的两个不等式做出以下解释:

    1. 第一个不等式来自 Jensen 不等式,原式是这样的:(ell( heta) geq sum_{i} sum_{z^{(i)}} Q_{i}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta ight)}{Q_{i}left(z^{(i)} ight)}) 对任意的$ Q_i $ 与( heta) 均成立,这里为什么对任意的 ( heta) 都成立, 是因为 (Q) 本质也是由 ( heta) 决定的.
    2. 对于第二个不等式, 就是我们前面 M 步骤的结果, 最大化下面的极大似然函数的时候得到$ heta^{(t+1)}$ 所以第二个式子显然:

    [sum_{i} sum_{z^{(i)}} Q_{i}^{(t)}left(z^{(i)} ight) log frac{pleft(x^{(i)}, z^{(i)} ; heta^{(t)} ight)}{Q_{i}^{(t)}left(z^{(i)} ight)} ]

    综上证明了 EM 算法的收敛性问题

    混合高斯分布

    在一般的情况下,对于所得到的样本集,(X=left{x_{1}, dots, x_{N} ight}),我们的目标是最大化似然函数,通过最大化似然函数来获取参数的值。这是似然函数往往取对数表示就是:

    [egin{aligned} L( heta | X) &=log left(prod_{i=1}^{N} pleft(x_{i} | heta ight) ight) \ &=sum_{i=1}^{N} log pleft(x_{i} | heta ight) end{aligned} ]

    这个过程的参数估计可以描述成:

    [hat{ heta}=arg max _{ heta} L( heta | X) ]

    这个结果是可以直接计算出来的,那么引入混合高斯分布会是什么样呢?

    混合高斯分布:

    简单的说就是,非混合的情况下,我们的数据集满足高斯分布,这样用极大似然就直接算出高斯分布中的参数就可以了。那么混合的情况下就是,数据集是由多个满足高斯分布的数据集线性组合起来,这个时候我们理解为:有 (k) 个满足不同参数的高斯分布的原数据集,并且 (sum_{j=1}^{k} pi_{j}=1).

    此时:

    [p(oldsymbol{x})=sum_{k=1}^{K} pi_{k} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{k}, oldsymbol{Sigma}_{k} ight) ]

    举个例子说明:

    我们假设男生的身高符合正态分布, (mathcal{N_1}left(oldsymbol{mu}_{1}, oldsymbol{Sigma}_{1} ight)) 女生的身高符合正态分布 (mathcal{N_2}left(oldsymbol{mu}_{2}, oldsymbol{Sigma}_{2} ight)) , 现在我们得到一个关于身高的样本数据集, 但是不知道男生与女生的比例, 要拟合这个身高数据集符合哪种分布, 这里就可以使用混合高斯分布,

    [p(oldsymbol{x})= pi_{1} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{1}, oldsymbol{Sigma}_{1} ight) + pi_{2} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{2}, oldsymbol{Sigma}_{2} ight) ]

    其中 (pi_1)(pi_2) 就表示男生与女生的比例. 也表示 (p(x)) 对不同高斯分布的可能性,

    在上面的问题中, (pi_1)(pi_2) 是未知参数 , 那么这个例子与一个问题很类似, 那就是贝叶斯分布的问题, 如果我们在已知后验概率的情况下, 就可以很简单的计算先验概率了, 也就是说, 如果我们已知男生满足的高斯分布的参数, (mu_1)({Sigma}_{1}) 和 女生高斯分布的参数 (mu_2)({Sigma}_{2}) , 这时就可以反过来使用贝叶斯计算出 (pi_1)(pi_2) 的值.

    混合高斯分布与贝叶斯

    对于一个混合高斯分布来说, (pi_1, pi_2, dots, pi_n) 对这个分布来说也是未知参数, 这里我们从贝叶斯的角度来改写这个未知参数的形式, 这样做的目的是为了更好的使用下面的 EM 算法以及 EM 算法的计算方式, 对于上面的 (N) 维混合高斯分布, 我们取一个向量(N) 维向量 ({z}) 表示取得的, 其中 $ z_k(1 leq k leq K) $, 只能取 0 与 1 两个值, 其中 (p(z_i == 1)) 表示第 i 中高斯分布被选中的概率, 而 (p(z_i = 0)) 表示未被选中的概率, 那么用向量 (z) 表示混合高斯分布就有两种方式, 一种是类似于

    [[1,1,1,0,0,dots, 1, 0, 1] ]

    的形式, 这样直接使用向量的乘法的到概率,

    而另一种是

    [[0,1,0, 0dots, dots, 0] ]

    的形式, 这里只有一个 (z_i) 为 1, 因为我们要使用贝叶斯分布, 反过来使用全概率公式, 所以我们使用第二种形式, 在这种形式下, (z_i) 之间是相互独立的, 所以 (z) 向量的概率可以写成:

    [p(oldsymbol{z})=pleft(z_{1} ight) pleft(z_{2} ight) ldots pleft(z_{K} ight)=prod_{k=1}^{K} pi_{k}^{z_{k}} ]

    然后对应于 (z_i) 的那一个分布是服从高斯分布的, 可以将对应的 (z_i) 分布的概率表示如下:

    [pleft(oldsymbol{x} | z_{i}=1 ight)=mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight) ]

    进而上式有可以写成如下形式:

    [p(oldsymbol{x} | oldsymbol{z})=prod_{i=1}^{N} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight)^{z_{i}} ]

    我们知道, 在向量 (z) 中, 假设 (z_j = 1) 那么对于所有的 (z_i (i eq j)) 都有 (z_i = 0), 所以上式成立很显然, 我们得到了 (p(z))(p(x|z)) 就可以得到 (p(x)) 的公式了:

    [egin{aligned} p(oldsymbol{x}) &=sum_{oldsymbol{z}} p(oldsymbol{z}) p(oldsymbol{x} | oldsymbol{z}) \ &=sum_{oldsymbol{z}}left(prod_{k=1}^{K} pi_{k}^{z_{k}} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{k}, oldsymbol{Sigma}_{k} ight)^{z_{k}} ight) \ &=sum_{k=1}^{K} pi_{k} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{k}, oldsymbol{Sigma}_{k} ight) end{aligned} ag{1} ]

    公式 (1) 本质上使用的是全概率公式, 因为向量 (z) 的表示方式我们取的是第二种, 上式中 (p(z)p(x|z)) 的结果是 (p(x,z)) 所以要用全概率公式相加, 从上面的式子我们惊奇的发现, 我们 (p(x)) 改写成了后验概率的形式, 因为我们划分了所有的 (z_i) , 在不同的向量 (z) 的基础上得到的概率 (p(x)), 这时反过来就可以求先验概率:

    [egin{aligned} gammaleft(z_{i} ight) &=pleft(z_{i}=1 | oldsymbol{x} ight) \ &=frac{pleft(z_{i}=1 ight) pleft(oldsymbol{x} | z_{i}=1 ight)}{pleft(oldsymbol{x}, z_{i}=1 ight)} \ &=frac{pleft(z_{i}=1 ight) pleft(oldsymbol{x} | z_{i}=1 ight)}{sum_{j=1}^{N} pleft(z_{j}=1 ight) pleft(oldsymbol{x} | z_{j}=1 ight)}( ext { 全概率公式 }) \ &=frac{pi_{i} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{oldsymbol{i}} ight)}{sum_{j=1}^{N} pi_{j} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{j}, oldsymbol{Sigma}_{j} ight)}( ext {结合}(1) . end{aligned} ag{2} ]

    EM 算法用于混合高斯分布

    前面的 EM 算法是一般情况下的 EM 算法求解以及证明的过程, 那么在实际的计算的过程中, EM 算法求解 GMM 问题是什么样的呢? 下面的公式 (3) 是基本的公式,

    [p(oldsymbol{x|pi, mu, Sigma})=sum_{k=1}^{K} pi_{i} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight) ag{3} ]

    这个公式与 EM 算法最终的最大化的公式是一样的(因为 Jensen 不等式的条件), 所以对这个公式求极大似然,

    [prod_{i}^{N} p(oldsymbol{x|pi, mu, Sigma}) = prod_{n}^{N} sum_{i}^{N} pi_{i} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight) ag{4} ]

    然后对这个公式求对数, 以及求偏导, 可以很容易得到下面的公式:

    [0=-sum_{n=1}^{N} frac{pi_{i} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight)}{sum_{j} pi_{j} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{j}, oldsymbol{Sigma}_{j} ight)} oldsymbol{Sigma}_{i}^{-1}left(oldsymbol{x}_{n}-oldsymbol{mu}_{i} ight) ag{5} ]

    对于公式(5), 前面部分是通常的对数求导, 将 (sum^{K} pi_{i} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight)) 放在分母, 然后后面的 (oldsymbol{Sigma}_{i}^{-1}left(oldsymbol{x}_{n}-oldsymbol{mu}_{i} ight)) 就是一维的正态分布的求导公式, 这样就可以求出 (mu_i) 的值了,

    [oldsymbol{mu}_{i}=frac{1}{N_{i}} sum_{n=1}^{N} gammaleft(z_{n i} ight) oldsymbol{x}_{n} ag{6.1} ]

    其中

    [N_i = sum_{n=1}^N gammaleft(z_{n i} ight) ag{6.2} ]

    对上面公式的解释:

    1. (gammaleft(z_{n i} ight)) 表示对于样本 (x_n) 属于第 (i) 类的概率, 在男生女生问题中就是对于某一个样本属于男生或者女生的概率
    2. (N_i) 就是对于所有样本, 属于第 (i) 类的概率的个数总和
    3. (mu_i) 就将这个值进行加权平均, 在公式 6.1 中, (x_n) 表示权值, (gammaleft(z_{n i} ight)) 表示概率

    同理, 使用极大似然的方法可以求得:

    [oldsymbol{Sigma}_{i}=frac{1}{N_{i}} sum_{n=1}^{N} gammaleft(z_{n i} ight)left(x_{n}-oldsymbol{mu}_{i} ight)left(x_{n}-oldsymbol{mu}_{i} ight)^{T} ]

    最后是 (pi_1, pi_2, dots, pi_n) , 这个的极大值是条件极大值, 条件限制是 (sum_{1}^{N} pi_i = 1) , 也就是在求偏导的时候需要加入拉格朗日算子, 使用拉格朗日乘子法:

    [ln p(oldsymbol{x} | oldsymbol{pi}, oldsymbol{mu}, oldsymbol{Sigma})+lambdaleft(sum_{i=1}^{N} pi_{i}-1 ight) ag{7.1} ]

    求上式关于 (pi_i) 的极大似然函数, 得到:

    [0=sum_{n=1}^{N} frac{mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight)}{sum_{j} pi_{j} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{j}, oldsymbol{Sigma}_{j} ight)}+lambda ag{7.2} ]

    公式 7.2 求导的过程就是, 第一部分与前面对 (mu) 求导类似, 这里分子不变很显然, 对于结果两边同乘以 (pi_i) 得到:

    [0=sum_{n=1}^{N} frac{pi_{i} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight)}{sum_{j} pi_{j} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{j}, oldsymbol{Sigma}_{j} ight)}+lambda pi_{i} ]

    结合公式 (2) 可以得到:

    [0 = N_i + lambda pi_i \ 0 = sum_{i = 1}^{N} N_i + lambda sum_{i=1}^{N} pi_i \ 0 = N+k ag{7.3} ]

    最终得到就是:

    [pi_i = frac{N_i}{N} ]

    EM 算法求解的迭代过程

    EM 算法的不断迭代收敛过程就是不断地计算 E 步骤与 M 步骤, 直到期望的似然函数收敛,

    E Step

    根据当前的 (pi_i) , (mu_i) , (Sigma_i) 计算:

    [y(z_nk) = frac{pi_{i} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight)}{sum_{j}^{N} pi_{j} mathcal{N}left(oldsymbol{x}_{n} | oldsymbol{mu}_{j}, oldsymbol{Sigma}_{j} ight)} ]

    同时需要记录的值为:

    [p(oldsymbol{x})=sum_{i=1}^{N} pi_{i} mathcal{N}left(oldsymbol{x} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight) ]

    M Step

    使用上面的公式更新参数:

    [oldsymbol{mu}_{i+1}=frac{1}{N_{i}} sum_{n=1}^{N} gammaleft(z_{n i} ight) oldsymbol{x}_{n} \ oldsymbol{Sigma}_{i+1}=frac{1}{N_{i}} sum_{n=1}^{N} gammaleft(z_{n i} ight)left(x_{n}-oldsymbol{mu}_{i} ight)left(x_{n}-oldsymbol{mu}_{i} ight)^{T} \\ pi_{i+1} = frac{N_i}{N} ag{8} ]

    其中:

    [N_i = sum_{n=1}^N gammaleft(z_{n i} ight) ag{6.2} ]

    最后再计算,

    [ln p(oldsymbol{x} | oldsymbol{pi}, oldsymbol{mu}, oldsymbol{Sigma})=sum_{n=1}^{N} ln left{sum_{i=1}^{N} pi_{i} mathcal{N}left(oldsymbol{x}_{i} | oldsymbol{mu}_{i}, oldsymbol{Sigma}_{i} ight) ight} ]

    判断是否收敛, 不收敛, 返回至 E Step, 然后在 M Step 继续更新, 直到收敛.

  • 相关阅读:
    PAT (Advanced Level) Practice 1055 The World's Richest (25 分) (结构体排序)
    PAT (Advanced Level) Practice 1036 Boys vs Girls (25 分)
    PAT (Advanced Level) Practice 1028 List Sorting (25 分) (自定义排序)
    PAT (Advanced Level) Practice 1035 Password (20 分)
    PAT (Advanced Level) Practice 1019 General Palindromic Number (20 分) (进制转换,回文数)
    PAT (Advanced Level) Practice 1120 Friend Numbers (20 分) (set)
    从零开始吧
    Python GUI编程(TKinter)(简易计算器)
    PAT 基础编程题目集 6-7 统计某类完全平方数 (20 分)
    PAT (Advanced Level) Practice 1152 Google Recruitment (20 分)
  • 原文地址:https://www.cnblogs.com/wevolf/p/11039101.html
Copyright © 2011-2022 走看看