0

我有这个用winBUGS编写的代码：

n <- 100

x1 <- rbinom(n,1,.7)

x2 <- rbinom(n,1,.5)

sum(x1)

sum(x2)

model{

x1 ~ dbin(p1, n) x2 ~ dbin(p2, n) p1 ~ dbeta(a1, b1) p2 ~ dbeta(a2,b2)

diff <- p1 - p2 p.value <- step(diff)

} list(n = 100, x1 = 70, x2 = 54, a1 = 1, b1 = 1, a2 = 1, b2 = 1)

我在如何在 R/JAGS 中执行此操作时遇到问题。事实上，我什至不完全确定这段代码试图做什么（我认为计算后验？）。之前没用过winBUGS，对R也是新手。这也是我第一次上贝叶斯课，一介绍代码就迷茫了。

此外，我将如何计算比例差异的后验均值和标准差？或者如何找到p1大于的后验概率p2，以及是否显着？

4

1 回答 1

2

如果您在设置rjags时需要帮助，我希望您的助教或同学可以帮助您。这个答案的重点是解释代码的作用。

出于教学目的，让我们对这些数据进行简单的叙述。假设我们有两枚硬币，我们想知道一枚硬币是否比另一枚更有可能翻转正面。代码分为两部分。

生成/正向模拟

第一个是前向模拟，我们已经知道两枚硬币正面朝上的概率是多少。

## how many times we will flip each coin
n <- 100

## coin 1 has 70% chance to land on heads
## simulate n flips
x1 <- rbinom(n, 1, .7)

## coin 2 has 50% chance to land on heads
## simulate n flips
x2 <- rbinom(n, 1, .5)

## count number of coin 1 heads
sum(x1) # 70

## count number of coin 2 heads
sum(x2) # 54

后验概率

现在我们采用这些模拟数据并尝试反转实验。也就是说，如果我们观察到一枚硬币有 70 个正面，另一枚硬币有 54 个正面，我们能说一下每枚硬币翻转正面的概率吗？具体来说，代码会问，“在关于每枚硬币的真实概率的统一假设下，硬币 1 比硬币 2 更有可能正面朝上的（后验）概率是多少？ ”

计算可能性

JAGS 将使用 MCMC 执行此操作，其中样本取自所有可能的参数配置（在本例中为p1和的值p2）的空间，并根据每个配置生成观察数据的可能性进行加权。因此，JAGS 代码需要完成的主要事情是定义如何生成参数值，以及如何在给定这些值的情况下计算数据的可能性。这是模型代码的第一部分。

计算感兴趣的值

模型部分代码的第二部分是测试我们提出的问题。这是通过计算一些对配置的可能性没有贡献的额外变量来完成的。相反，它们会告诉我们有关正在采样的配置的信息。具体来说，diff将跟踪两个硬币概率之间差异的分布情况，p.value并将跟踪p1大于 的频率p2。

model{

## COMPUTE LIKELIHOODS
## compute likelihood that heads resulted from coins
## with given probabilities after `n` coin flips
x1 ~ dbin(p1, n) 
x2 ~ dbin(p2, n) 

## SAMPLE PARAMETERS
## randomly select coin probabilities from Beta distributions
## Note that since these are all 1, it's really just a Uniform[0,1]
p1 ~ dbeta(a1, b1)
p2 ~ dbeta(a2, b2)

## HYPOTHESIS TEST
## compute how often coin1's probability is greater than coin2's
diff <- p1 - p2 
p.value <- step(diff)
} 

## INPUT VALUES
list(n = 100, x1 = 70, x2 = 54, a1 = 1, b1 = 1, a2 = 1, b2 = 1)
于 2019-03-06T16:47:46.360 回答