机器学习---算法基础(十二)LDA主题模型

参考文章:
一文详解LDA主题模型
小白都能看懂的蒙特卡洛方法以及python实现
机器之心——马尔可夫模型
马尔可夫链蒙特卡罗算法(MCMC)
MCMC(一)蒙特卡罗方法
MCMC(二)马尔科夫链
MCMC(三)MCMC采样和M-H采样
通俗理解LDA模型

数学基础

Gamma分布

Beta分布与Dirichlet分布

随机样本采样

蒙特卡洛方法

基本思想

蒙特卡罗方法也称统计模拟方法,是1940年代中期由于科学技术的发展和电子计算机的发明,而提出的一种以概率统计理论为指导的数值计算方法。是指使用随机数(或更常见的伪随机数)来解决很多计算问题的方法。

蒙特卡洛思想可以粗略的分成两类。

  • 求解问题的内在随机性。借助计算机模拟这种随机过程
  • 求解问题可以转化为某种随机分布的特征数,比如随机事件出现的概率,或者随机变量的期望值。通过随机抽样的方法,以随机事件出现的频率估计其概率,或者以抽样的数字特征估算随机变量的数字特征,并将其作为问题的解。这种方法多用于求解复杂的多维积分问题。
蒙特卡洛求积分

如果需要求一个分布在某个区间内的可能性,一般我们会采用如下的方式进行计算:
I = ∫ a b f ( x ) d x I=\int_a^b{f(x)dx} I=abf(x)dx
对于一个未知的分布,我们在不知道概率密度的情况系就很难计算积分。对于这种情况我们可以使用一种方式来计算某一种分布。
在这里插入图片描述在这里插入图片描述当我们在[a,b]之间随机取一点x时,它对应的函数值就是f(x)。接下来我们就可以用f(x) * (b - a)来粗略估计曲线下方的面积,也就是我们需要求的积分值,当然这种估计(或近似)是非常粗略的。

通过计算我们可以估计出,这个区间的概率分布。
S = 1 4 ( b − a ) ∑ 1 4 f ( x i ) S=\frac14(b-a)\sum_1^4{f(x_i)} S=41(ba)14f(xi)

当我们将这个概率上的分布取无限多个时,我们假设假设概率分布为 g ( x ) g(x) g(x),抽样的概率密度函数为f(x)且满足概率密度为 ∫ a b f ( x ) d x = 1 \int_a^bf(x)dx=1 abf(x)dx=1是在之前的抽样中其实是有一个假设的——g(x)在区间[a,b]上是连续分布的。但是很多情况下分布并不是连续的,所以取样采用的是概率分布的形式进行计算)。那么令 g ∗ ( x ) = g ( x ) f ( x ) g^\ast(x)=\frac{g(x)}{f(x)} g(x)=f(x)g(x),原有的积分可以写成如下形式:
I = ∫ a b g ∗ ( X ) f ( x ) d x I=\int_a^b g^\ast(X)f(x)dx I=abg(X)f(x)dx

由于采样时为离散的采样,所以计算公式可以修改为:
I = 1 N ∑ N g ∗ ( X ) f ( x ) I=\frac{1}{N}\sum_N g^\ast(X)f(x) I=N1Ng(X)f(x)
一般取f(x)为均匀分布,当区间在[a,b]时所得到的积分为
I = b − a N ∑ N g ∗ ( X ) I=\frac{b-a}{N}\sum_N g^\ast(X) I=NbaNg(X)

拒绝接受采样(Acceptance-Rejection Sampling)

根据上一节我们可以知道,当我们了解样本的分布为常用分布时,我们可以通过蒙特卡洛方法来求解区间的概率值。

但是当我们无法了解样本的分布的时候我们就难以使用蒙克卡罗方法了。对于概率分布不是常见的分布,一个可行的办法是采用接受-拒绝采样来得到该分布的样本。既然 p ( x ) p(x) p(x) 太复杂在程序中没法直接采样,那么我设定一个程序可采样的分布 q ( x ) q(x) q(x)比如高斯分布,然后按照一定的方法拒绝某些样本,以达到接近 p ( x ) p(x) p(x) 分布的目的,其中 q ( x ) q(x) q(x)叫做 proposal distribution(拒绝接收采样)

