PageRank, in depth: the random surfer, dangling nodes, a sparse implementation, personalised PageRank and running it in production

The random surfer and the update rule

Picture a surfer who, at each step, follows a random out-link from the current page with probability d, the damping factor, conventionally 0.85, and with probability 1 - d jumps to a page chosen from a teleport distribution, uniform by default. A page's PageRank is the long-run fraction of time the surfer spends there. Written as an update over all n pages:

r_new[v] = (1 - d) * p[v]
         + d * sum(r[u] / outdeg(u) for u with an edge u -> v)
         + d * (sum of r[u] over dangling u) * p[v]

The first term is the teleport. The second is the link-following. The third handles dangling nodes, pages with no out-links: without it their rank would leak out of the system each step and the scores would not sum to one. The convention here returns that mass to the teleport distribution p, which is what most libraries do by default.

Why the teleport? Without it, two problems arise. A group of pages that link only to each other, a rank sink, would absorb all rank. And a graph that is not strongly connected, which real graphs never are, would not have a unique answer. With teleportation, the chain can reach every page from every other, the stationary distribution exists and is unique, and power iteration converges to it. For the structure of real graphs see strongly connected components.

Worked example: five pages by hand

Five pages: A links to B and C, B links to C, C links to A, D links to C and E, and E links nowhere. Start with every page at 0.2 and do one step by hand with d = 0.85. E is dangling and holds 0.2, so d x 0.2 / 5 = 0.034 goes to every page, on top of the teleport floor (1 - 0.85) / 5 = 0.03. Page C receives half of A (0.1), all of B (0.2) and half of D (0.1), so its new value is 0.03 + 0.034 + 0.85 x 0.4 = 0.404. D receives no links, so it gets only 0.064. After that one step the vector is A 0.234, B 0.149, C 0.404, D 0.064, E 0.149, and it still sums to one.

Iterating to convergence gives A 0.3502, B 0.1884, C 0.3654, D 0.0396 and E 0.0564. Read the result. C is first because three pages feed it. A is a close second with a single in-link, because that link is C's only out-link, so A receives everything C has. Importance comes from who links to you and how many other places they link to, not from counting links. D, with no in-links, sits at the floor plus its share of E's spread rank: 0.03 + 0.85 x 0.0564 / 5 = 0.0396, which you can check exactly.

A sparse implementation

The same computation as sparse matrix-vector products. Build a matrix whose column u holds 1 / outdeg(u) in each row v that u links to; then each iteration is one sparse multiply plus two vector operations. This runs as written with NumPy and SciPy and prints the values above.

import numpy as np
import scipy.sparse as sp

def pagerank(src, dst, n, d=0.85, p=None, tol=1e-10, max_iter=200):
    """Power iteration. src[i] -> dst[i] are edges between node ids 0..n-1."""
    out_deg = np.bincount(src, minlength=n).astype(float)
    w = 1.0 / out_deg[src]                        # each edge carries 1/outdeg(source)
    M = sp.csr_matrix((w, (dst, src)), shape=(n, n))
    dangling = out_deg == 0
    p = np.full(n, 1.0 / n) if p is None else p / p.sum()
    r = p.copy()
    for it in range(1, max_iter + 1):
        leaked = r[dangling].sum()                # mass on nodes with no out-links
        r_new = d * (M @ r) + (d * leaked + (1 - d)) * p
        err = np.abs(r_new - r).sum()             # L1 change
        r = r_new
        if err < tol:
            return r, it
    raise RuntimeError(f"no convergence after {max_iter} iterations, err={err:.2e}")

names = ["A", "B", "C", "D", "E"]
edges = [("A", "B"), ("A", "C"), ("B", "C"), ("C", "A"), ("D", "C"), ("D", "E")]
ix = {k: i for i, k in enumerate(names)}
src = np.array([ix[u] for u, _ in edges])
dst = np.array([ix[v] for _, v in edges])
r, it = pagerank(src, dst, len(names))
print(it, {k: round(float(x), 4) for k, x in zip(names, r)})
# 45 {'A': 0.3502, 'B': 0.1884, 'C': 0.3654, 'D': 0.0396, 'E': 0.0564}

Three details are deliberate. Duplicate edges are summed by the CSR constructor, so dedupe first if a repeated link should not count twice. The tolerance is on the L1 change, which is the right norm for a probability vector. And the function raises instead of returning a half-converged vector that looks like an answer.

