PyMC 中的自定义先验

本文介绍了PyMC 中的自定义先验的处理方法，对大家解决问题具有一定的参考价值，需要的朋友们下面随着小编来一起学习吧！ 问题描述

假设我想在 PyMC 中的两个变量 a 和 b 上放置一个自定义先验，例如:

Say I want to put a custom prior on two variables a and b in PyMC, e.g.:

p(a,b)∝(a+b)^(−5/2)

(关于这种先验选择背后的动机，请参阅这个答案)

(for the motivation behind this choice of prior, see this answer)

这可以在 PyMC 中完成吗?如果是这样怎么办?

Can this be done in PyMC? If so how?

举个例子，我想在下面的模型中在 a 和 b 上定义这样的先验.

As an example, I would like to define such prior on a and b in the model below.

import pymc as pm # ... # Code that defines the prior: p(a,b)∝(a+b)^(−5/2) # ... theta = pm.Beta("prior", alpha=a, beta=b) # Binomials that share a common prior bins = dict() for i in xrange(N_cities): bins[i] = pm.Binomial('bin_{}'.format(i), p=theta,n=N_trials[i], value=N_yes[i], observed=True) mcmc = pm.MCMC([bins, ps])

更新

按照 John Salvatier 的建议，我尝试了以下操作(请注意，我在 PyMC2，虽然我很乐意切换到 PyMC3)，但我的问题是:

Update

Following John Salvatier's advice, I try the following (note that I am in PyMC2, although I would be happy to switch to PyMC3), but my questions are:

我应该导入什么才能正确继承Continuous?

在 PyMC2 中，我还需要坚持使用 Theano 表示法吗?

最后，我以后如何告诉我的 Beta 分布 alpha 和 beta 有一个来自这个多变量分布的先验?

What should I import so that I can properly inherit from Continuous?

In PyMC2, do I still need to stick to Theano notation?

Finally, how can I later tell my Beta distribution that alpha and beta have a prior from this multivariate distribution?

导入 pymc.Multivariate.Continuous

import pymc.Multivariate.Continuous

类 CustomPrior(Continuous):"""p(a,b)∝(a+b)^(−5/2)

class CustomPrior(Continuous): """ p(a,b)∝(a+b)^(−5/2)

:Parameters: None :Support: 2 positive floats (parameters to a Beta distribution) """ def __init__(self, mu, tau, *args, **kwargs): super(CustomPrior, self).__init__(*args, **kwargs) def logp(self, a,b): return np.log(math.power(a+b),-5./2)

推荐答案

在 PyMC2 中，诀窍是将 a 和 b 参数放在一起:

In PyMC2, the trick is to put the a and b parameters together:

# Code that defines the prior: p(a,b)∝(a+b)^(−5/2) @pm.stochastic def ab(power=-2.5, value=[1,1]): if np.any(value <= 0): return -np.Inf return power * np.log(value[0]+value[1]) a = ab[0] b = ab[1]

这个笔记本有一个完整的例子.