问题标签 [mcmc]

我目前正在尝试将一些 JAGS 模型移植到 STAN。我收到一些奇怪的错误“stan::prob::exponential_log(N4stan5agrad3varE): Random variable is nan:0, but must not be nan!” 并调试那些我想知道一些本地参数的值。

在 JAGS 中，我可以为任何变量设置监视器。STAN 仅监控参数。但是参数不能有赋值（如果我理解正确的话）。

那么如何监控中间变量呢？

我还粘贴了模型代码，以防有人看到我犯的愚蠢错误。但是请注意，我知道可以将相同的模型表示为双指数（具有两个速率）的 CDF。这是我计划的简化形式。

这是一些虚拟数据：

我试图在 Stan 中运行这个模型。我有一个正在运行的 JAGS 版本（返回高度自相关的参数），并且我知道如何将其公式化为双指数的 CDF（具有两个速率），这可能会毫无问题地运行。但是，我想将此版本用作类似但更复杂模型的起点。

到目前为止，我怀疑这样的模型在 Stan 中是不可能的。可能由于采用布尔值之和引入的离散性，Stan 可能无法计算梯度。

有谁知道是这种情况，还是我在这个模型中以错误的方式做了其他事情？我将得到的错误粘贴到模型代码下方。

非常感谢提前一月

这是一些虚拟数据：

以下是错误：

在运行时：

这是此模型的工作 JAGS 版本：

关于 min() 和 max()：见这篇文章https://stats.stackexchange.com/questions/130978/observed-node-inconsistent-when-binomial-success-rate-exactly-one?noredirect=1 #comment250046_130978。

即使我使用pickle后端，在Python中遵循明显简单的MCMC代码也会导致大量内存使用（> 15GB）。每当我在 pymc 中使用观察到的变量数组时，都会发生这种情况。知道为什么会这样吗？

我不确定如何在 R 中保存 coda (mcmc.list) 对象。其他人也问过类似的问题，但我发现给出的答案并不是特别清楚。理想情况下，我想将 coda 对象保存为 R.data 文件或文本文件（例如 csv），这样我就可以重新导入它并分析 JAGS 链，而无需重新运行模型（大约需要 30 分钟）我的电脑）。现在我的尾声对象“coda.samples”看起来像这样：

如您所见，它是三个矩阵的列表，其中包含 1094 个参数的 3334 个估计值（即 3 个长度为 3334 的链）。我想存储这个 coda 对象，这样我就可以将它调用回 R 中，而不必每次都重新运行模型。我还想保留三个独特链的事实。

我有一个关于 Metropolis-Hastings 算法的简单问题。假设分布只有一个变量x，x的取值范围为s=[-2^31,2^31]。

在抽样过程中，我需要提出一个新的 x 值，然后决定是否接受。

如果我想自己实现它，如何确定\epsilon的值。

基本的解决方案是从 Uniform[-2^31,2^31] 中选择一个值并将其设置为 \epsilon。如果值范围像 [-inf, inf] 那样无界怎么办？

当前的 MCMC 库（例如 pymc）如何解决这个问题？

下面是定义两个随机伯努利随机变量的一种方法，一个依赖于另一个带有装饰器的变量。该模型旨在：

在 PyMC 中使用装饰器是：

这是最正确和最紧凑的方法吗？是否可以通过将节点显式定义为pymc.Bernoulli变量而不是使用stochastic装饰器来保存一些代码？

上面的代码似乎是多余的，必须始终定义如何从伯努利 RV 中采样，但这也许是不可避免的？

最后是否有必要从numpylike调用采样函数，np.random.binomial或者 PyMC 中是否有采样函数？

鉴于以下模型，我的问题是如何S知道关于alpha、beta和的任何信息theta？我已经看到MCMC了在单独的文件（即作为 Python 模块）中指定模型的示例，这对我来说很有意义。但是在这里我没有S明确地传递任何数据。只是想了解这是如何工作的。

我有一些启用 MPI 的 python MCMC 采样代码，可以触发对单独内核的并行可能性调用。因为它（必然 - 不要问）拒绝采样，所以我只需要一个 np 样本就可以成功开始下一次迭代，并且过去通过这种方法非常高兴地实现了 ~ np x 加速。

我已将此应用于一个新问题，其中可能性调用 f2py 包装的 fortran 子例程。在这种情况下，在每次迭代中，其他 np-1 进程等待最慢（有时非常慢）的结果返回，即使其中一个 np-1 已经可以接受。

所以我怀疑我需要将一条消息传递给所有非获胜（在速度方面）的进程以终止，以便下一次迭代可以开始，并且我需要明确执行此操作的最佳方法的一些细节，如下所示。

python代码是这样的。采样器是 PyMultiNEST。

广播应该通过主进程吗？

现在，棘手的部分是如何在 F90 代码中接收终止信号。大概如果代码总是在监听（while循环？）它会减慢很多 - 但我无论如何应该使用类似的东西：

然后如何在收到消息后最好地终止该进程？

最后，我是否需要在 F 代码中做任何事情以使下一次迭代重新启动 OK/产生新进程？

谢谢！

我目前正在用 R 编写一个 MCMC 程序来估计 Rasch 模型参数。为此，我在 Gibbs 采样器中使用了 Metropolis-hastings 算法。

在下面的代码中，给出了项目参数的提案函数的一部分。

我的问题是使用if(runif(1)<= exp(d.l.p[i]-d.l.c[i]))来确定是否应该接受提案值是否正确？我知道在非对数似然的情况下，您可以使用if(runif(1)<= d.l.p[i]/d.l.c[i])来确定新值的接受度。

由于这更像是一个概念问题而不是编码问题，因此我省略了其余代码。但是，如果需要所有代码，我很乐意提供。

提前Tnx！约斯特

想象一个简单的增长模型。我如何获得领先一步的预测？