Как рассчитать обратную интегральную функцию бета-распределения в Java

Я ищу библиотеку Java / реализацию, которая поддерживает вычисление обратной кумулятивной функции распределения для бета-распределения (или квантили оценки)with reasonable precision.

Конечно я пробовалapache commons математика, но в версии 3 все еще есть некоторыепроблемы с точностью, Ниже проблема, которая приводит к этому вопросу, подробно описана.

Предположим, я хочу рассчитать вероятный интервал бета-распределения с большим количеством испытаний. Вapache commons math ...

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));

который доставляет

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

Проблема в том, что 2,5 процентиль и медиана в то же время и выше среднего.

Для сравненияR-packagebinom обеспечивает

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

иR-packagestats

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

Чтобы подкрепить результаты R, вот чтоWolfram Alpha сказал мне

InverseBetaRegularized[0.025,10007+1,161750-10007+1] => 0.06070354631... InverseBetaRegularized[0.975,10007+1,161750-10007+1] => 0.06305170794...

Final note on the requirements:

I need to run a lot of these calculations. Hence any solution should not take longer than 1s (which is still a lot compared to the 41ms of (albeit wrong) apache commons math). I am aware that one can use R within java. For reasons I won't elaborate here, this is the last option if anything else (pure java) fails.

Update 21.08.12

Похоже на то что проблема была исправлена или, по крайней мере, улучшена в 3.1-SNAPSHOT из apache-commons-math. Для случая использования выше

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

Update 23.02.13

Хотя на первый взгляд этот вопрос и его ответы могут быть слишком локализованы, я думаю, что он очень хорошо иллюстрирует, что некоторые численные проблемы не могут быть решены (эффективно) с помощью подхода «что первым приходит на ум». Поэтому я надеюсь, что он остается открытым.

Ответы на вопрос(3)

Ваш ответ на вопрос