问题标签 [mcmc]

我正在使用 R2jags 和 CODA 为我的 MCMC 链生成一些诊断统计信息，但我遇到了麻烦。我想按如下方式运行 MCMC：

错误如下：

我正在使用 RStudio 0.97.551 和 R 版本 3.0.0 (2013-04-03)。

我很感激任何帮助！

这是我的 R 脚本：

假设我们有一个关于X的先验（例如 X ~ Gaussian）和一个前向算子y = f(x)。进一步假设我们通过一个实验观察了y，并且这个实验可以无限重复。假设输出Y为高斯 (Y ~ Gaussian) 或无噪声 (Y ~ Delta(observation))。

给定观察结果，如何持续更新我们对X的主观知识程度？我用 PyMC 尝试了以下模型，但似乎我遗漏了一些东西：

后验没有收敛到xtrue。

我有以下问题：

有 12 个样本，每个样本大约 20000 个元素来自未知分布（有时分布不是单峰的，因此很难自动估计分布的分析族）。基于这些分布，我计算出不同的数量。如何以最有效（和最简单）的方式探索目标数量的分布？

为了绝对清楚，这里有一个简单的例子：数量 A 等于 B*C/D

B，C，D是根据未知规律分布的，但是我有来自它们分布的样本，并且基于这些样本我想计算A的分布。所以实际上我想要的是一个工具来探索基于目标数量的分布关于变量的样本。

我知道有 MCMC 算法可以做到这一点。但是有人知道在 Python 或 C 中一个好的 MCMC 采样器实现吗？或者还有其他解决问题的方法吗？

格言

首先，让我说我缺乏科学数学或统计学方面的经验——所以这可能是一个非常有名的问题，但我不知道从哪里开始。

我有一个函数f(x1, x2, ..., xn)，我需要猜测 x'ses 并找到f. 该函数具有以下属性：

• 总数或参数通常在 40 到 60 左右，因此不可能采用蛮力方法。

• 每个 x 的可能值范围从 0.01 到 2.99

• 该函数是稳定的，这意味着更高的 f 值意味着对参数的猜测更好，反之亦然。

到目前为止，我在 python 中实现了一个非常基本的方法。它最初将所有参数设置为 1，随机猜测新值并检查 f 是否高于以前。如果不是，则回滚到以前的值。在具有 10,000 次迭代的循环中，这似乎以某种方式起作用，但结果可能远非完美。

我们将不胜感激有关如何改进对最佳参数的搜索的任何建议。当谷歌搜索这个问题时，linke MCMC 出现了，但这似乎是一种非常先进的方法，我什至需要很多时间才能理解该方法。基本的提示或概念比详细的方法和算法更能帮助我。

是否有可能将 pymc 用于以下任务：信号是我要估计的两个变量的非线性函数，噪声是非平稳的，在不同的时刻相关？噪声的分布是具有零均值和随时间变化的协变函数的高斯分布。信号是确定性的。

我是 pymc 的新手，不知道在这种情况下如何构造适当的似然函数。

将不胜感激任何建议。

我正在编写我的 R mcmc 代码的 C 版本。抛出错误的部分代码如下：

blaio.h 中的函数已经使用网站中的 makefile 构建，我是从那里得到的。

我正在使用以下内容进行编译：

它抛出以下错误：

我在 C 中没有嵌套函数。为什么会抛出这样的错误？如果是链接问题，它不应该抛出未找到函数的错误吗？

编辑：我被要求发布解决方案。在我的初始代码中，我没有编译上面的代码（而是另一个目录中的旧版本）。（是的，这是可能的最愚蠢的错误）。我正在重新发布代码的工作版本。除了查看 mcmc 采样步骤如何编码以及如何从 R 调用它的示例之外，我认为它没有太大价值。但是，就分组贝叶斯套索而言，此代码是不完整的。它还有代码（带有用于打印 blas 中使用的矩阵的来源）。我在大约一周内就开始使用 C 进行编码，所以对所有内容（打印代码除外）都持保留态度。如果您有更漂亮的版本（例如放置单独的文件以进行打印和编译命令），请继续进行更改。

我在这个练习中学到了 3 件事：

• 从 R 调用 C

• 使用 blas

• 使用 gsl

代码将在稍后演变为解决第 4 点： - 从 C 调用 R

