An independent Python implementation of the spatial ecotype algorithm published as SpatialEcoTyper (Zhang et al., Nature 2026, doi:10.1038/s41586-026-10452-4).
Spatial EcoTyper is an unsupervised framework for discovering spatial ecotypes — recurrent multicellular communities — from single-cell spatial transcriptomics. It builds one similarity network over spatial neighbourhoods per cell type, fuses them with similarity network fusion, clusters the fused graph, and then integrates spatial clusters across samples with consensus NMF to recover the ecotypes that recur.
Upstream is R + Seurat only. This is a complete Python re-implementation,
validated function-by-function against R 4.4.3 / Seurat 5.4.0 / NMF 0.28 /
SpatialEcoTyper 1.0.4, with the parity gate committed before any algorithmic
Python was written (data/manifest.yaml).
Every number below is measured, not asserted. Reproduce with
pytest tests/test_exact_match.py.
| Output | R function | Parity class | Gate | Measured |
|---|---|---|---|---|
| Filtered matrix | PreprocessST |
deterministic-strict | < 1e-13 | 0 (bit-exact) |
| Spatial metacells | GetSpatialMetacells |
deterministic | < 1e-8 | 3.6e-15 |
| Similarity network | getSN / GetSNList |
deterministic | < 1e-8 | 1.0e-15 (4/9 cell types); see Known divergences |
| Fused network | SNF2 |
deterministic | < 1e-8 | 4.2e-17 |
| Rank transform | rankSparse |
deterministic-strict | < 1e-13 | 0 (bit-exact) |
| Blocked matmul | matrixMultiply |
deterministic | < 1e-8 | 0 (bit-exact) |
| Weighted z-score | Znorm |
deterministic | < 1e-8 | 1.4e-14 |
| Cell-type PCs | GetPCList |
embedding (Procrustes) | ≥ 0.95 | 1.000 |
| SE labels, 1 sample | SpatialEcoTyper |
clustering (ARI) | ≥ 0.85 | 0.983 (6k cells) / 0.980 (27.9k cells) |
| Cell-level SE labels | AnnotateCells |
classification | ≥ 0.99 | 1.000 |
NMF basis W |
NMFGenerateW |
factorization (matched Pearson) | ≥ 0.95 | 1.000, max abs diff 7.5e-15 |
NMF loadings H |
NMFpredict |
factorization | ≥ 0.95 | 1.000, max abs diff 8.9e-16 |
| Integrated matrix | Integrate |
deterministic | < 1e-8 | 5.6e-16 (real 2-sample data) |
| Conserved SE labels | nmfClustering |
clustering (ARI) | ≥ 0.85 | 1.000; consensus matrix bit-identical |
| Coassociation index | Coassociation |
deterministic | < 1e-8 | 7.8e-16 |
| Coassociation p-values | CoassociationTest |
inference | ρ ≥ 0.90 | ρ = 1.000, Z max abs diff 4.4e-15 |
| Colocalization Z | Colocalization |
ordinal (Pearson) | ≥ 0.99 | 1.000, max abs diff 1.8e-13 |
| Colocalization meta | ColocalizationMetaAnalysis |
ordinal | ≥ 0.99 | 1.000 |
| Normalised Moran's I | ComputeNormalizedMoranI |
ordinal | ≥ 0.95 | 1.000, max abs diff 1.7e-13 |
| Recovery metrics | ComputeMetrics |
deterministic | < 1e-8 | 5.6e-16 |
| SE abundance by SN | ComputeSEAbundanceBySN |
deterministic | < 1e-8 | 5.5e-16 (payload columns) |
| Palettes | getColors |
byte equality | 100% | 323/323 cases identical to pals 1.10 |
R-function coverage: 35/35 exported functions (100%), 9/10 reachable internal helpers.
Full audit in AUDIT.md; full evidence in
RECONSTRUCTION_REPORT.md.
git clone https://github.com/omicverse/py-spatialecotyper
cd py-spatialecotyper
pip install -e ".[dev]"A wheel and sdist are built and twine check-clean
(dist/pyspatialecotyper-0.1.0-py3-none-any.whl).
Dependencies: numpy, scipy, pandas, scikit-learn, matplotlib, and
nmf-rs (the omicverse Rust port of
R's NMF, bit-equivalent to the brunet multiplicative updates within 1e-12).
A pure-NumPy fallback for the NMF kernels is built in, so the package still
works if the Rust wheel is unavailable for your platform.
import pandas as pd, scipy.sparse as sp
from pyspatialecotyper import SpatialEcoTyper, normalize_data
counts = pd.read_csv("Melanoma1_subset_counts.tsv.gz", sep="\t", index_col=0)
meta = pd.read_csv("Melanoma1_subset_scmeta.tsv", sep="\t", index_col=0)
meta = meta.loc[counts.columns] # X, Y, CellType columns required
normdata = normalize_data(sp.csc_matrix(counts.to_numpy(float)))
se = SpatialEcoTyper(normdata, meta, gene_names=list(counts.index),
radius=50, resolution=0.5, nfeatures=300).run()
print(se) # SpatialEcoTyper(1400 spatial neighborhoods, 8 SEs)
se.result.spot_metadata # neighbourhood-level metadata with an 'SE' column
se.result.metadata # single-cell metadata with an 'SE' column
se.result.fused # the fused, rank-transformed similarity graphFunction names map from R to snake_case; argument names and defaults are
verbatim. pyspatialecotyper.R_FUNCTION_MAP is the machine-readable
dictionary, and the mechanical transliterations (spatial_eco_typer,
nm_fpredict, …) are bound as aliases so code translated line-by-line from R
keeps working.
import pyspatialecotyper as se
expdat, meta, genes, cells = se.preprocess_st(normdata, meta, min_cells=5, min_features=10)
mc, mc_cols = se.get_spatial_metacells(expdat, meta, k=20, radius=50, gene_names=genes)
emb = se.get_pc_list(mc, mc_cols, genes, nfeatures=300)
nets, spots = se.get_sn_list(emb, npcs=20, k=50, min_cts_per_region=2)
fused = se.snf2([w for w, _ in nets.values()], t=10)
fused = se.rank_sparse(fused)Multi-sample integration (upstream Tutorial 2):
res = se.multi_spatial_ecotyper(
data_list={"SKCM": norm1, "CRC": norm2},
metadata_list={"SKCM": meta1, "CRC": meta2},
gene_names=genes, nmf_ranks=range(4, 13), nrun_per_rank=10, seed=1)
res["cluster_SE"] # conserved SE per spatial cluster
res["metadata"] # single-cell metadata with the conserved SE labelSE recovery in a new dataset, and bulk deconvolution:
Ws = se.nmf_generate_w_list(scdata, scmeta, gene_names=genes, Sample="Sample")
calls = se.recover_se(newdata, genes, cells, celltypes, Ws, min_score=0.6)
fracs = se.deconvolute_se(bulk, genes, samples, W, w_names, se_names)| R function | Python | Notes |
|---|---|---|
SpatialEcoTyper |
spatial_ecotyper / SpatialEcoTyper class |
single-sample discovery |
MultiSpatialEcoTyper |
multi_spatial_ecotyper |
end-to-end multi-sample |
IntegrateSpatialEcoTyper |
integrate_spatial_ecotyper |
integration from existing results |
PreprocessST, Znorm |
preprocess_st, znorm |
|
GetSpatialMetacells |
get_spatial_metacells |
|
SNF2, matrixMultiply, rankSparse |
snf2, matrix_multiply, rank_sparse |
|
NMFGenerateW, NMFGenerateWList, NMFpredict, nmfClustering |
nmf_generate_w, nmf_generate_w_list, nmf_predict, nmf_clustering |
|
RecoverSE, DeconvoluteSE, AggregateRecoverModels, LoocvPredict |
recover_se, deconvolute_se, aggregate_recover_models, loocv_predict |
|
Coassociation, CoassociationTest |
coassociation, coassociation_test |
|
Colocalization, ColocalizationMetaAnalysis |
colocalization, colocalization_meta_analysis |
|
ComputeMetrics, ComputeNormalizedMoranI |
compute_metrics, compute_normalized_moran_i |
|
ComputeSEAbundanceBySN, SmoothSEAbundances |
compute_se_abundance_by_sn, smooth_se_abundances |
|
CreatePseudobulks, InferNCells, AnnotateCells |
create_pseudobulks, infer_ncells, annotate_cells |
|
AverageMarkerExpression |
average_marker_expression |
gene sets must be supplied (upstream's are in the R package's inst/extdata) |
SpatialView, HeatmapView, CooccurrenceHeatmapView, drawRectangleAnnotation |
spatial_view, heatmap_view, cooccurrence_heatmap_view, draw_rectangle_annotation |
matplotlib |
getColors, mostFrequent |
get_colors, most_frequent |
Two substrate modules have no upstream Python equivalent and were written from the R/C++ sources:
pyspatialecotyper.rrandom— R's Mersenne-Twisterunif_rand(bit-identical), the post-R-3.6 rejection sampler behindsample()(identical), and the inversionrnorm(2.2e-16). This is what makes the permutation tests, the NMF restart seeds and the down-sampling reproducible rather than merely similar.pyspatialecotyper._modularity— Seurat'sComputeSNNand a port ofModularityOptimizer.cppincluding itsJavaRandomLCG, soFindClustersreproduces Seurat's Louvain exactly (element-wise identical labels given the same graph, verified on 31/31 parameter combinations).
# 1. reference (needs R 4.4.3 + Seurat 5.4.0 + NMF 0.28 + SpatialEcoTyper 1.0.4)
Rscript tests/r_reference_driver.R data/Melanoma1_subset_counts.tsv.gz \
data/Melanoma1_subset_scmeta.tsv \
reference_out/ci 6000
Rscript tests/r_nmf_driver.R reference_out/nmf
Rscript tests/r_nmfclust_driver.R reference_out/nmf
Rscript tests/r_stats_driver.R reference_out/stats
# 2. the pre-registered gate
pytest tests/test_exact_match.py -vEvery threshold is read out of data/manifest.yaml; no test hard-codes a
number, so the gate cannot be quietly widened.
getSNon tied nearest neighbours — the one gate this port fails. Spatial neighbourhoods can have exactly identical cell-type-specific metacell profiles (adjacent grid spots drawing the samek = 20nearest cells of a rare type). Their distances tie at the k-th neighbour boundary, where only one of them fits, andRANN's ANN kd-tree andscipy.spatial.cKDTreebreak the tie differently. Both are correct k-NN. Measured, feeding R's own PC embeddings to both sides: on the 6k-cell fixture, 5 of 9 cell types diverge,max abs0.108 on a 0-0.5 scale, and 72 of 88,436 stored entries differ (99.919% identical); on the full 27.9k-cell fixture, 314 of 483,550 stored entries differ — 99.935% identical, withnnzequal in both implementations for all 9 cell types. Downstream,SNF2on R's own networks still matches at 4.2e-17 and the end-to-end ARI is 0.98. Reported as a failure rather than papered over; the gate was not widened.FindNeighborsuses Annoy in Seurat 5, exact k-NN here. Annoy is approximate; measured on the full fixture, 95.5% of cells get an identical neighbour set and 99.7% of (cell, neighbour) pairs are shared, and the exact search's summed neighbour distance is ≤ Annoy's for 2133/2133 cells. This is the sole source of the end-to-end ARI being 0.98 rather than 1.00 — given R's own SNN graph, the Louvain labels are element-wise identical.NMFGenerateW's feature selection is not stable in R either. R sets no seed beforeNMF::rnmf, soW_0differs run to run. The fittedWis unaffected (KL is convex inWwhenHis fixed — seeMATH.md; measured R-vs-Rmax abs8.9e-16, Pearson 1.000), but the downstreamdelta > 0feature filter can flip near-threshold genes: two R runs on a smaller fixture kept 212 vs 218 features (Jaccard 0.92).Colocalization(ncores > 1)is not reproducible in R.mclapplyforks and each child reseeds from its PID. Two R runs with the sameset.seed(1)differ bymax abs11.29 (Pearson 0.991). Atncores = 1R is bit-reproducible, and that is what the port matches;colocalization()acceptsncoresfor signature parity and consumes the stream sequentially.normalization.method = "SCT"raisesNotImplementedError— it would require a port ofSCTransform, which is out of scope.- The bundled MERSCOPE recovery models (
inst/extdata/*.rds) are not redistributed.recover_seanddeconvolute_serequire explicitWs/W; the default-model branch, including the publishedSE01..SE11 → SE1..SE9relabelling, is not reachable without those files. ComputeMetricsdrops single-SE samples in R; the port keeps them. When a sample contains exactly one distinct true SE, R'sPrecision[, -1]degrades the data.frame to an unnamed vector, the subsequentmatch(ses, rownames(.))is allNA, and that sample silently vanishes from the cross-sample average. Measured effect on a deliberately degenerate partition containing one single-SE tile:max abs0.0878. The port does not reproduce this; on any input without a single-SE sample the two agree at 5.6e-16.nmfClusteringcophenetic coefficient differs from R in the 5th decimal (0.960832 vs 0.960812) on identical consensus matrices —hclustandscipy.cluster.hierarchy.averageorder equal merges differently. Cluster assignments are unaffected (ARI 1.000).
Full melanoma fixture (500 genes x 27,907 cells, 2133 spatial neighbourhoods), single node, 17 cores:
| R 4.4.3 | py-spatialecotyper | speed-up | |
|---|---|---|---|
SpatialEcoTyper end-to-end |
251.8 s | 84.1 s | 3.0x |
SNF2 alone |
122.9 s | (see ITERATION_LOG.md) |
|
Colocalization (same seed, same nperm) |
18.4 s | 4.5 s | 4.1x |
The Acceleration loop is logged in ITERATION_LOG.md with
per-iteration admissibility evidence; one rewrite was rejected and rolled back.
Perturbation bounds for the single bounded-ε rewrite are derived in
MATH.md.
examples/compare_R_vs_Python.ipynb— pipeline-level parity, one visualisation per manifest output.examples/tutorial_melanoma.ipynb— Python-only walkthrough, one subsection per public function.examples/function_by_function_R_parity.ipynb— R⇄Python dictionary with a parameter table and a numerical comparison per function.examples/evolution.ipynb— the Acceleration iterations, one panel each.
Built with the omicverse-rebuildr protocol and
depends on omicverse/rust-NMF for the
KL multiplicative-update kernels. Motivated by omicverse issue #760.
If you use this package, cite the original method:
Zhang, W. et al. Spatial ecotypes of the tumour microenvironment. Nature (2026). doi:10.1038/s41586-026-10452-4
and, optionally, this port:
py-spatialecotyper: a pure-Python port of SpatialEcoTyper. https://github.com/omicverse/py-spatialecotyper
GPL-3.0-or-later, the licence omicverse itself carries (see LICENSE).
This package is an independent implementation of a published algorithm, written from its mathematics rather than by translating the R sources, so it is not a derivative of the upstream distribution and does not carry its terms. The upstream R package is separately distributed by Stanford under a non-commercial agreement; that agreement governs their software, not this one.
Nothing here changes how the method should be cited — see Citation.