MCMC基本原理与应用(一)

背景

贝叶斯统计推断是一个非常艺术的东西,他的先验结果(prior knowledge)是一种专家的经验,而贝叶斯公式给出的后验分布是贝叶斯推断的基本工具。每个大学生都应该学过数理统计,其中提及的贝叶斯统计推断都是一些基于简单先验和简单后验的结果。

但当我们一旦碰到复杂的先验(比如高维参数、复杂先验分布),我们用贝叶斯公式得出的后验分布将变得非常复杂,在计算上会非常困难。为此,先人提出了MCMC算法方便我们可以对任何后验分布进行计算或推断。其思想就是其名字:两个MC。

第一个MC: Monte Carlo(蒙特卡洛)。这个简单来说是让我们使用随机数(随机抽样)来解决计算问题。在MCMC中意味着:后验分布作为一个随机样本生成器,我们利用它来生成样本(simulation),然后通过这些样本对一些感兴趣的计算问题(特征数,预测)进行估计。

第二个MC:Markov Chain(马尔科夫链)。第二个MC是这个方法的关键,因为我们在第一个MC中看到,我们需要利用后验分布生成随机样本,但后验分布太复杂,一些package中根本没有相应的随机数生成函数(如rnorm(),rbinom())怎么办?答案是我们可以利用Markov Chain的平稳分布这个概念实现对复杂后验分布的抽样。

Markov Chain

为了能顺利阐述MCMC算法,这里就简单地讲一下所涉及到的马尔科夫链概念。因为都是我个人的理解,这里所讲不敢说都是正确的,希望各位童鞋能够独立思考。

1.什么是Markov Chain

随机过程${X_n},X_n$的状态空间为有限集或可列集,比如当$X_n=i$,就称过程在n处于状态i。
定义:如果对任何一列状态$i_0,i_1,...,i_{n-1},i,j$,及对任何n≥0,随机过程${X_n}$满足Markov性质:
$$P{X_{n+1}=j|X_0=i_0,...X_{n-1}=i_{n-1},X_n=i}\=P{X_{n+1}=j|X_n=i}$$

转移概率
$P{x_{n+1}=j|X_n=i}$成为Markov链的一步转移概率$P_{ij}^{n,n+1}$,当这个概率与n无关,则记之为$P_{ij}$

转移概率矩阵P
这个矩阵的元素就是一步转移概率$P_{ij}$

结论
一个Markov链可以由它的初始状态以及转移概率矩阵P完全确定。(证明略,自行百度或翻书)

2.n步转移概率

所谓n步转移概率就是从状态i走n步正好到状态j的概率,我们记为$P_{ij}^{(n)}$。
利用概率分割的思想,由基础概率论中全概率公式可以得到$$P_{ij}{(n)}=\Sigma_{k=0}{\infty}P_{ik}P_{kj}^{(n-1)}$$ 写成矩阵形式就是$$P^{(n)}=P\times P^{(n-1)}$$
进一步推广,我们就推出了著名的Chapman-Kolmogorov方程:$$P_{ij}{(n+m)}=\Sigma_{k=0}{\infty}P_{ik}{(n)}P_{kj}{(m)}$$写成矩阵形式就是$$P{(m+n)}=P{(m)}\times P^{(n)}$$

3.Markov链的极限定理和平稳分布

现在我们已经有了n步转移概率的概念,一个很简单的想法就是如果n趋向无穷会怎么样?这个问题也是后面极限分布以及平稳分布的来源,更是MCMC算法中第二个MC的关键。
要回答这个问题首先要掌握几个关键概念,我先列出来,如果不熟悉的可以自行百度或翻书:

  1. 互达性($i \leftrightarrow j$)
  2. 周期性(d(i))
  3. 常返与瞬过($f_{ii}=1$)
  4. 常返时($\T_i$)

3.1互达性、周期性

几个重要结论(证明略,自行百度,或者call me)

  1. $i \leftrightarrow j,\Rightarrow d(i)=d(j)$
  2. $若存在d(i)<\infty,则存在N,对所有n>N,有P_{ii}^{(nd(i))}>0$
  3. $P_{ji}^{(m)}>0,\Rightarrow 存在N,对所有n>N,有P_{ii}^{(m+nd(i))}>0$
  4. 对于非周期不可约Markov链的转移概率矩阵P,存在N,当n≥N时,$P^{(n)}>0$

