Implementing a more sophisticated pymc3 model

[ahc_ahc]

Hi, I have worked through the basic Pymc3 tutorials, but now I need some help building something a little bit more sophisticated. I’d like to produce posterior distributions on three variables that are combined to produce a single observable parameter in the following way:

1. Using the value of one of the variables, I calculate the upper limit of integration by solving a polynomial equation via the bisection method.

2. I numerically integrate a function that is independent of the variables from 0 to the calculated upper limit.

3. I multiply the result by powers of each of the variables to produce a rate of observations.

4. I make 0, 1, 2, … draws from a Poisson distribution with shape parameter of rate*time

In this way I hope to constrain the three variables that make up the rate. How do I implement this second step with pymc3 distributions?

[junpenglao]

For numerical integration, you should check out the following posts:

Custom theano Op to do numerical integration

I think I worked out how to do that. I implemented a more general integration op: from scipy.integrate import quad import theano import theano.tensor as tt import numpy as np class Integrate(theano.Op): def __init__(self, expr, var, *extra_vars): super().__init__() self._expr = expr self._var = var self._extra_vars = extra_vars self._func = theano.function( [var] + list(extra_vars), self._expr, on_unused_input='ig…

PyMC3 : Using a parameter as a limit of integral

Try rewrite the function f into a theano function, remove the numpy.vectorize, and use the integration theano.Ops from below.

[ahc_ahc]

Thanks, I’ll try to implement that and see what problems I run into. Any clue for computing the limit of integration? It amounts to a root-finding problem where the input variable changes the root (and so changes the limit of integration).

[ahc_ahc]

Unfortunately, I’m already running into issues - I get AttributeError: module 'theano' has no attribute 'Op' due to the class definition line.

[ahc_ahc]

It looks like this is a persistent and fiddly problem - I will ping back here once I’ve resolved it, thanks.

[jessegrabowski]

Blockquote
Any clue for computing the limit of integration? It amounts to a root-finding problem where the input variable changes the root (and so changes the limit of integration).

I wrote an implementation of newton’s method in Aesara for a project, you might find it helpful. Basically you need to write a function for one step of your optimizer, then use a scan together with an until condition.

[ahc_ahc]

@junpenglao i’m having some serious problems getting theano working - is there an installation guide somwhere? Following the steps on the theano sphinx page isn’t working and others’ advice is contradictory. Sorry that this isn’t really on-topic.

[junpenglao]

If following GitHub - pymc-devs/pymc: Probabilistic Programming in Python: Bayesian Modeling and Probabilistic Machine Learning with Aesara does not work, could you please open a new topic?

[ahc_ahc]

Thanks @junpenglao. I have posted a new topc here where I dump the output from trying to run one of the simple polynomial integrals in the linked thread.