具体采用过程如下,设定一个方便采样的常用概率分布函数 q ( x ) q(x) q(x),以及一个常量 k k k,使得 p ( x ) p(x) p(x) 总在 k q ( x ) k_q(x) kq(x) 的下方。下图所示。
在这里插入图片描述首先,采样得到 q ( x ) q(x) q(x)的一个样本 z 0 z_0 z0,采样方法如第三节。然后,从均匀分布 ( 0 , k q ( z 0 ) ) (0,k_q(z_0)) (0,kq(z0))中采样得到一个值 u u u。如果 u u u落在了上图中的灰色区域,则拒绝这次抽样,否则接受这个样本 z 0 z_0 z0。重复以上过程得到n个接受的样本 z 0 , z 1 , . . . z n − 1 z_0,z_1,...z_{n−1} z0,z1,...zn1,则最后的蒙特卡罗方法求解结果为:

1 n ∑ i = 0 n − 1 f ( z i ) p ( z i ) \frac1n\sum_{i=0}^{n-1}\frac{f(z_i)}{p(z_i)} n1i=0n1p(zi)f(zi)

通过拒绝采样和蒙特卡洛方法我们依然不容易求出所需要的样本点。

  1. 对于一些二维分布p(x,y),有时候我们只能得到条件分布p(x|y)和p(y|x)和,却很难得到二维分布p(x,y)一般形式,这时我们无法用接受-拒绝采样得到其样本集。
  2. 对于一些高维的复杂非常见分布p(x1,x2,…,xn),我们要找到一个合适的q(x)和k非常困难。

通过上述可以看出,对于复杂概率模型,蒙特卡洛与拒绝采样方法难以应付,此时需要请出马尔科夫模型来对复杂概率模型进行样本的采样。

马尔科夫模型

马尔科夫性质

马尔可夫模型是一类随机过程,该过程有如下的特点:已知有几种状态,各种状态之间可以互相转化。 状态空间中经过从一个状态到另一个状态的转换的随机过程。该过程要求具备“无记忆”的性质:下一状态的概率分布只能由当前状态决定,在时间序列中它前面的事件均与之无关。 这种特定类型的“无记忆性”称作马尔可夫性质。马尔科夫链作为实际过程的统计模型具有许多应用。其公式如下:
在这里插入图片描述
模型在t+1时刻的状态条件概率仅仅依赖于t时刻。

马尔科夫模型

「马尔可夫模型」是指基于马尔可夫性质的模型,其假设一个给定过程的未来状态仅取决于当前状态。根据系统状态是否完全可被观测以及系统是自动的还是受控的,可以将常见的马尔可夫模型分成四种:马尔可夫链、隐马尔可夫模型(HMM)、马尔可夫决策过程(MDP)和部分可观测马尔可夫决策过程(POMDP),具体可见下表。
在这里插入图片描述最简单的马尔可夫模型是马尔可夫链。它用一个随时间变化的随机变量来模拟系统的状态。在这种情况下,马尔可夫性质表明,这个变量的分布只取决于之前状态的分布。当马尔可夫链的状态只能部分观察到,这就是隐马尔可夫模型。隐马尔可夫模型常用的用途是语音识别,它是大多数现代自动语音识别系统的基础。
建立一个马尔科夫模型如下所示:
在这里插入图片描述 其状态之间的转移我们用转移概率矩阵P表示:
在这里插入图片描述
  此时,我们给定一个初始状态,然后经过该状态转移矩阵的转换,最终会收敛到一个稳定的状态,具体如马尔科夫链定理所示
在这里插入图片描述当马尔科夫链达到平稳分布的情况下,我们得到的性质为:从某个时刻开始,每个状态之间发生的概率不再变化,我们成这个最后不在变化的概率为马尔科夫链的平稳分布。

由于马尔科夫链能收敛到平稳分布, 于是有了一个想法:如果我们能构造一个转移矩阵为P的马氏链,使得该马氏链的平稳分布恰好是p(x), 那么我们从任何一个初始状态x0出发沿着马氏链转移, 得到一个转移序列 x0, x1, x2,⋯xn, xn+1⋯, 如果马氏链在第n步已经收敛了,于是我们就得到了 π(x) 的样本xn, xn+1⋯(也就是从第n步收敛时开始,之后的x都服从同一个平稳分布,我们可以将这个分布设定为我们的目标采样分布)。

从上面可以看出马尔科夫链的平稳分布收敛主要依赖于状态转移矩阵,所以关键是如何构建状态转移矩阵,使得最终的平稳分布是我们所要的分布。想做到这一点主要依赖于细致平稳定理
在这里插入图片描述

马尔可夫采样

如果我们得到了某个平稳分布所对应的马尔科夫链状态转移矩阵,我们就很容易采用出这个平稳分布的样本集。

