R语言马尔可夫MCMC中的METROPOLIS HASTINGS,MH算法抽样(采样)法可视化实例

简介: R语言马尔可夫MCMC中的METROPOLIS HASTINGS,MH算法抽样(采样)法可视化实例

原文链接:http://tecdat.cn/?p=26324 

介绍

Metropolis Hastings 算法是一种非常简单的算法,用于从难以采样的分布中生成样本。

假设我们要从分布 π 中进行采样,我们将其称为“目标”分布。为简单起见,我们假设 π是实线上的一维分布,尽管它很容易扩展到一维以上(见下文)。

MH 算法通过模拟马尔可夫链来工作,其平稳分布为 π。这意味着,从长远来看,来自马尔可夫链的样本看起来像来自 π的样本。正如我们将看到的,该算法非常简单和灵活。


MH算法

转移核

要实现 MH 算法,用户必须提供一个“转移核”Q。转移核只是一种在 给定当前位置(例如 x)的情况下随机移动到空间中新位置(例如 y)的方式。也就是说,Q 是给定 x 在 y 上的分布,我们将其写成 Q(y|x)。在许多应用中,Q将是一个连续分布,在这种情况下 Q(y|x) 将是 y 上的密度,因此∫Q(y|x)dy=1(对于所有 x)。

例如,从当前位置 x 生成新位置 y 的一种非常简单的方法是向 x添加一个 N(0,1) 随机数。即设置y=x+N(0,1),或者转移y|x∼N(x,1)。所以

image.png

这种在当前位置x加上一些随机数得到y的核,在实际中经常使用,被称为“随机游走”核。

MH算法

使用转移核 Q 从目标分布 π 中采样的 MH 算法包括以下步骤:

  • 初始化,X1=x1 。
  • 对于 t=1,2,…
  • 从 Q(y|xt)中采样 y。将 y 视为 xt+1 的“建议”值。
  • 计算
  • image.png
  • A通常被称为“接受概率”。
  • 以概率 A“接受”提议的值,并设置 xt+1=y。否则设置 xt+1=xt。
  • Metropolis 算法
    请注意,上面给出的示例随机游走建议 Q 对于所有 x,y 满足 Q(y|x)=Q(x|y) 任何满足这一点的建议都称为“对称”。当 Q 是对称时,MH 算法中 A 的公式 简化为:
    image.png

该算法的这种特殊情况,具有 Q 对称,首先由 Metropolis 等人在 1953 年提出,因此它有时被称为“Metropolis 算法”。

示例

为了帮助理解 MH 算法,我们现在做一个简单的例子:我们实现算法以从指数分布中采样:

image.png

当然,以其他方式从指数分布中采样会容易得多;我们只是用它来说明算法。

请记住,π 被称为“目标”分布,因此我们调用函数来计算 π  target

现在我们实现 MH 算法,使用上面提到的简单正态随机游走转移核 Q。

这是代码:

x = rep(0,10000)
x1
1 = 3     #初始化;我任意地将其设置为3
for(i in 2:10000){
  if(){
    xi
i = proposed_x       # 以最小(1,A)的概率接受移动。
  } else {
    xii
i = current_x        # 否则就 "拒绝 "移动,并留在原地。
  }
}

运行此代码后,我们可以绘制马尔可夫链 x 访问的位置(有时称为轨迹图)。

image.png



请记住,我们设计此算法是为了从指数分布中采样。这意味着(只要我们运行算法足够长的时间!)x 的直方图应该看起来像一个指数分布。在这里我们检查一下:

hist(x)
lines

image.png

x 中的值的直方图确实提供了与指数分布的紧密拟合。

结束语

MH 算法的一个特别有用的特性是,即使 只知道π 是一个常数,它也可以实现:也就是说,对于一些已知的 f,π(x)=cf(x) , 但未知常数 c。这是因为该算法仅通过比率 image.png 依赖于π 。

