MCMC_关键概念理解

  • 关于概念和符号,参考前文:MCMC_方法示例
  • 以下是个人在理解MCMC方法中遇到的疑惑,以及个人的见解。

问题一

回顾如下公式,既然对于 α ( i , j ) \alpha(i,j) α(i,j)我们将之理解为从一条马氏链过渡到另一条马氏链的跳转率(也作接受率),并且实际操作中是通过物理过程:从均匀分布中抽样这个动作来实现的。那么式中的 q ( i , j ) q(i,j) q(i,j),作为转移矩阵中的一个点,它也是一个取值于(0,1)的概率值,它的作用形式是怎样的?
p ( i ) q ( i , j ) α ( i , j ) = p ( j ) q ( j , i ) α ( j , i ) p(i)q(i,j)\alpha(i,j)=p(j)q(j,i)\alpha(j,i) p(i)q(i,j)α(i,j)=p(j)q(j,i)α(j,i)

解答

回顾实际采样过程:

  • a 用X表示由MCMC过程产生的状态向量, Q 0 Q_0 Q0表示状态转移矩阵( D 0 D_0 D0表示抽样分布)。初始化 X 0 = x 0 X_0=x_0 X0=x0,其中 x 0 x_0 x0是从状态集(分布)中获取的随机样本。
  • b 对 t = 0 , 1 , 2 , ⋯   , n t = 0, 1, 2, \cdots,n t=0,1,2,,n重复如下步骤进行采样:
    • b-1 第 t t t时刻有 X t = x t X_t=x_t Xt=xt, 采样 y ∼ q ( x ∣ x t ) y\sim q(x|x_t) yq(xxt)(此时状态转移矩阵为 Q t = q ( x ∣ x t ) Q_t=q(x|x_t) Qt=q(xxt))。
    • b-2 从 U ( 0 , 1 ) U(0,1) U(0,1)中抽取一个样本点 u u u
    • b-3 如果 u &lt; α ( x t , y ) = m i n ( p ( j ) q ( j , i ) p ( i ) q ( i , j ) , 1 ) u&lt;\alpha(x_t, y)=min(\frac{p(j)q(j,i)}{p(i)q(i,j)},1) u<α(xt,y)=min(p(i)q(i,j)p(j)q(j,i),1)则接受转移 X t + 1 = y X_{t+1}=y Xt+1=y
    • b-4 否则拒绝转移, X t + 1 = x t X_{t+1}=x_t Xt+1=xt, Q t + 1 = Q t Q_{t+1} = Q_t Qt+1=Qt

其中b-1步中获得的样本 y y y是从条件分布 q ( x ∣ x t ) q(x|x_t) q(xxt)中随机抽样的。因此实现 q ( i , j ) q(i,j) q(i,j)发生作用的物理过程就是在这一步中所进行的随机抽样。

问题二

为什么在解读 α ( i , j ) \alpha(i, j) α(i,j)的时候要把它理解为一个物理过程?

解答

这实际上是要理解概率之于样本或者抽样过程的关系。

  • 概率是对所描述事物随机性的度量,是“描述”,是写在纸上的概率规则。
  • 抽样过程是从分布中依照概率获取真实样本的过程,其结果是确定的样本,是“实体”,因此可看做物理过程。

由于抽样过程的结果是实实在在的样本,它并不是一个概率值,例如:

  • 对于抛硬币的实验而言,每一次抛掷得到的结果是确定的正或反,即使从数学上每一面向上的概率不同,但每一次实验的结果也不可能是一个概率值。
  • 同样在通过马氏链进行的抽样中,即便能够直接给出稳态分布,但是要获取哪怕一个样本都要通过随机的方式进行物理抽样的。对于前文中马氏链部分的细节,即使知道了媳妇儿有0.613的概率不在“愤怒”状态,但正真见着她时,她必定是处于一个确定的心情状态的,这个时候如果要做一个预判的话,那就进行一次随机抽样把抽样结果当做“上帝的指示”。
  • “于是每天携带一块秒表,进门前按下并读取毫秒数。如果数值小于0.613放心大胆进门,否则等15分钟再按一次。”

问题三

如何理解拒绝跳转这一现象?我们知道如果发生一次拒绝跳转,就把上一次抽样结果当做本次抽样结果放到抽样链里去。相当于该样本重复了一次,那么在最后得到的样本中是否要去除这些重复的样本(或者说只取连续游程中的第一个)?