3.2常返、瞬过、常返时

  1. 引入重要概率$f_{ij}{(n)}$表示从i出发在n步**首次**到达j的概率,我们约定$f_{ii}{(0)}=0$

  2. $定义f_{ij}=\Sigma_{n=1}{\infty}f_{ij}{(n)}$表示从i出发最终到j的概率。

  3. $定义:如果f_{ii}=1,则状态i是常返的,否则是瞬过的$

  4. $我们记T_i为首次返回i的步长,且\mu_i=ET_i$

  5. $若\mu_i=\infty,则称状态i是零常返,否则正常返$

3.3 极限定理

是时候回答上面那个问题了!(摘自网络共享PPT)

另外对于一个正常返、非周期的状态(也称为遍历ergodic),我们有结论:$$对所有i\leftrightarrow j,有limP_{ji}{(n)}=limP_{ii}{(n)}=\frac{1}{\mu_i}$$

到这里,我们可以想象了,n步转移概率矩阵最终的极限形式应该是由相同的行向量组成的!!!我们把那个行向量称为极限概率分布。 但是,问题还没完,要计算极限概率,就要根据极限定理,求出很多东西(如$f_{ii}^{(n)}$),实际中并不方便。所以,我们要引入平稳分布这个概念。最关键的部分来了!!!!

3.3 平稳分布

什么是平稳分布?它和求极限概率分布有什么关系呢?

定义:Markov链有转移概率矩阵P,如果有一个概率分布${\pi_i }满足\pi_j =\Sigma_{i=0}^{\infty}\pi_i P_{ij}$,则称为这个Markov链的平稳分布。这个定义用矩阵形式写出来就是π*P=π.

这个定义内容很丰富:如果一个过程的初始状态$X_0$有平稳分布π,我们可以知道对所有n,$X_n$有相同的分布π。再根据Markov性质可以得到,对任何k,有$X_n,X_{n+1},...,X_{n+k}$的联合分布不依赖于n,显然这个过程是严格平稳的,平稳分布也由此得名!!

重要定理

若一个markov链中所有的状态都是遍历的,则对所有i,j有$limP_{ij}{(n)}=\pi_j存在,且\pi={\pi_j,j≥0}就是平稳分布$。反之,拖一个不可约Markov链只存在一个平稳分布,且这个Markov链的所有状态都是遍历的,则这个平稳分布就是该Markov链的极限概率分布。$$limP_{ij}{(n)}=\pi_j$$

上面这条定理就给出了通过求平稳分布来求Markov链极限概率分布的简单方法。到这里我对第二个MC的解读也就结束了。下一期我们就要开始具体讲讲怎么应用MCMC算法对后验分布进行抽样和推断。

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 211,042评论 6 490
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 89,996评论 2 384
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 156,674评论 0 345
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,340评论 1 283
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,404评论 5 384
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,749评论 1 289
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,902评论 3 405
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,662评论 0 266
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,110评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,451评论 2 325
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,577评论 1 340
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,258评论 4 328
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,848评论 3 312
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,726评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,952评论 1 264
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,271评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,452评论 2 348

推荐阅读更多精彩内容

  • 本文主要摘抄整理自Rickjin的《LDA数学八卦》 统计模拟中一个很重要的问题就是给定一个概率分布$p(x)$,...
    iV0id阅读 419评论 0 1
  • # 背景 贝叶斯统计推断是一个非常*艺术*的东西,他的先验结果(prior knowledge)是一种专家的经验,...
    大仙Dylon阅读 667评论 0 1
  • 原文传送门:也谈MCMC方法与Gibbs抽样 MCMC,即传说中的Markov Chain Mento Carlo...
    willheng阅读 10,977评论 2 14
  • 本系列第三篇,承接前面的《浅谈机器学习基础》和《浅谈深度学习基础》。 自然语言处理绪论 什么是自然语言处理? 自然...
    我偏笑_NSNirvana阅读 17,545评论 2 68
  • 我从来不曾忘记过,只是不愿想起。
    天命之光阅读 175评论 0 0