随机采样和随机模拟：吉布斯采样Gibbs Sampling实现高斯分布参数推断

怼烎@ 2022-08-22 14:20 293阅读 0赞

http://[blog.csdn.net/pipisorry/article/details/51539739][blog.csdn.net_pipisorry_article_details_51539739]

## 吉布斯采样的实现问题    ##

本文主要说明如何通过吉布斯采样来采样截断多维高斯分布的参数（已知一堆截断高斯分布的数据，推断其参数( μ , Σ )）。

关于吉布斯采样的介绍文章都停止在吉布斯采样的详细描述上，如[随机采样和随机模拟：吉布斯采样Gibbs Sampling][Gibbs Sampling](why)但并没有说明吉布斯采样到底如何实现的(how)？

也就是具体怎么实现从下面这个公式采样？

![Center][]

下面介绍如何为多维正态分布构建一个吉布斯采样器，来采样截断多维高斯分布的参数。

[皮皮blog][blog]

# 模型的表示 #

**思路**：要求吉布斯采样的两个维度mu和sigma的条件概率，就要先求出在数据x下mu和sigma的联合概率密度（似然函数），而这个就先要知道在mu和sigma下x的概率密度（就是高斯分布呀）以及mu和sigma的先验概率分布。

y为均值为μ协方差矩阵为 Σ的( N × 1)正态随机向量，其概率密度函数pdf为：

![Center 1][]

\[[概率论：高斯分布][Link 1]\]

假设x是y的截断truncated版本，有同样的参数location and scale parameters μ and Σ , 只是截断到区域 R = \{ ( a i < x i < b i ), i = 1,2,... N \} ，(a 可以为 −∞ or b 可以为 +∞ )。

**x的pdf为：**

![Center 2][]

I R ( x )是指示器函数：当x在区域R中时值为1，否则为0；

P ( μ , Σ ) \]− 1是正则化常数（由于截断），为

![Center 3][]

让 x = ( x 1 , x 2 ,..., x T ) 表示截断多维高斯分布pdf f ( x | μ , Σ ) 的随机采样

**则随机采样x的pdf（μ 和 Σ的似然函数）为**

![Center 4][]

## 高斯分布的先验和数据产生后验分布    ##

**( μ , Σ )的先验分布**

设为noninformative diffuse prior(扩散先验密度)，As a prior pdf for ( μ , Σ ) , we will use the conventional noninformative diffuse prior (see, for example, Zellner 1971, p.225)

故**( μ , Σ )的后验分布**

Combining this prior with the likelihood function yields the posterior pdf for ( μ , Σ )

![Center 5][]

得到后验分布后，我们可以从这些后验pdfs中采样，并利用这些采样值估计边缘后验分布和矩（moments）。

然而其中的和式没法消去，后面将介绍如何利用隐含变量消去它。

[皮皮blog][blog]

# 单变量高斯分布The Univariate Case #

单变量高斯分布中，正则化因子式3变为：

![Center 6][]

where Φ (.) is the standard normal cumulative distribution function (cdf).

(μ , σ 2)的后验pdf (式6)为

![Center 7][]

然而从上式没法得到f ( σ 2 | x ) and f ( μ | x )：It is not possible to analytically integrate out μ or σ 2 from this pdf to obtain the marginal posterior pdfs f ( σ 2 | x ) and f ( μ | x ) . Also, because the conditional posterior pdfs f ( σ 2 | μ , x ) and f ( μ | σ 2 , x ) are not recognisable, it is not possible to set up a Gibbs’ sampling algorithm that draws from these conditional pdfs.

所以我们添加对应x的隐含变量y：令y = ( y 1 , y 2 ,..., y T ) ′ ， 可看作非截断正态分布N ( μ , σ 2 )的采样，并且和截断观测值x有直接的（确定性的deterministic）对应。

考虑采样观测值y t from N ( μ , σ 2 )和x t from N ( μ , σ 2 ) × I ( a , b ) ( x t )的逆累积分布函数（inverse cdf）方法：

给定均匀分布抽样Ufrom (0,1)，y t and x t的抽样分别为：

![Center 8][]

![Center 9][]

。。。一堆推导（略）

得到x和y之间的确定deterministic关系

![Center 10][]

根据贝叶斯理论，我们可以**重写μ , σ 2 and y的联合后验分布（joint posterior pdf）为**

![Center 11][]

由于f ( x | y , μ , σ 2 ) = 1 when (12) holds, and is zero otherwise.     The remaining terms on the right side of (13) involve y not x.则

![Center 12][]

where y is the sample mean of the y t and the relationship between y and x is given by (12).

式14中，**吉布斯采样需要的条件后验pdfs分别为**：

