9.3.14 Random variables: random_variable randvar
The randvar command produces an object representing a random
variable. The value(s) can be generated subsequently by calling
sample (see Section 9.3.1), rand (see
Section 9.3.1), randvector (see Section 9.3.15)
or randmatrix (see Section 9.3.16).
random_variable is a synonym for randvar.
-
randvar takes a sequence of arguments:
distspec, which can specify a probability distribution.
The distributions and their arguments are:
-
Uniform distribution (see Section 9.4.2):
Arguments:
-
uniform or uniformd.
- Either:
-
Optionally, mean=µ, to specify a mean of µ.
- Optionally, stddev=σ, to specify a standard
deviation of σ.
- Optionally, variance=σ2, to specify a
variance of σ2.
or:
-
a and b, two numbers specifying the end points of a
range.
The range can also be specified by a..b or
range=a..b.
- Binomial distribution (see Section 9.4.3):
Arguments:
-
binomial.
- Either:
-
n, a positive integer.
- p, a probability (a number between 0 and 1).
or:
-
Optionally, mean=µ, to specify a mean of µ.
- Optionally, stddev=σ, to specify a standard
deviation of σ.
- Optionally, variance=σ2, to specify a
variance of σ2.
- Negative binomial distribution (see Section 9.4.4):
Arguments:
-
negbinomial.
- n, a positive integer.
- p, a probability (a number between 0 and 1).
- Multinomial distribution (see Section 9.3.4):
Arguments:
-
randmultinomial.
- [p0,p1,…,pj], a list of probabilities with
p0+⋯ + pj=1.
- Optionally, [a0,a1,…,aj], a list of
possible return values.
- Poisson distribution (see Section 9.4.6):
Arguments:
-
poisson.
- λ, a real number
or
mean=µ, the mean value.
- Standard normal distribution (see Section 9.3.6):
Argument:
- Normal distribution (see Section 9.4.7):
Arguments:
-
normald.
- Either no arguments (for the standard normal distribution)
or one of:
-
Optionally, mean=µ, to specify a mean of µ.
- Optionally, stddev=σ, to specify a standard
deviation of σ.
- Optionally, variance=σ2, to specify a
variance of σ2.
or:
-
a and b, two numbers specifying the end points of a
range.
The range can also be specified by a..b or
range=a..b.
- Student’s distribution (see Section 9.4.8):
Arguments:
-
student.
- n, an integer (the degrees of freedom).
- χ2 distribution (see Section 9.4.9):
Arguments:
-
chisquare.
- n, an integer (the degrees of freedom).
- Fisher-Snédécor distribution (see
Section 9.4.10):
Arguments:
-
fisher, fisherd, or snedecor.
- n1 and n2, integers (the degrees of freedom).
- Gamma distribution (see Section 9.4.11):
Arguments:
-
gammad.
- a and b, real numbers.
- Beta distribution (see Section 9.4.12):
Arguments:
-
betaa.
- a and b, real numbers.
- Geometric distribution (see Section 9.4.13):
Arguments:
-
geometric.
- p, a probability (number between 0 and 1).
- Cauchy distribution (see Section 9.4.14):
Arguments:
-
cauchy or cauchyd.
- a and b, real numbers.
- Exponential distribution (see Section 9.4.15):
Arguments:
-
exponential or exponentiald.
- λ, a real number.
- Weibull’s distribution (see Section 9.4.16):
Arguments:
-
k, an integer.
- λ, a real number.
- Discrete distribution:
Arguments:
-
W=[w1,w2,…,wn], a list of nonnegative weights.
- Optionally, V=[v1,v2,…,vn], a list of values.
or:
-
[[v1,w1],[v2,w2],…,[vn,wn]], a list of of
object-weight pairs.
or:
-
f, a nonnegative function.
- a..b or range=a..b with real numbers a and b, a range specification.
- Optionally, N, a positive integer or
V=[v0,v1,v2,…,vn], a list of values with n=b−a (here a
and b have to be integers).
The weights are automatically scaled by the inverse of their sum to
obtain the values of the probability mass function. If a function f
is given instead of a list of weights, then wk=f(a+k) for
k=0,1,…,b−a unless N is given, in which case wk=f(xk)
where xk=a+(k−1) b−a/N and k=1,2,…,N. The resulting
random variable X has values in {0,1,…,n−1} for 0-based
modes (e.g. Xcas) resp. in {1,2…,n} for 1-based modes
(e.g. Maple). If the list V of custom objects is given, then V[X]
is returned instead of X. If N is given, then vk=xk for
k=1,2,…,N.
- randvar(distspec) returns an object
representing a random variable.
Examples.
-
Define a random variable with a
Fisher-Snedecor distribution (two degrees of freedom):
Input:
X:=random_variable(fisher,2,3)
Output:
To generate some values of X:
Input:
rand(X)
or:
sample(X)
Output:
Input:
randvector(5,X)
or:
sample(X,5)
Output:
⎡
⎣ | 2.2652,0.1397,6.3320,1.0556,0.2995 | ⎤
⎦ |
- Define a random variable with multinomial distribution:
Input:
M:=randvar(multinomial,[1/2,1/3,1/6],[a,b,c])
Output:
multinomial, | ⎡
⎢
⎢
⎣ | | , | | , | | ⎤
⎥
⎥
⎦ | , | ⎡
⎣ | a,b,c | ⎤
⎦ |
Input:
randvector(10,M)
Output:
⎡
⎣ | b,b,b,b,b,b,a,a,b,b | ⎤
⎦ |
Some continuous distributions can be defined by specifying their first
and/or second moment.
Examples.
-
Input:
randvector(10,randvar(poisson,mean=5))
Output:
⎡
⎣ | 7,2,5,6,7,9,8,4,3,4 | ⎤
⎦ |
- Input:
randvector(5,randvar(weibull,mean=5.0,stddev=1.5))
Output:
⎡
⎣ | 1.6124,3.2720,7.02627,5.5360,3.1929 | ⎤
⎦ |
- Input:
X:=randvar(binomial,mean=18,stddev=4)
Output:
- Input:
X:=randvar(weibull,mean=12.5,variance=1)
Output:
weibulld | ⎛
⎝ | 3.08574940721,13.9803128143 | ⎞
⎠ |
Input:
mean(randvector(1000,X))
Output:
- Input:
G:=randvar(geometric,stddev=2.5)
Output:
geometric | ⎛
⎝ | 0.327921561087 | ⎞
⎠ |
Input:
evalf(stddev(randvector(1000,G)))
Output:
- Input:
randvar(gammad,mean=12,variance=4)
Output:
Uniformly distributed random variables can be defined by specifying
the support as an interval.
Examples:
-
Input:
randvector(5,randvar(uniform,range=15..81))
Output:
⎡
⎣ | 77.0025,77.7644,63.2414,52.0707,66.3837 | ⎤
⎦ |
- Input:
rand(randvar(uniform,e..pi))
Output:
The following examples demonstrate various ways to define a discrete
random variable.
Examples:
Discrete random variables can be used to approximate custom continuous
random variables. For example, consider a probability density function
f as a mixture of two normal distributions on the support
S=[−10,10]. You can sample f in N=10000 points in S.
Input:
F:=normald(3,2,x)+normald(-5,1,x):; |
c:=integrate(F,x=-10..10):; |
f:=unapply(1/c*F,x):; |
X:=randvar(f,range=-10..10,10000):;
|
Now generate 25000 values of X and plot a histogram:
Input:
R:=sample(X,25000):; |
hist:=histogram(R,-10,0.1):; |
PDF:=plot(f(x),display=red+line_width_2):; |
hist,PDF
|
Output:
Sampling from discrete distributions is fast: for instance, generating 25 million
samples from the distribution of X which has about 10000 outcomes
takes only couple of seconds. In fact, the sampling complexity is
constant. Also, observe that the process isn’t slowed down by spreading
it across multiple calls of randvector.
Input:
for k from 1 to 1000 do randvector(25000,X); od:;
Evaluation time: 2.12
Independent random variables can be combined in an expression,
yielding a new random variable. In the example below, you define a
log-normally distributed variable Y from a variable X with standard
normal distribution.
Input:
X:=randvar(normal):; mu,sigma:=1.0,0.5:; |
Y:=exp(mu+sigma*X):; |
L:=randvector(10000,Y):; |
histogram(L,0,0.33)
|
Output:
It is known that E[Y]=eµ+σ2/2. The mean of L
should be close to that number.
Input:
mean(L); exp(mu+sigma^2/2)
Output:
In case a compound random variable is defined as an expression
containing several independent random variables X,Y,… of the
same type, you sometimes need to prevent its evaluation when
passing it to randvector and similar functions.
Example.
Input:
X:=randvar(normal):; Y:=randvar(normal):;
If you want to generate, for example, the random variable X/Y, you
would have to forbid automatic evaluation of the latter expression;
otherwise it would reduce to 1 since X and Y are both
normald(0,1).
Input:
randvector(5,eval(X/Y,0))
Output (for example):
⎡
⎣ | −0.358479277895,5.03004946974,−5.5414073892,−0.885656967277,−2.63689662108 | ⎤
⎦ |
To save typing, you can define Z with eval(∗,0) and
pass eval(Z,1) to randvector or
randmatrix.
Input:
Z:=eval(X/Y,0):; randvector(5,eval(Z,1))
Output (for example):
⎡
⎣ | 0.404123429613,−4.06194898981,0.00356038536404,1.61619003525,−2.85682173195 | ⎤
⎦ |
Parameters of a distribution can be entered as symbols to allow
(re)assigning them at any time.
Input:
purge(lambda):; |
X:=randvar(exp,lambda):; |
lambda:=1:;
|
Now execute the following command line several times in a row. The
parameter λ is updated in each iteration.
Input:
r:=rand(X); lambda:=sqrt(r)
Output (by executing the above command line three times):
8.5682,2.9272 |
1.5702,1.2531 |
0.53244,0.72968
|