Inverse poisson
10-04-2015, 06:52 PM
Post: #1
 lrdheat Senior Member Posts: 512 Joined: Feb 2014
Inverse poisson
Accuracy of calculation (inclined to go with the Prime)...

Inverse Poisson with mean of 1E9 on Prime 1,000,073,566
Casio 9860gii 1,000,011,502
(Casio can't calculate for mean of 5E8, Prime hangs)

WP 34s. Way off!

TI nSpire Mean of 1E9 out of domain for calculation

Can anyone find an accurate value to use as a benchmark?
10-04-2015, 07:43 PM (This post was last modified: 10-04-2015 09:09 PM by Dieter.)
Post: #2
 Dieter Senior Member Posts: 2,397 Joined: Dec 2013
RE: Inverse poisson
(10-04-2015 06:52 PM)lrdheat Wrote:  Inverse Poisson with mean of 1E9 on Prime 1,000,073,566
Casio 9860gii 1,000,011,502
(Casio can't calculate for mean of 5E8, Prime hangs)

Fine - but for which probability ?-)
I assume it must be something like 0.01 resp. 0.99.

(10-04-2015 06:52 PM)lrdheat Wrote:  WP 34s. Way off!

I agree the Poisson (and Binomial) functions, especially the quantile functions, should get an update. Pauli and I had a discussion on this some months ago. The current, not yet completely finished code I proposed back then returns for p=0.01 a quantile of 999.926.434,656. This result may be flawed by inaccuracies in the regularized Gamma function used internally, so beware. Edit: result is correct – see below.

(10-04-2015 06:52 PM)lrdheat Wrote:  TI nSpire Mean of 1E9 out of domain for calculation

The problem seems to be the Poisson mass function where e–µ · µx has to be evaluated. For µ=1E9 the first factor evaluates to zero in most environments, but there is a workaround by working with logs. ;-)

I checked the 34s result (using my improved, but currently not implemented algorithm) with Wolfram Alpha. Et voilà: the double precision result for p=0.01 and µ=1E9 is exact in 32 out of 34 digits. Which is about as good as it gets without additional guard digits:

x = 999926434.65620471615794098329370...

You may check this yourself by entering Q(1+x, 1E9) at Wolfram Alpha where x is the quantile you look for, i.e. here 999926434.6562...

N.B.: Q(a, x) is the regularized incomplete Gamma function, which can be seen as a continuous (real) extension of the Poisson CDF. For integer x the Poisson CDF(x, µ) equals Q(x+1, µ). If you want an integer quantile, round up the result to get 999926435.

Dieter
10-05-2015, 01:05 AM
Post: #3
 lrdheat Senior Member Posts: 512 Joined: Feb 2014
RE: Inverse poisson
Sorry about not including probability...Using average of 1E9, at p=.99, what is correct value for inverse poisson? I'm guessing the Prime is closer than the Casio 9860gii.

Your new implementation sounds like a winner!
10-05-2015, 06:11 AM (This post was last modified: 10-05-2015 09:44 PM by Dieter.)
Post: #4
 Dieter Senior Member Posts: 2,397 Joined: Dec 2013
RE: Inverse poisson
(10-05-2015 01:05 AM)lrdheat Wrote:  Sorry about not including probability...Using average of 1E9, at p=.99, what is correct value for inverse poisson? I'm guessing the Prime is closer than the Casio 9860gii.

This depends on whether you look at the left hand or right hand sum, i.e. the probability P(x≤X) or P(x>X). ;-) The results I posted refer to the former.

For p=0.99 there is a problem with the 34s regularized Gamma function (which returns values > 1 – sounds like a bug). Despite this my program returns an integer result of 1000073565 < x < 1000073566. With a corrected Gamma function the program returns...

x = 1000073565.814426761131158060038507

...which, according to Wolfram Alpha, seems to be the correct 34-digit value for the 99%-quantile.

So the Prime is correct.

Dieter

Edit: sorry, the originally posted result was off by one... #-)
10-07-2015, 12:45 PM (This post was last modified: 10-07-2015 12:46 PM by Dieter.)
Post: #5
 Dieter Senior Member Posts: 2,397 Joined: Dec 2013
RE: Inverse poisson
And BTW...

(10-04-2015 06:52 PM)lrdheat Wrote:  (Casio can't calculate for mean of 5E8, Prime hangs)

...the quantile for µ=5E8 and p=0.99 is
x = 500052018.9551617032468798129971636.

Or 500052019 if you want an integer result.

Dieter
 « Next Oldest | Next Newest »

User(s) browsing this thread: 1 Guest(s)