randomStable

Generates pseudorandom numbers from a stable distribution.

Synopsis

randomStable (nRandom, alpha, bprime)

Required Arguments

int nRandom (Input)
Number of random numbers to generate.
float alpha (Input)
Characteristic exponent of the stable distribution. This parameter must be positive and less than or equal to 2.
float bprime (Input)
Skewness parameter of the stable distribution. When bprime = 0, the distribution is symmetric. Unless alpha = 1, bprime is not the usual skewness parameter of the stable distribution. bprime must be greater than or equal to - 1 and less than or equal to 1.

Return Value

An integer array of length nRandom containing the random deviates.

Description

Function randomStable generates pseudorandom numbers from a stable distribution with parameters alpha and bprime. alpha is the usual characteristic exponent parameter α and bprime is related to the usual skewness parameter β of the stable distribution. With the restrictions 0 < α ≤ 2 and - 1 ≤ β ≤ 1, the characteristic function of the distribution is

\[φ(t) = \exp\left[-|t|α \exp(πiβ(1 - |1 - α|)\mathrm{sign}(t)/2)\right] \phantom{...} \text{ for } α ≠ 1\]

and

\[φ(t) = \exp\left[-|t|(1 + 2iβ \ln|t|)\mathrm{sign}(t)/π)\right] \phantom{...} \text{ for } α = 1\]

When \(\beta=0\), the distribution is symmetric. In this case, if \(\alpha =2\), the distribution is normal with mean 0 and variance 2; and if \(\alpha =1\), the distribution is Cauchy.

The parameterization using bprime and the algorithm used here are due to Chambers, Mallows, and Stuck (1976). The relationship between bprime = β’ and the standard β is

\[β' = -\tan(π(1 - α)/2) \tan(-πβ(1-|1-α|)/2) \phantom{...} \text{ for } α ≠ 1\]

and

\[β' = β \phantom{...} \text{ for } α = 1\]

The algorithm involves formation of the ratio of a uniform and an exponential random variate.

Example

In this example, randomStable is used to generate five pseudorandom symmetric stable variates with characteristic exponent 1.5. The tails of this distribution are heavier than those of a normal distribution, but not so heavy as those of a Cauchy distribution. The variance of this distribution does not exist, however. (This is the case for any stable distribution with characteristic exponent less than 2.)

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

n_random = 5
alpha = 1.5
bprime = 0.0
randomSeedSet(123457)
r = randomStable(n_random, alpha, bprime)
writeMatrix("Stable random deviates", r, noColLabels=True)

Output

 
                    Stable random deviates
      4.409        1.056        2.546        5.672        2.166