假设我们任意初始的概率分布是 π 0 ( x ) π_0(x) π0(x), 经过第一轮马尔科夫链状态转移后的概率分布是 π 1 ( x ) π_1(x) π1(x),。。。第i轮的概率分布是 π i ( x ) π_i(x) πi(x)。假设从第i轮开始再经过n轮后马尔科夫链收敛到我们的平稳分布 π ( x ) π(x) π(x),即:
π n ( x ) = π n + 1 ( x ) = π n + 2 ( x ) = . . . = π ( x ) π_n(x)=π_{n+1}(x)=π_{n+2}(x)=...=π(x) πn(x)=πn+1(x)=πn+2(x)=...=π(x)
对于每个分布 π i ( x ) π_i(x) πi(x),我们有:
π i ( x ) = π i − 1 ( x ) P = π i − 2 ( x ) P 2 = π 0 ( x ) P i π_i(x)=π_{i−1}(x)P=π_{i−2}(x)P^2=π_0(x)P^i πi(x)=πi1(x)P=πi2(x)P2=π0(x)Pi

现在我们可以开始采样了,首先,基于初始任意简单概率分布比如高斯分布 π 0 ( x ) π_0(x) π0(x)采样得到状态值 x 0 x_0 x0,基于条件概率分布 P ( x ∣ x 0 ) P(x|x_0) P(xx0)采样状态值 x 1 x_1 x1,一直进行下去,当状态转移进行到一定的次数时,比如到n次时,我们认为此时的采样集 ( x n , x n + 1 , x n + 2 , . . . ) (x_n,x_{n+1},x_{n+2},...) (xn,xn+1,xn+2,...)即是符合我们的平稳分布的对应样本集,可以用来做蒙特卡罗模拟求和了。

总结下基于马尔科夫链的采样过程:

1)输入马尔科夫链状态转移矩阵P,设定状态转移次数阈值n1,需要的样本个数n2

2)从任意简单概率分布采样得到初始状态值x0

3)for t=0 to n 1 + n 2 − 1 n_1+n_2−1 n1+n21: 从条件概率分布 P ( x ∣ x t ) P(x|x_t) P(xxt)中采样得到样本 x t + 1 x_{t+1} xt+1

样本集 ( x n 1 , x n 1 + 1 , . . . , x n 1 + n 2 − 1 ) (x_{n_1},x_{n_1+1},...,x_{n_1+n_2−1}) (xn1,xn1+1,...,xn1+n21)即为我们需要的平稳分布对应的样本集。

小结:
如果假定我们可以得到我们需要采样样本的平稳分布所对应的马尔科夫链状态转移矩阵,那么我们就可以用马尔科夫链采样得到我们需要的样本集,进而进行蒙特卡罗模拟。但是一个重要的问题是,随意给定一个平稳分布π,如何得到它所对应的马尔科夫链状态转移矩阵P呢? 这是个大问题。我们绕了一圈似乎还是没有解决任意概率分布采样样本集的问题。

MCMC采样

通过上一节的描述,给定一个概率平稳分布π, 很难直接找到对应的马尔科夫链状态转移矩阵P。而只要解决这个问题,我们就可以找到一种通用的概率分布采样方法,进而用于蒙特卡罗模拟。解决这个问题的办法:MCMC采样和它的易用版M-H采样,以及Gibbs采样

对于无法直接获得的状态转移矩阵P,可以借用拒绝采样的思想进行计算。

由于一般情况下,目标平稳分布π(x)和某一个马尔科夫链状态转移矩阵Q不满足细致平稳条件,即:
π ( i ) Q ( i , j ) ≠ π ( j ) Q ( j , i ) \pi(i)Q(i,j)\neq\pi(j)Q(j,i) π(i)Q(i,j)=π(j)Q(j,i)
我们可以对上式做一个改造,使细致平稳条件成立。方法是引入一个 α ( i , j ) α(i,j) α(i,j),使上式可以取等号,即:
π ( i ) Q ( i , j ) α ( i , j ) = π ( j ) Q ( j , i ) α ( j , i ) \pi(i)Q(i,j)\alpha(i,j)=\pi(j)Q(j,i)\alpha(j,i) π(i)Q(i,j)α(i,j)=π(j)Q(j,i)α(j,i)

