# Translating user-defined joint-distribution from PyMC to PyMC3

I'm attempting to set up a simple beta binomial hierarchical model with an uninformative prior in PyMC3. I've read that the uninformative prior for this model should have alpha and beta hyper-parameters sampled from the following function: $$f(a,b) \sim (a+b)^{-5/2}$$ where a, b > 0. I stumbled across a way to represent this in PyMC (first block of code below), but am not quite sure how to express it in PyMC3 (attempt in second block of code below).

Given that I don't think my PyMC3 attempt is correct, I'd like to know the right way to set this up.

#This is the function I'm strugling to translate to PyMC3
@pymc.stochastic(dtype=np.float64)
def hyperpriors(value=[1.0, 1.0]):
a, b = value[0], value[1]
if a <= 0 or b <= 0:
return -np.inf
else:
return np.log(np.power((a + b), -2.5))

#hyperpriors
a = hyperpriors[0]
b = hyperpriors[1]

#prior
true_rates = pymc.Beta('true_rates', a, b, size=10)

#likelihood
observed_values = pymc.Binomial('observed_values', trials, true_rates, observed=True, value=successes)

# This is what we observed
trials = np.array([100, 100, 100, 100, 100, 100, 100, 100, 100, 100])
successes = np.array([40, 44, 47, 54, 63, 46, 44, 49, 58, 50])

model = pymc.Model([a, b, true_rates, observed_values])
mcmc = pymc.MCMC(model)


Here is my attempt in PyMC3. I end up with the following error:

'ValueError: length not known: Elemwise{le,no_inplace} [id A] ''
|ab [id B] |DimShuffle{x} [id C] ''
|TensorConstant{0} [id D]'.

I'm new to Theano and am struggling to debug this one, even after reading posts with similar error messages.

with pm.Model() as test:

#joint distribution for uninformative prior
def ab_dist(value=[1.0,1.0]):
return T.switch(any(T.le(value, 0)), -np.Inf, T.log(np.power((value[0] + value[1]), -2.5)))

#joint distribution represented in PyMC3
ab = pm.DensityDist('ab', ab_dist, shape=2, testval = [1,1])

#hyperpriors
a = ab[0]
b = ab[1]

#prior
p = pm.Beta('p', a, b, shape=trials.shape[0])

#likelihood
y = pm.Binomial('y', n = trials, p = p, observed = successes)

# This is what we observed
trials = np.array([100, 100, 100, 100, 100, 100, 100, 100, 100, 100])
successes = np.array([40, 44, 47, 54, 63, 46, 44, 49, 58, 50])

trace = pm.sample(2000)

• I think you need to convert the if in hyperpriors to theano IfElse or Switch as the whole dist needs to be a theano expression. – twiecki Nov 12 '16 at 21:07
• Thanks for your response @twiecki! I've converted the function to a T.switch as you described, but am now getting the following error message: 'ValueError: length not known: Elemwise{le,no_inplace} [id A] '' |ab [id B] |DimShuffle{x} [id C] '' |TensorConstant{0} [id D]'. I'm new to Theano and struggling to debug despite some research. – Eddie Nov 13 '16 at 19:26