Test it with properties rather than only with the example. Scores must be non-negative and sum to one on any graph. On a directed cycle every node must score exactly 1 / n, by symmetry. Relabelling the nodes must permute the scores and change nothing else. A node with no in-links must score exactly the teleport floor plus its share of dangling mass. And for small random graphs the result must match a dense eigenvector solve of the full Google matrix to within the tolerance:

def dense_reference(src, dst, n, d=0.85):
    G = np.zeros((n, n))
    out_deg = np.bincount(src, minlength=n)
    for u, v in zip(src, dst):
        G[v, u] += 1.0 / out_deg[u]
    G[:, out_deg == 0] = 1.0 / n                  # dangling columns teleport uniformly
    G = d * G + (1 - d) / n                       # every column now sums to one
    vals, vecs = np.linalg.eig(G)
    v = np.real(vecs[:, np.argmax(np.real(vals))])
    return v / v.sum()

The dense version costs n squared memory, so it is a test oracle for graphs of a few hundred nodes, never a production path.

Personalised PageRank

Change the teleport distribution and you get personalised PageRank (PPR): the surfer jumps back to a chosen set of seed nodes instead of anywhere. Scores then measure proximity to the seeds rather than global importance. In the example, putting all teleport mass on D gives A 0.2558, B 0.1087, C 0.3009, D 0.2348 and E 0.0998. D jumps from last to third, E nearly doubles from 0.056 to 0.100 though it still finishes last, and A and B, reachable from D only through C, both fall.

PPR is how recommendation systems answer "what is near this user" and how graph-based retrieval for LLMs spreads relevance from entities matched in a query to related passages; see graph RAG versus vector search. For one seed set you can run the same code with a different p. For many, use local push algorithms, which touch only nodes with non-negligible score, or Monte Carlo random walks that sample the surfer directly; both trade exactness for running time proportional to the neighbourhood rather than the whole graph.

Convergence and the damping factor

Each iteration shrinks the L1 error by at least a factor of d, so the worst case is about log(tol / 2) / log(d) iterations: roughly 145 for a tolerance of 1e-10 at d = 0.85, and many more as d approaches one. Real graphs usually converge faster, as the example's 45 shows. Higher damping gives scores that follow the link structure more closely but converge more slowly and become more sensitive to small changes; 0.85 remains the common default for good reasons.

Choose the tolerance from how the scores are used. If consumers only read the top 1,000 ranks, stop when that ranking is stable between iterations; there is no value in the eleventh decimal place. Warm-starting from yesterday's vector usually cuts iterations sharply because graphs change slowly, but always run to the tolerance rather than a fixed count, so a large change in the graph still converges.

Scaling and the production pipeline

Memory is the first constraint. CSR stores one integer index and one float weight per edge plus n + 1 row pointers, so a billion edges with 32-bit indices and weights is about 8 GB, which fits on one large machine. You can avoid storing weights altogether by multiplying by r / outdeg first and using an unweighted adjacency. Below a few billion edges, a single machine with a good sparse library is usually faster and simpler than a cluster.

Beyond that, PageRank is the canonical vertex-centric program: each vertex sends r / outdeg along its edges, then sums its messages, in supersteps separated by a barrier. That is the model described in bulk synchronous parallel computing and implemented by Pregel-style systems and Spark GraphX. Partitioning to keep edges local matters more than raw cores, since the work per superstep is trivial and the communication is not.

Failure modes

• Dropping dangling mass. Scores sum to less than one and decay over iterations; rankings may look plausible while being wrong. Assert the sum each run.
• Inconsistent conventions. Libraries differ on dangling handling and on whether scores sum to one or to n. Comparing scores across tools without checking is a common source of false regressions.
• Edge direction reversed. Building the matrix as (src, dst) instead of (dst, src) computes something else entirely. Test on a graph where the answer is known.
• Link spam and farms. Rank is manipulable by creating many nodes that link to a target. Mitigations include dedupe by owner, down-weighting new or untrusted nodes, and biasing the teleport toward trusted seeds, the idea behind TrustRank.
• Silent input breakage. A truncated edge dump moves every score. Compare the edge count and the top-k overlap with the previous run before publishing.

Trade-offs

What to do next

1. Run the code above and reproduce the five-page values, then compute the first step by hand once.
2. Load one real graph you own, dedupe edges, and assert that scores sum to one after each run.
3. Choose and write down your dangling convention, damping and tolerance, and keep them fixed across tools.
4. Add validation before publishing: edge count drift, score sum, and top-k overlap with the last run.
5. Try personalised PageRank from a few seed nodes and inspect whether the neighbourhoods make sense.
6. Only when one machine is too small, move to a vertex-centric engine and compare its output with the single-machine version on a sample.