t-distribution
The t-distribution (Student’s t-distribution) is a continuous probability distribution that arises in estimating the mean of a normally distributed population when the sample size is small and the population standard deviation is unknown. The larger the sample, the more the t-distribution resembles a normal distribution.
Like the F-distribution the t-distribution can be obtained using the regularized incomplete beta function with the following inputs: X = d.f. * t / (d.f. * t * t), a = d.f. / 2 and b = 0.5, where t is the value obtained in the t-test, and d.f. is the degrees of freedom. The source code below applies that method to estimate the probability of the null hypothesis for a given t-value. The code includes the gamma function and the incomplete beta function, which are both used.
[code language=”vb”]
#COMPILE EXE
#REGISTER NONE
#DIM ALL
‘
FUNCTION PfromT(BYVAL T AS DOUBLE,BYVAL df AS DOUBLE) AS DOUBLE
LOCAL Beta, betain AS EXT
IF ABS(T) > 0.0 AND df > 0.0 THEN
CALL IncomplBeta(df / (df + T * T), df * 0.5#, 0.5#, Beta, betain)
FUNCTION = betain
ELSE
FUNCTION = 1
END IF
END FUNCTION
‘
FUNCTION PBMAIN
LOCAL Result AS STRING, P, t, df AS EXT
LOCAL i1 AS LONG
Result = INPUTBOX$("Enter t-value, d.f.", "Get p-value")
i1 = INSTR(Result, ",")
t = VAL(LEFT$(Result, i1 – 1))
df = VAL(MID$(Result,i1 + 1))
P = PfromT(t, df)
MSGBOX "P=" + FORMAT$(P," ##.##################"), ,"Result:"
END FUNCTION
‘
SUB IncomplBeta(BYVAL X AS EXT, BYVAL P AS EXT, BYVAL Q AS EXT, _
BYREF beta AS EXT, BYREF betain AS EXT)
‘
‘ Derived from FORTRAN code based on:
‘ algorithm AS 63 Appl. Statist. (1973), vol.22, no.3
‘ Computes incomplete beta function ratio for arguments
‘ X between zero and one, p (=a) and q (=b) positive.
‘ Returns log beta and Regularized Incomplete Beta
‘
LOCAL ns, indx AS LONG
LOCAL zero, one, acu AS EXT
LOCAL psq, cx, xx, pp, qq, term, ai, rx, temp AS EXT
‘ define accuracy and initialise
zero = 0.0 : one = 1.0 : acu = 1.0E-18
betain = 0.0 : beta = 0.0
‘ test for admissibility of arguments
IF(p <= zero OR q <= zero) THEN EXIT SUB
IF(x < zero OR x > one) THEN EXIT SUB
IF(x = zero OR x = one) THEN EXIT SUB
‘ calculate log of beta by using function GammLn
beta = gammln(p) + gammln(q) – gammln(p + q)
betain = x
‘ change tail if necessary
psq = p + q
cx = one – x
IF (p < psq * x) THEN
xx = cx
cx = x
pp = q
qq = p
indx = 1
ELSE
xx = x
pp = p
qq = q
indx = 0
END IF
term = one
ai = one
betain = one
ns = qq + cx * psq
‘ use Soper’s reduction formulae.
rx = xx / cx
temp = qq – ai
IF (ns = 0) THEN rx = xx
DO
term = term * temp * rx / (pp + ai)
betain = betain + term
temp = ABS(term)
IF(temp <= acu AND temp <= acu * betain) THEN EXIT DO
ai = ai + one
ns = ns – 1
IF (ns >= 0) THEN
temp = qq – ai
IF (ns = 0) THEN rx = xx
ELSE
temp = psq
psq = psq + one
END IF
LOOP
‘ calculate Regularized Incomplete Beta
betain = betain * EXP(pp * LOG(xx) + (qq – one) * LOG(cx) – beta) / pp
IF indx = 1 THEN betain = one – betain
END SUB
‘
FUNCTION GammLn(BYVAL x AS EXT) AS EXT
‘ Returns Ln(Gamma()) or 0 on error
‘ Based on Numerical Recipes gamma.h
‘ Lanczos, C. 1964, "A Precision Approximation
‘ of the Gamma Function," SIAM Journal on Numerical
‘ Analysis, ser. B, vol. 1, pp. 86-96.
LOCAL j AS LONG, tmp, y, ser AS EXT
DIM cof(0 TO 13) AS LOCAL EXT
cof(0) = 57.1562356658629235
cof(1) = -59.5979603554754912
cof(2) = 14.1360979747417471
cof(3) = -0.491913816097620199
cof(4) = 0.339946499848118887e-4
cof(5) = 0.465236289270485756e-4
cof(6) = -0.983744753048795646e-4
cof(7) = 0.158088703224912494e-3
cof(8) = -0.210264441724104883e-3
cof(9) = 0.217439618115212643e-3
cof(10) = -0.164318106536763890e-3
cof(11) = 0.844182239838527433e-4
cof(12) = -0.261908384015814087e-4
cof(13) = 0.368991826595316234e-5
IF x <= 0.0 THEN FUNCTION = 0.0 : EXIT FUNCTION ‘ Bad argument
y = x
tmp = x + 5.2421875
tmp = (x + 0.5) * LOG(tmp) – tmp
ser = 0.999999999999997092
FOR j = 0 TO 13
y = y + 1
ser = ser + cof(j)/y
NEXT j
FUNCTION = tmp + LOG(2.5066282746310005 * ser / x)
END FUNCTION
[/code]