The Gamma function, Γ(x)

Γ(x) is a generalisation of the factorial function to real and complex numbers. Like x!, Γ(x+1) = x * Γ(x).

Mathematically, if z.re > 0 then Γ(z) = ∫₀^∞ t^z-1e^-t dt

Natural logarithm of the gamma function, Γ(x)

Returns the base e (2.718...) logarithm of the absolute value of the gamma function of the argument.

For reals,  logGamma is equivalent to log(fabs(gamma(x))).

The sign of Γ(x).

Returns -1 if Γ(x) < 0, +1 if Γ(x) > 0, NAN if sign is indeterminate.

Note that this function can be used in conjunction with logGamma(x) to evaluate gamma for very large values of x.

Beta function

The  beta function is defined as

 beta(x, y) = (Γ(x) * Γ(y)) / Γ(x + y)

Digamma function

The  digamma function is the logarithmic derivative of the gamma function.

 digamma(x) = d/dx logGamma(x)

Incomplete beta integral

Returns incomplete beta integral of the arguments, evaluated from zero to x. The regularized incomplete beta function is defined as

 betaIncomplete(a, b, x) = Γ(a + b) / ( Γ(a) Γ(b) ) * ∫₀^x t^a-1(1-t)^b-1 dt

and is the same as the the cumulative distribution function.

The domain of definition is 0 <= x <= 1. In this implementation a and b are restricted to positive values. The integral from x to 1 may be obtained by the symmetry relation

betaIncompleteCompl(a, b, x ) =  betaIncomplete( b, a, 1-x )

The integral is evaluated by a continued fraction expansion or, when b * x is small, by a power series.

Inverse of incomplete beta integral

Given y, the function finds x such that

betaIncomplete(a, b, x) == y

Newton iterations or interval halving is used.

Incomplete gamma integral and its complement

These functions are defined by

 gammaIncomplete = ( ∫₀^x e^-t t^a-1 dt )/ Γ(a)

gammaIncompleteCompl(a,x) = 1 -  gammaIncomplete(a,x) = (∫_x^∞ e^-t t^a-1 dt )/ Γ(a)

In this implementation both arguments must be positive. The integral is evaluated by either a power series or continued fraction expansion, depending on the relative values of a and x.

Inverse of complemented incomplete gamma integral

Given a and p, the function finds x such that

gammaIncompleteCompl( a, x ) = p.

Error function

The integral is

 erf(x) = 2/ √(π) ∫₀^x exp( - t²) dt

The magnitude of x is limited to about 106.56 for IEEE 80-bit arithmetic; 1 or -1 is returned outside this range.

Complementary error function

 erfc(x) = 1 - erf(x) = 2/ √(π) ∫_x^∞ exp( - t²) dt

This function has high relative accuracy for values of x far from zero. (For values near zero, use erf(x)).

Normal distribution function.

The normal (or Gaussian, or bell-shaped) distribution is defined as:

normalDist(x) = 1/√ π ∫_-∞^x exp( - t²/2) dt = 0.5 + 0.5 * erf(x/sqrt(2)) = 0.5 * erfc(- x/sqrt(2))

To maintain accuracy at values of x near 1.0, use  normalDistribution(x) = 1.0 -  normalDistribution(-x).

Inverse of Normal distribution function

Returns the argument, x, for which the area under the Normal probability density function (integrated from minus infinity to x) is equal to p.