Markov chains: the transition matrix, the stationary distribution and PageRank

Let (Ω, 𝓕, P) be a probability space, S a finite non-empty set (the set of states) and X₀, X₁, X₂, … : Ω → S a sequence of random variables (X_t is the state at time t ∈ ℕ₀). The sequence (X_t)_{t ≥ 0} is a Markov chain if and only if for every t ∈ ℕ₀ and all states i₀, …, i_t−1, i, j ∈ S such that P(X_t = i, X_t−1 = i_t−1, …, X₀ = i₀) > 0 the Markov property holds:

P(X_t+1 = j | X_t = i, X_t−1 = i_t−1, …, X₀ = i₀) = P(X_t+1 = j | X_t = i).

A Markov chain is homogeneous if and only if there are numbers p_ij ∈ [0, 1] (i, j ∈ S) such that

P(X_t+1 = j | X_t = i) = p_ij   for all t ∈ ℕ₀ and i, j ∈ S with P(X_t = i) > 0.

In words: the distribution of the next state depends only on the current state, not on the path by which the chain got there; for a homogeneous chain it does not depend on the time t either. The distribution π₀(i) = P(X₀ = i) is called the initial distribution.

Let n ≥ 1. A matrix A = (a_ij) ∈ ℝ^n×n is stochastic if and only if

1. a_ij ≥ 0 for all i, j ∈ {1, …, n};
2. ∑_j=1ⁿ a_ij = 1 for every i ∈ {1, …, n} (every row sums to 1).

Let (X_t)_{t ≥ 0} be a homogeneous Markov chain with set of states S = {1, …, n} and the numbers p_ij from the definition above. Its transition matrix is

P = (p_ij)_i,j=1ⁿ ∈ ℝ^n×n,   p_ij = P(X_t+1 = j | X_t = i).

The transition matrix is stochastic: for a fixed i the row (p_i1, …, p_in) is the probability distribution of the state X_t+1 given X_t = i. Conversely, every stochastic matrix together with an initial distribution π₀ determines a homogeneous Markov chain. (Using the same letter P for probability and for the matrix is standard; the context always makes clear which is which.)

Today is sunny (π₀ = (1, 0, 0)). What is the probability that the next three days will be sunny, cloudy, rain? By the formula for a trajectory

P = π₀(sunny) · p_{sunny→sunny} · p_{sunny→cloudy} · p_{cloudy→rain} = 1 · 0.6 · 0.3 · 0.3 = 0.054

(colours as in the formula for a trajectory above).

That is 5.4%. The sum over all 3³ = 27 possible triples of days is 1 — precisely because every row of the matrix sums to 1.

Let (X_t)_{t ≥ 0} be a homogeneous Markov chain with set of states S = {1, …, n} and transition matrix P ∈ ℝ^n×n. For t ∈ ℕ₀ let π_t = (π_t(1), …, π_t(n)) ∈ ℝ^1×n be the row vector of the distribution of X_t, i.e. π_t(i) = P(X_t = i). Then for every t ∈ ℕ₀

1. π_t+1 = π_t P, i.e. π_t+1(j) = ∑_i=1ⁿ π_t(i) p_ij for every j ∈ S;
2. π_t = π₀ P^t;
3. for every m ∈ ℕ and i, j ∈ S with P(X_t = i) > 0, P(X_t+m = j | X_t = i) = (P^m)_ij.

So the entry (i, j) of the matrix P^m is the probability of moving from i to j in m steps. Point 1 is the law of total probability with the partition by the value of X_t; point 2 follows from it by induction.

Tomorrow (after a sunny today) π₁ = (0.6, 0.3, 0.1). For the day after tomorrow, component by component — a column of the matrix multiplied by the vector π₁:

• sunny: 0.6 · 0.6 + 0.3 · 0.3 + 0.1 · 0.2 = 0.36 + 0.09 + 0.02 = 0.47,
• cloudy: 0.6 · 0.3 + 0.3 · 0.4 + 0.1 · 0.4 = 0.18 + 0.12 + 0.04 = 0.34,
• rain: 0.6 · 0.1 + 0.3 · 0.3 + 0.1 · 0.4 = 0.06 + 0.09 + 0.04 = 0.19.

