Help testing modification to distribution

Development
1

It bothered me (more than it probably should!) that the Gamma supports parameterization with mu, sd as well as the standard alpha, beta, but the InverseGamma does not.
So I modified the InverseGamma by analogy to the Gamma function, adding InverseGamma.get_alpha_beta(). I am trying to test this, but I’m afraid I really don’t get how pymc3_matches_scipy works. This is what I am trying to do:

    def test_inverse_gamma(self):
        self.pymc3_matches_scipy(
            InverseGamma, Rplus, {'alpha': Rplus, 'beta': Rplus},
            lambda value, alpha, beta: sp.invgamma.logpdf(value, alpha, scale=beta))
        # here's the new stuff...
        self.pymc3_matches_scipy(
            InverseGamma, Rplus, {'mu': Rplus, 'sd': Rplus},
            lambda value, mu, sd: sp.invgamma.logpdf(value, InverseGamma.get_alpha_beta(mu, sd)[0],
                                                     scale=InverseGamma.get_alpha_beta(mu, sd)[1]))

This gives me an error that I’m calling InverseGamma.get_alpha_beta() with missing arguments, mu and sd.
I think this means that the lambda function is being called in a way that I don’t expect or understand. Could someone please clue me in about this? I’d like to submit this enhancement as a pull request, but I would like to feel confident that it works by testing it.

P.S. I note that get_alpha_beta for the Gamma function is a normal method, but as far as I can tell, it should really be either a function local to the initializer, or a static method, because it doesn’t need access to any information local to the object. Right? Should this be changed, or is it something needed for backwards compatibility?

2

Can you share the implementation and stack trace for InverseGamma.get_alpha_beta? The test suite gets pretty convoluted in there, but it sounds from what you’re saying like your implementation is throwing the exception.

3

Sure. Here’s get_alpha_beta():

    @staticmethod
    def get_alpha_beta(alpha, beta, mu, sd):
        if (alpha is not None):
            if (beta is not None):
                pass
            else:
                beta = 1
        elif (mu is not None) and (sd is not None):
            alpha = (2 * sd**2 + mu**2)/sd**2
            beta = mu * (mu**2 + sd**2) / sd**2
        else:
            raise ValueError('Incompatible parameterization. Either use '
                             'alpha and (optionally) beta, or mu and sd to specify '
                             'distribution.')

        return alpha, beta

…and this is the message I get from pytest:

______________________________________________________ TestMatchesScipy.test_inverse_gamma _______________________________________________________

self = <pymc3.tests.test_distributions.TestMatchesScipy object at 0x1261176a0>

    def test_inverse_gamma(self):
        self.pymc3_matches_scipy(
            InverseGamma, Rplus, {'alpha': Rplus, 'beta': Rplus},
            lambda value, alpha, beta: sp.invgamma.logpdf(value, alpha, scale=beta))
        self.pymc3_matches_scipy(
            InverseGamma, Rplus, {'mu': Rplus, 'sd': Rplus},
            lambda value, mu, sd: sp.invgamma.logpdf(value, InverseGamma.get_alpha_beta(mu, sd)[0],
                                                     scale=InverseGamma.get_alpha_beta(mu, sd)[1]))
    
>   def test_pareto(self):

pymc3/tests/test_distributions.py:709: 
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 
pymc3/tests/test_distributions.py:438: in pymc3_matches_scipy
    self.check_logp(model, value, domain, paramdomains, logp, decimal=decimal)
pymc3/tests/test_distributions.py:448: in check_logp
    assert_almost_equal(logp(pt), logp_reference(pt), decimal=decimal, err_msg=str(pt))
pymc3/tests/test_distributions.py:437: in logp
    return scipy_dist(**args)
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 

value = array(1.5), mu = array(100.), sd = array(0.9)

>   def test_pareto(self):
E   TypeError: get_alpha_beta() missing 2 required positional arguments: 'mu' and 'sd'

pymc3/tests/test_distributions.py:709: TypeError

when I look at the test code, it looks like I should be passing pymc3_matches_scipy a function of three arguments, value, mu, and sd (i.e., value and the parameters of the distribution), but this is what’s giving me the error. Looking at the code for TestMatchesScipy it looks almost like this is used in a curried way – like the parameters to the distribution are passed first, and then the points to evaluate the logp are plugged in separately, but TBH I really don’t understand the testing code at all.

Actually, a quick follow up: the type error makes it look like the function was passed a pair of singleton arrays for mu and sd which my code definitely wouldn’t expect. But the error message says that it’s missing two arguments. So this error message is pretty confusing. It shows three arguments, value, mu, and sd, getting arrays of 1.5, 100, and 0.9, respectively, but then it says that two of the arguments are missing.

Final follow-up: the call to the scipy function is made by scipy_dist(**args) (shouldn’t that be **kwargs to be more pythonic?).
So my guess is that this is going awry because of something involving my simple lambda getting an argument dictionary.
But if that’s the case, since my lambda looks like others in distribution_tests.py, I’m still puzzled by the error. That is, the following test form in the same test definition works correctly:

        self.pymc3_matches_scipy(
            InverseGamma, Rplus, {'alpha': Rplus, 'beta': Rplus},
            lambda value, alpha, beta: sp.invgamma.logpdf(value, alpha, scale=beta))

… so I’m still a little lost here.

6

It looks like your code expects 4 arguments all the time, and the test will pass only 2. Do you need to set default arguments of None to alpha, beta, mu, and sd?

7

D’oh, you are right. I didn’t really expect that static method to be called except by the initializer, which ensures that all four arguments are passed. And… I was looking so hard for something complicated in the test framework that I missed the obvious! Thank you so much. I’ll get this pushed to my repo and set up the pull request.