On the consistent and scalable detection of spatial patterns

TL;DR

The study addresses the inconsistent and computationally bottlenecked landscape of spatial-pattern detection in spatial omics. It unifies major approaches into a single quadratic-form framework $Q_n = \mathbf{z}^\top \mathbf{K} \mathbf{z}$ and analyzes when such tests are consistent, showing that they are capable only of detecting mean shifts unless the kernel is strictly definite. By adopting a CAR-based kernel $\mathbf{K} = (\mathbf{I} - \rho \tilde{\mathbf{W}})^{-1}$, the authors ensure a positive spectrum and universal consistency, while introducing scalable implementations (sparse solvers with Hutchinson trace estimation and FFT-based kernels) to handle millions of locations. They extend the framework to the bivariate case for co-expression and demonstrate, through simulations and real data (Visium HD and tumor lineage datasets), that Moran's I suffers from spectral cancellation and Moran-based methods can miss true patterns that CAR-based Q-tests detect. The work provides a principled, scalable foundation for reliable spatial pattern detection with broad applicability in spatial omics and single-cell lineage tracing.

Abstract

Detecting spatial patterns is fundamental to scientific discovery, yet current methods lack statistical consensus and face computational barriers when applied to large-scale spatial omics datasets. We unify major approaches through a single quadratic form and derive general consistency conditions. We reveal that several widely used methods, including Moran's I, are inconsistent, and propose scalable corrections. The resulting test enables robust pattern detection across millions of spatial locations and single-cell lineage-tracing datasets.

Figures (4)

• Figure 1: Theoretical analysis of spatial variability tests.(a) Summary of main theoretical results. (b) Q-test can only detect mean-shift patterns (left); variance shifts (right) remain invisible. (c) Spatial patterns can be modeled as deterministic functions in a Hilbert space decomposable into a frequency basis. The discrete kernel matrix converges to an integral operator with a corresponding spectral decomposition. (d) A hypothetical kernel spectrum illustrating two sources of false negatives: "blind spots" (signal aligns with the kernel's null space) and "cancellation" (signal energy is split between positive and negative eigenvalues summing to zero). (e) Spectral cancellation of Moran's I. The adjacency kernel with a fixed distance cutoff has an indefinite and oscillatory spectrum. (f) Graph Laplacian and its inverse (CAR) enforce strict definiteness, resolving cancellation but restricting sensitivity to high or low frequency regimes.
• Figure 2: Detecting heritability and plasticity in tumor lineage tracing data.(a) Schematic of intrinsic ($Q$) and relational ($R$) pattern detection showing four hypothetical genes: (1) heritable (phylogeny-aligned); (2) plastic (unaligned); (3) co-expressed with (1); and (4) anti-coexpressed with (1). (b) Cluster-level heritability concordance (Spearman $\rho$ between Moran’s I and CAR-$Q$) across 16 tumors. (c) Cell-state heritability in tumor 3730_NT_T2 showing $Q$-statistic comparison (left) and tumor phylogeny (right). (d) Distribution of cell-state heritability Z-scores (CAR-$Q$) across tumors. Red dashed lines indicate significance thresholds ($Z=\pm1.96$). (e) Gene-level performance heritability concordance. Left: correlation of test statistics for genes stratified by Moran’s I significance (p-val $< 0.05$). Right: cumulative recovery curves showing fraction of CAR-significant genes (p-val $< 0.05$) recovered when ranked by Moran’s I in decreasing order. (f) Gene expression heritability (left) and sparsity (right) in tumor 3430_NT_T2. Red: genes significant by both methods; blue: significant only by CAR (p-val $< 0.05$). (g) Log-normalized expression of discordant genes from (f). Left: phylogeny with tracks showing cell-state annotation and log-normalized expression (inner to outer: Rgs4, Col4a1, Tmem37, Cd200). Right: expression distributions by cell state. (h) Bivariate co-expression analysis. Left: Z-scores (Moran-$R$ vs. CAR-$R$) of gene pairs with regression outliers highlighted. Right: heatmap of gene modules clustered by CAR-$R$; green boxes indicate outlier pairs. (i) Spectral decomposition of a discordant gene pair (Spink5, Gdpd2). Left: graph spectral coefficients. Middle: eigenmode contribution to Moran's R showing cancellation between Mode 95 (red, positive) and Mode 605 (blue, negative). Right: phylogeny showing log-normalized expression of both genes and the two eigenmodes. (j) Proportion of method-specific outlier gene pairs across varying outlier thresholds ($\sigma$). Boxplots in b, d, e, g, j show median (center line), interquartile range (box), and 1.5× interquartile range (whiskers).
• Figure S1: Accuracy and efficiency of spectral approximation.(a) Runtime breakdown of Q-test operations across different kernels. Exact eigen-decomposition for $p$-value computation (Liu method) represents the primary bottleneck. (b) Comparison of $p$-values derived from the fast Welch-Satterthwaite approximation (y-axis) and the Liu method (x-axis) under both null and spatially variable scenarios, each with 50 replicates. (c) FFT-based computation of kernel spectra on hexagonal grids. Open: open boundary condition; Torus: periodic boundary condition.
• Figure S2: Scaling spatial variability tests to millions of spots.(a) Spectral analysis on rectangular grids. Normalized eigenvalues computed via FFT (torus boundary) closely approximate the exact matrix eigenvalues (open boundary) and theoretical values across all kernels. (b) Runtime and peak memory usage averaged over five replicates, performed on an M1 Macbook Air with 16GB RAM. By design, the CAR kernel switches to implicit mode for $n>5,000$. (c) Gene-ranking experiment results. A ground truth pattern was corrupted with varying noise levels to generate 100 genes, which were ranked using either the Q-statistic per kernel or its $p$-value (Welch approximation for matrix kernel, Liu method for FFT kernel). Boxplots show median (center line), interquartile range (box), and 1.5× interquartile range (whiskers). (d) Test consistency in a Visium HD mouse small intestine sample. Significance was called at p-adj $< 0.01$. For pairwise comparison, we extracted the top 2,000 genes ranked by each method and computed the Spearman correlation of the Q-statistic on the union gene set. (e) Detailed comparison between CAR and Moran’s I. (Left) Scatter plot of Z-scores (Q-statistic normalized by its mean and variance from the null) with a fitted linear regression model. (Right) Per-expression-level test agreement. Top bar shows the fraction of significant SV genes called by CAR (p-adj $< 0.01$) in each gene group. (f) Cross-resolution analysis (8 µm vs. 16 µm). Comparison of Z-scores highlights outliers like Reg3b (red), which exhibits a "Swiss roll" artifact distinct from immune-cluster signals like Coro1a (blue). (g) Radial power density profiles comparing genes with different spatial patterns highlighted in (f). (h) Spatial co-expression analysis using the R-statistic $\mathbf{z}^\top \mathbf{K} \mathbf{y}$. (Left) Scatter plot of Z-scores (R-statistic normalized by its mean and variance from the null) with a fitted linear regression model. (Right) Fourier coefficient analysis of two example genes, confirming that strong negative associations (high-frequency patterns, $\lambda \approx -1$) are rare in biological data, mitigating the theoretical limitations of Moran’s R.

Theorems & Definitions (20)

• Definition 2.1: Spatial Variability
• Proposition 2.1
• Definition 2.2: Q-statistic
• Theorem 1: Q-statistic measures mean dependence
• proof
• Corollary 1.1: Q-test consistency
• proof
• Definition 2.3: Conditional Spatial Variability
• Definition 3.1: Spatial Pattern
• Definition 3.2: Spatial Kernel Operators