问题是什么样的 α ( i , j ) α(i,j) α(i,j)可以使等式成立呢?其实很简单,只要满足下两式即可:
α ( i , j ) = π ( j ) Q ( j , i ) \alpha(i,j)=\pi(j)Q(j,i) α(i,j)=π(j)Q(j,i)
α ( j , i ) = π ( i ) Q ( i , j ) \alpha(j,i)=\pi(i)Q(i,j) α(j,i)=π(i)Q(i,j)
这样,我们就得到了我们的分布 π ( x ) π(x) π(x)对应的马尔科夫链状态转移矩阵P,满足:
P ( i , j ) = Q ( i , j ) α ( i , j ) P(i,j)=Q(i,j)\alpha(i,j) P(i,j)=Q(i,j)α(i,j)
也就是说,我们的目标矩阵P可以通过任意一个马尔科夫链状态转移矩阵Q乘以 α ( i , j ) α(i,j) α(i,j)得到。 α ( i , j ) α(i,j) α(i,j)我们有一般称之为接受率。取值在[0,1]之间,可以理解为一个概率值。即目标矩阵P可以通过任意一个马尔科夫链状态转移矩阵Q以一定的接受率获得。这个很像我们在蒙特卡罗方法中讲到的接受-拒绝采样,那里是以一个常用分布通过一定的接受-拒绝概率得到一个非常见分布,这里是以一个常见的马尔科夫链状态转移矩阵Q通过一定的接受-拒绝概率得到目标转移矩阵P,两者的解决问题思路是类似的。

总结下MCMC的采样过程。

1)输入我们任意选定的马尔科夫链状态转移矩阵Q,平稳分布π(x),设定状态转移次数阈值 n 1 n_1 n1,需要的样本个数 n 2 n_2 n2

2)从任意简单概率分布采样得到初始状态值 x 0 x_0 x0

3)for t=0 to n 1 + n 2 − 1 : n_{1}+n_2−1: n1+n21:
  a) 从条件概率分布 Q ( x ∣ x t ) Q(x|x_t) Q(xxt)中采样得到样本 x ∗ x_∗ x
  
  b) 从均匀分布采样 u ∼ u n i f o r m [ 0 , 1 ] u∼uniform[0,1] uuniform[0,1]
  
  c) 如果 u < α ( x t , x ∗ ) = π ( x ∗ ) Q ( x ∗ , x t ) u<α(x_t,x_∗)=π(x_∗)Q(x_∗,x_t) u<α(xt,x)=π(x)Q(x,xt), 则接受转移 x t → x ∗ xt→x∗ xtx,即 x t + 1 = x ∗ x_t+1=x_∗ xt+1=x
  
  d) 否则不接受转移,即xt+1=xt即为我们需要的平稳分布对应的样本集。

上面这个过程基本上就是MCMC采样的完整采样理论了,但是这个采样算法还是比较难在实际中应用,为什么呢?问题在上面第三步的c步骤,接受率这儿。由于 α ( x t , x ∗ ) α(x_t,x_∗) α(xt,x)可能非常的小,比如0.1,导致我们大部分的采样值都被拒绝转移,采样效率很低。有可能我们采样了上百万次马尔可夫链还没有收敛,也就是上面这个 n 1 n_1 n1要非常非常的大,这让人难以接受,怎么办呢?这时就轮到我们的M-H采样出场了。

MCMC(马尔科夫—蒙特卡洛算法)中的代表——MH算法

M-H采样是Metropolis-Hastings采样的简称,这个算法首先由Metropolis提出,被Hastings改进,因此被称之为Metropolis-Hastings采样或M-H采样。

M-H采样解决了我们上一节MCMC采样接受率过低的问题。
我们回到MCMC采样的细致平稳条件:
π ( i ) Q ( i , j ) α ( i , j ) = π ( j ) Q ( j , i ) α ( j , i ) π(i)Q(i,j)α(i,j)=π(j)Q(j,i)α(j,i) π(i)Q(i,j)α(i,j)=π(j)Q(j,i)α(j,i)
我们采样效率低的原因是 α ( i , j ) α(i,j) α(i,j)太小了,比如为0.1,而 α ( j , i ) α(j,i) α(j,i)为0.2。即:
π ( i ) Q ( i , j ) × 0.1 = π ( j ) Q ( j , i ) × 0.2 π(i)Q(i,j)×0.1=π(j)Q(j,i)×0.2 π(i)Q(i,j)×0.1=π(j)Q(j,i)×0.2
这时我们可以看到,如果两边同时扩大五倍,接受率提高到了0.5,但是细致平稳条件却仍然是满足的,即:
π ( i ) Q ( i , j ) × 0.5 = π ( j ) Q ( j , i ) × 1 π(i)Q(i,j)×0.5=π(j)Q(j,i)×1 π(i)Q(i,j)×0.5=π(j)Q(j,i)×1
这样我们的接受率可以做如下改进,即:
α ( i , j ) = m i n { π ( j ) Q ( j , i ) π ( i ) Q ( i , j ) , 1 } \alpha(i,j)=min\{\frac{\pi(j)Q(j,i)}{\pi(i)Q(i,j)},1\} α(i,j)=min{π(i)Q(i,j)π(j)Q(j,i),1}

