randomGamma

../../_images/OpenMp_27.png

Generates pseudorandom numbers from a standard gamma distribution.

Synopsis

randomGamma (nRandom, a)

Required Arguments

int nRandom (Input)
Number of random numbers to generate.
float a (Input)
Shape parameter of the gamma distribution. This parameter must be positive.

Return Value

An array of length nRandom containing the random standard gamma deviates.

Description

Function randomGamma generates pseudorandom numbers from a gamma distribution with shape parameter a and unit scale parameter. The probability density function is

\[f(x) = \frac{1}{\mathit{\Gamma}(a)}x^{a-1}e^{-x}\text{ for } x \geq 0\]

Various computational algorithms are used depending on the value of the shape parameter a. For the special case of \(a=0.5\), squared and halved normal deviates are used; for the special case of \(a=1.0\), exponential deviates are generated. Otherwise, if a is less than 1.0, an acceptance-rejection method due to Ahrens, described in Ahrens and Dieter (1974), is used. If a is greater than 1.0, a ten-region rejection procedure developed by Schmeiser and Lal (1980) is used.

Deviates from the two-parameter gamma distribution with shape parameter a and scale parameter b can be generated by using randomGamma and then multiplying each entry in r by b. The following statements (in single precision) would yield random deviates from a gamma \((a,b)\) distribution.

The Erlang distribution is a standard gamma distribution with the shape parameter having a value equal to a positive integer; hence, randomGamma generates pseudorandom deviates from an Erlang distribution with no modifications required.

Function randomSeedSet can be used to initialize the seed of the random number generator; function randomOption can be used to select the form of the generator.

Example

In this example, randomGamma generates five pseudorandom deviates from a gamma (Erlang) distribution with shape parameter equal to 3.0.

from numpy import *
from pyimsl.stat.randomGamma import randomGamma
from pyimsl.stat.randomSeedSet import randomSeedSet
from pyimsl.stat.writeMatrix import writeMatrix

n_random = 5
a = 3.0
randomSeedSet(123457)
r = randomGamma(n_random, a)
writeMatrix("Gamma(3) random deviates", r)

Output

 
                   Gamma(3) random deviates
          1            2            3            4            5
      6.843        3.445        1.853        3.999        0.779