Inconsistency in graphs used by `sc.tl.leiden` and `sc.metrics.modularity` leads to modularity mismatch

[jpintar]

Please make sure these conditions are met

• I have checked that this issue has not already been reported.
• I have confirmed this bug exists on the latest version of scanpy.
• (optional) I have confirmed this bug exists on the main branch of scanpy.

What happened?

The igraph modularity recorded by sc.tl.leiden and the one computed by sc.metrics.modularity differ, when the reasonable expectation is for them to be equal.

This is because the graphs over which they are computed are not the same. When the graph is constructed during sc.tl.leiden, connectivities are passed to sc._utils.get_igraph_from_adjacency, resulting in a graph in which each pair of connected vertices is connected by two edges (possibly directed when flavor is not igraph). In sc.metrics.modularity, on the other hand, connectivities are passed to ig.Graph.Weighted_Adjacency, and the in resulting graph, each pair of connected vertices is connected by a single edge.

Minimal code sample

# /// script
# requires-python = ">=3.12"
# dependencies = [
#   "scanpy@git+https://github.com/scverse/scanpy.git@main",
#   "igraph",
# ]
# ///
#
# This script automatically imports the development branch of scanpy to check for issues
import igraph as ig
import scanpy as sc


adata = sc.datasets.pbmc3k_processed()
sc.tl.leiden(adata, flavor="igraph", n_iterations=10, directed=False, random_state=1234)

modularity_recorded = adata.uns["leiden"]["modularity"]
modularity_metrics = sc.metrics.modularity(adata, labels="leiden", mode="calculate")

# graph creation equivalent to sc.metrics.modularity
graph_1 = ig.Graph.Weighted_Adjacency(adata.obsp["connectivities"], mode=ig.ADJ_UNDIRECTED)
modularity_1 = graph_1.modularity(adata.obs["leiden"].array.codes, "weight")

# graph creation equivalent to sc.tl.leiden
graph_2 = sc._utils.get_igraph_from_adjacency(adata.obsp["connectivities"], directed=False)
modularity_2 = graph_2.modularity(adata.obs["leiden"].array.codes, "weight")

# creation methods agree for directed graphs
graph_3 = ig.Graph.Weighted_Adjacency(adata.obsp["connectivities"], mode=ig.ADJ_DIRECTED)
graph_4 = sc._utils.get_igraph_from_adjacency(adata.obsp["connectivities"], directed=True)

# tests
assert graph_3.isomorphic(graph_4)

assert not graph_2.isomorphic(graph_1)  # passes, should fail
assert graph_2.subisomorphic_vf2(graph_1)
assert graph_2.ecount() == 2 * graph_1.ecount()  # passes, should fail

assert modularity_1 == modularity_metrics
assert modularity_2 == modularity_recorded
assert modularity_metrics == modularity_recorded, "modularities unequal"  # fails

Error output

AssertionError: modularities unequal

Versions

Details

| Package | Version |
| ------- | ------- |
| igraph  | 1.0.0   |
| scanpy  | 1.12.1  |
| anndata | 0.12.11 |

