Lecture 3

Generating Random Variables

1/1
Where does “randomness” come from ?

RANDOM
NUMBER
GENERATOR

True random Pseudo-random

generator generator

Human sources Natural sources

2/1
True random generator

Pros:
They are truly random numbers.

Cons:
• Random numbers are generated slowly (compared to pseudo-random
generator).
• Need a huge storage to record.
• Systematic error.

3/1
Pseudo-random number generator (Computer
algorithms)
▶ Pseudo-Random Number Generator (PRNG) are computer algorithms to
generate sequences of Unif[0,1] random variables u1 , u2 , · · ·, un .
PRNG is defined by the tuple {S, s1 , f , g} where
• S is the state space (set of all possible values during the process).
• s1 ∈ S is the seed.
• f : S → S.
• g : S → [0, 1].

▶ Uniform random variables u1 , u2 , · · · , un are simulated as follows:

• sk = f (sk−1 ). The common algo: sk = a ∗ sk−1 + b (mod m)
• uk = g(sk ) = sk /(m − 1) ∈ [0, 1]

4/1
Pseudo-random number generator

SEED

f f f f f
s1 s2 s3 ··· sd s1

g g g g g

u1 u2 u3 ··· ud u1

▶ If seed is not specified then R initialises it using a value taken from the
system clock.

5/1
Pseudo-random number generator

• Since sk ∈ {0, . . . , m − 1} it follows that after some finite number (of

at most m), the whole sequence will begin to repeat.

6/1
Pseudo-random number generator

• Since sk ∈ {0, . . . , m − 1} it follows that after some finite number (of

at most m), the whole sequence will begin to repeat.
• In general m should be chosen to satisfy three criteria:
◦ For any initial seed, the resultant sequence has the "appearance" of
being a sequence of independent Uni[0, 1] random variables.
◦ For any initial seed, the number of variables that can be generated
before repetition begins is large.
◦ The values can be computed efficiently on a digital computer.

6/1
Pseudo-random number generator

• Since sk ∈ {0, . . . , m − 1} it follows that after some finite number (of

at most m), the whole sequence will begin to repeat.
• In general m should be chosen to satisfy three criteria:
◦ For any initial seed, the resultant sequence has the "appearance" of
being a sequence of independent Uni[0, 1] random variables.
◦ For any initial seed, the number of variables that can be generated
before repetition begins is large.
◦ The values can be computed efficiently on a digital computer.

• One can choose m to be a sufficiently large number that can be fitted

to the computer word size. For example, the ideal choice for a 32-bit
word machine is m = 231 − 1.

6/1
Generating continuous random variables

Proposition. Let U be a uniform [0, 1] random variable. For any

continuous cumulative distribution function F , the random variable X
defined by
X = F −1 (U)
has distribution F .

The inverse transform algorithm

1. Simulate a random number U in [0, 1].
2. Return X = F −1 (U) where F −1 is the inverse function of F .
(Using uniroot function in R if there is no close form of F −1 ).
3. Repeat the above steps n times to have n realizations of X.

7/1
Generate continuous random variables
Example 1. A continuous random variable X ≥ 0, has distribution function:

1 2
FX (x ) = P(X ≤ x ) = 1 − − , x ≥ 0.
3(x + 1) 3(x + 1)2

• Is there a close form solution for X = F −1 (U) ?

• Calculate P(X ≤ 1) from the sample and compare to the true
probabilities.

▶ If F −1 does not have a closed-form, we can use uniroot function in R.

Command uniroot(f,c(a,b))$root returns the root of equation
f (x ) = 0 which lies in the interval (a, b).

Note: f (a) and f (b) must have opposite signs.

8/1
Generate continuous random variables
Solution.
1 t 2t 2
Closed-form solution: Let t = then 1 − − = u.
x +1 3 3
p
2 1 + 24(1 − u) − 1
→ 2t + t − 3(1 − u) = 0 → t =
4
4
→ X=p −1
1 + 24(1 − U) − 1

F_inv <- function(u){

return(4/(sqrt(1+24*(1-u))-1)-1)
}
U <- runif(1e5)
X <- F_inv(U)
mean(X<=1)

[1] 0.6672
9/1
Generate continuous random variables

Example 2. The total expenses of company ABC is calculated by the sales

amount X as follow:

Total expense = x 0.1 + x 0.2 + x 0.3

where X is a random variable with distribution function

P(X ≤ x ) = 1 − 0.3e −x − 0.7e −1.5x

• What is the expectation of total expense?

• What is VaR90% (Value at Risk) of the total expense?

10 / 1
Generate continuous random variables
n=1e4
u <- runif(n)
x <- rep(0,n)
for (i in 1:n){
F1 <- function(x){
return(1-0.3*exp(-0.5*x)-0.7*exp(-1.5*x)-u[i])}
x[i]<-uniroot(F1,c(0,100))$root
}
totalexpense <- xˆ0.1+xˆ0.2+xˆ0.3
mean(totalexpense)

[1] 2.738066
quantile(totalexpense,0.9) ### 90% Value at Risk

90%
3.609704

11 / 1
Generate continuous random variables

Variable Density Generate n obsers

Uniform(a, b) 1/(b − a) runif (n, a, b)

Exponential(λ) λe −λx rexp(n, λ)

