Spreadshit addiction [General Sta­tis­tics]

posted by yicaoting  – NanKing, China, 2011-10-01 18:39 (4964 d 09:48 ago) – Posting: # 7405
Views: 15,983

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,424 posts in 4,927 threads, 1,668 registered users;
39 visitors (0 registered, 39 guests [including 7 identified bots]).
Forum time: 04:28 CEST (Europe/Vienna)

My doctor gave me six months to live,
but when I couldn’t pay the bill
he gave me six months more.    Walter Matthau

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