与 stan 一样,pymc3 中的直接对数概率增量(目标 += ...)

direct log probability increment (target += ...) in pymc3 as in stan

在概率编程语言 stan 中,给定 data/params 个块:

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))