概率图模型课堂笔记：2.4 取样方法
从前一章可知，用精确推理来获得系统的整体或局部概率分布（比如某些选定变量的边缘分布，比如某种条件分布），往往计算量过于庞大，而采样往往是一种代价更低的近似方法。
1.1 简单取样
1. 用独立等概率取样的方式取M个样本D，$E_p[f]\approx \frac1M \sum_{m=1}^{M}f(x[m])$
2. 对离散分布生成采样的方法，设$X\in {x_1,x_2,\cdots,x_M}, P(X=X_i)=\theta_i$, 把$[0,1]$分成M个区间, $I_1=[0,\theta_1), I_2=[\theta_1, \theta_1+\theta_2),\cdots$，生成$[0,1]$随机数，落在$I_k$表示抽样$x_k$
3. 大数定律两个Bound
(1) 取样落在$[p-\epsilon, p+\epsilon]$之外的概率小于等于$2e^{-2M\epsilon^2}$
若要取样落在$[p-\epsilon, p+\epsilon]$之外的概率小于$\delta$，则$M\ge=\frac{ln(2/\delta)}{2\epsilon^2}$
(2) 取样落在$[p(1-\epsilon), p(1+\epsilon)]$之外的概率小于等于$2e^{-2Mp\epsilon^2/3}$
若要取样落在$[p(1-\epsilon), p(1+\epsilon)]$之外的概率小于$\delta$，则$M\ge=\frac{3ln(2/\delta)}{p\epsilon^2}$
所以对于p很小的情况，都会造成需要大量样本。
4. 采样生成方法：对于贝叶斯图来说，生成样本自顶向下就可以了。
(1)对root结点$X$，直接根据$X$的概率分布$P(X)$来生成
(2)对于非root结点$X$，因为其父亲结点已经生成，所以在条件概率$P(X|Pa_X)$中选择相应的行来生成。
(3)所有变量生成完成后，一个样本就生成出来了。
(4)重复(1)-(3), 生成M个样本
(5)要估算$P(Y=y)$，用所有样本中$Y=y$样本个数除以M。注意M要满足大数定律的Bound
(6)要估算$P(Y=y|E=e)$，方法和(5)一样，遇到$E\ne e$的样本，直接抛弃。这样生成的样本都是$E=e$的。如果$E=e$的概率很小，那么生成会很费时间。因为最后能保留的样本只有$MP(e)$个。
5. 对马尔科夫图来说，4这个方法不能用。因为没有根节点。
1.2 马尔可夫链蒙特卡洛
思想是：$\mathbf X$每个可能的取值$\mathbf x$都是一个状态，状态之间以一定的概率转移。一开始每个状态的分布率（可以认为是一个概率向量）目标概率P很可能是不相等的，但是粒子在马尔科夫网络中进行不断地转移，概率会逼近目标概率，达到稳态(stationary)，在此之后，再取样，就是符合P的分布了。
1. 马尔可夫链定义了：
(1)每个变量的初始概率$P^0$
(2)状态之间互相转换的概率$T(x\to x')$。且$\sum_{x'}T(x\to x')=1$
2. 转移等式(Temporal Dynamics)：$P^{(t+1)}(X^{(t+1)}=x')=\sum_P^{(t)}(X^{(t)}=x)T(x\to x')$
3. 马尔可夫链收敛的条件：对任意的$x'$，有：$P^{(t+1)}(x')=P^{(t)}(x')=\sum_xP^{(t)}(x)T(x\to x')$，达到这个状态后的概率我们称为收敛概率。
所以，如果要证明构造出来的某个马尔科夫链的收敛概率为P，只要把P代入上面的式子，证明$P(x')=\sum_xP(x)T(x\to x')$
4. 正则马尔可夫链
(1)如果存在一个$k$（先选定$k$），通过$k$步，能从任意状态到另一任意状态，我们称这个马尔可夫链是正则(regular)的
(2)一个正则的马尔可夫链一定可以收敛
(3)正则的一个充分条件是：任意两个状态之间都有路径，且任意一个状态都有回路。
1.3马尔科夫链的使用
1. 假设已经针对某分布$P$（如某马尔科夫网表示的分布）生成了一个收敛概率为$P$的马尔科夫链
2. 通过一个初始分布$P^{(0)}$生成C个样本$x^{(c,0)}$
3. 对每个样本$x^{(c,0)}$，同步生成各自的样本序列。直到认为已经收敛。判断收敛的方法：
(1)针对同一个样本序列，在不同时间段的窗口中，得出来的分布接近一致的
(2)针对不同的样本序列，在同一个窗口中，它们的分布已经接近一致。
这个通常更好，因为有可能状态集分为两个区域，跨区的转移可能性很小。
用多图来描述，C个序列在同一个评价概率（比如选定一个$P(x_i)$）上趋向相等，可以认为收敛
用单图描述，建立C维坐标，如果发现所有序列的评价概率收拢在对角线上，可以认为收敛
4. 收敛之后，继续各个样本序列的生成，再生成M个样本（M的大小参照前一节）。拿这M个样本计算概率
1.4 吉布斯采样
1. 目标分布：$P_{\phi}(X_1,X_2,\cdots,X_n)$。步骤：
(1)初始状态：$\mathbf{x}=\{x_1,x_2,...,x_n\}$
(2)单次迭代过程：for i=1 to n: Sample $x_i$ by $P_{\phi}(X_i|x_{-i})$
2. 这个过程实际上是运行在一个马尔可夫链上，链中各个状态就是目标分布中的$\mathbf{X}_i$各个可能的取值，可转移的两个状态之间只有一个变量$x_i$不同，转移概率是$P_{\phi}(X_i|x_{-i})$。当然，每个状态向自己也是可以转移的。这个就符合正则马尔科夫链的两个充分条件。
同时，可以证明这个网络的收敛概率是$P_{\phi}$（？？？）
3. 这个过程中取样的代价也很低，为$P_{\phi}(X_i|x_{-i})=\frac{P_{\phi}(X_i,x_{-i})}{P_{\phi}(x_{-i})}\propto \prod_{j:X_i\in Scope[C_j]}\phi_j(X_i,x_{j,-i})$，也就是不用计算整体的联合分布，只要把和$X_i$有关的factor乘起来再代入其他变量。
4. 对于ＸＯＲ吉布斯很难收敛，因为给定$Y$的情况下，${x_1,x_2}$从一个状态跳到另外一个状态，不可能像迭代中描述的那样，只改变其中一个变量的值。
5. 对于图分块来说，也很难，因为很多块都是和相邻块有强相关性，比如某片区域全是水块或者全是天空块的可能性都不小，但如果只改变其中一个，就会导致一块水被天包围，这个联合概率很小，很难取样到。因此一次只改变其中一个变量的值，来把一整片水块变成天空块，也是接近不可能。
1.5
1. 可逆链：对于任何一对状态来说，如果互相迁移的“量”是相等的$P(x)T(x\to x')=P(x')T(x'\to x)$，则称为可逆链。因为“往来”等量，所以可逆链一定已经收敛，用数学描述：两边对$x$进行求和，就得到了收敛式。
但是收敛链不一定可逆。考虑一个$A\to B\to C\to A$的顺向转移概率为0.8，反向为0.2的环。
2. 令$\pi(\mathbf X)$表示收敛概率，可逆状态用$\pi(x)T(x\to x')=\pi(x')T(x'\to x)$表示。
我们把$T(x\to x')$拆为两部分：$Q(x\to x')$和$A(x\to x')$。$Q$表示提议概率，$A$表示接受度。系统以$Q(x\to x')$的概率提出转移，被接受的概率有$A(x\to x')$，如果被接受，则转移，否则，留在原地。
(1)若$x\ne x'$，则$T(x\to x')=Q(x\to x')A(x\to x')$
(2)若$x=x'$，则$T(x\to x)=Q(x\to x)+\sum_{y\ne x}Q(x\to y)(1-A(x\to y))$
3. 现在考虑$x\ne x'$的情形。对可逆链来说，平衡时我们会达到$\pi(x)Q(x\to x')A(x\to x')=\pi(x')Q(x'\to x)A(x'\to x)$，也就是
$\frac{A(x\to x')}{A(x'\to x)}=\frac{\pi(x')Q(x'\to x)}{\pi(x)Q(x\to x')}$
现在假设$A(x\to x')<A(x'\to x)$，$\frac{A(x\to x')}{A(x'\to x)}=\rho<1$。
那么对$A(x\to x')$和$A(x'\to x)$来说，只要满足比例$\rho$的任意取值对都是合法的。我们令$A(x\to x')=\rho,A(x'\to x)=1$
对于这种取值方式，我们可以扩展到一般的情况，定义$A(x\to x')=\min[1, \frac{\pi(x')Q(x'\to x)}{\pi(x)Q(x\to x')}]$
4. Q的取值准则
(1) Q不能为0，否则计算A会除0
(2) Q必须广泛分布，这样才能远程跳转。但是对一些不太可能的跳转，A通常会低
5. 考虑一个最小距离和的二部图匹配问题，每个左侧点$X_i$都要指派唯一的右侧点$v_i=Y_j$，
(1)如果指派合法，系统概率是$e^(-\sum_i dist(i,v_i))$，即总距离越大，概率越小。
(2)如果指派不合法，概率为0
一个状态转移算法是：
(1)为$X_i$随机指派一个$v_i$
(2)如果$v_i$原来指派给了$X_j$，则原$X_j$执行(1)，否则完成。
对比：
如果用上一节的吉布斯采样，很难收敛。
2.1 时间模型的推导
1. 预备知识：
(1)贝叶斯定理：$P(A|B)=\frac{P(B|A)P(A)}{P(B)}$
(2)条件贝叶斯定理：$P(A|B,C)=\frac{P(B|A,C)P(A|C)}{P(B|C)}$
条件贝叶斯定理就是在(1)的每个式子的竖线右边都加上C。
2. 对时间模型来说，我们希望做推断的时候不需要追溯太久之前的状态和数据，否则要维护巨大的状态空间。所幸的是这是可以做到的。考虑一个通用的时间状态模型。
(1) 标准状态方程：$\sigma^{(t)}(S^{(t)})=P(S^{(t)}|o^{(1:t)})$
(2) 现在要基于$\sigma^{(t)}$计算$\sigma^{(t+1)}$
(3) $\sigma^{(t+1)}=P(S^{(t+1)}|o^{(1:t)},o^{(t+1)})$
$=\frac{P(o^{(t+1)}|S^{(t+1)})P(S^{(t+1)}|o^{(1:t)})}{P(o^{(t+1)}|o^{(1:t)})}$ // 条件贝叶斯定理
$=\frac{P(o^{(t+1)}|S^{(t+1)})\sum_{s^{(t)}}P(S^{(t+1)}|s^{(t)})P(s^{(t)}|o^{(1:t)})}{P(o^{(t+1)}|o^{(1:t)})}$
$=\frac{P(o^{(t+1)}|S^{(t+1)})\sum_{s^{(t)}}P(S^{(t+1)}|s^{(t)})\sigma^{(t)}(s^{(t)})}{P(o^{(t+1)}|o^{(1:t)})}$
其中分母只是用来normalize的。
看到所有运算步骤只是用到前一个状态，说明可以用较少的状态空间来实现推断。
3. 观察机器人的状态模型。以及这个思想的应用
4. 对时间模型来说，在同一个时间t里，不同的状态变量并不是独立的，通过追溯到t,t-1或更早，往往能找到它们之间的active trail，所以往往只能对时间t所有的状态变量做一个完全的联合分布。状态空间随着状态数量指数增长。因此精确推导的计算量还是很大。
5. 对碟子模型来说，也是一样，同一个碟子里不同个体之间也是相关的，往往只能把它们当作一个整体作联合分布。
6. 考虑一个网页分类的模型，传统机器学习算法一般是根据网页本身内容来，但是用概率图来说，可以定义相互连接的网页之间类别的factor。比如学生和院系之间的连接更多，而学生之间的连接则很少。运用这个可以显著提高性能。
