0
$\begingroup$

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)
$\endgroup$
  • $\begingroup$ 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. $\endgroup$ – twiecki Nov 12 '16 at 21:07
  • $\begingroup$ 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. $\endgroup$ – Eddie Nov 13 '16 at 19:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.