Γ(a + b) a−1
Beta(a, b) x (1 − x )b−1 rbeta(n, a, b)
Γ(a) Γ(b)
β α α−1 −βx
Gamma(α, β) x e rgamma(n, α, β)
Γ(α)
(x − µ)2
!
1
Normal (µ, σ) √ exp − rnorm(n, µ, σ)
2πσ 2σ 2
Log-normal (µ, σ) Ln(X ) ∼ Normal(µ, σ) rlnorm(n, µ, σ)
α βα Install Pareto package
Pareto(α, β)
x α+1 rPareto(n, β, α)

12 / 1
Generate discrete random variables
Example 3. Let X be a Bernoulli random variable where

P(X = 1) = p
P(X = 0) = 1 − p

To simulate X
1. Generate U ∼ Uniform(0,1).
2. If U ≤ p return 1; else return 0.
p=0.3
u <- runif(1)
x <- ifelse(u<=p, 1, 0)
print(x)

[1] 0

13 / 1
Generate discrete random variables

Generate a vector contains n Bernoulli random variables

p=0.3; n=1e4
u <- runif(n)
x <- ifelse(u<=p, 1, 0)
print(x[1:20])

[1] 0 0 1 0 0 1 0 0 0 0 1 0 0 0 0 0 1 1 0 0
mean(x) #should be close to p

[1] 0.2986

14 / 1
Generate discrete random variables
The general method:

Discrete random variable X has probability mass function (pmf)

p0 = P(X = 0)
p1 = P(X = 1)
···
pk = P(X = k)
···
k
Let qk = P(X ≤ k) = pi with k = 0, 1, 2, · · ·
P
i=0

To simulate X :
1. Generate an uniform random variable U in [0, 1].
2. If qk−1 ≤ U < qk , we assign X ← k.
15 / 1
Generate discrete random variables
∞
•
P
pi = 1
i=0
k
• qk = P(X ≤ k) =
P
pi
i=0

▶ In practice, notice that we can check the condition qk−1 ≤ U < qk by

running a loop over k = 0, 1, 2, . . . and finding the first index k such that
U < qk .

16 / 1
Generate discrete random variables

17 / 1
Generate discrete random variables

Example 4. A discrete random variable X (X ≥ 0) has pmf:

1 2 1
P(X = k) = p(1 − p)k + e −1
3 3 k!
where p = 0.2.
a) Simulate a sample of N = 106 of value of X .
b) Calculate from the sample: P(X = 0), P(X = 1), · · · , and compare
them with the true probabilities.
c) If you notice, X is a mixed distribution of a Geometric distribution and
a Poisson distribution. Using this comment, can you suggest another
way to simulate X ?

18 / 1
Generate discrete random variables
Code 1:
m=100; p=0.2; N=10ˆ6
calculate_probability <- function(k){
1/3*p*(1-p)ˆk + 2/3*exp(-1)/factorial(k)}
pmf <- sapply(0:m, calculate_probability)
cdf <- cumsum(pmf) #calculate the cumulative sum qk
#Notice that qk=cdf[k+1] since R start indexing from 1

U <- runif(N,0,1)
X <- rep(0,N)
for (i in 1:N) {
k <- 0 #find the first index k such that U < qk
while (U[i] >= cdf[k+1]) {
k <- k+1
X[i] <- k }
}
mean(X==0) #should be close to P(X=0)=0.312
19 / 1
Generate discrete random variables

Code 2:
Once we have the pmf P(X = k), k = 0, . . . , m, we can simply using the
sample function to generate the discrete values corresponding to its given
probability.
m=100; p=0.2; N=10ˆ6
calculate_probability <- function(k){
1/3*p*(1-p)ˆk + 2/3*exp(-1)/factorial(k)}
pmf <- sapply(0:m, calculate_probability)

X <- sample(0:m, N, replace=TRUE, prob=pmf)

mean(X==0)

[1] 0.312381

20 / 1
Generate discrete random variables

Variables pmf P(X = k) Generate N obsers

Bernoulli(p) pk (1 − p)(1−k) , k ∈ {0, 1} rbinom(N, 1, p)
n
k
Binomial(n, p) k
p (1 − p)(n−k) rbinom(N, n, p)
Geometric(p) p (1 − p)k rgeom(N, p)
r +k−1
r
Negbinomial(r , p) k
p (1 − p)k rnbinom(N, r , p)
λk
Poisson (λ) e −λ rpois(N, λ)
k!

21 / 1
Exercises

Exercise 1. A random variable X has cumulative distribution function:




0 x <0
 x2

0≤x <1

F (x ) = 4
x


 2 1≤x <2

1 2≤x



2
Calculate E(e X + e X ).

22 / 1
Exercises

Exercise 2. Given that X1 , X2 , · · · , X5 are 5 independent Pareto random

variables with parameters:

α1 = 1.2; β1 = 2.0
α2 = 1.3; β2 = 1.8
α3 = 1.4; β3 = 1.6
α4 = 1.5; β4 = 1.4
α5 = 1.6; β5 = 1.2

Let S5 = X1 + X2 + X3 + X4 + X5 . Calculate
a) P(S5 > 10); P(S5 > 103 ); P(S5 > 105 ).
b) VaR99% (S5 ).

23 / 1
Exercises

Exercise 3. (Collective model) Total claim amount of a car driver in a

year is modeled by 
0
 N=0
SN = N
Xi N ≥ 1
P


i=1

where N is number of claims in the year and Xi is the claim amount of the
i th accident. N and Xi are independent variables. The following assumptions
are given about the number of claim N and the claim amounts Xi :
• N is a Poisson random variable with parameter λ = 2.
• Xi are i.i.d Pareto random variables with parameters α = 1.5 and
β = 2.5.
Calculate P(SN > 20) and VaR99% (SN )

24 / 1