通过这个微小的改造,我们就得到了可以在实际应用中使用的M-H采样算法过程如下:
1)输入我们任意选定的马尔科夫链状态转移矩阵Q,平稳分布π(x),设定状态转移次数阈值 n 1 n_1 n1,需要的样本个数 n 2 n_2 n2
2)从任意简单概率分布采样得到初始状态值 x 0 x_0 x0
3)for t=0 to n 1 + n 2 − 1 n_1+n_2−1 n1+n21:
  a) 从条件概率分布 Q ( x ∣ x t ) Q(x|x_t) Q(xxt)中采样得到样本 x ∗ x_∗ x
  b) 从均匀分布采样 u ∼ u n i f o r m [ 0 , 1 ] u∼uniform[0,1] uuniform[0,1]
  c) 如果 u < α ( x t , x ∗ ) = m i n π ( j ) Q ( j , i ) π ( i ) Q ( i , j ) , 1 u<α(x_t,x_∗)=min{π(j)Q(j,i)π(i)Q(i,j),1} u<α(xt,x)=minπ(j)Q(j,i)π(i)Q(i,j),1, 则接受转移 x t → x ∗ x_t→x_∗ xtx,即 x t + 1 = x ∗ x_{t+1}=x_∗ xt+1=x
  d) 否则不接受转移,即 x t + 1 = x t x_{t+1}=x_t xt+1=xt 样本集 ( x n 1 , x n 1 + 1 , . . . , x n 1 + n 2 − 1 ) (x_{n1},x_{n_1+1},...,x_{n_1+n_2−1}) (xn1,xn1+1,...,xn1+n21)

即为我们需要的平稳分布对应的样本集。很多时候,我们选择的马尔科夫链状态转移矩阵Q。如果是对称的,即满足 Q ( i , j ) = Q ( j , i ) Q(i,j)=Q(j,i) Q(i,j)=Q(j,i),这时我们的接受率可以进一步简化为:
α ( i , j ) = m i n { π ( j ) π ( i ) , 1 } \alpha(i,j)=min\{\frac{\pi(j)}{\pi(i)},1\} α(i,j)=min{π(i)π(j),1}

M-H采样完整解决了使用蒙特卡罗方法需要的任意概率分布样本集的问题,因此在实际生产环境得到了广泛的应用。但是在大数据时代,M-H采样面临着两大难题:

  1. 我们的数据特征非常的多,M-H采样由于接受率计算式 π ( j ) Q ( j , i ) π ( i ) Q ( i , j ) \frac{\pi(j)Q(j,i)}{\pi(i)Q(i,j)} π(i)Q(i,j)π(j)Q(j,i)的存在,在高维时需要的计算时间非常的可观,算法效率很低。同时α(i,j)一般小于1,有时候辛苦计算出来却被拒绝了。能不能做到不拒绝转移呢?

  2. 由于特征维度大,很多时候我们甚至很难求出目标的各特征维度联合分布,但是可以方便求出各个特征之间的条件概率分布。这时候我们能不能只有各维度之间条件概率分布的情况下方便的采样呢?

Gibbs采样解决了上面两个问题,因此在大数据时代,MCMC采样基本是Gibbs采样的天下,接下来就来讨论Gibbs采样。

Gibbs算法

M-H采样已经可以很好的解决蒙特卡罗方法需要的任意概率分布的样本集的问题。但是M-H采样有两个缺点:一是需要计算接受率,在高维时计算量大。并且由于接受率的原因导致算法收敛时间变长。二是有些高维数据,特征的条件概率分布好求,但是特征的联合分布不好求。因此需要一个好的方法来改进M-H采样,这就是我们下面讲到的Gibbs采样。

我们讲到了细致平稳条件:如果非周期马尔科夫链的状态转移矩阵P和概率分布π(x)对于所有的i,j满足:
π ( i ) P ( i , j ) = π ( j ) P ( j , i ) \pi(i)P(i,j)=\pi(j)P(j,i) π(i)P(i,j)=π(j)P(j,i)
则称概率分布π(x)是状态转移矩阵P的平稳分布。在M-H采样中我们通过引入接受率使细致平稳条件满足。 现在我们换一个思路。

从二维的数据分布开始,假设π(x1,x2)是一个二维联合数据分布,观察第一个特征维度相同的两个点 A ( x , y 1 ) A(x,y_1) A(x,y1) B ( x , y 2 ) B(x,y_2) B(x,y2),容易发现下面两式成立:
π ( x 1 , y 1 ) π ( y 2 ∣ x 1 ) = π ( x 1 ) π ( y 1 ∣ x 1 ) π ( y 2 ∣ x 1 ) \pi(x_1,y_1)\pi(y_2\vert x_1)=\pi(x_1)\pi(y_1\vert x_1)\pi(y_2\vert x_1) π(x1,y1)π(y2x1)=π(x1)π(y1x1)π(y2x1)