In each line the first factor is π₁(i) and the second p_ij (colours as in the formula above). The sum is 1 again. In three days the same procedure gives (0.422, 0.353, 0.225).

Let P ∈ ℝ^n×n (n ≥ 2) be the transition matrix of a homogeneous Markov chain and assume that

1. P has n mutually distinct (in general complex) eigenvalues λ₁, …, λ_n, numbered so that λ₁ = 1 and |λ₂| ≥ |λ₃| ≥ … ≥ |λ_n|;
2. |λ_k| < 1 for all k = 2, …, n.

For k = 2, …, n let u_k ∈ ℂ^1×n be a left eigenvector belonging to λ_k (u_kP = λ_ku_k), and let π ∈ ℝ^1×n be the left eigenvector belonging to λ₁ = 1 whose components sum to 1 — the stationary distribution (definition below). Then for every initial distribution π₀ there are numbers a₂, …, a_n ∈ ℂ such that for all t ∈ ℕ₀

π_t = π₀P^t = π + a₂λ₂^tu₂ + … + a_nλ_n^tu_n.

Reason. Eigenvectors belonging to distinct eigenvalues form a basis, so π₀ = a₁π + ∑_k≥2 a_ku_k, and multiplying by P^t multiplies the k-th component by λ_k^t. For k ≥ 2 the components of u_k sum to 0 (from u_kP1 = u_k1, where 1 is the column of ones, and λ_k ≠ 1), so comparing the sums gives a₁ = 1.

The deviation π_t − π therefore shrinks to zero roughly like |λ₂|^t, where λ₂ is the second largest eigenvalue in absolute value.

The sum of the eigenvalues is the trace of the matrix, 0.6 + 0.4 + 0.4 = 1.4, and their product is the determinant, which comes out as 0.01. One eigenvalue is 1 (the rows sum to 1), so the other two satisfy λ₂ + λ₃ = 0.4 and λ₂λ₃ = 0.01. They are the roots of λ² − 0.4λ + 0.01 = 0:

λ_2,3 = 0.2 ± √0.03 = 2 ± √3/10,   so λ₂ ≈ 0.373 and λ₃ ≈ 0.027.

The deviation of the probability of sunshine from its stationary value 0.393, for a sunny start, is 0.0286, 0.0107 and 0.0040 after 3, 4 and 5 days — each day about 0.373 times smaller (colours as in the formula above). After a week it is 0.0006, and −0.0005 for a rainy start: the two beginnings differ only in the fourth decimal place. The deviation shrinks to one hundredth in t days, where 0.373^t = 0.01, so

t = ln 0.01/ln 0.373 ≈ 4.7 days.

Let P = (p_ij) ∈ ℝ^n×n be the transition matrix of a homogeneous Markov chain with set of states S = {1, …, n}. A row vector π = (π₁, …, π_n) ∈ ℝ^1×n is a stationary distribution of the chain if and only if

1. π_j ≥ 0 for every j ∈ S and ∑_j=1ⁿ π_j = 1 (π is a distribution on S);
2. π P = π, i.e. π_j = ∑_i=1ⁿ π_i p_ij for every j ∈ S.

Condition 2 says that π is a left eigenvector of the matrix P with eigenvalue 1. In words: if X₀ has the distribution π, then by the theorem on the distribution after t steps every X_t has the distribution π too — a step of the chain does not change the distribution.

Let (X_t)_{t ≥ 0} be a homogeneous Markov chain with a finite set of states S and transition matrix P. For i, j ∈ S and m ∈ ℕ let p_ij^(m) = (P^m)_ij denote the probability of moving from i to j in m steps.

