Copyright © 2023 by Karl Sigman

1 Simulating Markov chains

Here we give a quick introduction on how to simulate Markov chains in discrete time with
emphasis on a discrete state space.

1.1 Review of the definition of a Markov Chain

We shall mainly assume that the state space S = ZZ = {. . . , −2, −1, 0, 1, 2, . . .}, the integers, or
a proper subset of the integers.
Definition 1.1 A stochastic process {Xn } is called a Markov chain if for all times n ≥ 0 and
all states i0 , . . . , i, j ∈ S,

P (Xn+1 = j|Xn = i, Xn−1 = in−1 , . . . , X0 = i0 ) = P (Xn+1 = j|Xn = i) (1)

= Pij .

(1) mathematically expresses the Markov property which in words states that “The future, given
(conditional on) the present state, is independent of the past”. Here, the present state is that
of Xn , the future is {Xn+1 , Xn+2 , . . .} and the past is {X0 , . . . , Xn−1 }.
Pij denotes the probability that the chain, whenever in state i (regardless of what time n it is),
moves next (one unit of time later) into state j independent of the past, and is referred to as a
one-step transition probability. Since we are assuming that the transition probabilities do not
depend on the time n, in particular, using n = 0 in (1) yields

Pij = P (X1 = j|X0 = i).

The square matrix P = (Pij ), i, j ∈ S, is called the one-step transition matrix, and since
when leaving state i the chain must move to one of the states j ∈ S, each row sums to one
(e.g., forms a probability distribution): For each i ∈ S
X
Pij = 1.
j∈S

For any time n ≥ 0, the ith row of the transition matrix P is the probability mass
function (pmf ) of Xn+1 if Xn = i.
More generally, a discrete-time Markov process is any stochastic process {Xn } that satisfies
The Markov Property: for any time n ≥ 0,
Conditional on the rv Xn , the future sequence of rvs {Xn+1 , Xn+2 , . . .} is indepen-
dent of the past sequence of rvs {X0 , . . . , Xn−1 }.
More general state spaces are then allowed such as the case when S = R (the real numbers), or
even S = Rd (Euclidean d space.)

1
1.2 General Simulation Algorithm for Discrete-Time Markov Chains
The key to simulating a Markov chain is understanding what we explained earlier: For any
time n ≥ 0, the ith row of the transition matrix P is the pmf of Xn+1 if Xn = i. By
the Markov property, we thus see that if Xn = i (no matter what time n this happens and what
state i it is), then independent of the past, Xn+1 can be simulated by independently simulating
a random variable Yi distributed as the ith row of P , and setting Xn+1 = Yi . Specifically the
pmf of Yi is given by
P (Yi = j) = Pi,j , j ∈ S. (2)
For example,
∞
X
E(Xn+1 | Xn = i) = E(Yi ) = jPi,j ,
j=−∞

and represents the expected value of the next state visited if the state is now i.
Using the discrete inverse transform method (on the CDF Fi (x) = P (Yi ≤ x), x ∈ R), we
generate a U (uniform over (0, 1)) and set Yi = Fi−1 (U ), which is given by setting
j−1
X j
X
Yi = j, if Pi,m < U ≤ Pi,m , j ∈ Z. (3)
m=−∞ m=−∞

To see that this works (i.e., that the Yi defined in (3) has the distribution of the ith row of
P ): Fix an i.
If we define a(j) = j−1
P Pj
m=−∞ Pi,m and b(j) = m=−∞ Pi,m , then b(j) − a(j) = Pi,j and thus

P (Yi = j) = P (a(j) < U ≤ b(j)) = b(j) − a(j) = Pi,j .

This proves that Yi indeed has the distribution of the ith row of P ; the algorithm works. (The
(a(j), b(j)] subintervals as j varies are disjoint and cover all of (0, 1); thus U must fall in the
j th one for some j.)
When the state space is non-negative (S = {0, 1, 2, . . .}), this becomes: Set Yi = 0 if U ≤ Pi,0
and otherwise set
j−1
X X j
Yi = j, if Pi,m < U ≤ Pi,m , j ≥ 1.
m=0 m=0

