2

我正在将一个非常简单的示例从 pymc 2.3 转换为 pymc 3.0,并且似乎无法弄清楚如何从预测后验分布中采样(或获取 MAP)。按照文档中的建议(7.3 模型检查和诊断:拟合优度),我可以使用 pymc 2.3 通过添加未观察到的随机变量从该分布中采样。这是笔记本的链接。一切似乎运作良好。

但是,当我尝试在pymc 3.0中执行此操作时,会发生两件奇怪的事情。

  1. MAP 值甚至不接近,好像未观察到的随机性正在影响最小化,而未观察到的随机性的 MAP 值是错误的?
  2. NUTS 采样器不会改变未观察到的随机变量的值,因此轨迹只是单个值 10。

当然,如果观察到的随机变量不存在,这个值是有意义的。在 pymc 3.0 中,为预测后验找到 MAP 值并从中采样有何变化?

更新:

这是一个最小的示例,说明了我做错了什么:

import pymc as mc
with mc.Model() as model:
    p = mc.Beta('p',2,2)
    surv_sim = mc.Binomial('surv_sim',n=20,p=p)
    surv = mc.Binomial('surv',n=20,p=p,observed=15)

with model:
    step = mc.step_methods.HamiltonianMC(vars=model.vars) #Again must specificy
                                                          #model.vars or else
                                                          #only continuous values will
                                                          #be sampled

with model:
    trace = mc.sample(100000,step)

mc.traceplot(trace);

with model:
    map_est = mc.find_MAP(vars=model.vars) #Must Specify model.vars, or else
                                           #only continuous stochastics will
                                           #be fit, however this fit will be 
                                           #horrible

我想我自己已经部分回答了这个问题。首先是抽样。据我所知,NUTS 采样器无法处理离散变量。然而,Metropolis 采样器肯定可以处理离散变量(已在 issue #235 中解决),HamiltonianMC 采样器似乎也可以处理离散变量。

然而,离散随机指标的 MAP 估计注定是糟糕的。即使您指定使用所有变量(不是默认值)。因为离散随机对数概率函数返回任何离散随机的底的对数概率,所以无论使用哪个 scipy 最小化函数,都会陷入局部最优。最小化只会优化连续随机变量,同时保持离散变量固定,因为离散变量中的任何小步骤都不会导致对数概率的提高。我不确定这是否是一个错误,或者只是当您同时具有离散变量和连续变量时查找 MAP 估计值的基本限制。

更新 2:

如果您使用多个步骤,采样效果会更好。对于上述模型,您可以使用。

with model:
    step1 = mc.step_methods.NUTS(vars=[p])
    step2 = mc.step_methods.Metropolis(vars=[surv_sim])

with model:
    trace = mc.sample(100000,[step1,step2])
4

1 回答 1

0

最近有很多关于 PyMC3 的开发活动,所以我想我会为你的问题发布一个小答案,即使你在一年前就已经知道了这一点。

首先,该pymc3.find_MAP函数现在尝试处理离散变量,因此要找到 MAP,您只需执行以下操作:

import pymc3 as pm

m = pm.Model()
with m:
    p        = pm.Beta('p',2,2)
    surv_sim = pm.Binomial('surv_sim', n=20, p=p)
    surv     = pm.Binomial('surv', n=20, p=p, observed=15)

    map_est = pm.find_MAP()

map_est

对我来说,这会导致以下结果

{'p': array(0.6190476210094416), 'surv_sim': 10.0} 

这似乎不正确。也许这值得一个错误报告。

从后验分布中采样现在需要明确选择步进方法,并且对于离散变量pymc3.step_methods.Metropolis是可接受的选择。它也适用于本例中的连续变量:

with m:
    step1 = pm.step_methods.Metropolis(vars=[p])
    step2 = pm.step_methods.Metropolis(vars=[surv_sim])
    trace = pm.sample(10000, [step1, step2])
于 2015-06-18T04:12:42.633 回答