Spreadshit addiction [General Sta­tis­tics]

posted by Helmut Homepage – Vienna, Austria, 2011-09-30 18:52 (4965 d 13:04 ago) – Posting: # 7397
Views: 16,330

Dear yicaoting & ElMaestro!

❝ FUNCTION NcTinv() […]


❝ BTW: I have tried it, the result was satisfactory.


Hhm, but how?
I must confess that I have only a somewhat outdated version of Excel (2000), but I had no problems in running FARTSIE.xls.
In FARTSSIE I added a new sheet and 5, 5.123, 6 to A1:A3 and B1:B3 =TINV(0.05,A1…A3).
I get (as expected) 2.57058, 2.57058, 2.4469. If I add to C1:C3 = NCtInv(0.05,A1…A3) the debugger throws an error in VBA Module 1Error in compiling: Sub or Function not defined pointing to
NCtInv = Solve("tDummy", p, x, y) in
Public Function NCtInv(p As Double, df As Double, Optional Delta) As Double
Solve() calls (?) the Solver-Add-in (Solver.xla), which I have installed (proof: FARTSSIE works fine…). The noncentrality parameter Delta is optional – don’t know what’s going on here (IMHO without noncentrality this section of the code should not be called at all).

BTW, NCt() works fine:
NCt(2.570582,5)     → 0.97500
NCt(2.552129,5.123) → 0.97500
NCt(2.446912,6)     → 0.97500


if compared to R:
t  <- c(2.570582,2.552129,2.446912)
df <- c(5,5.123,6)
pt(t,df)
[1] 0.975 0.975 0.975


@ElMaestro:
A crude workaround based on the first section on the macro…
df in cells A1:A3 (5, 5.123, 6).
In cells B1…B3 =SQRT(A1…A3*BETAINV(1-0.05,0.5,A1…A3/2)/(1-BETAINV(1-0.05,0.5,A1…A3/2))). I got:

df      Excel 2000  OO Calc 3.0.10    R 2.13
5        2.570578      2.570582      2.570582
5.123    2.552128      2.552129      2.552129
6        2.446912      2.446912      2.446912


The VBA-section for fractional df without noncentrality shift is:
If p = 0.5 Then p = 0.5000000001
  If p < 0.5 Then
    x = Application.BetaInv(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)
    NCtInv = Sqr(df * x / (1 - x))
End If

You should be able to adapt the lines above into a series of nested if()s into OO calc (like I did with the crowbar above).


P.S.: Don’t know what I’m doing here. Just fiddlin’ around. :smoke:

Edit: Toys for boys. The nested Ifs in the VBA-code as a one-liner gives a pretty long formula. Use the simple one from above only for α 0.05. Results are identical to the 8th digit; It tried dfs between 1 and 6. OO calc performed slightly better than Excel 2000 (maybe this issue disappeared in later versions - see here). The largest deviation I got was for t(0.05,<2) ~0.00002. Funny that the formula based on the inverse β-distribution outperformed the one on the inverse t in Excel:
0.05, df=1 (diff to R and OO calc)
based on BETAINV() 12.70618521 (-0.00001953)
based on TINV()    12.70615030 (-0.00005444)

The full story: in A1 df, in B2 p=1-α/2, in C2 =IF(OR(B1<=0,B1>=1,A1<=0),NA(), IF(B1=0.5, SQRT(A1*BETAINV(1-2*(1-0.5000000001),0.5,A1/2)/(1-BETAINV(1-2*(1-0.5000000001),0.5,A1/2))), IF(B1<0.5,-SQRT(A1*BETAINV(1-2*B1,0.5,A1/2)/(1-BETAINV(1-2*B1,0.5,A1/2))), SQRT(A1*BETAINV(1-2*(1-B1),0.5,A1/2)/(1-BETAINV(1-2*(1-B1),0.5,A1/2))))))
  1. p ≤0 | p ≥1 | df ≤0 throw a #NA
  2. workaround for p=0.5 (preventing error in BETAINV(0;0.5;A1/2)
  3. p <0.5
  4. p >0.5 (else)
Test case for your famous df=5.123
df       p      Excel 2000    OO Calc      R
5.123   0.000          #NA          #NA         -Inf
5.123   0.025  -2.55212813  -2.55212890  -2.55212890
5.123   0.500   0.00000038   0.00022634   0
5.123   0.975   2.55212813   2.55212890   2.55212890
5.123   1.000          #NA          #NA          Inf
0       0.975          #NA          #NA          NaN


I hate Microsoft’s partly translations – which were carried over to OO Calc. On a German system NA()NV(), IF()WENN(), SQRT()WURZEL().

Dif-tor heh smusma 🖖🏼 Довге життя Україна! [image]
Helmut Schütz
[image]

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

Complete thread:

UA Flag
Activity
 Admin contact
23,424 posts in 4,927 threads, 1,669 registered users;
44 visitors (0 registered, 44 guests [including 4 identified bots]).
Forum time: 07:56 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