| Dependency        | Version      |
| ----------------- | ------------ |
| parso             | 0.8.7        |
| numpy             | 2.4.3        |
| asttokens         | 3.0.1        |
| jedi              | 0.19.2       |
| cycler            | 0.12.1       |
| google-crc32c     | 1.8.0        |
| threadpoolctl     | 3.6.0        |
| traitlets         | 5.14.3       |
| fast-array-utils  | 1.4.1        |
| prompt_toolkit    | 3.0.52       |
| jupyter_client    | 8.8.0        |
| python-dateutil   | 2.9.0.post0  |
| Pygments          | 2.20.0       |
| kiwisolver        | 1.5.0        |
| packaging         | 26.2         |
| executing         | 2.2.1        |
| ipykernel         | 7.2.0        |
| donfig            | 0.8.1.post1  |
| debugpy           | 1.8.20       |
| platformdirs      | 4.9.6        |
| h5py              | 3.16.0       |
| six               | 1.17.0       |
| pure_eval         | 0.2.3        |
| PyYAML            | 6.0.3        |
| ipython           | 9.13.0       |
| pyparsing         | 3.3.2        |
| decorator         | 5.2.1        |
| scipy             | 1.17.1       |
| scikit-learn      | 1.8.0        |
| tornado           | 6.5.5        |
| psutil            | 7.2.2        |
| zarr              | 3.2.0        |
| comm              | 0.2.3        |
| numcodecs         | 0.16.5       |
| llvmlite          | 0.47.0       |
| typing_extensions | 4.15.0       |
| pillow            | 12.2.0       |
| natsort           | 8.4.0        |
| stack_data        | 0.6.3        |
| pyzmq             | 27.1.0       |
| setuptools        | 82.0.1       |
| session-info2     | 0.4.1        |
| wcwidth           | 0.7.0        |
| pytz              | 2026.1.post1 |
| msgpack           | 1.1.2        |
| legacy-api-wrap   | 1.5          |
| numba             | 0.65.1       |
| matplotlib        | 3.10.9       |
| scverse-misc      | 0.0.5        |
| appnope           | 0.1.4        |
| joblib            | 1.5.3        |
| pandas            | 2.3.3        |
| texttable         | 1.7.0        |
| jupyter_core      | 5.9.1        |

| Component | Info                                                                              |
| --------- | --------------------------------------------------------------------------------- |
| Python    | 3.13.13 | packaged by conda-forge | (main, Apr  8 2026, 02:29:07) [Clang 19.1.7 ] |
| OS        | macOS-26.4.1-arm64-arm-64bit-Mach-O                                               |
| CPU       | 10/10 logical CPU cores, arm                                                      |
| GPU       | No GPU found                                                                      |
| Updated   | 2026-05-03 15:31                                                                  |

[ilan-gold]

@jpintar I'm not sure these are really meant to produce identical results. It does seem reasonable that they would, but not everyone calculating modularity will have used leiden on the nieghborhood graph.

I do think it could be reasonable to check uns for leiden results and warn about this fact though. Or return them if they are present/calculated.

Maybe you could highlight the use-case a bit - why run modularity again?

[jpintar]

Definitely, people may want to calculate modularity of the kNN graph with respect to any partition, not just Leiden, and even Leiden can be computed in different ways. So it's not going to be a general expectation that the scores recorded in uns match the ones from sc.metrics.modularity, even for Leiden. And it's true – I don't have a good use-case for running sc.metrics.modularity on a Leiden partition now that modularity is already being recorded in uns (which was added in the same PR, I believe), except, perhaps, to compare different flavors against a common yardstick.

What actually puzzles me more is the implied lack of agreement about what the relevant kNN graph even is. For the purposes of sc.tl.leiden with directed=False it's the undirected multigraph with two edges between each pair of connected vertices (weighted as in connectivities), but for the purposes of sc.metrics.modularity, it's the undirected simple graph.

It seems the leidenalg implementation of Leiden (always?) returns the same partition for both those graphs, but the igraph implementation does not, which has extra implications for when it becomes the default in Scanpy v2.

My inclination would be to:

• use sc._utils.get_igraph_from_adjacency in sc.metrics.modularity and wherever else a kNN graph is being built from connectivities, for consistency
• change sc._utils.get_igraph_from_adjacency to create a simple graph when directed=False and the input matrix is symmetric (a breaking change for flavor=igraph, so perhaps in Scanpy v2?)

[ilan-gold]

Very interesting @jpintar, so I think this problem is independent of leiden as the flavor. If I change your flavor to leidenalg, I have the same fundamental problem - the metrics from modularity and leiden don't agree. This makes sense to me from the description of the issue/the code. The underlying problem is sc._utils.get_igraph_from_adjacency.

But you raise an interesting question - do we need to use the multigraph for leiden? This seems like a distinct question from the above.

I'll have a look into things algorithmically. I'm still not sure above the get_igraph_from_adjacency change. Pining @flying-sheep for visibility.

[flying-sheep]

I see you’re using scanpy 1.12.1, which contains the modularity fix in #4045. Good! That bug actually caused modularities to differ significantly.

