首页 > 解决方案 > 使用 NUTS 初始化的 PyMC3 贝叶斯推理

问题描述

我正在尝试使用 ODE 模型实现简单的贝叶斯推理。我想使用 NUTS 算法进行采样,但它给了我一个初始化错误。我对 PyMC3 不太了解,因为我是新手。请看一下并告诉我出了什么问题。

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
import seaborn
import pymc3 as pm
import theano.tensor as T
from theano.compile.ops import as_op

#Actual Solution of the Differential Equation(Used to generate data)
def actual(a,b,x):
    Y =  np.exp(-b*x)*(a*np.exp(b*x)*(b*x-1)+a+b**2)/b**2
    return Y

#Method For Solving the ODE
def lv(xdata, a=5.0, b=0.2):
    def dy_dx(y, x):
        return a*x - b*y
    y0 = 1.0
    Y, dict  = odeint(dy_dx,y0,xdata,full_output=True)
    return Y

#Generating Data for Bayesian Inference
a0, b0 = 5, 0.2
xdata = np.linspace(0, 21, 100)
ydata = actual(a0,b0,xdata)

# Adding some error to the ydata points
yerror = 10*np.random.rand(len(xdata))
ydata += np.random.normal(0.0, np.sqrt(yerror))
ydata = np.ravel(ydata)

@as_op(itypes=[T.dscalar, T.dscalar], otypes=[T.dvector])

def func(al,be):
    Q = lv(xdata, a=al, b=be)
    return np.ravel(Q)

# Number of Samples and Initial Conditions
nsample = 5000
y0 = 1.0

# Model for Bayesian Inference
model = pm.Model()
with model:
    # Priors for unknown model parameters
    alpha = pm.Uniform('alpha', lower=a0/2, upper=a0+a0/2)
    beta = pm.Uniform('beta', lower=b0/2, upper=b0+b0/2)

    # Expected value of outcome
    mu = func(alpha,beta)


    # Likelihood (sampling distribution) of observations
    Y_obs = pm.Normal('Y_obs', mu=mu, sd=yerror, observed=ydata)

    trace = pm.sample(nsample, nchains=1)


pm.traceplot(trace)
plt.show()

我得到的错误是

Auto-assigning NUTS sampler...
Initializing NUTS using jitter+adapt_diag...
Initializing NUTS failed. Falling back to elementwise auto-assignment.

任何帮助将非常感激

标签: python-2.7bayesianpymc3mcmcstatistical-sampling

解决方案


推荐阅读