1

我试图让一个简单的 PyMC2 模型在 PyMC3 中工作。我已经让模型运行,但模型对变量给出了非常不同的 MAP 估计。这是我的 PyMC2 模型:

import pymc
theta = pymc.Normal('theta', 0, .88)

X1 = pymc.Bernoulli('X2', p=pymc.Lambda('a', lambda theta=theta:1./(1+np.exp(-(theta-(-0.75))))), value=[1],observed=True)
X2 = pymc.Bernoulli('X3', p=pymc.Lambda('b', lambda theta=theta:1./(1+np.exp(-(theta-0)))), value=[1],observed=True)

model = pymc.Model([theta, X1, X2])
mcmc = pymc.MCMC(model)
mcmc.sample(iter=25000, burn=5000)
trace = (mcmc.trace('theta')[:])
print "\nThe MAP value for theta is", trace.sum()/len(trace)

这似乎按预期工作。我在弄清楚如何在 PyMC3 中使用等效的 pymc.Lambda 对象时遇到了各种各样的麻烦。我最终遇到了确定性对象。以下是我的代码:

import pymc3

with pymc3.Model() as model:

    theta = pymc3.Normal('theta', 0, 0.88)
    X1 = pymc3.Bernoulli('X1', p=pymc3.Deterministic('b', 1./(1+np.exp(-(theta-(-0.75))))), observed=[1])
    X2 = pymc3.Bernoulli('X2', p=pymc3.Deterministic('c', 1./(1+np.exp(-(theta-(0))))), observed=[1])

    start=pymc3.find_MAP()
    step=pymc3.NUTS(state=start)
    trace = pymc3.sample(20000, step, njobs=1, progressbar=True)

pymc3.traceplot(trace)

我遇到的问题是我使用 PyMC2 对 theta 的 MAP 估计值约为 0.68(正确),而 PyMC3 给出的估计值约为 0.26(不正确)。我怀疑这与我定义确定性函数的方式有关。PyMC3 不允许我使用 lambda 函数,所以我只需要内联编写表达式。当我尝试使用 lambda theta=theta:... 我得到这个错误:

AsTensorError: ('Cannot convert <function <lambda> at 0x157323e60> to TensorType', <type 'function'>)

和Theano有什么关系??任何建议将不胜感激!

4

4 回答 4

2

当您在 Deterministic 中使用 theanotensor而不是 numpy 函数时,它会起作用。

import pymc3
import theano.tensor as tt

with pymc3.Model() as model:

    theta = pymc3.Normal('theta', 0, 0.88)
    X1 = pymc3.Bernoulli('X1', p=pymc3.Deterministic('b', 1./(1+tt.exp(-(theta-(-0.75))))), observed=[1])
    X2 = pymc3.Bernoulli('X2', p=pymc3.Deterministic('c', 1./(1+tt.exp(-(theta-(0))))), observed=[1])

    start=pymc3.find_MAP()
    step=pymc3.NUTS(state=start)
    trace = pymc3.sample(20000, step, njobs=1, progressbar=True)

print "\nThe MAP value for theta is", np.median(trace['theta'])

pymc3.traceplot(trace);

这是输出:

在此处输入图像描述

于 2015-08-31T12:39:17.540 回答
0

MAP 值不定义为分布的平均值,而是定义为最大值。你可以通过pymc2以下方式找到它:

M = pymc.MAP(model)
M.fit()
theta.value

返回array(0.6253614422469552)

find_MAP这与您在 中找到的 MAP 一致,pymc3您称之为start

{'theta': array(0.6253614811102668)}

哪个是更好的采样器的问题是不同的,并且不依赖于 MAP 的计算。MAP 计算是一种优化。请参阅:https ://pymc-devs.github.io/pymc/modelfitting.html#maximum-a-posteriori-estimates为pymc2.

于 2015-08-31T21:41:37.363 回答
0

此外,NUTS 不处理离散变量。如果要使用 NUTS,则必须拆分采样器:

step1 = pymc3.NUTS([theta])
step2 = pymc3.BinaryMetropolis([X1,X2])

trace = pymc3.sample(10000, [step1, step2], start)

编辑:错过了 'b' 和 'c' 是内联定义的。从 NUTS 函数调用中删除它们

于 2015-08-31T16:21:57.600 回答
0

以防其他人有同样的问题,我想我找到了答案。在尝试了不同的采样算法后,我发现:

  • find_MAP 给出了错误的答案
  • NUTS 采样器给出了错误的答案
  • Metropolis 采样器给出了正确答案,耶!

我在其他地方读到 NUTS 采样器不适用于 Deterministic。我不知道为什么。也许 find_MAP 也是这种情况?但现在我会坚持使用 Metropolis。

于 2015-08-31T04:54:45.807 回答