Spreadshit addiction [General Sta­tis­tics]

posted by yicaoting  – NanKing, China, 2011-10-01 18:39 (5384 d 02:43 ago) – Posting: # 7405
Views: 19,802

Dear HS and ElMaestro, thanks for ElMaestro's sharing of that great paper and HS's useful links for R.
To improve the precision of Dr. Russel Lenth's NctInv() in Excel, I have to find a solution to improve the Excel's BetaInv() function instead of direct using of BetaInv() in MS Excel.

I have tested the BetaInv(0.001,2,4) in MS Excel 2007, it returns 0.010101318359375 which is identical to that in Excel 2003. How disappointing with MS's so called improvment. Actually, they did nothing! I will test it in MS Excel 2010 if possible.

I found a solution to improve the precision of Excel's BetaInv() function to a level of MS's so called 10^(-15). It is based on Excel's built-in BetaDist() function which has obviousely higher precision than BetaInv().

'-----------------------------------------------------------------
Public Function NCtInv(p As Double, df As Double) As Double
    Dim X As Double, y As Double, Incr As Double
    If (p <= 0) Or (p >= 1) Or (df <= 0) Then
        MsgBox "Illegal argument", vbCritical, "NCtInv"
        Exit Function
    End If
        If p = 0.5 Then p = 0.5000000001
        If p < 0.5 Then
            'X = Application.BetaInv(1 - 2 * p, 0.5, df / 2)
            X = BetaInverse(1 - 2 * p, 0.5, df / 2)
            NCtInv = -Sqr(df * X / (1 - X))
        Else
            'X = Application.BetaInv(1 - 2 * (1 - p), 0.5, df / 2)
            X = BetaInverse(1 - 2 * (1 - p), 0.5, df / 2)
            NCtInv = Sqr(df * X / (1 - X))
        End If
End Function
'-----------------------------------------------------------------
'-----------------------------------------------------------------
Public Function BetaInverse(p As Variant, Alpha As Variant, Beta As Variant) As Variant
Dim X As Variant
Dim A As Variant
Dim B As Variant
Dim Precision As Variant
X = 0
A = 0
B = 1
Precison = 10 ^ (-15)
Do While B - A > Precision
    X = (A + B) / 2
        If Application.BetaDist(X, Alpha, Beta) > p Then
            B = X
        Else
            A = X
        End If
Loop
BetaInverse = X
End Function
'-----------------------------------------------------------------

How to quickly leave a blank at the head of each row in BEBAC? My code is with blank at several rows, but after I submit the post, the blanks disappeared. *) See note. [HS]

Using the improved BetaInverse(), I get:
Excel's BetaInv(0.001,2,4)=     0.010101318359375
My     BetaInverse(0.001,2,4)=  0.0101017878834719
It can be seen that my BetaInverse(0.001,2,4) is much closer to R's qbeta(0.001,2,4) 0.01010178788373775407572 than Excel's BetaInv(0.001,2,4).

Thus, the new results for ElMaestro's df 5.123 example were refreshed as:
df      p    Excel 2003   OO Calc     R
5.123  0.000      #NA          #NA      -Inf
5.123  0.025  -2.55212890  -2.55212890  -2.55212890
5.123  0.500   0.00000006   0.00022634   0
5.123  0.975   2.55212890   2.55212890   2.55212890
5.123  1.000      #NA          #NA       Inf
0.000  0.975      #NA          #NA      NaN

Excel's calculation is performed with improved NctInv() containing improved BetaInverse()
Now, the bias has been greatly decreased. May be acceptable in most cases.



Complete thread:

UA Flag
Activity
 Admin contact
23,655 posts in 4,993 threads, 1,571 registered users;
355 visitors (0 registered, 355 guests [including 12 identified bots]).
Forum time: 21:23 CEST (Europe/Vienna)

The real struggle is not between the right and the left
but between the party of the thoughtful
and the party of the jerks.    Jimmy Wales

The Bioequivalence and Bioavailability Forum is hosted by
BEBAC Ing. Helmut Schütz
HTML5