www.digitalmars.com         C & C++   DMDScript  

digitalmars.D.bugs - [Issue 22502] New: Potential error where


          Issue ID: 22502
           Summary: Potential error where
                    gives different results from Wolfram Alpha for
                    parameters aa and bb being half integer and yy0  being
           Product: D
           Version: D2
          Hardware: x86_64
                OS: Linux
            Status: NEW
          Severity: minor
          Priority: P1
         Component: phobos
          Assignee: nobody puremagic.com
          Reporter: alex dengenius.com

So I was trying to match a table in this article
(https://www.jstor.org/stable/2676784), in particular Table 5. which basically
corresponds to betaIncompleteInverse(x+1/2,N-x+1/2,0.05/2).  Using Wolfram
Alpha I can match the table, but I am getting numbers which seem surprisingly
off when using betaIncompleteInverse.

For example: 

InverseBetaRegularized[0.05/2, 7+1/2, 1+1/2] == 0.546280678967585... and
betaIncompleteInv(7+1/2,1+1/2,0.05/2) ~= 0.590384

So I did a few sanity checks using the test cases in the unittests in the file
For example:

InverseBetaRegularized[0.0109343152340992, 8, 10] == 0.2 just like
betaIncompleteInv(8, 10, 0.010_934_315_234_099_2L) ~= 0.2

Based on that I think I am using the function correctly, it is still possible I
am doing something wrong though.

Using the program:
import std.mathspecial;
import std.stdio;

real jefferys(real x, real N, real alpha=0.05){
 return betaIncompleteInverse(x+1/2,N-x+1/2,alpha/2);

int main(){
    for(real N =7; N< 10; N++){
        for(real x=1; x < N;x++){
         writef("(%f,%f) %f\n", N,x,jefferys(x,N));
    return 0;
and Wolfram Alpha like
just varying the parameters and Table 5 from the article I linked I get:

 (N,x)             Wolfram               D             Table 5
 (7,1)             0.0158765...         0.004211..      0.016
 (7,2)             0.0647282...         0.043272..      0.065
 (7,3)             0.1388642...         0.118117..      0.139
 (7,4)             0.2345012...         0.222778..      0.234

And so on.

Sorry if this is just some kind of user error on my part.

Nov 10 2021