π ( x 1 , y 2 ) π ( y 1 ∣ x 1 ) = π ( x 1 ) π ( y 2 ∣ x 1 ) π ( y 1 ∣ x 1 ) \pi(x_1,y_2)\pi(y_1\vert x_1)=\pi(x_1)\pi(y_2\vert x_1)\pi(y_1\vert x_1) π(x1,y2)π(y1x1)=π(x1)π(y2x1)π(y1x1)

因此从上式可以推断出模型如下:
π ( x 1 , y 1 ) π ( y 2 ∣ x 1 ) = π ( x 1 , y 2 ) π ( y 1 ∣ x 1 ) \pi(x_1,y_1)\pi(y_2\vert x_1)=\pi(x_1,y_2)\pi(y_1\vert x_1) π(x1,y1)π(y2x1)=π(x1,y2)π(y1x1)

π ( A ) π ( y 2 ∣ x ) = π ( B ) π ( y 1 ∣ x ) \pi(A)\pi(y_2\vert x)=\pi(B)\pi(y_1\vert x) π(A)π(y2x)=π(B)π(y1x)

观察上式再观察细致平稳条件的公式,我们发现在 x = x 1 x=x_1 x=x1这条直线上,如果用条件概率分布 π ( y ∣ x 1 ) π(y|x_1) π(yx1)作为马尔科夫链的状态转移概率,则任意两个点之间的转移满足细致平稳条件!这真是一个开心的发现,同样的道理,在 y = y 2 y=y_2 y=y2这条直线上,如果用条件概率分布 π ( x 1 ∣ y 2 ) π(x_1|y_2) π(x1y2)作为马尔科夫链的状态转移概率,则任意两个点之间的转移也满足细致平稳条件。那是因为假如有一点 C ( x 2 , y 1 ) C(x_2,y_1) C(x2,y1),我们可以得到:

π ( A ) π ( x 2 ∣ y 1 ) = π ( C ) π ( x 1 ∣ y 1 ) \pi(A)\pi(x_2\vert y_1)=\pi(C)\pi(x_1\vert y_1) π(A)π(x2y1)=π(C)π(x1y1)

基于上面的发现,我们可以这样构造分布π(x1,x2)的马尔可夫链对应的状态转移矩阵P:
P ( A → B ) = π ( x 2 ( B ) ∣ x 1 ( 1 ) )    i f   x 1 ( A ) = x 1 ( B ) = x 1 ( 1 ) P(A\rightarrow B)=\pi(x_2^{(B)}\vert x_1^{(1)})\text{ }\text{ }if\text{ }x_1^{(A)}=x_1^{(B)}=x_1^{(1)} P(AB)=π(x2(B)x1(1))ifx1(A)=x1(B)=x1(1)

P ( A → C ) = π ( x 1 ( C ) ∣ x 2 ( 1 ) )    i f   x 2 ( A ) = x 2 ( C ) = x 2 ( 1 ) P(A\rightarrow C)=\pi(x_1^{(C)}\vert x_2^{(1)})\text{ }\text{ }if\text{ }x_2^{(A)}=x_2^{(C)}=x_2^{(1)} P(AC)=π(x1(C)x2(1))ifx2(A)=x2(C)=x2(1)

二维Gibbs采样

利用上一节找到的状态转移矩阵,我们就得到了二维Gibbs采样,这个采样需要两个维度之间的条件概率。具体过程如下:
1)输入平稳分布 π ( x 1 , x 2 ) π(x_1,x_2) π(x1,x2),设定状态转移次数阈值 n 1 n_1 n1,需要的样本个数 n 2 n_2 n2
2)随机初始化初始状态值 x 1 ( 0 ) x_1^{(0)} x1(0) x 2 ( 0 ) x_2^{(0)} x2(0)
3)for t=0 to n 1 + n 2 − 1 n_1+n_2−1 n1+n21 :
  a) 从条件概率分布 P ( x 2 ∣ x 1 ( t ) ) P(x_2|x_1^{(t)}) P(x2x1(t))中采样得到样本 x 2 ( t + 1 ) x_2^{(t+1)} x2(t+1)
  b) 从条件概率分布 P ( x 1 ∣ x 2 ( t + 1 ) ) P(x_1|x_2^{(t+1)}) P(x1x2(t+1))中采样得到样本 x 1 ( t + 1 ) x_1^{(t+1)} x1(t+1)
