RNHYP
Generates pseudorandom numbers from a hypergeometric distribution.
Required Arguments
N — Number of items in the sample. (Input)
N must be positive.
M — Number of special items in the population, or lot. (Input)
M must be positive.
L — Number of items in the lot. (Input)
L must be greater than both N and M.
IR — Vector of length NR containing the random hypergeometric deviates. (Output)
Each element of IR can be considered to be the number of special items in a sample of size N drawn without replacement from a population of size L that contains M such special items.
Optional Arguments
NR — Number of random numbers to generate. (Input)
Default: NR = size (IR,1).
FORTRAN 90 Interface
Generic: CALL RNHYP (N, M, L, IR [, …])
Specific: The specific interface name is S_RNHYP.
FORTRAN 77 Interface
Single: CALL RNHYP (NR, N, M, L, IR)
Description
Routine RNHYP generates pseudorandom numbers from a hypergeometric distribution with parameters N, M, and L. The hypergeometric random variable X can be thought of as the number of items of a given type in a random sample of size N that is drawn without replacement from a population of size L containing M items of this type. The probability function is
for x = max(0, N ‑ L + M), 1, 2, …, min(N, M)
If the hypergeometric probability function with parameters N, M, and L evaluated at N ‑ L + M (or at 0 if this is negative) is greater than the machine epsilon (AMACH(4) (Reference Material)), and less than 1.0 minus the machine epsilon, then RNHYP uses the inverse CDF technique. The routine recursively computes the hypergeometric probabilities, starting at x = max(0, N ‑ L + M) and using the ratio f(X = x + 1)/f(X = x) (see Fishman 1978, page 457).
If the hypergeometric probability function is too small or too close to 1.0, then RNHYP generates integer deviates uniformly in the interval [1, L ‑ i], for i = 0, 1, …; and at the i‑th step, if the generated deviate is less than or equal to the number of special items remaining in the lot, the occurrence of one special item is tallied and the number of remaining special items is decreased by one. This process continues until the sample size or the number of special items in the lot is reached, whichever comes first. This method can be much slower than the inverse CDF technique. The timing depends on N. If N is more than half of L (which in practical examples is rarely the case), the user may wish to modify the problem, replacing N by L ‑ N, and to consider the deviates in IR to be the number of special items not included in the sample.
Comments
The routine
RNSET can be used to initialize the seed of the random number generator. The routine
RNOPT can be used to select the form of the generator.
Example
In this example, RNHYP is used to generate five pseudorandom deviates from a hypergeometric distribution to simulate taking random samples of size 4 from a lot containing 20 items of which 12 are defective. The resulting hypergeometric deviates represent the numbers of defectives in each of the five samples of size 4.
USE RNHYP_INT
USE UMACH_INT
USE RNSET_INT
IMPLICIT NONE
INTEGER NR
PARAMETER (NR=5)
!
INTEGER IR(NR), ISEED, L, M, N, NOUT
!
CALL UMACH (2, NOUT)
N = 4
M = 12
L = 20
ISEED = 123457
CALL RNSET (ISEED)
CALL RNHYP (N, M, L, IR)
WRITE (NOUT,99999) IR
99999 FORMAT (' Hypergeometric random deviates: ', 5I8)
END
Output
Hypergeometric random deviates: 4 2 3 3 3