Adding debug output to your test case above, we see modularity_metrics=0.6667003742937355 ≠ modularity_recorded=0.6667003742937329, which is pretty close!

What actually puzzles me more is the implied lack of agreement about what the relevant kNN graph even is. For the purposes of sc.tl.leiden with directed=False it's the undirected multigraph with two edges between each pair of connected vertices (weighted as in connectivities), but for the purposes of sc.metrics.modularity, it's the undirected simple graph.

Sure, there are probably inconsistencies in what we do, but depending on what you want to do, there isn’t just the kNN graph.

• the basic kNN distance graph, as it comes out of a kNN search, is a directed graph: for a given cell pair AB, A might be among the kNN of B while B has >k closer neighbors than A. It comes with distances that can be used as weights.
• from this graph, several others can be inferred based on the algorithm – discard direction, distances, or replace the distances with other weights (often connectivities).

[flying-sheep]

I’m pretty sure having a multigraph is not intentional. #4112 would fix things in a consistent way, but I checked, and fully switching to Weighted_Adjacency causes no qualitative difference between the results. The following tests fail, but all edges seem to be identical, only the graph layout and the dendrogram layout changes, in other words only nondeterministic randomness changes.

FAILED tests/test_paga.py::test_paga_plots[dense--paga] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[dense-continuous-func1] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[dense-continuous_obs-func2] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[dense-continuous_multiple-func3] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_pie[dense] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[sparse--paga] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[sparse-continuous-func1] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[sparse-continuous_obs-func2] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_plots[sparse-continuous_multiple-func3] - AssertionError: Image files did not match.
FAILED tests/test_paga.py::test_paga_pie[sparse] - AssertionError: Image files did not match.
FAILED tests/test_plotting.py::test_heatmap_var_as_dict - AssertionError: Image files did not match.

and then there’s tests/test_clustering.py::test_leiden_random_state[igraph-random_state], where the “different seeds can create different clusterings” part fails (as it also did when igraph 1.0 came out, which we fixed by using a different seed)

[jpintar]

Sure, there are probably inconsistencies in what we do, but depending on what you want to do, there isn’t just the kNN graph.

Yes, definitely. I should have said: given an adjacency matrix, there should be a single agreed-upon graph. I agree that #4113 is the right fix.

fully switching to Weighted_Adjacency causes no qualitative difference between the results

I'm not sure how generally that will hold. Consider the following:

import igraph as ig
import numpy as np
import scanpy as sc
import scipy as sp
from sklearn.metrics import adjusted_rand_score


ARIs = []
rng = np.random.default_rng(1234)
for _ in range(100):
    adjacency = sp.sparse.random(3000, 3000, density=0.01, format="csr", dtype=np.float64, rng=rng)
    adjacency = (adjacency + adjacency.T) / 2  # symmetrise
    adjacency.setdiag(0)
    adjacency.eliminate_zeros()

    graph_1 = ig.Graph.Weighted_Adjacency(adjacency, mode=ig.ADJ_UNDIRECTED)
    with sc._utils.random._set_igraph_rng(1234):
        part_1 = graph_1.community_leiden(
            weights="weight",
            objective_function="modularity",
            n_iterations=10,
        )

    graph_2 = sc._utils.get_igraph_from_adjacency(adjacency, directed=False)  # before fix in #4113
    with sc._utils.random._set_igraph_rng(1234):
        part_2 = graph_2.community_leiden(
            weights="weight",
            objective_function="modularity",
            n_iterations=10,
        )

    ARIs.append(adjusted_rand_score(part_1.membership, part_2.membership))

res = sp.stats.describe(ARIs)
for k, v in res._asdict().items():
    print(f"{k}: {v}")

On my Macbook, this gives

nobs: 100
minmax: (np.float64(0.4416185660532235), np.float64(0.8114675312566055))
mean: 0.6586289702883189
variance: 0.005895196302277514
skewness: -0.2208622838603807
kurtosis: -0.19022826811685212

so it does look like igraph finds quite different partitions for the multigraph and simple graph built from the same (symmetric) adjacency matrix.