解答

假设状态集合为 S = [ s 1 , s 2 , s 3 ] S=[s_1, s_2, s_3] S=[s1,s2,s3],满足细致平稳条件的 Q ′ Q\prime Q的矩阵结构如下式所示,可知由于乘以 α ( i , j ) \alpha(i,j) α(i,j)而丢失的概率质量全部转移到了对角线上,相当于放大了状态转移到自己的概率。而这会降低马氏链收敛到稳态的速度。
Q ′ = Q ⋅ A = { 1 − ∑ j = 1 3 p ( 1 , j ) α ( 1 , j ) p ( 1 , 2 ) α ( 1 , 2 ) p ( 1 , 3 ) α ( 1 , 3 ) p ( 2 , 1 ) α ( 2 , 1 ) 1 − ∑ j = 1 3 p ( 2 , j ) α ( 2 , j ) p ( 2 , 3 ) α ( 2 , 3 ) p ( 3 , 1 ) α ( 3 , 1 ) p ( 3 , 2 ) α ( 3 , 2 ) 1 − ∑ j = 1 3 p ( 3 , j ) α ( 3 , j ) } Q\prime=Q\cdot A= \\ \left\{ \begin{matrix} 1-\sum_{j=1}^3 p(1,j)\alpha(1,j) &amp; p(1,2)\alpha(1,2) &amp; p(1,3)\alpha(1,3) \\ p(2,1)\alpha(2,1) &amp; 1-\sum_{j=1}^3 p(2,j)\alpha(2,j) &amp; p(2,3)\alpha(2,3) \\ p(3,1)\alpha(3,1) &amp; p(3,2)\alpha(3,2) &amp; 1-\sum_{j=1}^3 p(3,j)\alpha(3,j) \end{matrix} \right\} Q=QA=1j=13p(1,j)α(1,j)p(2,1)α(2,1)p(3,1)α(3,1)p(1,2)α(1,2)1j=13p(2,j)α(2,j)p(3,2)α(3,2)p(1,3)α(1,3)p(2,3)α(2,3)1j=13p(3,j)α(3,j)

假设t-1步链上的样本为 s ( t − 1 ) = s 1 s^{(t-1)}=s_1 s(t1)=s1,同时t步已抽取样本 s t = s 2 s^t=s_2 st=s2,则此时 Q ′ ∣ s t = s 2 Q\prime|s^t=s_2 Qst=s2具有如下形式,由下式可知这就是跳转的真实过程。

Q ′ ∣ ( s ( t − 1 ) = s 1 , s t = s 2 ) = ( Q ⋅ A ) ∣ ( s ( t − 1 ) = s 1 , s t = s 2 ) Q\prime|(s^{(t-1)} = s_1,s^t=s_2)=(Q\cdot A)|(s^{(t-1)} = s_1,s^t=s_2) Q(s(t1)=s1,st=s2)=(QA)(s(t1)=s1,st=s2)

= { ( 1 − ∑ j = 1 3 p ( 1 , j ) α ( 1 , j ) ) ∣ t ( p ( 1 , 2 ) α ( 1 , 2 ) ) ∣ t ( p ( 1 , 3 ) α ( 1 , 3 ) ) ∣ t … … … … … … } = \left\{ \begin{matrix} (1-\sum_{j=1}^3 p(1,j)\alpha(1,j))|_t &amp; (p(1,2)\alpha(1,2))|_t &amp; (p(1,3)\alpha(1,3))|_t \\ \dots &amp; \dots &amp; \dots \\ \dots &amp; \dots &amp; \dots \end{matrix} \right\} \\ =(1j=13p(1,j)α(1,j))t(p(1,2)α(1,2))t(p(1,3)α(1,3))t

= { 1 − p ( 1 , 2 ) α ( 1 , 2 ) p ( 1 , 2 ) α ( 1 , 2 ) 0 … … … … … … } =\left\{ \begin{matrix} 1-p(1,2)\alpha(1,2) &amp; p(1,2)\alpha(1,2) &amp; 0 \\ \dots &amp; \dots &amp; \dots \\ \dots &amp; \dots &amp; \dots \end{matrix} \right\} =1p(1,2)α(1,2)p(1,2)α(1,2)0

问题四

