Random Walks, in depth: gambler's ruin, the square-root law, walks on graphs, hitting and mixing times, sampling and randomized 2-SAT

The walk on a line: gambler's ruin

Start on the integers. A walker at position i moves to i+1 with probability p and to i-1 with probability q = 1-p. With absorbing barriers at 0 and N, this is the gambler's ruin: you have i dollars, bet one dollar per round, and stop when broke or when you reach N.

For the fair walk (p = 1/2), call P(i) the probability of reaching N before 0. Conditioning on the first step gives P(i) = (P(i-1) + P(i+1)) / 2, so P is linear in i, and the boundary values P(0) = 0, P(N) = 1 force P(i) = i/N. The same first-step trick for the expected duration D(i) = 1 + (D(i-1) + D(i+1)) / 2 with D(0) = D(N) = 0 gives D(i) = i(N - i).

For a biased walk, with r = q/p not equal to 1, the success probability is (1 - r^i) / (1 - r^N). The bias matters enormously. With p = 0.49, i = 50 and N = 100, the fair answer of 0.5 drops to about 0.12. That single formula is why a small house edge wins in the long run.

import random

def ruin(i, N, p=0.5, trials=100_000):
    wins = steps = 0
    for _ in range(trials):
        x = i
        while 0 < x < N:
            x += 1 if random.random() < p else -1
            steps += 1
        wins += x == N
    return wins / trials, steps / trials

print(ruin(3, 10))         # about (0.30, 21): i/N and i(N-i)
print(ruin(50, 100, 0.49)) # about 0.12, not 0.5

Worked example: starting with 3 and stopping at 0 or 10, you win with probability 0.3 and the game lasts 21 rounds on average. Doubling everything to 6 and 20 keeps the win probability at 0.3 but makes the expected length 84 rounds. Duration grows with the square of the scale.

The walk on a line: gambler's ruin

Start on the integers. A walker at position i moves to i+1 with probability p and to i-1 with probability q = 1-p. With absorbing barriers at 0 and N, this is the gambler's ruin: you have i dollars, bet one dollar per round, and stop when broke or when you reach N.

For the fair walk (p = 1/2), call P(i) the probability of reaching N before 0. Conditioning on the first step gives P(i) = (P(i-1) + P(i+1)) / 2, so P is linear in i, and the boundary values P(0) = 0, P(N) = 1 force P(i) = i/N. The same first-step trick for the expected duration D(i) = 1 + (D(i-1) + D(i+1)) / 2 with D(0) = D(N) = 0 gives D(i) = i(N - i).

For a biased walk, with r = q/p not equal to 1, the success probability is (1 - r^i) / (1 - r^N). The bias matters enormously. With p = 0.49, i = 50 and N = 100, the fair answer of 0.5 drops to about 0.12. That single formula is why a small house edge wins in the long run.

import random

def ruin(i, N, p=0.5, trials=100_000):
    wins = steps = 0
    for _ in range(trials):
        x = i
        while 0 < x < N:
            x += 1 if random.random() < p else -1
            steps += 1
        wins += x == N
    return wins / trials, steps / trials

print(ruin(3, 10))         # about (0.30, 21): i/N and i(N-i)
print(ruin(50, 100, 0.49)) # about 0.12, not 0.5

Worked example: starting with 3 and stopping at 0 or 10, you win with probability 0.3 and the game lasts 21 rounds on average. Doubling everything to 6 and 20 keeps the win probability at 0.3 but makes the expected length 84 rounds. Duration grows with the square of the scale.

Spread and recurrence

Remove the barriers. After n fair steps of size one, the position S_n is a sum of n independent plus-or-minus ones, so E[S_n] = 0 and E[S_n^2] = n. The typical distance from the start is sqrt(n), not n. To travel distance d, a random walk needs about d^2 steps. That square law shows up again and again: it is why a local search that wanders randomly is slow to cross a large space, and why mixing times on path-like graphs grow quadratically.

