Modularity, in depth: the configuration null model, a worked score, the modularity matrix, tested code, resolution limits and how to use Q honestly

Definition from the null model

Take an undirected graph with m edges, adjacency matrix A and node degrees k_i. Imagine cutting every edge into two half-edges, or stubs, and reconnecting all 2m stubs uniformly at random. This is the configuration model. A stub of node i lands on one of node j's k_j stubs with probability about k_j/2m, and i has k_i stubs, so the expected number of edges between i and j is k_ik_j/2m.

Modularity sums observed minus expected over all ordered pairs in the same community and normalises by 2m:

Q = (1 / 2m) * sum over i, j of [ A_ij - k_i * k_j / (2m) ] * delta(c_i, c_j)

Grouping the terms by community gives the form you compute with. Let L_c be the number of edges with both ends in community c and d_c the sum of degrees of its nodes. Then Q = sum over c of [ L_c/m - (d_c/2m)² ]. The first term is the fraction of edges inside c; the second is the fraction expected inside c under the null model. Only community-level totals matter, which is why algorithms can update Q in constant time per move.

The range follows. Q is below 1 because the expected term is positive whenever a community has edges. Putting every node in one community gives exactly 0. Every node alone gives -sum (k_i/2m)², which is negative. Brandes and coauthors proved that Q is never below -1/2, and with c communities it cannot exceed 1 - 1/c.

Worked example

The graph above has m = 7 and degrees 2, 2, 3, 3, 2, 2, so 2m = 14. Score four partitions.

Read the singleton row carefully: a partition with no internal edges is penalised, not scored zero, because the null model still expects some edges inside each group. The all-in-one row is zero by construction, since observed and expected internal fractions are both 1. Good partitions therefore have to beat two trivial baselines, and Q measures by how much they beat the random expectation rather than how dense the groups are in absolute terms. A clique of three nodes in a graph of a million edges contributes almost its full internal fraction, while the same clique in a seven-edge graph pays a large expected term.

The natural split wins clearly. Moving node 2 across the bridge loses 0.235: it gives up two internal edges for one, and community B's degree total grows, which raises its expected term faster than its edge count.

The modularity matrix and spectral splits

Define the modularity matrix B = A - k k^T/2m. Each entry says how much more connected a pair is than chance. Every row of B sums to zero, so the all-ones vector is an eigenvector with eigenvalue 0. For a split into two groups, encode membership as s_i = +1 or -1; then Q = s^T B s / 4m.

Newman's leading-eigenvector method relaxes s to real values: the s that maximises s^TBs under a norm constraint is the eigenvector of the largest eigenvalue, and the signs of its entries give the split. If the largest eigenvalue is zero or below, no split improves Q and the group is left whole. Recursive bisection with a generalised matrix for subgroups, followed by node-moving refinement, gives good partitions on graphs of tens of thousands of nodes. It connects modularity to spectral methods and to eigenvector centrality, where the same power iteration machinery applies.

Computing Q and the move gain

The function below computes Q for a weighted edge list, and the second computes the change from moving one node, the quantity every local optimiser evaluates millions of times. Both accept the resolution parameter gamma discussed later.

from collections import defaultdict

def modularity(edges, part, gamma=1.0):
    """edges: (u, v, w) undirected; part: node -> community."""
    m, deg = 0.0, defaultdict(float)
    inside, tot = defaultdict(float), defaultdict(float)
    for u, v, w in edges:
        m += w
        deg[u] += w
        deg[v] += w
        if part[u] == part[v]:
            inside[part[u]] += w
    for node, k in deg.items():
        tot[part[node]] += k
    return sum(inside[c] / m - gamma * (tot[c] / (2 * m)) ** 2 for c in tot)

def move_gain(k_i, k_ia, k_ib, d_a, d_b, m, gamma=1.0):
    """Change in Q when node i leaves community a (d_a includes i) for b.
    k_ia, k_ib: edge weight from i to other members of a and to b."""
    return (k_ib - k_ia) / m - gamma * k_i * (d_b - d_a + k_i) / (2 * m * m)

edges = [(0, 1, 1), (0, 2, 1), (1, 2, 1), (3, 4, 1), (3, 5, 1), (4, 5, 1), (2, 3, 1)]
part = {0: "A", 1: "A", 2: "A", 3: "B", 4: "B", 5: "B"}
assert abs(modularity(edges, part) - 5 / 14) < 1e-12
gain = move_gain(k_i=3, k_ia=2, k_ib=1, d_a=7, d_b=7, m=7)
moved = {**part, 2: "B"}
assert abs(modularity(edges, moved) - modularity(edges, part) - gain) < 1e-12