如何确定Burn-In过程的长度,即:该舍弃多少次抽样?

解答

视问题复杂度,由于 α ( i , j ) \alpha(i,j) α(i,j)的存在延长了抽样过程中马氏链收敛到稳态时间。

<think>好的,我需要帮助用户了解MCMC算法的概念和编程实现。首先,我得回忆一下MCMC的基础知识,确保自己理解正确。MCMC是马尔可夫链蒙特卡罗的缩写,主要用于从复杂概率分布中抽样,特别是在高维空间或者解析解难以求得的情况下。 接下来,我需要整理MCMC的核心概念。比如,马尔可夫链的平稳分布,这是保证抽样结果正确性的关键。还有Metropolis-Hastings算法,这是最经典的MCMC方法之一。另外,Gibbs抽样也是常用的一种方法,适用于条件分布容易采样的情形。 然后,编程实现部分。用户可能希望看到具体的代码示例,Python是比较常用的语言,所以我应该用Python来写。Metropolis-Hastings算法需要先定义目标分布,比如正态分布,然后选择提议分布,比如均匀分布或者正态分布。代码的结构包括初始化参数、循环迭代、接受或拒绝样本等步骤。需要注意的是调整步长参数,以确保合适的接受率,通常推荐在20%-50%之间。 验证算法的正确性也很重要。可以通过绘制抽样结果的直方图,并与目标分布的概率密度函数对比,观察是否匹配。此外,可能需要提到收敛诊断的方法,比如Gelman-Rubin统计量,但代码示例中可以简单用图形对比来展示。 用户可能的问题包括如何选择提议分布,如何处理收敛问题,MCMC的应用场景,以及不同算法之间的比较,比如MH和Gibbs。这些可以作为相关问题列在最后,帮助用户进一步探索。 需要检查LaTeX格式是否正确,所有数学公式都用$...$或$$...$$包裹。例如,目标分布的定义,接受概率的计算等。确保代码块语法正确,注释清晰,方便用户理解。 最后,确保回答结构清晰,分点讲解概念和实现步骤,代码示例后附上相关问题和引用。引用部分要根据提供的参考信息添加,但用户提供的引用似乎与MCMC无关,可能需要调整或者确认是否有误。这里可能用户提供的引用有误,但按照要求还是在末尾加上引用标识。</think>### MCMC算法基础与实现 马尔可夫链蒙特卡罗(MCMC)是一种通过构建马尔可夫链从复杂概率分布中抽样的方法[^1]。其核心思想是通过迭代生成样本,使得马尔可夫链的平稳分布等于目标分布。常见算法包括**Metropolis-Hastings**和**Gibbs抽样**。 #### 关键步骤 1. **目标分布**:定义需抽样的分布$p(x)$。 2. **提议分布**:选择过渡分布$q(x'|x)$生成候选样本。 3. **接受概率**:计算候选样本的接受概率$\alpha = \min\left(1, \frac{p(x')q(x|x')}{p(x)q(x'|x)}\right)$。 4. **迭代更新**:以概率$\alpha$接受候选样本,否则保留当前样本。 #### Python实现示例(Metropolis-Hastings) 以下代码从正态分布$N(0,1)$中抽样: ```python import numpy as np import matplotlib.pyplot as plt def target_distribution(x): return np.exp(-x**2 / 2) / np.sqrt(2 * np.pi) # 标准正态分布 def metropolis_hastings(n_samples, initial=0, step_size=1): samples = [initial] current = initial for _ in range(n_samples): proposal = current + np.random.normal(0, step_size) accept_ratio = target_distribution(proposal) / target_distribution(current) if np.random.rand() < accept_ratio: current = proposal samples.append(current) return samples # 生成样本并验证 samples = metropolis_hastings(10000) plt.hist(samples, bins=50, density=True) x = np.linspace(-4, 4, 100) plt.plot(x, target_distribution(x), 'r') plt.show() ``` **代码说明**: - `step_size`控制提议分布的步长,影响接受率[^2]。 - 直方图应与红色曲线(目标分布)吻合,验证算法正确性。 #### 优化与注意事项 - **步长调整**:接受率通常在20%-50%之间最优。 - **预烧期**:丢弃前$N$个样本以避免初始值影响。 - **相关性**:可通过稀疏采样(如每隔10个样本取一次)降低样本相关性。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值