Szukam biblioteki / implementacji java, która obsługuje obliczanie odwrotnej funkcji rozkładu skumulowanego dla dystrybucji beta (aka oszacowanie kwantyli)z rozsądną precyzją.

Oczywiście, że próbowałemapache commons matematyka, ale w wersji 3 nadal wydaje się, że jest ich trochęproblemy z precyzją. Poniżej opisano problem, który prowadzi do tego pytania.

Załóżmy, że chcę obliczyć wiarygodny przedział dystrybucji beta z wieloma próbami. Wapache commons matematyka ...

final int trials = 161750;
final int successes = 10007;
final double alpha = 0.05d;

// the supplied precision is the default precision according to the source code
BetaDistribution betaDist = new BetaDistribution(successes + 1, trials - successes + 1, 1e-9);

System.out.println("2.5 percentile :" + betaDist.inverseCumulativeProbability(alpha / 2d));
System.out.println("mean: " + betaDist.getNumericalMean());
System.out.println("median: " + betaDist.inverseCumulativeProbability(0.5));
System.out.println("97.5 percentile :" + betaDist.inverseCumulativeProbability(1 - alpha / 2d));

który dostarcza

2.5 percentile :0.062030402074808505
mean: 0.06187249616697166
median: 0.062030258659508855
97.5 percentile :0.06305170793994147

Problem polega na tym, że 2.5 percentyl i mediana są w tym samym czasie większe niż średnia.

W porównaniu zR-pakietbinom dostarcza

binom.confint(10007+1,161750+2,methods=c("agresti-coull","exact","wilson"))
         method     x      n      mean      lower      upper
1 agresti-coull 10008 161752 0.0618725 0.06070873 0.06305707
2         exact 10008 161752 0.0618725 0.06070317 0.06305756
3        wilson 10008 161752 0.0618725 0.06070877 0.06305703

iR-pakietstatystyki

qbeta(c(0.025,0.975),10007+1,161750-10007+1)
[1] 0.06070355 0.06305171

Po drugie wyniki R, oto coWolfram Alpha powiedział mi

InverseBetaRegularized [0,025,10007 + 1,161750-10007 + 1] => 0,06070354631 ...InverseBetaRegularized [0.975,10007 + 1,161750-10007 + 1] => 0,06305170794 ...

Ostatnia uwaga dotycząca wymagań:

Muszę wykonać wiele tych obliczeń. Stąd każde rozwiązanie nie powinno trwać dłużej niż 1s (co jest wciąż dużo w porównaniu z matematyką apache commons 41ms (aczkolwiek błędną)).Zdaję sobie sprawę, że można używać R w java. Z powodów, których nie omówię tutaj, jest to ostatnia opcja, jeśli coś innego (czysta java) zawiedzie.

Aktualizacja 21.08.12

Wydaje się że problem został naprawiony lub przynajmniej poprawiony w 3.1-SNAPSHOT apache-commons-math. Do powyższej listy zastosowań

2.5 percentile :0.06070354581340706
mean: 0.06187249616697166
median: 0.06187069085946604
97.5 percentile :0.06305170793994147

Aktualizacja 23.02.13

Chociaż na pierwszy rzut oka to pytanie i jego odpowiedzi mogą być zbyt zlokalizowane, myślę, że bardzo dobrze ilustruje to, że niektórych problemów numerycznych nie można rozwiązać (efektywnie) za pomocą podejścia „co pierwszy przychodzi do umysłu-hakera”. Mam nadzieję, że pozostanie otwarty.