使用R语言进行Metroplis-in-Gibbs采样和MCMC运行分析
2021-02-12 11:17
标签:基本算法 process 简单 有趣 title 慢慢 tin name 频繁 原文链接:http://tecdat.cn/?p=12200 对于许多模型,例如物流模型,没有共轭先验。因此,吉布斯采样不适用。 这篇文章展示了我们如何使用Metropolis-Hastings(MH)从每次Gibbs迭代中的非共轭条件后验对象中进行采样–比网格方法更好的替代方法。 我将说明该算法,给出一些R代码结果,然后分析R代码以识别MH算法中的瓶颈。 此示例的模拟数据是包含?患者的横截面数据集。有一个二元结果,?一个二元治疗变量,?一个因子age。年龄是具有3个等级的分类变量。我用贝叶斯逻辑回归建模: ? 对于大都市吉布斯来说,这是一个相当现实的示例: 让我们看一下该模型的(非标准化)条件后验。我不会进行推导,但是它遵循我以前的帖子中使用的相同过程。 此条件分布不是已知分布,因此我们不能简单地使用Gibbs从中进行采样。相反,在每个gibbs迭代中,我们需要另一个采样步骤来从该条件后验中提取。第二个采样器将是MH采样器。 目标是从中取样?。请注意,这是4维密度。 MH采样器的工作方式如下: 因此,总是接受产生更高条件的后验评估的提议。但是,有时仅接受具有较低密度评估的提案-提案的相对密度评估越低,其接受的可能性就越低。 经过多次迭代,从后验的高密度区域开始的抽样被接受,并且被接受的序列“爬升”到高密度区域。一旦序列到达此高密度区域,它将趋于保持在那里。因此,这也类似于模拟退火。 这种表示法很容易扩展到我们的4维示例:提案分布现在是4维多元高斯模型。代替标量方差参数,我们有一个协方差矩阵。因此,我们的建议是系数的向量。从这个意义上讲,我们运行的是Gibbs –使用MH每次迭代绘制整个系数块。 这是我们感兴趣的4个参数的MCMC链。红线表示真实值。 ? 有一些改进的空间: 关于贝叶斯范式的好处是,所有推断都是使用后验分布完成的。现在,系数估计值是对数刻度,但是如果我们需要比值比,则只需对后验取幂。如果我们想要对比值比进行区间估计,那么我们就可以获取指数后验平局的2.5%和97.5%。 下面是使用R分析,显示了这一点。for循环运行Gibbs迭代。在每个Gibbs迭代中,我都调用函数rcond_post_beta_mh(),该函数使用MH从参数向量的条件后验中得出图形。 ? 深入研究rcond_post_beta_mh(),我们看到子例程log_cond_post_beta()是MH运行中的瓶颈。此函数是beta载体的对数条件后验密度,将其评估两次。 ? 使用R语言进行Metroplis-in-Gibbs采样和MCMC运行分析 标签:基本算法 process 简单 有趣 title 慢慢 tin name 频繁 原文地址:https://www.cnblogs.com/tecdat/p/12731381.html
模型
非规范条件后验
?Metroplis-in-Gibbs采样
结果
profvis(expr = {
for(i in 2:gibbs_iter){
# sample from posterior of phi
gibbs_res[i,p+1]
文章标题:使用R语言进行Metroplis-in-Gibbs采样和MCMC运行分析
文章链接:http://soscw.com/index.php/essay/54413.html