For example, if S = {0, 1, 2, } and

 
1/2 1/3 1/6
P =  0 1/5 4/5  ,
1/3 0 2/3

then we can generate Y0 via 

0 if U ≤ 1/2,

Y0 = 1 1/2 < U ≤ 5/6,

2 if U > 5/6.


We can generate Y1 via (

1 if U ≤ 1/5,
Y1 =
2 if U > 1/5.

2
We can generate Y2 via (
0 if U ≤ 1/3,
Y2 =
2 if U > 1/3.

In the following algorithm, whenever we say “generate a Yi ”, we mean doing so

using the inverse transform method (3) using a new independent uniform U . Our
objective is to simulate the MC out to a desired time N to obtain (X0 , X1 , X2 , . . . , XN ) as the
output of our algorithm.

General Algorithm for simulating a Markov chain with given transition matrix P = (Pi,j ) up to
the first N steps:
1. Choose an initial value, X0 = i. Set n = 1, enter N ≥ 1.
2. Generate a Yi , and set Xn = Yi .
3. If n < N , then set i = Xn , set n = n + 1 and go back to (2).
4. Stop. Output (X0 , X1 , . . . , XN ).

Remark 1.1 It is the algorithm above that gave us the proof of the nice (theoretical) result
that every Markov chain can be represented as a recursion; Xn+1 = f (Xn , Un ), n ≥ 0, for some
function f (i, u) and an iid sequence of continuous uniform rvs {Un : n ≥ 0}. (Contained in our
Lecture Notes on Markov chains.) The function is

f (i, u) = Fi−1 (u), i ∈ S.

1.3 Examples in special settings

The general algorithm we gave is useful when one explicitly has the transition matrix P at
hand. In some cases, however, a simulation algorithm can be carried out easier by taking
advantage of the particular recursive structure of the given Markov chain, and the chain need
not be discrete-valued; instead it can be a discrete-time Markov process with a general state
space such as R. We give some examples here.
1. Simple random walk: With X0 = i chosen, we need to only generate iid copies of the
increments {∆n : n ≥ 1} distributed as P (∆ = 1) = p, P (∆ = −1) = 1 − p for
a given fixed p ∈ (0, 1); for example set ∆ = 1, if U ≤ p; set ∆ = −1 if U > p.
Then we just use Xn+1 = Xn + ∆n+1 , n ≥ 0 out to how ever many steps N we desire.
We could P even first generate all N iid ∆i : 1 ≤ i ≤ N , store them, and then define
Xn = i + nm=1 ∆m , 1 ≤ n ≤ N ; but then of course this has the disadvantage of forcing
us to first store all N iid copies of the ∆n as opposed to only one at a time.
2. General iid random walk: In this case, we allow the iid {∆n : n ≥ 1} to have any desired
probability distribution F (x) = P (∆ ≤ x), x ∈ R, and of course we must already have
in advance a simulation algorithm for generating iid copies distributed as this F . Such a
random walk still has the Markov property, even though it might no longer be discrete-
valued. For example, we might want F to be the normal distribution with a given mean
µ and variance σ 2 , the N (µ, σ 2 ) distribution. In this case, our initial condition could be
X0 = x ∈ R, and once again we just use Xn+1 = Xn + ∆n+1 , n ≥ 0 out to however many
steps N we desire by sequentially generating iid copies of ∆n distributed as F .

3
3. Product of iid rvs: If {Lk : k ≥ 1} is an iid sequence of rvs, then, with X0 the initial value
(assumed independent of the sequence {Lk },

Xn = X0 L1 × L2 × · · · Ln , n ≥ 0,

defines a Markov chain as we know since we have the recursion

Xn+1 = Xn Ln+1 , n ≥ 0.

Thus we need only sequentially simulate the Lk and use the recursion out N steps analo-
gous to Examples (1) and (2) above to obtain (X0 , X1 , X2 , . . . , XN ). This covers, for exam-
ple, the Binomial Lattice Model for risky assets, in which P (L = u) = p, P (L = d) = 1−p.