1

在我的模型中，我需要使用复杂的 python 函数从一组父变量中获取确定性变量的值。

有可能这样做吗？

以下是一个 pyMC3 代码，它显示了我在简化情况下尝试做的事情。

import numpy as np
import pymc as pm


#Predefine values on two parameter Grid (x,w) for a set of i values (1,2,3)
idata = np.array([1,2,3])
size= 20
gridlength = size*size
Grid = np.empty((gridlength,2+len(idata)))
for x in range(size):
    for w in range(size):
        # A silly version of my real model evaluated on grid.
        Grid[x*size+w,:]= np.array([x,w]+[(x**i + w**i) for i in idata])

# A function to find the nearest value in Grid and return its product with third variable z
def FindFromGrid(x,w,z):
    return Grid[int(x)*size+int(w),2:] * z

#Generate fake Y data with error
yerror = np.random.normal(loc=0.0, scale=9.0, size=len(idata))
ydata = Grid[16*size+12,2:]*3.6 + yerror  # ie. True x= 16, w= 12 and z= 3.6


with pm.Model() as model:

    #Priors
    x = pm.Uniform('x',lower=0,upper= size)
    w = pm.Uniform('w',lower=0,upper =size)
    z = pm.Uniform('z',lower=-5,upper =10)

    #Expected value
    y_hat = pm.Deterministic('y_hat',FindFromGrid(x,w,z))

    #Data likelihood
    ysigmas = np.ones(len(idata))*9.0 
    y_like = pm.Normal('y_like',mu= y_hat, sd=ysigmas, observed=ydata)

    # Inference...
    start = pm.find_MAP() # Find starting value by optimization
    step = pm.NUTS(state=start) # Instantiate MCMC sampling algorithm
    trace = pm.sample(1000, step, start=start, progressbar=False) # draw 1000 posterior samples using NUTS sampling


print('The trace plot')
fig = pm.traceplot(trace, lines={'x': 16, 'w': 12, 'z':3.6})
fig.show()

当我运行此代码时，在 y_hat 阶段出现错误，因为int()函数内部的FindFromGrid(x,w,z)函数需要整数而不是 FreeRV。

从预先计算的网格中查找y_hat很重要，因为我的 y_hat 真实模型没有要表达的分析形式。

我之前曾尝试使用 OpenBUGS，但在这里我发现在 OpenBUGS 中无法做到这一点。PyMC 有可能吗？

更新

基于 pyMC github 页面中的示例，我发现我需要将以下装饰器添加到我的FindFromGrid(x,w,z)函数中。

@pm.theano.compile.ops.as_op(itypes=[t.dscalar, t.dscalar, t.dscalar],otypes=[t.dvector])

这似乎解决了上述问题。但我不能再使用 NUTS 采样器了，因为它需要渐变。

大都会似乎没有融合。

在这种情况下，我应该使用哪种步骤方法？

4

1 回答 1

1

您找到了正确的解决方案as_op。

关于收敛：您是否正在使用pm.Metropolis()而不是pm.NUTS()偶然？这不能收敛的一个原因是，Metropolis()默认情况下，联合空间中的样本，而 Metropolis 中的 Gibbs 通常更有效（这是 pymc2 中的默认设置）。话虽如此，我刚刚合并了这个：httpsMetropolis ://github.com/pymc-devs/pymc/pull/587 它将默认情况下和采样器的默认行为更改为Slice非阻塞（因此在 Gibbs 内）。其他类似的采样器NUTS主要设计用于对关节空间进行采样，但仍默认为阻塞。您始终可以使用 kwarg 显式设置它blocked=True。

无论如何，用最新的 master 更新 pymc，看看收敛性是否有所改善。如果没有，请尝试Slice采样器。

于 2014-10-18T09:23:08.787 回答