Does the walker come back to its start? George Pólya proved in 1921 that on the integer lattice in one or two dimensions the answer is yes, with probability one, infinitely often: the walk is recurrent. In three or more dimensions it is transient. In three dimensions the probability of ever returning is about 0.34, so about two thirds of walkers wander off for good. The intuition is a count: the chance of being back at the origin after 2n steps is about n^(-d/2), and the sum over n diverges only when d is at most 2.

Walks on graphs

On a graph, a simple random walk moves from vertex u to a uniformly chosen neighbour. It is a Markov chain with transition matrix P = D^-1 A, where A is the adjacency matrix and D the diagonal matrix of degrees.

For a connected undirected graph with m edges, the stationary distribution is π(v) = deg(v) / 2m. Checking it takes one line: the flow into v is the sum over neighbours u of π(u) / deg(u) = 1 / 2m, and there are deg(v) such neighbours. So in the long run the walk spends time in proportion to degree. The expected return time to v is 1 / π(v) = 2m / deg(v).

Two conditions make the walk actually converge to π from any start. Connectivity makes π unique. Aperiodicity is needed too: on a bipartite graph the walk alternates sides forever and never settles. The standard fix is the lazy walk, which stays put with probability 1/2 and otherwise steps; it has the same stationary distribution and always converges.

import random
from collections import Counter

G = {0: [1, 2, 3], 1: [0, 2], 2: [0, 1, 3], 3: [0, 2, 4], 4: [3]}
m = sum(len(n) for n in G.values()) // 2           # 6 edges

def walk(G, start, steps, lazy=True):
    v, seen = start, Counter()
    for _ in range(steps):
        if not lazy or random.random() < 0.5:
            v = random.choice(G[v])
        seen[v] += 1
    return seen

seen = walk(G, 4, 1_000_000)
for v in G:
    print(v, round(seen[v] / 1_000_000, 3), round(len(G[v]) / (2 * m), 3))
# vertex 0: 0.25 vs 3/12; vertex 4: about 0.083 vs 1/12

That degree bias is the first thing to remember when you use walks on real data. A crawler that follows random links visits popular pages far more often than obscure ones, so averages computed over its visits are averages over edges, not vertices.

Hitting, commute and cover times

Three time quantities describe how fast a walk explores.

• Hitting time H(u, v): expected steps from u to first reach v. It is not symmetric. On a path of n vertices, going from one end to the other takes about n^2 steps, the square law again.
• Commute time C(u, v) = H(u, v) + H(v, u). It has an exact formula: C(u, v) = 2m * R_eff(u, v), where R_eff is the effective resistance between u and v when every edge is a one-ohm resistor. Many parallel paths mean low resistance and a short commute.
• Cover time: expected steps to visit every vertex. For any connected graph it is at most 2m(n - 1), by walking around a spanning tree, where each tree edge's commute time is at most 2m.

The extreme case is the lollipop: a clique on about 2n/3 vertices joined to a path of n/3. Starting in the clique, reaching the far end of the path takes on the order n^3 steps, because the walk keeps getting pulled back into the dense part. Cover-time bounds are the reason a randomized algorithm can decide s-t connectivity in logarithmic space: walk for 4mn steps and, by Markov's inequality, you reach t with probability at least 1/2 if it is reachable. The electrical view also powers uniform spanning-tree sampling, as shown in Random Spanning Tree, in depth.

Mixing time

Mixing time asks how many steps until the walk's distribution is close to π regardless of where it started, usually measured as total variation distance below 1/4. For a lazy walk it is controlled by the spectral gap γ = 1 - λ2, where λ2 is the second-largest eigenvalue of the transition matrix:

t_mix(eps) <= (1 / gap) * ln(1 / (eps * pi_min))

Graphs with a bottleneck, such as two dense communities joined by one edge, have a tiny gap and mix slowly: the walk stays in one community for a long time before crossing. Expanders, where every set has many edges leaving it, have a constant gap and mix in O(log n) steps. That is why a lazy walk on a well-connected social graph forgets its start in a few dozen steps while a walk on a long road network needs thousands.

In practice you rarely compute λ2 for a large graph. Instead, run several walks from very different starting points and watch a statistic you care about, such as the average degree of visited vertices; when the runs agree, you are past the mixing time for that statistic. That is the same diagnostic used for MCMC samplers, and the Metropolis acceptance rule behind them is developed in Simulated Annealing, in depth.