样本集 { ( x 1 ( n 1 ) , x 2 ( n 1 ) ) , ( x 1 ( n 1 + 1 ) , x 2 ( n 1 + 1 ) ) , . . . , ( x 1 ( n 1 + n 2 − 1 ) , x 2 ( n 1 + n 2 − 1 ) ) } \{(x_1^{(n_1)},x_2^{(n_1)}),(x_1^{(n_1+1)},x_2^{(n_1+1)}),...,(x_1^{(n_1+n_2-1)},x_2^{(n_1+n_2-1)})\} {(x1(n1),x2(n1)),(x1(n1+1),x2(n1+1)),...,(x1(n1+n21),x2(n1+n21))}即为我们需要的平稳分布对应的样本集。
整个采样过程中,我们通过轮换坐标轴,采样的过程为:
( x 1 ( 1 ) , x 2 ( 1 ) ) → ( x 1 ( 1 ) , x 2 ( 2 ) ) → ( x 1 ( 2 ) , x 2 ( 2 ) ) → . . . → ( x 1 ( n 1 + n 2 − 1 ) , x 2 ( n 1 + n 2 − 1 ) ) (x_1^{(1)},x_2^{(1)})\rightarrow(x_1^{(1)},x_2^{(2)})\rightarrow(x_1^{(2)},x_2^{(2)})\rightarrow...\rightarrow(x_1^{(n_1+n_2-1)},x_2^{(n_1+n_2-1)}) (x1(1),x2(1))(x1(1),x2(2))(x1(2),x2(2))...(x1(n1+n21),x2(n1+n21))

用下图可以很直观的看出,采样是在两个坐标轴上不停的轮换的。当然,坐标轴轮换不是必须的,我们也可以每次随机选择一个坐标轴进行采样。不过常用的Gibbs采样的实现都是基于坐标轴轮换的。

在这里插入图片描述

多维Gibbs采样

上面的这个算法推广到多维的时候也是成立的。比如一个n维的概率分布 π ( x 1 , x 2 , . . . x n ) π(x_1,x_2,...x_n) π(x1,x2,...xn),我们可以通过在n个坐标轴上轮换采样,来得到新的样本。对于轮换到的任意一个坐标轴 x i x_i xi上的转移,马尔科夫链的状态转移概率为 P ( x i ∣ x 1 , x 2 , . . . , x i − 1 , x i + 1 , . . . , x n ) P(x_i|x_1,x_2,...,x_{i−1},x_{i+1},...,x_n) P(xix1,x2,...,xi1,xi+1,...,xn)即固定n−1个坐标轴,在某一个坐标轴上移动。

文本建模

一篇文档,可以看成是一组有序的词的序列 d = d= d= . 从统计学角度来看,文档的生成可以看成是上帝抛掷骰子生成的结果,每一次抛掷骰子都生成一个词汇,抛掷N词生成一篇文档。在统计文本建模中,我们希望猜测出上帝是如何玩这个游戏的,这会涉及到两个最核心的问题:

1. 上帝都有什么样的骰子;
2. 上帝是如何抛掷这些骰子的;

第一个问题就是表示模型中都有哪些参数,骰子的每一个面的概率都对应于模型中的参数;第二个问题就表示游戏规则是什么,上帝可能有各种不同类型的骰子,上帝可以按照一定的规则抛掷这些骰子从而产生词序列。

PLSA模型

PLSA模型在单词—>文章的骰子之前增加了一个隐变量,文章分类。
在PLSA看来,在写一篇文章的时候分几个部分。一共投两个骰子,分别决定从doc—topic和topic—word。
其文章生成的过程如下:

  1. 现有两种类型的骰子,一种是doc-topic骰子,每个doc-topic骰子有K个面,每个面一个topic的编号;一种是topic-word骰子,每个topic-word骰子有V个面,每个面对应一个词;
  2. 现有K个topic-word骰子,每个骰子有一个编号,编号从1到K;
  3. 生成每篇文档之前,先为这篇文章制造一个特定的doc-topic骰子,重复如下过程生成文档中的词:
    3.1. 投掷这个doc-topic骰子,得到一个topic编号z;
    3.2 选择K个topic-word骰子中编号为z的那个,投掷这个骰子,得到一个词;

由于文档之间相互独立,很容易写出整个语料的生成概率。求解PLSA 可以使用著名的 EM 算法进行求得局部最优解。

在这里插入图片描述

LDA模型

对于一篇文档d中的每一个单词,LDA根据先验知识 确定某篇文档的主题分布θ,然后从该文档所对应的多项分布(主题分布)θ中抽取一个主题z,接着根据先验知识 确定当前主题的词语分布ϕ,然后从主题z所对应的多项分布(词分布)ϕ中抽取一个单词w。然后将这个过程重复N次,就产生了文档d。

