1

我目前正在尝试使用 PyMC 进行模型检查，其中我的模型是伯努利模型，并且我之前有 Beta。我想做一个（i）gof图以及（ii）计算后验预测p值。

我已经让我的代码使用二项式模型运行，但我很难找到使伯努利模型工作的正确方法。不幸的是，在任何地方都没有我可以使用的示例。我的代码如下所示：

import pymc as mc
import numpy as np
alpha = 2
beta = 2
n = 13
yes = np.array([1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0,0,0])

p = mc.Beta('p',alpha,beta)
surv = mc.Bernoulli('surv',p=p,observed=True,value=yes)
surv_sim = mc.Bernoulli('surv_sim',p=p)

mc_est = mc.MCMC({'p':p,'surv':surv,'surv_sim':surv_sim})
mc_est.sample(10000,5000,2)

import matplotlib.pylab as plt
plt.hist(mc_est.surv_sim.trace(),bins=range(0,3),normed=True)
plt.figure()
plt.hist(mc_est.p.trace(),bins=100,normed=True)

mc.Matplot.gof_plot(mc_est.surv_sim.trace(), 10/13., name='surv')

#here I have issues
D = mc.discrepancy(yes, surv_sim, p.trace())
mc.Matplot.discrepancy_plot(D)

我遇到的主要问题是确定discrepancy函数的预期值。仅使用p.trace()在这里不起作用，因为这些是概率。不知何故，我需要在这里合并样本量，但我正在努力以与二项式模型类似的方式来做到这一点。我也不太确定，如果我做得gof_plot正确。

希望有人可以在这里帮助我！谢谢！

4

1 回答 1

0

根据discrepancy函数doc string，参数是：

observed : Iterable of observed values (size=(n,))
simulated : Iterable of simulated values (size=(r,n))
expected : Iterable of expected values (size=(r,) or (r,n))

所以你需要纠正两件事：

1）修改您的模拟结果以具有大小n（例如，在您的示例中为 13）：

surv_sim = mc.Bernoulli('surv_sim', p=p, size=n)

p.trace()2）用方法封装你bernoulli_expval：

D = mc.discrepancy(yes, surv_sim.trace(), mc.bernoulli_expval(p.trace()))

（bernoulli_expval只是吐回去p。）

通过这两个更改，我得到以下信息：

在此处输入图像描述

于 2015-08-17T21:13:51.093 回答