Sampling a graph you cannot download

Suppose you can only explore a huge graph by asking for a vertex's neighbours, for example a social network behind an API, and you want a uniform sample of vertices. A simple walk is biased toward high degree. Two fixes work.

• Reweighting. Keep the simple walk and weight each visited vertex by 1 / deg(v) when computing averages. Cheap and unbiased in the limit.
• Metropolis-Hastings walk. Propose a random neighbour w of v and accept the move with probability min(1, deg(v) / deg(w)); otherwise stay at v. The stationary distribution becomes uniform.

def mh_uniform_walk(G, start, steps):
    v, seen = start, Counter()
    for _ in range(steps):
        w = random.choice(G[v])
        if random.random() < len(G[v]) / len(G[w]):
            v = w                        # accept; otherwise stay put
        seen[v] += 1
    return seen

seen = mh_uniform_walk(G, 4, 1_000_000)
print({v: round(seen[v] / 1_000_000, 3) for v in G})   # each about 0.2

Discard the first stretch of the walk as burn-in, and remember that consecutive samples are correlated: a million steps is worth far fewer than a million independent samples, and the exchange rate is set by the mixing time.

An algorithm that is a random walk: 2-SAT

Random walks also prove algorithms fast. Christos Papadimitriou's algorithm for 2-SAT starts with any assignment and, while some clause is unsatisfied, picks one and flips one of its two variables at random.

def random_walk_2sat(n, clauses, max_rounds=None):
    """clauses: list of (a, b) literals; +k means x_k, -k means not x_k (1-based)."""
    x = [random.random() < 0.5 for _ in range(n + 1)]
    sat = lambda lit: x[abs(lit)] == (lit > 0)
    for _ in range(max_rounds or 2 * n * n):
        bad = [c for c in clauses if not (sat(c[0]) or sat(c[1]))]
        if not bad:
            return x[1:]
        lit = random.choice(random.choice(bad))
        x[abs(lit)] = not x[abs(lit)]
    return None                          # probably unsatisfiable

The analysis is gambler's ruin. Fix a satisfying assignment S and let k be the number of variables where the current assignment agrees with S. An unsatisfied clause has both literals false under the current assignment, but S makes at least one true, so at least one of the two variables disagrees with S. Flipping a random one increases k with probability at least 1/2. The walk on k from 0 to n is at least as favourable as a fair walk with a reflecting barrier at 0, whose expected time to reach n is at most n^2. Run 2n^2 flips and, by Markov's inequality, you fail to find a solution with probability at most 1/2; repeating drives the failure probability down exponentially. The same argument for 3-SAT gives a walk biased the wrong way, which is why Schöning's 3-SAT algorithm restarts frequently rather than walking long.

Failure modes

• Forgetting degree bias. Averaging over a simple walk's visits estimates an edge-weighted quantity. Reweight or use Metropolis-Hastings.
• Periodicity. A non-lazy walk on a bipartite graph never converges; counts oscillate between the two sides.
• Stopping before mixing. Short walks report the neighbourhood of the start. Use multiple starts and compare.
• Disconnected graphs and sinks. On a directed graph the walk can be trapped in a component or at a vertex with no out-edges; PageRank's teleport exists to fix exactly this.
• Trusting independence. Correlated samples make naive confidence intervals far too narrow; use batch means or effective sample size.
• Weak random number generators. Walks consume billions of draws; use a modern generator and seed each parallel walker independently.

What to do next

1. Run the gambler's ruin simulation and confirm i/N and i(N-i) yourself; then try p = 0.49 to see how strongly bias compounds.
2. Run the graph walk on a small graph of your own and check visits against deg(v)/2m.
3. Before using a walk to sample real data, decide whether you want vertex or edge averages, and add reweighting or a Metropolis-Hastings step accordingly.
4. Measure mixing empirically with several starting points before trusting any estimate.
5. Implement the randomized 2-SAT solver and compare its flip counts with the n^2 bound.
6. Continue with HITS, in depth and PageRank to see walk-based ranking, and with random spanning trees for the electrical view.