Conversation
ColCarroll
left a comment
There was a problem hiding this comment.
Sorry this sat for so long -- looks good to me, other than the one comment!
pymc3/distributions/dist_math.py
Outdated
There was a problem hiding this comment.
these two are actually implemented in theano now: https://github.com/Theano/Theano/blob/master/theano/scalar/basic_scipy.py
|
Okay that was a great comment! I used the same implementation of i0 in theano to implement the i0e function which is used in scipy to improve the numerical stability of the pdf for large values of value and nu. That means the tests can use a much larger domain now for testing! |
|
Looks like there is merging trouble. My fork is named Notice the |
…alues of nu and value Use ((x-b)**2)/2 + xb instead of (x**2 + b**2)/2 in the pdf for the rice distribution and include the np.exp(-xb) in the i0e to match the scipy implementation
|
I can't seem to reproduce the failing test case on my local machine. I'm using python 2.7 and set the FLOATX env var to float32. Any idea what's happening? |
|
That looks strange - like the output is -inf when you expect it to be like, -109. Replicating the tests is a little tough (you need to set |
|
Tests are passing now! |
|
Looks good to me -- thanks for the contribution, @MBrouns! |
Continuation from #2963. I rebased to get rid of the merge commit and got the tests to pass by changing the tested domain to
Rplusbig.The problem with large values of
valueandnuis that the regulari0function will overflow toinf. Scipy seems to go around this by slightly modifying the original function from:x * exp(-(x**2+b**2)/2) * I[0](x*b)to:x * np.exp(-(x-b)*(x-b)/2.0) * sc.i0e(x*b)The i0e function is then partitioned into two intervals: [0, 8] and (8, infinity) which are implemented slightly different:
If we want the Rice distribution to be able to handle larger values for
valueandnuwe would need to implement a similar strategy but I am unsure on how exactly to implement this in PyMC3. Would this PR be acceptable without that?