Details

Bug

Status: Closed

Major

Resolution: Fixed

None

None

None

None

Operating System: other
Platform: PC

36105
Description
The inverseCumulativeProbability for the Poisson distribution is calculated
through the standard implementation in AbstractIntegerDistribution. This
implementation fails with a continue fraction stack overflow for large lambda.
Instead, I suggest the following override in PoissonDistributionImpl
/**
 Calculate the inverse cumulative probability for the Poisson distribution
 based on a CornishFisher expansion approximation followed up by a line
 search.
 For the CornishFisher expansion see Abramawitz and Stegun
 'Handbook of Mathmatical Functions' pages 935 and 928
 @param prob the target probability
 @return the desired quantile
 @throws MathException when things go wrong
*/
private int inverseCumulativeProbability(double prob) throws MathException{
if (prob < 0.0  prob >= 1.0)
throw new MathException("Probability must be in rightside open interval
[0.0, 1.0)");
if (prob == 0) return 0; // there is nothing to calculate
// Use the CornishFisher Expansion with two terms to get a very good
approximation
// see Abramawitz and Stegun 'Handbook of Mathmatical Functions'
// pages 935 and 928
double mu = this.getMean(); // mean
double sigma = Math.sqrt(mu); // standard deviation
double gamma = 1.0/Math.sqrt(mu); // skewness
double z = new NormalDistributionImpl(0.0,
1.0).inverseCumulativeProbability(prob);
// this is the actual expansion
// the floor(... + 0.5) operation is the continuity correction
int y = (int) Math.floor(mu + sigma*(z + gamma*(z*z  1.0)/6.0) + 0.5);
// Given this starting point, line search to the right or left.
// Bisection search is not necessary, since the approximation is rarely
// off by more than 1 or 2
z = this.cumulativeProbability;
if ( z > prob) { // missed it to the right, search to the left
while(true)
} else { // missed it to the left, search to the right
while(true)
}
}