haskell / math-functions

Special mathematical functions
http://hackage.haskell.org/package/math-functions
BSD 2-Clause "Simplified" License
40 stars 28 forks source link

Loss of precision in incomplete beta when x=p/(p+q) #36

Closed Shimuuar closed 8 years ago

Shimuuar commented 8 years ago

There's discontinuity when switching between two areas for x, Here is plot for incompleteBeta 4.5 4.5

download

This of course adversely affects precision of inverse incomplete beta in area.

Shimuuar commented 8 years ago

Plot for \a → incompleteBeta a a (0.5+δ) - incompleteBeta a a (0.5-δ). E.g. gap size

download

Note sudden drop to 0 at a=10

Shimuuar commented 8 years ago

Turns out is due to loss of precision in logBeta for p < 10 && q < 10. It uses logGamma which isn't very precise but itself. Replacing logGamma with logGammaL in otherwise branch of logBeta fixed problem.

Shimuuar commented 8 years ago

Partially fixed in b4819949e5016f3cbcdc88d25556531b5c839a64 by switching to Lanczos approximation everywhere.

Shimuuar commented 8 years ago

b4819949e5016f3cbcdc88d25556531b5c839a64 markedly improved precision of gamma function for z<1, so I think it's safe to mark this as fixed.