BrevFeed
Sign in

47. Markov chains and stochastic matrices (bridge to probability course)

📖 Reading · 10 min
💡 Every code box below is live — edit it and hit Run.

This closing module trades theory for payoff — five applications, each built from a specific tool developed somewhere in Modules 16–22. This first one needs §22.3's power iteration and §22.4's Perron-Frobenius directly: a Markov chain, a system that moves between states with fixed probabilities, and whose long-run behavior turns out to be nothing but an eigenvector.

Stochastic matrices

A Markov chain on n states has a transition matrix P where P_{ij} is the probability of moving from state j to state i in one step. Every column sums to 1 (from state j, the system goes somewhere with total probability 1) and every entry is \ge0 — this is called a (column-)stochastic matrix.

If \vec x_k is the probability distribution over states at step k (a vector of nonnegative entries summing to 1), then

\vec x_{k+1}=P\vec x_k\ \Longrightarrow\ \vec x_k=P^k\vec x_0

— exactly §19.1's power formula, now applied to tracking probability over time instead of an abstract computation.

The stationary distribution

A stationary distribution \vec\pi satisfies P\vec\pi=\vec\pi — an eigenvector for eigenvalue 1. Every stochastic matrix has 1 as an eigenvalue: P^T\mathbb 1=\mathbb 1 (the all-ones vector, since each row of P^T — i.e. column of P — sums to 1), so 1 is an eigenvalue of P^T, and P^T,P share every eigenvalue (their characteristic polynomials are identical, since \det(P^T-\lambda I)=\det((P-\lambda I)^T)=\det(P-\lambda I), §16.4).

If P's entries are all strictly positive (or the chain is irreducible, §22.4), Perron-Frobenius guarantees \lambda=1 is the strictly dominant eigenvalue (every stochastic matrix has |\lambda_i|\le1 for all eigenvalues — not shown here — combined with Perron-Frobenius's guarantee of one dominant positive eigenvalue for positive matrices, that dominant eigenvalue must be exactly 1), with an eigenvector that can be chosen all-positive — exactly the probability-distribution shape a stationary distribution needs.

By §22.3's power iteration, this also proves: \vec x_k=P^k\vec x_0\to\vec\pi as k\to\infty, regardless of the starting distribution — the chain "forgets" where it started and converges to the same long-run behavior. This is the single most useful practical fact about Markov chains, and it's a direct corollary of two theorems this course already built, not a new one needing separate proof.

Doing it in Python

import numpy as np

# A 3-state weather model: columns = today, rows = tomorrow
# today:      sunny  cloudy  rainy
P = np.array([[0.7,  0.3,    0.2],   # tomorrow sunny
              [0.2,  0.4,    0.3],   # tomorrow cloudy
              [0.1,  0.3,    0.5]])  # tomorrow rainy

print("columns sum to 1:", np.allclose(P.sum(axis=0), 1))

x = np.array([1., 0., 0.])   # start: sunny for certain
for _ in range(30):
    x = P @ x

print("long-run distribution:", [round(v, 4) for v in x.tolist()])

# Confirm it's the eigenvector for eigenvalue 1
eigvals, eigvecs = np.linalg.eig(P)
idx = np.argmin(np.abs(eigvals - 1))
pi = eigvecs[:, idx].real
pi = pi / pi.sum()   # normalize to sum to 1 (a probability distribution)
print("eigenvector (normalized):", [round(v, 4) for v in pi.tolist()])

columns sum to 1: True
long-run distribution: [0.4565, 0.2826, 0.2609]
eigenvector (normalized): [0.4565, 0.2826, 0.2609]

Confirming convergence is independent of the starting distribution:

import numpy as np

P = np.array([[0.7, 0.3, 0.2], [0.2, 0.4, 0.3], [0.1, 0.3, 0.5]])

for name, x0 in [("all sunny", [1., 0., 0.]), ("all rainy", [0., 0., 1.]), ("uniform", [1/3, 1/3, 1/3])]:
    x = np.array(x0)
    for _ in range(30):
        x = P @ x
    print(f"{name}: converges to {[round(v, 4) for v in x.tolist()]}")

