Gamma probability distribution function.
Gamma(g,alpha,beta)
Generates a random variate following the Gamma probability distribution function, with parameters alpha and beta.
The pseudo-code algorithm used for this distribution is:
I. If alpha == 1, return y = exponential(1)
II. If alpha == 0.5, return y = normal(0,1)
III. If a < 1
b = (e + alpha)/e
1. Generate u1 as U(0,1). Let P =
beta*u1. If P > 1, go to step 3. Otherwise, go to step 2.
2. Let Y = P^(1/alpha), and generate u2
as U(0,1). If u2 &le e^(-Y), return y = Y. Otherwise, go to
step 1.
3. Let Y = -ln[(b - P)/alpha] and
generate u2 as U(0,1). If u2 &le Y^(alpha-1), return y = Y.
Otherwise, go to step 1.
IV. If a > 1 a = 1/sqrt(2*alpha
-1)
b = alpha -ln(4)
q = alpha + 1/a
theta = 4.5
d = 1 + ln(theta)
1. Generate u1 and u2 as IID
U(0,1).
2. Let V = a*ln[u1/(1-u1)], Y =
alpha*e^V, Z = u1^2*u2, and W = b + q*V - Y.
3. If W + d - theta*Z &ge 0, return
y = Y. Otherwise, go to step 4.
3. If W &ge ln(Z), return y = Y.
Otherwise, go to step 1.
V. x = beta*y;
(u,g) := Gamma(g,0.4,0.7);
function Gamma extends Var; end Gamma;