问题标签 [mcmc]

我试图在 PyMC3 中组合一个动态系统模型，以推断两个参数。该模型是基本的 SIR，常用于流行病学：

dS/dt = - r0 * g * S * I

dI/dt = g * I (r * S - 1 )

其中 r0 和 g 是要推断的参数。到目前为止，我根本无法走得很远。我见过的将这样的马尔可夫链组合在一起的唯一示例会产生关于递归太深的错误。这是我的示例代码。

任何帮助将非常感激。谢谢 ！

我正在更新一些我使用 pymc2 到 pymc3 的计算，当我的模型上有一些离散随机变量时，我遇到了采样器行为的一些问题。例如，考虑使用 pymc2 的以下模型：

它并不能真正代表任何东西，它只是一个模型，其中一个未观察到的变量是离散的。当我使用 pymc2 对该模型进行采样时，我得到以下结果：

但是当我对 PYMC3 进行同样的尝试时，我得到了这个：

看起来变量 A 根本没有被采样。我没有阅读太多关于 pymc3 中使用的采样方法的内容，但我注意到它似乎特别针对连续模型。这是否意味着它排除了模型上离散的未观察到的变量，或者有什么方法可以做我想做的事情？

假设我有（X，Y）点的随机集合：

并且我适合简单的线性回归：

如何获得y_est[0], y_est[1], y_est[2].. 的分布或统计数据（请注意，这些变量对应于每个输入y值的估计x值。

我想从数据中估计一个简单线性函数的参数和一个伽马分布的噪声项。（注意：这是https://stats.stackexchange.com/questions/88676/regression-with-unidirectional-noise的后续问题，但经过简化且更具体化）。假设我的观察数据生成如下：

如下所示：

我尝试使用 pymc 估计这些参数，如下所示：

但是，这给了我与真实值相去甚远的估计：

• 拦截：真：1.000，估计：3.281
• 斜率：真实：2.000，估计：-3.400

我做错了什么吗？

output = RBugsfit(..., coda=T, ...)输出一个 mcmc.list 对象，其中包含四个参数的后验分布样本及其样本后验均值。使用summary()我可以看到样本后验均值，但我想知道如何将样本后验均值从output程序中的变量中检索出来？谢谢！

我在执行模型的 MCMC 分析时遇到困难。我相信这与我在模型中的伽玛函数不完整有关。

我试图最小化高斯对数似然，但看起来步行者被困在他们的井里，而不是试图最小化似然函数。如下图所示，其中 y 轴是模型的参数，x 轴是步数。该图显示了步行者如何不探索参数空间。我添加了另一张图片来演示对参数空间的正确探索是什么样的。

参数空间的错误探索和参数空间的正确探索

我在下面添加了一些代码来演示我在做什么，其中 x、y 和 yerr 是大约 4000 个点的数组。该代码适用于其他模型，但仅在此模型上中断，因此它必须是其他模型所固有的。其他模型最明显的变化是增加了不完全伽马函数，否则它的函数形式与其他模型非常相似。

我正在拟合的模型具有以下形式：

注意我使用的是python包emcee（我会发布链接但据说我没有足够的声誉......）。我真的不明白为什么步行者拒绝为这个模型“走路”，而他们为其他模型做。非常感谢任何帮助，但我知道这是合理的利基领域。

我想知道为什么采样器在逐步采样时速度非常慢。例如，如果我运行：

采样速度很快。但是，如果我运行：

采样速度较慢（采样越多，速度越慢）。

如果您想知道我为什么要问这个.. 好吧，我需要逐步采样，因为我想在采样器的每一步之后对变量的值执行一些操作。

有没有办法加快速度？

先感谢您！

- - - - - - - - - 编辑 - - - - - - - - - - - - - - - - ------------------------------

在这里，我更详细地介绍了具体问题：

我有两个竞争模型，它们是更大模型的一部分，该模型具有一个分类变量，作为两者之间的“开关”。

在这个玩具示例中，我有观察到的向量“Y”，可以用泊松或几何分布来解释。分类变量“switch_model”在 = 0 时选择几何模型，在 =1 时选择泊松模型。

在每个样本之后，如果 switch_model 选择几何模型，我希望泊松模型的变量不要更新，因为它们不会影响可能性，因此它们只是逐渐消失。如果 switch_model 选择 Poisson 模型，则相反。

基本上，我在每一步所做的就是通过手动将其向后退一步来“更改”未选择模型的值。

我希望我的解释和注释代码足够清楚。如果您需要更多详细信息，请告诉我。

我正在使用多项式模型（我有 4 种模态），我想估计模型的所有参数（6 个参数，因为我有 2 个参数用于前三种模态）。我设法做到了。我现在的问题是我想对参数添加一个约束，以便我的第一种模态的概率随着我的协变量而降低x。

因此，我计算出我需要：

我正在使用runjags使用 JAGS 的 R 包。我希望每个样本都遵守这个约束。有可能这样做吗？

通常我使用 C/C++ 库，在每次迭代中，如果条件不满足，我会检查条件并重新采样此迭代的所有参数。但我找不到如何使用 JAGS 来做到这一点（出于协作原因，我不能使用 C）。

这是当前代码：

提前感谢您的帮助，

最好的

我正在尝试基于 Haskell 中的 Metropolis 算法编写一个 MCMC 程序，但我在从概率分布（生成伪随机数）采样和构建程序时遇到问题。现在我很高兴使用带有硬编码种子的生成器，而不是处理处理 IO 的复杂性。

看来我应该使用 state monad 来跟踪随机生成器状态、之前的马尔可夫链状态、卡方值和算法每一步之间的接受计数，然后最后收集所有马尔可夫链状态和最后一次接受数数。这是最好的/惯用的方法吗？如果是这样，程序的布局应该是什么（即提案功能和大都会阶梯功能的类型签名等）。

我见过一些处理随机数的示例程序，其中特定长度的随机数列表是从某种概率单子生成的，然后通过一些简单的函数进行线程化以执行计算。如果可能的话，我真的很想避免这种由内而外的程序形式。

编辑：暂时删除 WIP 代码。

我正在尝试使用 PyMC 学习简单离散 HMM 的参数。我正在对 HMM 上的Wiki页面中的雨晴模型进行建模。该模型如下所示：

我正在使用以下先验。

最初，我使用此设置创建一个训练集，如下所示。

但是，我无法理解如何使用 PyMC 学习这些参数。我开始如下。

可以在此处找到包含此代码的完整笔记本。

旁白：包含高斯 HMM 代码的要点真的很难理解！（未记录）

更新

根据以下答案，我尝试按如下方式更改我的代码：

下一步将是创建一个模型，然后运行 ​​MCMC 算法。但是，上面编辑的代码也不起作用。它给出了一个ZeroProbability error.

我不确定我是否正确解释了答案。