Simulations of stochastic phenomena often use random variables with a normal distribution. The Central Limit Theorem seemingly provides a justification for such a practice. But the sum of a large number of independently random variables does not always approach a normal distribution. The more general version of the Central Limit Theorem says that the limit will a stable distribution. The family of stable distribution includes the normal distribution as a special case. The more general family, identified by the French mathematician Paul Lévy, is characterized by four parameters.

• Alpha: usually α is called the stability parameter. For the normal distribution α = 2. Alpha has to be in the interval 0<α≤2.
• Beta: usually β is called the skewness parameter. For the normal distribution and any other symmetric distribution β = 0. Beta can have any real number value.
• Nu: ν is called the scale parameter or the dispersion parameter. For the normal distribution ν is equal to the standard deviation. For non-normal distribution ν has a value but it is not the same as the standard deviation, which for non-normal stable distributions is infinite. Nu can have any positive real number value.
• Delta: δ is called the mean or the measure of centrality. Delta can have any real number value.

For a normal distribution α=2, β=0, ν is equal to the standard deviation and δ is equal to the mean.

Three mathematical statisticians at Bell Laboratories in New Jersey developed an efficient algorithm for generating random variables with stable distributions with specified values for the above parameters. There names are J.M. Chambers, C.L. Mallows and B.W. Stuck. The algorithm will be referred to later as the CMS algorithm. It operates by taking two random variables which are uniformly distributed on the interval [0, 1] and combines them into a single variable which has the desired distribution.

The CMS takes as input the values of alpha and beta and generates a variable with μ (mean) zero and ν (dispersion) equal to unity. The output of the algorithm can be multiplied by any desired value of the dispersion parameter ν and an amount μ added to that result to give a random variable with the mean and dispersion.

The CMS Algorithm

• Define k(α) = 1 − |1−α|. Thus if α≤1 then k(α)=α and if α≥1 then k(α)=2−α.
• Compute φ₀ = −½β(k(α)/α).
• Transform β to β' by

  β' = β if α=1
  and
  β'= −tan(½π(1−α))tan(αφ₀) otherwise.

• Generate a random variable u uniformly distributed on the interval [0, 1] and compute φ = π(u−½).
• Compute ε=1−α and then τ = −εtan(αφ₀)
• Compute tan(½φ), tan(½εφ) and tan(½εφ)/(½εφ).
• Generate a random variable v which has a uniform distribution on the interval [0,1] and then compute w=−ln(v).
• Compute z = (cos(εφ)−tan(αφ₀)sin(εφ)/(wcos(φ)).
• Compute d = z^ε/α/ε
• Compute s = tan(αφ₀) + z^ε/α(sin(αφ)−tan(αφ₀)cos(αφ))/cos(φ)

The variable s has the desired distribution.

In the displays below samples of 2000 are generated using the CMS algorithm and the histogram compiled. Each time the refresh button is clicked new samples of 2000 are generated. The values that are beyond the lower and upper limits are tabulated at those limits.

First there is the normal distribution.

In this case the horizontal axis in terms of standard deviation units. Note the nearly complete absence of any values beyond four.

Note how much the distribution looks like that of a normal distribution except for the fact that there are many cases beyond four and even beyond eight.

This distribution is skewed to the left. Note the much higher frequency of the extreme cases of values of minus eight or beyond.