![Center 13][]

式16：[Inverse-gamma distribution][]， i\_gamma(alpha, beta)采样数据实际上就是gamma分布(alpha, beta^-1)采样数据取倒数。

IG的概率密度函数pdf为

![Center 14][]

综上，

## 吉布斯采样截断高斯分布参数算法    ##

从后验pdf中生成( μ , σ 2 )的步骤为

![Center 10][]

Note: 这一步对于非截断正态分布参数的推断当然是不需要的，这时12式就是y = x.(ab都取相应无穷大时)。

![Center 15][]

![Center 16][]

Note:

第3步就是

![Center 17][]

第4步就是

![Center 18][]

[皮皮blog][blog]

## 同步图计算并行化处理 ##

lz的不一定正确。。。

### 一维高斯参数推断 ###

高斯分布参数均值mu和方差sigma作为吉布斯采样的两个坐标维度  
吉布斯采样作为同步图节点，节点间是全连接的  
将高斯分布数据划分为节点数目个文件提供给对应节点采样新维度

使用已有方法产生截断高斯分布的数据点，将数据分成M份，分布在M个worker节点中。

**超步0**  
随机初始值作为mu和sigma的局部值，利用数据更新mu和sigma，传递mu和sigma的局部值给其它节点。（局部模型）  
**超步1+**

每个节点收到其它节点发来的( μ , σ 2 )信息，重新计算( μ , σ 2 )，记录到文件中，并作为新一次迭代的初始值。（全局模型，替换局部模型）

所有节点重新运行吉布斯采样截断高斯分布参数算法，得到部分训练数据的参数采样结果( μ , σ 2 )。

将新采样的mu和sigma的局部值传递给其它节点。

或者可以每隔几个Gibbs sampling迭代结束超步进行同步，并准备进入下一个超步。而不是一次迭代就进入下一下超步。

**判断程序是否终止**  
累加量acc，最大迭代步数。

\[Distributed Inference for Latent Dirichlet Allocation\]

最后lz有一个问题，吉布斯采样能用在连续的n维高斯分布采样中吗？如果可以如何实现，马尔可夫毯？

from: [http://blog.csdn.net/pipisorry/article/details/51539739][blog.csdn.net_pipisorry_article_details_51539739]

ref: W. Griffiths. A Gibbs Sampler for the Parameters of a TruncatedMultivariate Normal Distribution.

[blog.csdn.net_pipisorry_article_details_51539739]: http://blog.csdn.net/pipisorry/article/details/51539739
[Gibbs Sampling]: http://blog.csdn.net/pipisorry/article/details/51373090
[Center]: /images/20220722/5098a645d369446299883e1de78397b7.png
[blog]: http://blog.csdn.net/pipisorry
[Center 1]: /images/20220722/8def7b9b8fb249308e6c18903e9025ec.png
[Link 1]: http://blog.csdn.net/pipisorry/article/details/49516209
[Center 2]: /images/20220722/2d7d37739dff448e9ee69309ba8f5b93.png
[Center 3]: /images/20220722/ba580211f9ea405c9ea8b311f56fd6a2.png
[Center 4]: /images/20220722/7c2c7efb5d044484a2b0315f42026e1f.png
[Center 5]: /images/20220722/6e2f1e5a9a7d4a9f8c33fa220dd5ec6d.png
[Center 6]: /images/20220722/03f6c67fa61a4b44b2c73099004414d7.png
[Center 7]: /images/20220722/6d1eea7d16cd4991b5b9e7c986c54b95.png
[Center 8]: /images/20220722/ed863cbf39a548a3989931fbc891a77a.png
[Center 9]: /images/20220722/2b90793f296143078592a93dab5ec52c.png
[Center 10]: /images/20220722/1e390179b73746fea01e3924fa81dd56.png
[Center 11]: /images/20220722/ac683050689f4c3c84841dddbaa5e9dd.png
[Center 12]: /images/20220722/22dde882d10f4879a3830366cf410a57.png
[Center 13]: /images/20220722/9bbaceb37e4a41c4948b60bddb66f905.png
[Inverse-gamma distribution]: https://en.wikipedia.org/wiki/Inverse-gamma_distribution
[Center 14]: /images/20220722/a879e90bc7414044b1d8adc1d4352029.png
[Center 15]: /images/20220722/723f3fdd93f9456ba2fc509b807b74cd.png
[Center 16]: /images/20220722/bff9b0cfcf544d52bb6fe3409f6cdc2b.png
[Center 17]: /images/20220722/26d0b525c05d47a5a5a20b8568e54f95.png
[Center 18]: /images/20220722/af43eb97e77e400c94175c8717d0a711.png