Derive the gain by expanding the per-community form: leaving a removes k_ia inside edges and k_i degree from a, joining b adds k_ib and k_i, and the squared terms simplify to the expression in the code. We checked both functions against each other on a few thousand random small weighted graphs and random moves, a cheap test worth keeping in any implementation.

Variants

• Weighted. Replace A with edge weights and degree with strength (sum of incident weights); m becomes total weight. The code above already does this.
• Directed. Leicht and Newman's form uses out- and in-degrees: Q = (1/m) sum [A_ij - k_i^outk_jⁱⁿ/m] delta(c_i, c_j). Symmetrising a directed graph and using the undirected formula throws away who links to whom.
• Resolution. Reichardt and Bornholdt multiply the null-model term by gamma: Q_gamma = sum [L_c/m - gamma (d_c/2m)²]. Gamma above 1 favours more, smaller communities; below 1, fewer and larger.
• Other null models. Bipartite graphs (Barber's formulation), signed graphs and spatial networks each need a null model that matches what is random for that data. The structure Q = observed - expected stays the same.

Known limits

Resolution limit. Fortunato and Barthelemy showed that modularity maximisation cannot see communities below a scale set by the whole graph. The expected term for two small groups is tiny when m is large, so merging them gains more from the single edge between them than it loses. The standard example is a ring of cliques joined by single edges. For a ring of n cliques with l_c internal edges each, a direct calculation from the per-community form shows that merging neighbours in pairs scores higher once n exceeds 2(l_c + 1). For 5-node cliques (10 edges) that is any ring longer than 22: at 30 cliques, single cliques score 0.876 and pairs score 0.888, although each clique is obviously a community. Roughly, modules with fewer internal edges than about the square root of 2m are at risk.

Degeneracy. Good, de Montjoye and Clauset showed that real graphs typically have an enormous number of structurally different partitions whose Q is within a small fraction of the maximum. The optimum is a plateau, not a peak, so two runs of the same algorithm with different seeds can return different answers with nearly equal scores. Reporting one partition and its Q hides this.

Random graphs score well. Guimera, Sales-Pardo and Amaral showed that Erdos-Renyi random graphs, which have no planted structure, admit partitions with clearly positive modularity, especially when sparse. A Q of 0.4 is not by itself evidence of communities.

Hardness. Brandes and coauthors proved that finding the maximum-modularity partition is NP-hard. Every practical method is a heuristic, and its output depends on order and seed.

Using modularity in practice

Use modularity as a relative score with a baseline, not as an absolute grade.

1. Compare against a null distribution: rewire the graph many times preserving degrees, optimise each copy with the same algorithm, and report how far the real Q sits above that distribution.
2. Run the optimiser with many seeds and measure partition stability, for example with normalised mutual information between runs, or build a consensus partition from co-assignment frequencies.
3. Scan gamma and look for partitions that stay stable across a range of values, instead of trusting gamma = 1.
4. Never compare Q across graphs of different size or density; the achievable maximum depends on both.
5. Inspect communities whose internal edge count is near the square root of 2m; they may be merges forced by the resolution limit.

Prefer an optimiser that guarantees connected communities: Leiden fixes Louvain's known defect of sometimes returning internally disconnected groups. For data that is not a graph to begin with, it is often cleaner to cluster directly with k-means or DBSCAN than to build a similarity graph just to run modularity on it.

Common bugs

• Counting each undirected edge twice. Summing over an adjacency list where every edge appears in both directions doubles m; Q comes out wrong by a constant factor. Count m from a deduplicated edge list.
• Self-loops. Conventions differ on whether a self-loop adds 1 or 2 to degree. Pick one and apply it in both degree and internal-edge counts.
• Isolated nodes. They have degree zero and do not change Q, but some code divides by degree. Drop or handle them explicitly.
• Using Q to choose the number of clusters for a different algorithm without considering the resolution limit or the null distribution.

What to do next

1. Implement the two functions above and reproduce the 5/14 example.
2. Build a ring of 30 five-node cliques and confirm that pairs beat single cliques.
3. Run a modularity optimiser on your graph with 20 seeds and measure how much the partitions agree.
4. Generate degree-preserving rewired copies and compare their best Q with yours.
5. Scan gamma from 0.5 to 2 and record the number of communities and their stability.
6. Move on to Louvain and Leiden to see how the move gain drives a fast optimiser.