1. The chain is irreducible if and only if ∀ i, j ∈ S ∃ m ∈ ℕ : p_ij^(m) > 0.
2. Let i ∈ S be a state for which the set R(i) = {m ∈ ℕ : p_ii^(m) > 0} is non-empty. The period of the state i is d(i) = gcd R(i), the greatest common divisor of the lengths of all paths that return to i.
3. A state i is aperiodic if and only if d(i) = 1. The chain is aperiodic if and only if all its states are aperiodic.

In an irreducible chain all states have the same period, so it is enough to check one state.

The weather is both straight away: all entries of P are positive (m = 1 suffices) and p_ii > 0, so a return is possible in 1 step and the period is 1. By contrast, a chain with two states that regularly alternates between them (p₁₂ = p₂₁ = 1) is irreducible but has period 2: π₀P^t for the start (1, 0) keeps jumping between (1, 0) and (0, 1) forever, even though the stationary distribution (½, ½) exists.

Let (X_t)_{t ≥ 0} be a homogeneous Markov chain with a finite set of states S = {1, …, n} and transition matrix P that is irreducible and aperiodic. Then

1. the chain has exactly one stationary distribution π = (π₁, …, π_n), and π_j > 0 for every j ∈ S;
2. for every initial distribution π₀, lim_t→∞ π₀P^t = π (componentwise);
3. for every initial distribution and every state j ∈ S, with probability 1 lim_T→∞ |{t ∈ {0, …, T − 1} : X_t = j}|/T = π_j, i.e. the fraction of time spent in state j tends to π_j;
4. for every state j ∈ S the mean return time to j is m_j = E[min{t ≥ 1 : X_t = j} | X₀ = j] = 1/π_j.

We are looking for π = (a, b, c) with πP = π. Column by column of the matrix:

• a = 0.6a + 0.3b + 0.2c, so 4a = 3b + 2c,
• b = 0.3a + 0.4b + 0.4c,
• c = 0.1a + 0.3b + 0.4c, so 6c = a + 3b.

The equations are dependent (adding all three gives a + b + c = a + b + c), so we drop the middle one. From the third, a = 6c − 3b; substituting into the first gives 24c − 12b = 3b + 2c, so b = 22c/15 and a = 6c − 66c/15 = 24c/15. The ratio is a : b : c = 24 : 22 : 15, and the condition a + b + c = 1 gives

π = (24/61, 22/61, 15/61) ≈ (0.393, 0.361, 0.246).

By the ergodic theorem, the mean return time to rain is m_rain = 1/π_rain = 61/15 ≈ 4.07 days.

311 letters give 310 neighbouring pairs. Vowel–vowel 6, vowel–consonant 111, consonant–vowel 111, consonant–consonant 82. Hence

p̂_vow→vow = 6/117 = 0.051,   p̂_vow→con = 111/117 = 0.949,

p̂_con→vow = 111/193 = 0.575,   p̂_con→con = 82/193 = 0.425

(colours as in the formula above; 117 = 6 + 111 pairs start with a vowel, 193 = 111 + 82 with a consonant).

Máj (colours as in the formula above):

π_vow = 0.575/0.949 + 0.575 = 0.377.

More precisely: there are as many vowel→consonant transitions as consonant→vowel ones (111), so

π_vow = 111/193/111/117 + 111/193 = 1/193/117 + 1 = 117/310.

The actual share of vowels is 117/311 — it differs only in the denominator, because the last letter of the stanza (“r”) starts no pair.

Onegin: from Markov’s numbers (a vowel after a vowel 12.8%, so vowel→consonant 87.2%; a vowel after a consonant 66.3%) we get

π_vow = 0.663/0.872 + 0.663 = 0.432

— exactly the 43.2% of vowels that Markov counted directly.

Let N ≥ 1 be the number of pages, numbered 1, …, N, and let L ⊆ {1, …, N} × {1, …, N} be the set of links: (i, j) ∈ L if and only if page i links to page j. For every page i let k_i = |{j : (i, j) ∈ L}| be the number of its outgoing links, and let d ∈ (0, 1) be the damping factor (usually d = 0.85). The random surfer chain is the homogeneous Markov chain with states {1, …, N} and transition matrix P = (p_ij), where

