pymc3 中的直接对数概率增量目标 += ...，如 stan

如何解决pymc3 中的直接对数概率增量目标 += ...，如 stan

在概率编程语言stan中，给定数据/参数块：

data {
    int N;
    real[N] data;
}

params {
    real mu;
}

以下模型块是等效的：

“示例符号”：

model { 
    data ~ normal(mu,1);
}

“直接对数概率增量表示法”：

model {
    for (n in 1:N)
        target += normal_lpdf(data[n] | mu,1);
}

其中 target 表示由 MCMC (NUTS) 采样器随机更新的总对数密度。使用后一种表示法的好处是增加了如何定义对数概率模型的灵活性，特别是可以提供模型通过计算生成的样本（在上面的示例中 data[n] 但这也可以在更多上下文中使用).

示例符号可以在 pymc3 (python) 中应用为：

with pm.Model() as model:
    mu = pm.Flat('mu')
    x = pm.normal('x',mu=mu,sd=1.0,observed=data)

问题：如何在 pymc3 中应用相同的直接对数概率增量表示法，在其中指定样本？

解决方法

您可以使用 pm.Potential 运行它，如下所示：

import numpy as np
import pymc3 as pm

data = 5 + 3 * np.random.randn(20)


with pm.Model() as model:
    x = pm.Flat('x')
    pm.Potential('y',pm.Normal.dist(x,1).logp(data))

您可以通过以下方式验证这样做是否正确：

import scipy.stats as st

print(st.norm(0,1).logpdf(data).sum(),model.logp({'x': 0}))

请注意，Potential 必须包含一个 theano 表达式——pm.Normal.dist(x,1) 恰好在 theano 中实现。您也可以明确地编写日志 pdf：

import theano.tensor as tt

with pm.Model() as model:
    x = pm.Flat('x')
    pm.Potential('y',-0.5 * tt.sum((x - data) ** 2) 
                     - 0.5 * len(data) * np.log(2 * np.pi))