我必须在此模型中从广义逆高斯采样（此代码中排除）。但是，我看不到在 C 中执行此操作的简单函数。那时我可能不得不使用 R 的函数。所以变得很奇怪。R 调用此 C 函数，并且在此 C 函数中将调用其他一些 R 函数。如果有任何替代方法，您的意见将不胜感激，如果是这样，怎么做 - RInside？有什么容易的吗？

时间对我来说是最重要的。我必须用大约一百万个数据点运行这个模型，至少有 10,000 次 mcmc 迭代。我想在 5 分钟内完成，最好是 2 分钟。

本练习的第三阶段是并行化 z 中的采样步骤。

在 R 中模拟随机数据并调用 C 函数：

所以我必须运行在 R 中不可行的大型 mcmc。因此我在 C 中写出循环并使用 .C 从 R 中调用它

然而，在这个循环中，我偶尔会有采样步骤或建模步骤，其中有 R 中的标准包。从头开始用 C 编写这些东西对我来说太过分了。

这就是为什么我想在这里调用 R 函数。我将在“编写 R 扩展”中使用“8.1 Embedding R under Unix-likes”来做到这一点。

所以它变成了，R 被嵌入到 C 中，这本身就是来自 R 的 .C 调用。就以下方面而言，这种方法是否可行

• 简单，还是有更优雅的框架？

• 稳定性（R调用C，后者调用R）

• 可扩展性（R 调用 C，它调用 R），我将在亚马逊上运行多个这样的 mcmc 采样器？

• 内存使用情况？

我只对嵌入 C（不是 c++）的 R 感兴趣。速度对我来说非常重要。

我遇到了一个技术问题，看来我自己无法解决。我使用mcmcglmm包进行了估算。通过 results$Sol 我可以访问估计的后验分布。应用class()告诉我该对象属于“mcmc”类。使用as.data.frame()会生成一个嵌套数据框，其中包含其他数据框（一个数据框包含许多其他数据框）。我想rbind()主数据帧中的所有数据帧，以便生成一个数据帧（或者更确切地说是一个向量），其中包含所有后验分布的所有值以及（辅助）数据帧的名称作为行名。，有任何想法吗？我将不胜感激每一个提示！

更新：我没有设法为 stackoverflow 生成有用的数据集，对于所有这些采样链，这些数据集总是太大。如果您想帮助我，请考虑运行以下（示例）模型

result$Sol（一个 mcmc 对象）是所有链的存储位置。我想 rbind 所有链，以便拥有一个向量，其中所有后验分布的所有值和变量名称作为行名（或者因为不允许重复的行名，作为附加字符向量）。

我对 PyMC 比较陌生，我有一个关于 MCMC 采样器输出的快速问题。我希望找到 MCMC 采样器发现的变量的最可能值（后验的最大值）。有没有快速的方法来做到这一点？据推测，MCMC 采样器发现的后验最大值处的变量值可能与 PyMC 的 MAP 方法发现的变量值大不相同。

感谢开发者提供 PyMC。它对我的工作非常有用。这个问题也已经发布到 PyMC google group。

我正在尝试为我目前使用 R 包在 R 中运行的 JAGS 分析计算 Gelman 和 Rubin 的收敛诊断rjags。

例如，我想评估我的参数的收敛诊断beta。为此，我使用库 coda 和命令：

具有out_2MCMC多个链的 MCMC 列表对象，从而导致正确的输出而没有错误消息或任何内容。但是，由于我使用大量迭代作为老化，我想仅计算迭代子集的收敛诊断（仅老化后的部分！）。

为此，我尝试了：

这导致了以下错误：

因此，我尝试了：

但是，令人惊讶的是，这导致了以下错误：

由于这与我从 JAGS 分析中获得的 MCMC 列表对象相同，并且与我在评估所有迭代的收敛诊断时使用的相同（它工作得非常完美），我在这里看不到问题。

该函数本身仅提供在收敛诊断计算中使用系列的后半部分（迭代）的选项。由于我的老化阶段比这更长，不幸的是这对我来说还不够。

我想这是很明显的事情，我只是想念。有没有人有任何建议或提示？

由于代码很多，我没有提供 R 代码来运行完整的 2MCMC-JAGS 分析。我希望上面的代码能够很好地说明问题，也许有人以前遇到过同样的问题，或者识别出我的语法中的任何错误。但是，如果您觉得需要完整的代码来理解我的问题，我仍然可以提供运行 2MCM JAGS 分析的示例代码。