1. p_ij = 1 − d/N + d/k_i,   if (i, j) ∈ L;
2. p_ij = 1 − d/N,   if k_i ≥ 1 and (i, j) ∉ L;
3. p_ij = 1/N for all j,   if k_i = 0 (a page without links — the surfer jumps anywhere).

Every row sums to k_i · d/k_i + N · 1 − d/N = 1 (or N · 1/N = 1), so the matrix is stochastic. PageRank is the stationary distribution PR = (PR(1), …, PR(N)) of this chain, i.e. PR · P = PR and ∑_j PR(j) = 1.

In words: with probability d the surfer clicks on one of the links of the current page (all equally likely), otherwise they jump to any of the N pages. Because all p_ij > 0, the chain is irreducible and aperiodic, and by the ergodic theorem PageRank is unique.

For N = 6 and d = 0.85, 1 − d/N = 1 − 0.85/6 = 0.025. Only the link from C leads to A, and C has a single link (k_C = 1); only the link from A leads to B, and A has two (colours as in the formula above):

• PR(A) = 0.025 + 0.85 · 0.37351/1 = 0.34248,
• PR(B) = 0.025 + 0.85 · 0.34248/2 ≈ 0.1706,
• PR(C) = 0.025 + 0.85 · (PR(A)/2 + PR(B)/1 + PR(D)/2 + PR(E)/3 + PR(F)/1) = 0.025 + 0.85 · 0.41001 = 0.37351.

The values from the power method satisfy the equation — that is how you can tell the iteration has finished.

We convert the text of the stanza to lower case and replace punctuation with a space. In it, the letter “l” is followed by another character 31 times: 13 times by a space, 7 times by “a”, 6 times by “á”, twice by “o” and once each by “i”, “h”, “u”. By the estimation formula, then,

P(space | l) = 13/31 ≈ 0.42,   P(a | l) = 7/31 ≈ 0.23,   P(á | l) = 6/31 ≈ 0.19,

P(o | l) = 2/31 ≈ 0.06   and the others   1/31 ≈ 0.03

(colours as in the estimation formula for p̂_ij in the section on Máj).

Generation then does the same as the weather simulation: from the row of the current character it draws the next character. After “l” the word is most often at its end — “byl”, “lhal”, “kol”, “bol”.

Let V be a finite non-empty set of symbols (the vocabulary) and n ≥ 2 an integer. An n-gram model of order n is a collection of probability distributions q(· | v) on V, one for every context v = (v₁, …, v_n−1) ∈ Vⁿ⁻¹. A sequence of symbols W₁, W₂, … with values in V is generated by the model so that for every t ≥ n and w ∈ V

P(W_t = w | W₁, …, W_t−1) = q(w | W_t−n+1, …, W_t−1),

i.e. the next symbol depends only on the n − 1 previous ones. The sequence Y_t = (W_t−n+2, …, W_t) ∈ Vⁿ⁻¹ (t ≥ n − 1) is then a homogeneous Markov chain with set of states Vⁿ⁻¹: its state is the last n − 1 symbols. For n = 2 (bigrams) the state is a single symbol and q(w | v) is directly the transition probability p_vw.

With a vocabulary of V symbols, an n-gram model has |V|ⁿ⁻¹ states and |V| probabilities in each row. The stanza of Máj uses 32 different characters (letters and the space):

• bigrams (n = 2): 32 states, a 32 × 32 matrix = 1,024 numbers — more than the stanza has characters, so most of them come out as zero;
• trigrams: 32² = 1,024 states and over 32 thousand numbers — the stanza does not contain most pairs of letters at all, and their rows stay empty;
• a language model with a vocabulary of 50,000 tokens and a window of 1,000 tokens: 50,000^1,000 ≈ 10^4,699 states.

That is why a longer context needs either enormous data or a model that does not take the probabilities from a table but estimates them with a function — a neural network.