all sunny: converges to [0.4565, 0.2826, 0.2609]
all rainy: converges to [0.4565, 0.2826, 0.2609]
uniform: converges to [0.4565, 0.2826, 0.2609]

Worked example

Find the stationary distribution of P=\begin{pmatrix}0.9&0.4\\0.1&0.6\end{pmatrix} by solving P\vec\pi=\vec\pi directly.

P\vec\pi=\vec\pi\Rightarrow(P-I)\vec\pi=\vec0: \begin{pmatrix}-0.1&0.4\\0.1&-0.4\end{pmatrix}\vec\pi=\vec0\Rightarrow -0.1\pi_1+0.4\pi_2=0\Rightarrow\pi_1=4\pi_2.

Normalizing so \pi_1+\pi_2=1: 4\pi_2+\pi_2=1\Rightarrow\pi_2=0.2,\ \pi_1=0.8.

\boxed{\vec\pi=(0.8,0.2)}

Sanity check. P\vec\pi=\begin{pmatrix}0.9(0.8)+0.4(0.2)\\0.1(0.8)+0.6(0.2)\end{pmatrix}=\begin{pmatrix}0.72+0.08\\0.08+0.12\end{pmatrix}=\begin{pmatrix}0.8\\0.2\end{pmatrix}=\vec\pi ✓ — applying P leaves \vec\pi completely unchanged, confirming it's genuinely a fixed point (and both entries are positive and sum to 1, confirming it's a valid probability distribution too).

Your turn

1. Why must a stochastic matrix's dominant eigenvalue be exactly 1, never larger?

2. For P=I (staying in whatever state you start in, forever), what happens to power iteration's "forgets where it started" claim?

3. True or false: a stationary distribution is always unique.

Solutions

1. If some eigenvalue had |\lambda|>1, then P^k\vec x_0 would have a component growing without bound (§19.1's power formula) — but P^k\vec x_0 is always a valid probability distribution (entries \ge0 summing to 1), which can never grow unboundedly. So no eigenvalue can exceed 1 in magnitude, and 1 itself is always achieved (shown above), making it exactly the maximum.

2. It fails — P=I has every vector as an eigenvector for \lambda=1 (an n-dimensional eigenspace, the opposite of Perron-Frobenius's simple-eigenvalue requirement), so P^k\vec x_0=\vec x_0 for every k: the chain never moves at all, and the long-run distribution is just wherever it started. This is exactly why Perron-Frobenius's strict positivity (or irreducibility) hypothesis matters: P=I is about as reducible as a matrix can be (each state only reaches itself).

3. False, without irreducibility. A chain that splits into two disconnected groups of states (never crossing between them) has one stationary distribution per group — any convex combination of the two is also stationary, so there are infinitely many. Irreducibility (§22.4) is exactly the condition that rules this out and forces a unique stationary distribution.

Check yourself in code

For P=\begin{pmatrix}0.5&0.2\\0.5&0.8\end{pmatrix}, run power iteration for 40 steps from \vec x_0=(1,0) to find the stationary distribution.

Print exactly this:

stationary distribution: [0.2857, 0.7143]

import numpy as np

P = np.array([[0.5, 0.2], [0.5, 0.8]])
x = np.array([1., 0.])
for _ in range(40):
    x = P @ x

# print "stationary distribution: [...]" with x rounded to 4 decimal places

import numpy as np

P = np.array([[0.5, 0.2], [0.5, 0.8]])
x = np.array([1., 0.])
for _ in range(40):
    x = P @ x

print("stationary distribution:", [round(v, 4) for v in x.tolist()])

A stochastic matrix's columns are probability distributions, and it always has 1 as an eigenvalue; when positivity/irreducibility holds, Perron-Frobenius makes 1 strictly dominant, and power iteration guarantees convergence to the unique stationary distribution from any starting point — two theorems from Module 22, combined into the entire long-run theory of Markov chains for free.

Next: another eigenvector-based application, this time built on §21.4's SVD — Principal Component Analysis.

← Previous Next →