"We used Seurat" isn't reproducible: how your scRNA-seq package silently changes your results

A collaborator re-runs your single-cell analysis to double-check a figure before submission. They live in Python, you handed them a Seurat workflow, so they rebuild it in Scanpy from your methods section. Same counts matrix, same nominal steps: normalize, find variable genes, PCA, cluster, find markers. They get a different number of clusters, and roughly a third of your “significant” marker genes don’t come back. Nobody made a mistake. You just discovered that “we used Seurat” was never a complete description of what you did.

This isn’t a hypothetical. A 2026 Cell Systems study from the Pachter lab put the two standard single-cell stacks, Seurat and Scanpy, through the same data and measured how far their outputs drift. The short version: far enough that the choice of package, and even the version of that package, is a parameter that moves which genes and clusters you report. If you build single-cell pipelines for other people to trust, this is worth internalizing, because it changes what “reproducible” has to mean.

What the study actually measured

The work started from an unexpected observation. A member of the team, Lambda Moses, noticed that Seurat and Scanpy produced appreciably different PCA plots after each library’s default preprocessing on identical input, which is strange given that PCA is a mathematical procedure whose output shouldn’t depend on the implementation language. From there, the authors traced it step by step. A few of the numbers, because the specifics matter more than the headline:

Highly variable gene selection overlapped with a Jaccard index of just 0.22 under default settings. The two tools largely disagree about which genes are “variable.”
The nearest-neighbor graph had a median Jaccard index of 0.11, meaning the neighborhoods that clustering operates on are mostly not shared between the two.
Clustering agreed at an Adjusted Rand Index of 0.53 with default parameters. Better than chance, well short of “the same answer.”
Differential expression overlapped at a Jaccard of 0.62 for significant markers, with Seurat reporting roughly 50% more significant genes than Scanpy on the same data.

That last gap is worth reading correctly. Seurat reports more significant genes not because it is better powered or Scanpy is more conservative in any statistical sense, but because of default choices in how the two tools filter and correct. Seurat pre-filters candidate markers before testing and applies tie correction in its Wilcoxon test; Scanpy does neither by default. That is what drives the gap, and it happens despite Seurat using Bonferroni, the stricter of the two multiple-testing corrections (Scanpy defaults to the more permissive Benjamini-Hochberg). The 50% is a filtering-and-correction artifact, not a power difference.

The published summary framing is the one to remember. The authors state that the magnitude of the Seurat-versus-Scanpy difference is comparable to the variability introduced by sequencing fewer than 5% of reads, or analyzing fewer than 20% of cells, in benchmarking datasets. Read one way, switching packages perturbs your result about as much as discarding most of your data. The authors are careful to leave the interpretation open, though: the same result also says you can recover most of the signal from a small fraction of the reads and cells. Either the package differences matter biologically, or a similar answer is reachable from far less data, and this study doesn’t settle which. Either way, it isn’t a rounding difference you can wave off in a methods section.

One honest caveat that cuts against the scariest reading: quality control and normalization were largely equivalent between the packages when handed identical inputs. The divergence isn’t everywhere. It begins at feature selection and compounds from there, which, unfortunately, is exactly the part of the pipeline that feeds everything downstream.

One scoping note before you treat any of these numbers as constants. They come from a standard 10x PBMC benchmarking dataset (PBMC 10k, with PBMC 5k as validation), and the authors are explicit that the magnitudes are dataset-dependent. On the smaller PBMC dataset the downsampling equivalence shifted substantially, landing nearer 10 to 25% of reads rather than under 5%. Treat the figures as an existence proof of the effect, not as constants that transfer unchanged to your tissue.

Why it compounds instead of cancels

It’s tempting to hope these differences are noise that averages out. They don’t, and the reason is structural. Single-cell analysis is a chain where each step consumes the previous step’s output: variable genes define the feature space, PCA runs on those genes, the neighbor graph is built in PCA space, clustering partitions that graph, and DE tests the clusters. A 0.22 overlap in variable genes doesn’t stay contained. It reshapes the PCA, which reshapes the graph (0.11 overlap), which reshapes the clusters (0.53 ARI). By the time you’re testing markers, you’re testing different cells grouped differently in a different feature space. Small upstream disagreements become different biology downstream.

That’s also why “just use the defaults” is a false comfort. The defaults are where the tools disagree most, because each library’s authors made reasonable but different choices: a different variance-stabilization for HVGs, a different pseudocount in the fold-change, a different default resolution. The defaults are opinions, and Seurat’s opinions aren’t Scanpy’s.

A fair reading has to include the other half of the result, because most of these gaps are reconcilable. The HVG difference is a choice of algorithm, and each package’s default has an exact equivalent available in the other (seurat_v3 flavor in Scanpy, or mean.var.plot in Seurat). Matching the differential-expression settings (filtering, tie-correction, and multiple-testing method) raises the marker overlap from 0.62 to 0.73, and unifying the cluster assignments on top of that takes it to 0.99. The two packages implement Leiden clustering identically given the same input; the 0.53 ARI is Louvain under defaults, and much of current practice has moved to Leiden. What the study documents is not that the tools are irreconcilable, but that their defaults are, and that nobody is writing down which defaults they used. That framing is the useful one: the work is deliberate alignment, not fatalism.

