fix(neighbors): match Seurat's stored KNN/SNN graphs; the clustering divergence is not a defect

[shanikawm]

Closes the last item left open by #54: shanuz and Seurat disagreeing on clusters over an ifnb RPCA embedding that now matches to 30/30 PCs.

The divergence is not a defect

Ran Seurat's own RPCA embedding through both tools in three stages, with nn.method = "rann" so the neighbour search is exact on both sides.

stage	result
KNN indices	identical, cell for cell
SNN graph	agree off-diagonal to 2.8e-08 (float32 rounding)
community detection	this is the whole gap

Seurat's FindClusters runs n.start = 10 restarts of its own modularity optimiser and keeps the best; shanuz runs a single igraph multilevel pass. Seurat's partition scores 0.899903 to shanuz's 0.898336 under igraph's own modularity at γ=0.5 — so shanuz genuinely searches less hard, and 20 seeds never reach Seurat's value.

But the deeper optimum is the worse answer. shanuz's partition is a strict coarsening: no Seurat cluster is split, and one shanuz cluster absorbs two of Seurat's. Both are CD14 Mono, split along the batch axis:

Seurat cluster	cells	CTRL	STIM
0	2,607	73.8%	26.2%
1	1,740	16.7%	83.3%
dataset	13,999	46.8%	53.2%

The extra 0.17% modularity is spent re-discovering the batch effect integration just removed. Against seurat_annotations, shanuz scores ARI 0.9195 vs Seurat's 0.7368; batch mixing 0.8937 vs 0.8328.

Not a lucky seed: across 20 seeds shanuz gives 15 clusters 17 times, 85% beat Seurat's cell-type ARI, and none reaches its modularity. The Louvain search is therefore left alone and no n.start equivalent added — its main effect here would be to converge harder onto the batch split. Documented caveat: shanuz's single pass is more seed-sensitive (cell-type ARI 0.72–0.93 over those seeds).

Four real defects found proving that

All were invisible to find_clusters, whose igraph conversion takes the strict upper triangle and discarded exactly the entries that were missing. The graphs are stored objects users read directly.

	was	Seurat	now
nn graph	symmetrised	directed, nnz = n·k	directed
snn diagonal	dropped	SNN[i,i] = 1, all 13,999	kept
snn precision	float32 (~3e-08)	double	float64 (4.4e-16)
singletons	left alone	GroupSingletons	ported

Symmetrising was the costliest: Seurat's nn is the raw ranked table, so every row sums to k while column sums run 21–68 — that spread is the in-degree, and mat + mat.T erased it. Both in-tree consumers already symmetrise at the point of use, as Seurat's do.

Restoring the SNN diagonal meant following Seurat to the consumer too: RunUMAP.Graph opens with diag(x = object) <- 0, so run_umap now strips it rather than feeding the layout n zero-length self-edges.

Both graphs now reproduce Seurat exactly on ifnb — nn nnz 279980 with column sums [68, 24, 21, 23, 32], snn nnz 1120457 and sum 156062.938571.

Added

find_clusters(group_singletons=True) ports Seurat's rule: a size-1 cluster joins whichever non-singleton cluster it is most connected to, scored by mean SNN weight, candidate list fixed before the loop so one singleton cannot absorb another. ifnb has none at this resolution; it fires on sparser graphs.

Tests

11 new in tests/test_neighbors_graph_parity.py, 798 passing (was 787). Every guard mutation-tested — revert the fix, confirm a named test fails, 7/7 caught.

Two drafts were decorative and were rebuilt rather than kept green:

• the mean-vs-sum guard used w + w.T, which made both rules pick the same cluster;
• a guard claimed float32 straddled the prune threshold. It doesn't — a weak Python prune_snn is cast down to float32, so both sides of the comparison round identically. The code comment asserting that consequence was corrected too; float32 only ever changed the stored weights, not which edges survive.

A third trap worth recording: float32→float64 is a same-length edit, so a mutate/restore inside one mtime second left Python reusing the mutated .pyc and reporting a fix as caught when it wasn't. The harness now runs with bytecode disabled.

🤖 Generated with Claude Code