这个问题出现在贝叶斯应用中,其中后验分布与先验概率成正比,但比例常数通常是未知的。因此,MH 算法对于从后验分布进行采样以执行难以解析的贝叶斯计算特别有用。

相关文章
|
6天前
|
移动开发 算法 数据可视化
【视频】马尔可夫链蒙特卡罗方法MCMC原理与R语言实现|数据分享(上)
【视频】马尔可夫链蒙特卡罗方法MCMC原理与R语言实现|数据分享
|
6天前
|
数据可视化
R语言离散时间马尔可夫链(Markov chain)模型分类案例可视化分析
R语言离散时间马尔可夫链(Markov chain)模型分类案例可视化分析
|
6天前
|
存储 机器学习/深度学习 算法
R语言贝叶斯Metropolis-Hastings采样 MCMC算法理解和应用可视化案例
R语言贝叶斯Metropolis-Hastings采样 MCMC算法理解和应用可视化案例
|
6天前
|
机器学习/深度学习 算法 数据可视化
Matlab决策树、模糊C-均值聚类算法分析高校教师职称学历评分可视化
Matlab决策树、模糊C-均值聚类算法分析高校教师职称学历评分可视化
|
6天前
|
算法 数据可视化
圆填充( CIRCLE PACKING)算法圆堆图圆形空间填充算法可视化
圆填充( CIRCLE PACKING)算法圆堆图圆形空间填充算法可视化
|
6天前
|
算法
【视频】马尔可夫链蒙特卡罗方法MCMC原理与R语言实现|数据分享(下)
【视频】马尔可夫链蒙特卡罗方法MCMC原理与R语言实现|数据分享
|
6天前
|
算法 数据安全/隐私保护 计算机视觉
基于二维CS-SCHT变换和LABS方法的水印嵌入和提取算法matlab仿真
该内容包括一个算法的运行展示和详细步骤,使用了MATLAB2022a。算法涉及水印嵌入和提取,利用LAB色彩空间可能用于隐藏水印。水印通过二维CS-SCHT变换、低频系数处理和特定解码策略来提取。代码段展示了水印置乱、图像处理(如噪声、旋转、剪切等攻击)以及水印的逆置乱和提取过程。最后,计算并保存了比特率,用于评估水印的稳健性。
|
2天前
|
算法
m基于BP译码算法的LDPC编译码matlab误码率仿真,对比不同的码长
MATLAB 2022a仿真实现了LDPC码的性能分析,展示了不同码长对纠错能力的影响。短码长LDPC码收敛快但纠错能力有限,长码长则提供更强纠错能力但易陷入局部最优。核心代码通过循环进行误码率仿真,根据EsN0计算误比特率,并保存不同码长(12-768)的结果数据。
21 9
m基于BP译码算法的LDPC编译码matlab误码率仿真,对比不同的码长
|
4天前
|
算法
MATLAB|【免费】融合正余弦和柯西变异的麻雀优化算法SCSSA-CNN-BiLSTM双向长短期记忆网络预测模型
这段内容介绍了一个使用改进的麻雀搜索算法优化CNN-BiLSTM模型进行多输入单输出预测的程序。程序通过融合正余弦和柯西变异提升算法性能,主要优化学习率、正则化参数及BiLSTM的隐层神经元数量。它利用一段简单的风速数据进行演示,对比了改进算法与粒子群、灰狼算法的优化效果。代码包括数据导入、预处理和模型构建部分,并展示了优化前后的效果。建议使用高版本MATLAB运行。
|
6天前
|
算法 计算机视觉
基于高斯混合模型的视频背景提取和人员跟踪算法matlab仿真
该内容是关于使用MATLAB2013B实现基于高斯混合模型(GMM)的视频背景提取和人员跟踪算法。算法通过GMM建立背景模型,新帧与模型比较,提取前景并进行人员跟踪。文章附有程序代码示例,展示从读取视频到结果显示的流程。最后,结果保存在Result.mat文件中。

热门文章

最新文章