Versions are unstated parameters too

Here’s the part that catches even careful teams. It isn’t only which package. It’s which release. The same study found that version bumps within a single tool move results, and the effect concentrates in differential expression, the output that most often becomes a figure:

Seurat v5 vs. v4 produced sizeable differences in the set of significant markers and in log-fold-change estimates, traced to a change in how the pseudocount is handled.
Scanpy v1.9 vs. v1.4 shifted DE results, which the authors attribute to a change in marker filtering behavior between those releases.
Cell Ranger v7 vs. v6 changed results across all steps, largely because the default flipped on whether to count intronic reads.

So “we used Seurat v5” is closer to reproducible than “we used Seurat,” but it’s still under-specified if you don’t say v5-point-what, and if you don’t pin the upstream quantifier, because whether intron counts are in your matrix changes the numbers before Seurat or Scanpy ever sees them. As of this writing the moving targets are real and current: Seurat 5.5.1 (CRAN, June 2026) and Scanpy 1.12.2 (PyPI, June 2026) are the latest releases, and both lines have shipped multiple point releases in the past year. “The latest version” is a date-stamped claim, not a stable one.

”We used Seurat” is not an analysis spec

None of this means Seurat is wrong and Scanpy is right, or vice versa. The study’s point, and ours, isn’t that one tool is more correct. It’s that the tool and its version are parameters you’re setting whether you record them or not, exactly like a p-value threshold or a clustering resolution. You wouldn’t write “we filtered by adjusted p-value” without the cutoff. Writing “we used Scanpy” without the version is the same omission.

Practically, a single-cell methods section that someone else can actually reproduce records:

The package and its exact version: Seurat 5.5.1, not “Seurat v5.”
The upstream quantifier and version, and its consequential defaults: Cell Ranger vs. kb/kallisto-bustools vs. STARsolo, and whether intronic reads were counted.
The parameters you changed from default and the ones you didn’t, because “default” is a version-dependent value, not a constant.
A pinned environment that can rebuild all of the above, not a prose list a reader has to reconstruct by hand.

That last point is where good intentions usually fail. A methods paragraph is a description of an environment; it is not the environment. The only thing that reliably reproduces a result is something that can rebuild the software stack, which is the same discipline that makes any reproducible pipeline work: pin the tools, pin the versions, capture them in a container or lockfile, and put that under version control alongside the code.

How we treat this on real projects

When we run single-cell analysis as a service, tool-and-version pinning isn’t a nice-to-have we add if there’s time. It’s load-bearing, because clients hand results to reviewers and sometimes to regulators. A few habits do most of the work:

Pin the whole stack, not just the top layer. The container records Seurat/Scanpy and the quantifier and the reference transcriptome build. A result that can’t name its own reference genome isn’t reproducible even if the package version is pinned.
Record which package produced each number, not just “single-cell analysis.” When a marker gene lands in a figure, the provenance trail should say it came from Seurat 5.5.1’s FindMarkers with these parameters, so that if a reviewer re-runs in Scanpy and gets a shorter list, the discrepancy is explainable, not alarming. This is exactly the kind of traceability our Atlas portal is being built to surface as its analysis suite rolls out, so that a number never floats free of the tool and version that made it.
Sanity-check that a result is biology, not an artifact of the stack. If a cluster or a marker only appears under one tool’s defaults, that’s a flag to investigate (cross-check with the other library, vary the resolution, confirm the marker holds up) before it becomes a claim. A finding that survives a tool switch is one you can defend; one that doesn’t is one you need to understand before you publish it.

Our own read, which goes beyond what the study set out to test, is that a great many published methods sections are underspecified relative to what would actually be needed to re-run them, not through fraud or sloppiness, but because the field treated “we used Seurat” as a complete specification when it never was. The fix isn’t to pick the “right” tool. It’s to treat the tool and its version as what they are: parameters that belong in your record, pinned in your environment, and traceable in your results.

If you’re standing up single-cell pipelines you’ll need to defend two years from now (to a reviewer, a collaborator re-running in the other language, or your future self), our pipeline development service builds exactly this discipline in from the start: pinned stacks, recorded runs, and provenance down to the package that produced each number.

Cytogence is the bioinformatics division of KeyQ, Inc. We build single-cell and multi-omics analysis our clients can trust, trace, and re-run. See how we approach reproducible pipelines.

References

Rich JM, Moses L, Einarsson PH, Jackson K, Luebbert L, Booeshaghi AS, Antonsson S, Sullivan DK, Bray N, Melsted P, Pachter L. “The impact of package selection and versioning on single-cell RNA-seq analysis.” Cell Systems 2026 (published online 27 Mar 2026). DOI 10.1016/j.cels.2026.101560. Preprint: bioRxiv 10.1101/2024.04.04.588111.
Tang L. “Package selection and versioning matter.” Nature Methods 23, 872 (2026). DOI 10.1038/s41592-026-03108-8.
Seurat 5.5.1: CRAN (June 2026). Scanpy 1.12.2: PyPI (June 2026).