与PLSA不相同的一点是,LDA模型采用的是贝叶斯概率进行分析。为主题分布和词分布分别加了两个 Dirichlet 先验分布。两个Dirichlet分布的参数分别为 α → \overrightarrow\alpha α β → \overrightarrow\beta β

  • 在pLSA中,我们使用EM算法去估计“主题-词项”矩阵Φ(由转换得到)和“文档-主题”矩阵Θ(由转换得到)这两个参数,而且这两参数都是个固定的值,只是未知,使用的思想其实就是极大似然估计MLE。
  • 而在LDA中,估计Φ、Θ这两未知参数可以用变分(Variational inference)-EM算法,也可以用gibbs采样,前者的思想是最大后验估计MAP(MAP与MLE类似,都把未知参数当作固定的值),后者的思想是贝叶斯估计。贝叶斯估计是对MAP的扩展,但它与MAP有着本质的不同,即贝叶斯估计把待估计的参数看作是服从某种先验分布的随机变量。

由于LDA把要估计的主题分布和词分布看作是其先验分布是Dirichlet分布的随机变量,所以, 在LDA这个估计主题分布、词分布的过程中,它们的先验分布(即Dirichlet分布)事先由人为给定,那么LDA就是要去求它们的后验分布(LDA中可用gibbs采样去求解它们的后验分布,得到期望)!
pLSA模型的示意图:
在这里插入图片描述
LDA模型的示意图: 其中Φ表示词分布,Θ表示主题分布, 是主题分布Θ的先验分布(即Dirichlet 分布)的参数, 是词分布Φ的先验分布(即Dirichlet 分布)的参数,N表示文档的单词总数,M表示文档的总数
在这里插入图片描述

LDA模型推导过程

给定一个文档集合,w是可以观察到的已知变量,和是根据经验给定的先验参数,其他的变量z,θ和φ都是未知的隐含变量。Φ表示词分布,Θ表示主题分布, 是主题分布Θ的先验分布(即Dirichlet 分布)的参数, 是词分布Φ的先验分布(即Dirichlet 分布)的参数。根据LDA的图模型,可以写出所有变量的联合分布:
在这里插入图片描述
因为产生主题分布θ,主题分布θ确定具体主题,且产生词分布φ、词分布φ确定具体词,所以上述式子等价于下述式子所表达的联合概率分布 p ( w → , z → ) p(\overrightarrow w,\overrightarrow z) p(w ,z )(等于说我们最后的结果是要求出文章中的词与主题的联合分布概率):

在这里插入图片描述

其中,第一项因子表示的是根据确定的主题和词分布的先验分布参数采样词的过程,第二项因子是根据主题分布的先验分布参数采样主题的过程,这两项因子是需要计算的两个未知参数。由于这两个过程是独立的,所以下面可以分别处理,各个击破。
最后得到的结果是:
在这里插入图片描述
有了联合分布 p ( w → , z → ) p(\overrightarrow w,\overrightarrow z) p(w ,z ),便可以通过联合分布来计算在给定可观测变量 w 下的隐变量 z 的条件分布(后验分布) p ( z → ∣ w → ) p(\overrightarrow z\vert\overrightarrow w) p(z w )来进行贝叶斯分析。
换言之,有了这个联合分布后,要求解第m篇文档中的第n个词(下标为 的词)的全部条件概率就好求了。
先定义几个变量。 ¬ i \neg i ¬i表示除去 i i i 的词, w → = { w i = t , w ¬ i → } \overrightarrow w=\left\{w_i=t,\overrightarrow{w_{\neg i}}\right\} w ={wi=t,w¬i } z → = { z i = k , z ¬ i → } \overrightarrow z=\left\{z_i=k,\overrightarrow{z_{\neg i}}\right\} z ={zi=k,z¬i }
然后,排除当前词的主题分配,即根据其他词的主题分配和观察到的单词来计算当前词主题的概率公式为:

$$$$

最终求解的Dirichlet 分布期望为:
在这里插入图片描述

也就是上式的解为:
在这里插入图片描述
仔细观察上述结果,可以发现,式子的右半部分便是 p ( t o p i c ∣ d o c ) ∗ p ( t o p i c ∣ d o c ) p(topic \vert doc) * p(topic \vert doc) p(topicdoc)p(topicdoc),这个概率的值对应着 d o c → t o p i c → w o r d doc\rightarrow topic\rightarrow word doctopicword的路径概率。如此,K 个topic 对应着K条路径,Gibbs Sampling 便在这K 条路径中进行采样,如下图所示:

在这里插入图片描述

  • 0
    点赞
  • 5
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

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

余额充值