Language Reference

RANDMULTINOMIAL Function

RANDMULTINOMIAL (N, NumTrials, Prob);

The RANDMULTINOMIAL function is part of the IMLMLIB library . The RANDMULTINOMIAL function generates a random sample from a multinomial distribution, which is a multivariate generalization of the binomial distribution.

The input parameters are as follows:

N

is the number of observations to sample.

NumTrials

is the number of trials. $NumTrials[j] \geq 0$, for $j = 1 \ldots p$.

Prob

is a $1 \times p$ vector of probabilities with $0 < \mi{Prob}[j] \leq 1$ and $\Sigma _{j=1}^ p \mbox{Prob}[j] = 1$.

For each trial, $\mbox{Prob}[j]$ is the probability of event $E_ j$, where the $E_ j$ are mutually exclusive and $\Sigma _{j=1}^ p \mbox{Prob}[j] = 1$.

The RANDMULTINOMIAL function returns an $N \times p$ matrix that contains N observations of $\mbox{NumTrials}$ random draws from the multinomial distribution. Each row of the resulting matrix is an integer vector $\{ X_1 \, X_2 \, \ldots \, X_ p\} $ with $\Sigma X_ j = \mbox{NumTrials}$. That is, for each row, $X_ j$ indicates how many times event $E_ j$ occurred in $\mbox{NumTrials}$ trials.

If $X=\{ X_1 \, X_2 \, \ldots \, X_ p\} $ follows a multinomial distribution with n trials and probabilities $\rho = \{ \rho _1 \, \rho _2 \, \ldots \, \rho _ p\} $, then

  • the probability density function for x is

    \[ f(x; n, \rho )= \frac{n!}{ \prod _{i=1}^ p {x_ i!}} \prod _{i=1}^ p {{\rho _ i}^{x_ i}} \]

  • the expected value of $X_ i$ is $n \rho _ i$.

  • the variance of $X_ i$ is $n \rho _ i (1-\rho _ i)$.

  • the covariance of $X_ i$ with $X_ j$ is $-n \rho _ i\rho _ j$.

  • if $p=1$ then X is constant.

  • if $p=2$ then $X_1$ is Binomial$(n, \rho _1)$ and $X_2$ is Binomial$(n, \rho _2)$.

The following example generates 1,000 samples from a multinomial distribution with three mutually exclusive events. For each sample, 10 events are generated. Each row of the returned matrix x represents the number of times each event is observed. The example also computes the sample mean and covariance and compares them with the expected values. 

call randseed(1);
prob = {0.3,0.6,0.1};
NumTrials = 10;
N = 1000;
x = RandMultinomial(N,NumTrials,prob);

/* population mean and covariance */
Mean = NumTrials * prob`;
Cov = -NumTrials*prob*prob`;
/* replace diagonal elements of Cov with Variance */
Variance = NumTrials*prob#(1-prob);
do i = 1 to nrow(prob);
   Cov[i,i] = Variance[i];
end;

SampleMean = mean(x);
SampleCov = cov(x);
print SampleMean Mean, SampleCov Cov;

Figure 25.305: Estimated Mean and Covariance Matrix

SampleMean   Mean  
2.891 6.059 1.05 3 6 1

SampleCov   Cov  
2.0051241 -1.69126 -0.313864 2.1 -1.8 -0.3
-1.69126 2.3198388 -0.628579 -1.8 2.4 -0.6
-0.313864 -0.628579 0.9424424 -0.3 -0.6 0.9 




            

            


            
For further details about sampling from the multinomial distribution, see Gentle (2003), or Fishman (1996).