Modeling Microenvironment Trajectories on Spatial Transcriptomics with NicheFlow
Authors: Kristiyan Sakalyan
Alessandro Palma
Filippo Guerranti
Fabian J Theis
Stephan Günnemann
Venue: NeurIPS 2025 poster ()
Keywords: flow matching
spatial transcriptomics
microenvironments
trajectory inference
generative modeling
point cloud
Primary area: machine_learning_for_sciences
TL;DR: We develop a Variational Flow Matching algortihm to learn microenvironment trajectories on time-resolved spatial transcriptomics.

Abstract: Understanding the evolution of cellular microenvironments in spatiotemporal data is essential for deciphering tissue development and disease progression. While experimental techniques like spatial transcriptomics now enable high-resolution mapping of tissue organization across space and time, current methods that model cellular evolution operate at the single-cell level, overlooking the coordinated development of cellular states in a tissue. We introduce NicheFlow, a flow-based generative model that infers the temporal trajectory of cellular microenvironments across sequential spatial slides. By representing local cell neighborhoods as point clouds, NicheFlow jointly models the evolution of cell states and spatial coordinates using optimal transport and Variational Flow Matching. Our approach successfully recovers both global spatial architecture and local microenvironment composition across diverse spatiotemporal datasets, from embryonic to brain development.

Modeling Microenvironment Trajectories on Spatial Transcriptomics with NicheFlow
Kristiyan Sakalyan" Alessandro Palma" 1,2233 Filippo Guerranti
Fabian J.
Stephan Günnemann' 1,2,4
School of Computation, Information and Technology. Technical University of Munich Munich Centre for Machine Learning (MCML)
Institute of Computational Biology, Helmholtz Munich
'Munich Data Science Institute (MDSI), Technical University of Munich
TTUM School of Life Sciences Weihenstephan, Technical University of Munich Correspondence to: {k. sakalyan, palma, f. guerraatt /utum...t.m..... de
Abstract
Understanding the evolution cellular microenviromments spatiotemporal data is essential for deciphering tissue development and disease progression. While experimental techniques like spatial transcriptomics now enable high-resolution mapping of tissue organization across space and time, current methods that model cellular evolution operate at the single-cell level, overlooking the coordinated development of cellular states in tissue. We introduce NicheFlow, a flwwbbsedd generative model that infers the temporal trajectory of cellular microenvironments across sequential spatial slides. By representing local cell neighborhoods as point clouds, NicheFlow jointly models the evolution of cell states and spatial coordinates using optimal transport and Variational Flow Matching. Our approach successfully recovers both global spatial architecture and local microenvironment composition across diverse spatiotemporal datasets, from embryonic brain development'.
Introduction
Uncovering the principles governing tissue organization across space and time remains one of the most fundamental challenges in biology, with profound implications for evolutionary and developmental studies [1, 2]. While individual cells form the basic units of biological systems, they operate not in isolation but integral parts of spatially organized microenvironments, functionally distinct neighborhoods, or niches, shaped by cell-to-cell interactions and extracellular components [3-5]. These spatial microenvironments influence crucial biological processes, from tumor progression to immune infiltration and tissue regeneration [6-8].
Spatial transcriptomics (ST) has transformed our ability to investigate these tissue architectures by providing single-cell resolution mapping of gene expression while preserving spatial context [9-12]. This technological breakthrough has enabled researchers to examine the molecular underpinnings of tissue organization with unprecedented detail. However, ST provides only static snapshots of inherently dynamic biological systems. Time-resolved spatial analysis extends ST by capturing how gene expression patterns and cellular arrangements evolve across developmental stages or experimental time [13-15]. This temporal dimension offers critical insights into the development of tissue organization health and disease [16] that cannot be inferred from static observations alone.
"Equal contribution
Project page: https: cs. cit. tum.
cit.tum. www.tu/..eeeeeeeeemmommmmmm
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

Despite these technological advances, modeling trajectories rime-resolved spatial slides is compli- cated, as no direct correspondence exists between cells across slides due to the destructive nature of acquisition practices. Moreover, current computational methods fall short in modeling the evolution of tissue organization at the level of cellular microenvironments. Most approaches infer trajectories by modeling single-cell dynamics using velocity- based models [17-19] or optimal transport between individual cells [20, 21]. While effective at capturing cell evolution, these cell-centric methods fundamentally miss the coordinated evolution of structured niches within tissues.
This limitation presents critical research gap that we address with the following question:
How can model the spatiotemporal evolution of cellular microenvironments preserving both local neighborhood relationships and cellular state transitions?
To address this question, we directly model the dynamics of cellular neighborhoods as cohesive units rather than focusing on iso- lated cell trajectories. This approach aligns naturally with tissue-scale biological pro- cesses and enables principled learning of dynamics in structured, -igh-iimensionl,, and variably sized spatial domains.
to
t
Noisy ME Source Target ME
M
D
D
M
We introduce Niche Flow Matching (NicheFlow) (Fig. 1), a generative model for learning spatiotemporal dynamics of cel- lular niches from time-resolved spatial tran- scriptomies data. NicheFlow builds on re- cent advances in Flow Matching (FM) and Optimal Transport (OT) to operate over dis- tributions of microenvironments, which we represent as point clouds. NicheFlow en- ables accurate modeling of global spatial ar- chitecture and local microenvironment com- position within evolving tissues.
M O.0)
M
Figure 1: Overview of NicheFlow. At time 11, we gen- erate target microenvironment M by transforming Gaussian noise .M using a Variational Flow Match- ing model with posterior μd conditioned microenvironment M at lo- Source-target pairs are identified via entropic OT over pooled microenviron- ment coordinates and gene expression profiles.
source
Our contributions include:
A microenvironment rewererrrrrmmettteneeeeeeedd trajectory inference paradigm that shifts from modeling individ- ual cells in time to modeling niches as point clouds, enabling simultaneous prediction of spatial coordinates and gene expression profiles while preserving local tissue context.
A factorized Variational Flow Matching (VFM) approach with distributional families (Laplace for spatial coordinates, Gaussian for gene expression) that jointly trains spatial and cell state dynamics using a factorized loss, modeling spatial reconstruction and biological fidelity with tailored distributional assumptions.
A spatially-awan sampling strategy using OT between niche representations, enabling scalable training on large tissue sections while ensuring comprehensive coverage of heterogeneous regions.
Our approach consistently outperforms baselines in recovering cell-type organization and spatial structure across embryonic, brain development, and ageing datasets. NicheFlow enables principled learning of dynamics in structured, high-dimensional, and variably sized spatial domains, a challenge with parallels in other spatiotemporal modeling domains beyond biology.
2 Related Work
NicheFlow at the interface between generative models and spatiotemporal transcriptomic data.
FM and single-cell transcriptomics. We propose a model based on FM, framework introduced by several seminal works [22-24]. Specifically, we adopt a variational view of the FM objective, following Fijkelboom et al. [25], but extend it to mixed-factoriaa generation. Our method, NicheFlow, combines FM with OT, a pairing that has proven effective in modeling cellular data [26-29]. Unlike these models, however, we focus point clouds of
distributions for point cloud

---

ressive transeriptomic profiles. Closest to our approach is Wasserstein FM for point cloud generation [30], applied reconstruct cellular niches. Yet, that work does not address joint generation of spatial coordinates and cellular states, nor OT-based temporal trajectory prediction, both of which are central to our contribution.
Generative models for spatial transcriptomics. Generative models have been key to spatial tasks such as gene expression prediction from histology slides [31, 32]. integration with dissociated single-cell data [32], spatial imputation [33, 34], and perturbation [35]. More recently, LUNA [36] demonstrated strong performance in predicting single-cell spatial coordinates using diffusion models [37] conditioned on transcription data. While related, our model addresses the distinct task of inferring niche trajectories, enabling the simultaneous generation of coordinates and cellular states.
Trajectory inference for spatial transcriptomics. Previous work has explored learning trajectories from spatial slides. Pham et al. [38] proposed graph-hased spatiotemporal algorithm for pseudotime inference, while others leveraged dissue-resolved transcriptomics to estimate cell velocity [17-19]. Closer to our approach, Klein al. [20] and Bryan al. [21] use discrete OT to link cells across time and infer the evolution of cell states from spatially-resolved gene expression. Similarly, DeST-OT [39] aligns spatial slides with semi-relaxed OT couplings, preserving transcriptomic and spatial proximity between ancestor and descendant cells, while SpaTraak [19] uses Fused balancing transcriptomic and spatial differences based on spatial autocorrelation of features. Unlike our bani-batch deep learning model, these methods do not operate on entire microenvironment and rely on exact OT at the single-cell level, resulting in limitations in scalability and generalization.
OT [40],
3 Background
3.1 Optimal Transport with FM
FM [23] is generative model that transforms source density Po into a target density p1. It operates by learning a time-dependent velocity field (x)) for t E [0, 1] and x E RR,,, which generates a probability path {pe)o_oo.). This path is constructed such that the marginals at time t = and 1 match the source and target distributions, i.e., Po and Pl, respectively. The velocity field induces an Ordinary Differential Equation (ODE), whose solution o (a+) defines a flow map that transports samples from the source to the target distribution.
In practice, FM approximates (x) with time-conditioned neural network (z)) While the exact marginal velocity field uy (x) is intractable, it can be expressed in terms of data-conditioned velocity fields and joint distribution ग ( 21) over the source and target samples vo po and æ1 Pr:
ve(x)= = f ut (x zo. x1) (r
Do. x (T 21) dxy dæ
(I)
(x(
where (x To, x1) is a pre-defined interpolating probability path. Here, we consider the tractable probability path pa æ1) = 6(x= 9o (3r,, 3)) where ge (30,, aa1) (1- + fæ, is linear interpolation and 8 denotes a Dirac delta function, representing a deterministic conditional path.
Lipman et al. [23] show that regressing the conditional field wz( æo, x1) is equivalent to learning the marginal field (x) in expectation. Hence, the FM objective becomes the task of learning the velocity along the conditional probability path between any pair of source and target data points. For oo x1) has a closed form, and the FM loss is:
linear conditional probability path, the velocity (x
2 I|.e (gt (2 x1)) - at
(0) = 12-2 (1)
It (oo. x1)
(2)
In practice, the coupling Tt (x=, x1) is instantiated using sample pairs drawn from mini-batch
estimate. When one chooses Tt" as the OT coupling under a squared Fuclidean cost between samples from Po and pi, FM approximates the dynamic OT map between source and target densities [41, 26). Thus, the solution samples (æo, Φ1) ) ~ x² from the joint distribution approximately follow:
Po: 3
8 [
(x0))
(3)
where of is the solution of the ODE with velocity field of
3

---

3.2 Generative OT on incomparable source and target spaces
Klein et al. [29] generalize the OT FM formulation to settings where the souree and target distributions are defined on incomparable spaces and propose an approach to generative entropic OT using FM. Given a standard normal noise distribution with samples ~N Ip). the authors show that the following sampling procedure:
(4) defines generative model that implicitly samples from an Entropic OT (EOT) coupling. where (2 Øo) is FM model that maps noise to target samples, conditioned on source points. To achieve this, &8 trained using source-target pairs (#g, drawn from the EOT coupling nee with c denoting the entropic regularization parameter, which the model aims to approximate.
Po. ~ 8 (21 1021-1((z
Φf(2 mo))
Crucially, this formulation enables OT between distinct source and target spaces, as æo does not flow directly into I₁ but instead conditions the generation of target samples from noise.
3.3 Source-conditioned VFM
Consider the conditoend FM formulation in Sec. 3.2. Given conditioning source æo and noise-based generation, the marginal field in Eq. (1) can be written as:
Ug u1 xo) Ep, (m).m) [ur(w: [x1)]
where we drop the conditioning on By in the velocity field, as u, is entirely determined by the target 21 when generating from noise under linear probability paths.
Since (a (e1) is tractable [23], one can recast FM as a variational inference problem, following Eijkelboom et al. [25], by introducing parameterized approximation qf (2x x,, 22) to the true posterior Pt (æ) x, Zn). Integrating the expected velocity in Eq. (5) over t E [o, 1] enables the
generation of target points Ij from noise, conditioned on Zo.
During training, the source-conditioned Variational Flow Matching (VFM) loss is:
LSC (0) = -E,-1(0011), (#0.2) a=ppo( [log & (at1 x zo)]
(6)
where x. is entropic OT coupling modeling the joint distribution over source and target samples, and pr xi) interpolates between target samples and noise. In the generation phase, one samples 20 po and noise <n(x) = z to generate a target sample 271.
(0,, then simulates the marginal field in Eq. (5) starting from
Crucially, under the assumption that (x w1) is linear in æ1, which holds when using straight- -ine interpolation paths, the marginal field in Eq. (5) only depends on the posterior's 's first moment on 21
(æ 1Ep, (ur (=, [x1])
This implies that the VFM objective reduces to matching the first moment of the approximate posterior df (x) 2;, Z.) to that of the true posterior pz( x, zn). As result, the approximate posterior can be chosen fully factorized under a mean field assumption, since each dimension can be matched independently if the mean of the posterior is preserved; see App. C.2 for more details.
4 NicheFlow
We introduce a flow-based generative OT model to infer the temporal evolution of spatially resolved cellular microenvironments. More specifically, given spatial microenvironment represented as point cloud of cell states with their coordinates, NicheFlow predicts the corresponding tissue structure later time point. To delineate our approach, we define list of desiderata.
Generative model on structured data. Similar to prior work [30, 36], we consider generative model over structured point cloud data representing cellular microenvironments. This approach implicitly accounts for spatial correlations between cells, in contrast to models that study the evolution of spatial trajectories at the single-cell level [20].
Sub-regions and variable location. Crucially, for better memory efficiency, we do not consider an entire spatial slide for trajectory inference, but instead learn the dynamics of variably located

---

sub-regions. This design choice enables scalability and flexibility in modeling functional regions across different parts of the screened tissue.
Changes the number of nodes. To model the temporal evolution and densification of microen- vironments, we allow the source and target regions to differ in the number of nodes. We adopt formulation similar to Sec. 3.2, implementing OT FM between non-comparable spaces.
Flexible generative models for features and coordinates. We allow flexibility in the choice of generative models for the features and coordinates. To this end, we implement the approach described in Sec. 3.3, factorizing features and coordinates into separate posteriors from different families.
4.1 Data description and problem statement
We are given sequence of time-resolved spatial transcriptomic measurements across biological processes such as development ageing. For simplicity, we formulate our problem in terms of two consecutive time points indexed by s, such that s E {0, 1). though the model can be extended to collections of more consecutive discrete temporal measurements, Each dataset at a time point is full tissue slide that can be represented as an attributed point cloud:
P. {('', ) i= ,.l..,
(8)
where of E 2 denotes the 2D spatial coordinate of cell time 8, and x! E RD denotes its associated feature vector, typically corresponding to its gene expression profile or a low-dimensional representation thereof. Thus, each dataset is attributed point set two-dimensional space. Note that each slide can have variable number of cells, and no direct correspondence exists between single cells across subsequent slides.
To capture spatial context beyond individual cells, we define local microenvinonments as fixed-radius neighborhoods. Specifically, for each cell (c xx)) time s, construct neighborhood M consisting of all neighboring cells within a spatial radius r:
= [(c] m3) ||e]
ctl sr]
Let {M No and {M! be collections of source and target microenvironments at consecutive time points. Our goal is to train parameterized flow model of with i E [0. 1]. that generates target microenvironments conditioned on source point clouds. Specifically, to sample a target microenvironment with l: cells conditioned on variably sized source M we define sampling as:
M {[c], 21) | of ~N(00, I2). z; Io(,I Vi ......k].
(10)
M = (M IMM)
(11)
where M is point cloud composed of noisy coordinates and features, and M ' is a generated prediction for the evolution of M at the next time point. As explained in Sec. 3.1 and Sec. 3.2, we want our generative model to parameterize some notion of optimal entropic coupling ". between microenvironmentt across slides (see Sec. 4.2).
4.2 OT formulation
To train the flow map of to perform conditional EOT, define cost function that induces an optimal entropic coupling न. between source and target point clouds. This coupling is used to sample pairs of source and target microenvironment mini-batches during training. While there is no established notion of optimal cost in this setting, we propose to compute OT using source and target microenvironment representations based on the weighted average of features and coordinates.
Specifically, we compute a pooled representation for each microenvironment in the source and target slides by averaging spatial coordinates and gene expression features, weighted by tunable hyperparameter E [0, that balances spatial versus cellular state information:
1 M
m
Σ
Σ
(12)
I
(c)a))))M;
where || denotes concatenation, m E 2+D and s E (0. We then apply EOT on the sets ( Jin and {m]) using squared Euclidean cost and regularization parameter yiclding the
5

---

coupling - During training. we sample matched pairs (M,,'), _ computed over mini- batches, and use them as supervision for learning the conditional generative model. higher value of prioritizes feature similarity, while lower values favor proximity coordinates (see Fig. 12).
4.3 Mixed-factorized VFM
Once we have established a strategy for performing mini-batch OT, we proceed to describe our approach for learning the flow model of and simulating Eq. (11). ). To this end, we adopt variant of VFM (Sec. 3.3), originally developed for graph generation, and adapt to our point cloud setting.
In line with Sec. 3.2 and 3.3, we delineate an objective to train a parameterized, source-conditione posterior over target point clouds qf (M' M,, M.). where M and M represent source and target niches, and MM denotes a noisy point cloud at interpolation time t. Importantly, our posterior comes with the following characteristics: (i) the posterior factorized across the single points in a point cloud; (ii) the posterior is factorized across cellular features and coordinate dimensions; and (iii) the family of posteriors can be chosen differently between cellular features and coordinates.
Following (i), we model the variational distribution over M by factorizing it across individual points (c1,21) .M.. Moreover, we tackle (ii) and (iif) using the mean-field VFM assumption, modeling cellular state and positions separately (see App. C.1 and C.2 for theoretical justifications):
(
D M I d=1
(M- MM, (EM
II
M,
(13)
Here, f and r denote distinct approximate posterior families for cellular states (x1) and spatial posi- tions (cf), respectively. We use a Laplace distribution for f. due its concentration around the mean. which supports precise modeling of coordinate features, while Gaussian distribution is used for rī.
As explained in Sec. 3.3, if df is required to simulate the generative field in Eq. (5). Therefore, the posterior is replaced by time-dependent predictor (ft ff) n (M, .M)) of the mean features Fl and coordinates f of the target point clouds, learned minimizing the following loss:
FM with straight probability paths, only the first moment of
CNicheFlow (0)
E t- (M..M² (M.M))>-m YouT
Σ (cL
(lei (la--F') fill1 - _11011 -r?13)
(14)
(MIMM))
We derive the objective in App. C.4 and provide algorithms in App. E. Here, oe inputs noisy point cloud M and a source M and provides mean prediction vector for coordinates and feature dimen- sions, respectively indicated as fft and F(. We implement it as a point cloud transformer (see Sec. 4.4).
We highlight crucial aspect about Eq. (14). While the single dimensions are fully factorized the predictions from n4, every feature's mean is function of the whole noisy point cloud M well as the target microenvironment. M In other words, the predictions exploit structural information the point cloud to predict the individual posterior mean of each dimension. Like most approaches, method has modeling limitations that we outline in App. B.
4.4 Backbone architecture: Microenvironmenn transformer
To parameterize the conditional posterior mean naد from Sec. 4.3, we use permutaaiion-equiiviinn transformer architecture designed for point clouds with variable size.
Encoder-decoder structure. The model follows an encoder-decoder layout, processing the source microenvironment M via the encoder and predicting the posterior mean from noisy target M ~ (- MM)) via the decoder conditioned on the encoder's output.
Input embeddings. Each point is represented by features R and spatial coordinates R embedded separately and concatenated. Time is encoded using sinusoidal embeddings, linearly projected and broadeast across points.
Cross-attention conditioning. The encoder processes the embedded source M microenvironment using self-attention. The decoder operates on the noisy target M using self-attention, followed by

---

Table 1: Performance comparison across three biological datasets for SPFlow, RPCFlow and NicheFlow. Models are trained using the Conditional Flow Matching (CFM). Gaussian VFM (GVFM), or Laupsian-Laplacian VFM (GLVFM) strategies. Results are reported as mean + standard deviation over five evaluation runs on mouse embryonic development (MED), axolotl brain develop- ment (ABD), and mouse brain aging (MBA). For all experiments, we use fixed value of 0.1, enabling spatial location preservation across time.
MED
ABD
MBA
Madel
Obj
INN-FI
PSD (109)
SPD (100)
INN-FI
PSD (10")
SPD (109)
INN-FI 1
PSD J (10')
SPD (10')
LUNA
2004
0.331 AwEs
0.222 LoRe
SPFlow
CFM
0,272
1681 +3081
0.602
0.190 c00
2.494 -1005
1,119
0.205 ecron
1886 1886-
0.834 -180
SPFlow
GVIM
0,259
1.383
0.582
0.175
3,375
1.104
0,181
2.585
0.834
GLVFM
0.175
RFCFloo
13304
1
RPCFlow
GVTM
0.503 2uans
1,155 18:14
0.578 2::000
0.477 1:000
2.260 :1:077
1,036 218251
0.249 21.0
1.753 20063)
0.784 c000
RFCHlow
GLVFM
0.596
40000
0.979
0.586
400012
0.554 800
2,053
2.0044
1,038 282222
0,265
6000
1,723 00015
0.779
0000
NicheFlow
CFM
0.609 s19es
0.979 23:22
0,402
0.604 1018
2.086 1:0e
0.568 8180
0.283
1100
1.557 =
0.556
NacheFlow
GVIM
0.596 40002
0.991
0,405
0.574
2010
2:220
2.800
0.594
2.8046
0.268
R.000
1.661
0.531 200000
NicheFlow
GLVPM
0,664 -1011
0,883
0,398 8.223
0.628 0013
2,079 -100
0.576
0.285 68066
1.554
1.554-8107
0.532 0.132-01030
cross-attention to condition on the encoder's outputs. The cross-attention mechanism allows each target point to attend to all source points.
Output projection. The decoder outputs are linearly projected to yield posterior mean estimates (fi', Ff) for expression and coordinates.
5 Experiments
We propose quantitative and qualitative evaluations of our algorithm. Quantitatively, we test whether source-conditooedd samples generated by model preserve the biological structure and shape of future tissue states. Qualitatively, we demonstrate that our approach accurately captures compositional shifts in substructural components and developmental trajectories across time.
5.1
Quantitative evaluation
Our first research question is assess the impact of two core modeling choices in NicheFlow: (i) learning trajectories over spatial microenvironments, rather than independently for each cell and (ii) restricting source and target point clouds to spatially co-localized neighborhoods of cells, instead of sampling them randomly across the slide.
We NicheFlow and baseline FM approaches that do not incorporate (i) and (ii) simulate spatial trajectories conditioned on early time point observations. Assuming that spatial arrangements and the biological composition at later stages evolve from earlier slides, the global correspondence between predicted and true slides indicates the quality of the generative trajectory.
5.1.1 Training setup
Datasets. We assess model performance across three spatiotemporal datasets: (i) Mouse embryo- genesis [20, 8] and (ii) the axoloti brain development [42], two Stereo-seq datasets profiling the spatially-resolved cellular development of a mouse embryo and axolotl brain across three (E9.5, 5,, E10.5 and E11.5) and five time points, respectively. We also consider the (iii) mouse brain ageing dataset [43], profiled with MERFISH [44] across twenty time points (see Apps. F.I and F.2).
Dataset construction. For each dataset and time point $ S, we construct a set of cellular microenvironments by applying the fixed-radius neighborhood definition introduced in Sec. 4.1. Each microenvironment Mi is centered cell i and contains all cells within fixed radius r. This results in:
M (M, i 1,...,
where N₂ is the number of cells in the tissue at time s, and each M, is attributed point cloud encoding both spatial and gene expression information. We standardize coordinates for cross time comparability and reduce the normalized gene expression its top 50 Principal Components (PC).
Batching. We train NicheFlow with mini-batches source and target cellular point clouds. To ensure spatial diversity during training, we sample individual batches uniformly from within discrete
7

---

SPF CFM
CVRM
Target (10.5)
SPFlow CFM
SPriow SVFM
GLVFM
MichFFlow CFM
NicheFlow SVFM
NicheFlew GLVFM
Qualitative comparison of generated samples the embryonic development dataset (9.5- 11.5 days). We show source and target samples alongside predictions from SPFlow and NicheFlow with different objectives.
Figure 2:
regions of the slides computed with K-MMens clustering over the 2D coordinates (see Fig. 13 for visualization with different K values). From these regions, we collect M source and target microen- viromments and resample N S M matching pairs from the entropic OT coupling (M .M' ~ T as described in Sec. 3.1, where M and M denote the sampled sets (see Sec. 4.2 for details on our OT coupling).
Evaluation data. For consistent and reproducible evaluation, we discretize each tissue into fixed 2D grid and define evaluation microenvironments fixed-radius neighborhoods around the nearest cells each grid point. This guarantees full spatial coverage and ensures deterministic comparison across methods. See App. F.5 for details.
Multiple time-point. NicheFlow predicts piecewise trajectories between subsequent time points. Instead of learning one flow for each couple of subsequent slides, we train single model with additional conditioning on source and target labels (see App. F.6).
5.1.2 Quantitative evaluation metrics
Spatial structure. We quantify coordinate generation accuracy using two asymmetric distance metrics. The point-to-shape distance (PSD) measures how far predicted coordinates deviate from the true structure by computing the mean squared distance from each generated point to its nearest ground truth counterpart. In contrast, the shape-lo-point distance (SPD) evaluates how well the generated points cover the target region by averaging the squared distance from each ground truth point to its nearest generated point (see App. F.3 for a mathematical formulation of the metries).
Cell-type organization. To assess how well the model reconstructs the spatial organization of different cell types, we use a t-nearest-neighhor ( (INN) classification setup. Since the model generates only gene expression profiles and spatial coordinates, we assign cell type labels to generated cells using a classifier trained on ground truth gene expression data (see App. F.4). Each predicted cell is then matched to its nearest real cell, and report the weighted F1 score (INN-F1).
5.1.3 Models and results
Baselines. We compare against what we call SPFlow (Single-Point Flow), standard FM-based model that predicts temporal trajectories across slides at a single-cell level using an MI.P-based velocity field. We also consider RPCFlow (Random Point Cloud Flow), which has the same backbone as NicheFlow, but conditions on randomly sampled point clouds instead radius neighborhoods. Additionally, we include LUNA [36], diffusion model for spatial reconstruction from dissociated cells. Note that LUNA does not model temporal dynamics and only generates coordinates from noise with their respective biological annotations. Therefore, we use a reference for spatial generation accuracy via the INN-F1 metric rather than a proper baseline.
Ablations. We assess different training objectives by comparing standard Conditional Flow Matching (CFM) [26] with two variational formulations modeling posteriors over coordinates and features: Gaussian-only (GVFM) and Coussian-Laplace (GLVFM). The former uses Gaussian posteriors for coordinates and features. The latter uses the factorized formulation in Sec. 4.3. For the point-cloud- based methods,
fixed value of À 0.1 and sampled batches of 64 regions chosen from

---

Mapping spinal cord E10.5-E11.5
Mapping neural crest E9.5-E10.5
Source
NicheFlow mapping probability 0.1
Source
NicheFlow mapping probability A 0.5
Spinal cord
*euural ORST
Expected target
moscot mapping probability
Expected target
moscot mapping probability
Offtangtt
0.9
Ofi-targtt
Thoud manaecchytte
60
Figure 3: Left and right panels show the mapping of spinal cord (E10.5 to E11.5) and head neural crest cells (E9.5 to E10.5). In each panel, the left column shows source cells and expected targets, and the right column shows density contours of the most likely mapped regions Bar plots display transition probabilities to the most likely descendant cell types. For NicheFlow, contours represent the proportion of samples in generated point clouds assigned to real cell coordinates across 10 samples.
the ablations in App. D.5 and D.6. We present results with a single value for À as this experiment compares different models on reconstructing fixed tissue structures over time, for which we enforce spatial preservation (see Sec. 4.2). However, we also include model comparisons across multiple values of in Tab. 3 for completion.
Results. Our quantitative evaluation (Tab. 1) demonstrates that NicheFlow trained with the GLVFM objective consistently achieves strong performance across both spatial and semantic metrics. It outperforms all baselines in reconstructing spatial structure (PSD, SPD) and cell-type organization (INN-F1) on developmental datasets, while remaining competitive on ageing data. These results highlight the importance of structured microenvironment modeling for capturing the spatiotemporal dynamics of complex tissues. We complement our quantitative results visually in Fig. 2 and Fig. and 10 in the Appendix, where we show that SPFlow fails to capture tissue-level organization, producing blurry and spatially incoherent samples. contrast, NicheFlow generates predictions that preserve spatial structure and cell-type organization, despite learning only from local microenvironments. In App. D.4 and App. D.7, we additionally demonstrate that NicheFFlow produces conditionally consistent outputs with the source, while RPCFlow generates very diffused mappings across the slide, failing to preserve spatial consistency across time.
5.2 Qualitative evaluation and biological analysis
We explore the capabilities of NicheFFlow on the spatial trajectory inference task through qualita- tive and biological assessments. Specifically, we focus on validating whether our model captures compositional changes within fixed spatial structures and developmental trajectories.
Experimental setup. We train the model as described in Sec. 5.1.1. Using the mouse embryonic dataset [8], we select specific microenvironments as source niches for which we want to study the trajectory over time. Specifically, we propose two scenarios for the application of NicheFlow depending on the choice of the OT parameter À (see Sec. 4.2):
Compositional changes in fixed structures across time. We choose the evolution of the spinal cord of the embryo from E10.5 to 115. as example and set low value of À to prioritize the preservation of the spatial location in the trajectory.
2. Spatial and cellular development of immature cells. Developing cells may displace to different areas of the embryo, requiring a higher value we inspect how neural crest cells in the head evolve into mesenchymal and cranial structures.
account for gene expression.
case study,
Baseline. We compare NicheFlow with the spatiotemporal framework moscot [20], which models spatial trajectories at the single-cell level. In contrast to NicheFlow, moscot integrates spatial

---

coordinates directly into the OT formulation using a Fused hyperparameter a controls the trade-off between spatial and feature- distances (see App. F.8). A high a increases the influence of spatial distances, whereas in our framework this role is played by the hyperparameter À (see Sec. 4.2). We provide more details on the selection of the hyperparameter for our experiments in App. D.2. Notably, moscot relies on exact OT and learns a transition matrix between source and target samples. As such, it is not a generative model over point clouds like NicheFlow. However, given the overlap in downstream tasks, we consider the comparison relevant.
cost [45], where the
Evaluation and results. For both scenarios (1) and (2), we have prior knowledge of the ground truth regions that the source microenvironments are expected to occupy at later developmental stages, as well as their corresponding biological compositions. For both moscot and NicheFlow, we assess whether the transported mass of source samples concentrates within the conrec anatomical region at the target time point, and whether the predicted descendant cell types are biologically consistent (see App. F.8). Results are summarized in Fig. 3. When modeling the evolution of the spinal cord, moscot assigns considerable mass to unrelated regions such as the urogenital ridge and branchial arches, whereas NicheFlow correctly maps source niches to the maturing spinal cord. Similarly, NicheFlow captures the differentiation of neural crest cells into mesenchymal and cranial tissues within the head region, while moscot exhibits substantial off-target leakage towards lower regions.
Conclusions and Discussion
We introduce NicheFlow, a point-cloud-based generative model designed to capture the spatiotem- poral dynamics of cellular niches in time-resolved spatial transcriptomics data. Unlike methods that model single-cell trajectories independently, NicheFlow implicitly captures spatial correlations between cells by learning trajectories on variably sized local neighborhoods. To this end, we combine OT with a new version of VFM that factorizes features and coordinates into distinct posteriors from different distribution families. We showed that NicheFlow outperforms standard FM approaches at reconstructing spatial context from previous time points and improves the mapping of biological structures in time over established exact OT approaches. With the expected increase in the volume and quality of spatial data, modeling coordinated cellular state translations with generative models is promising avenue for generalization beyond spatiotemporal inference. We envision that models like NicheFlow will enable spatial perturbation prediction and modality translation tasks requiring principled parameterized maps beyond discrete OT to extrapolate and drive biological hypotheses.
Acknowledgments
A.P. is supported by the Helmholtz Association under the joint research school Munich School for Data Science (MUDS). A.P., F.G., S.G. and F.J.T. also acknowledge support from the German Federal Ministry of Education and Research (BMFTR) through grant numbers L01292 and LO0299 F.J.T. acknowledges support from the Helmholtz Association's 's Initiative and Networking Fund via the CausalCellDynamies project (grant number Interlabs-0029 and the European Union (ERC, DeepCell, grant number 101054957). The authors of this work take full responsibility for its content.
References
[1] Urs Mayr, Denise Serra, Prisca Liberali, Allon Klein, and Barbara Treutlein. Exploring single cells in space and time during tissue development, homeostasis and regeneration. Development, 146(12), 2019.
[2] Marietta Zinner, Ilya Lukonin, and Prisca Liberali. Design principles of tissue organisation: How single cells coordinate across scales. Current opinion in cell biology, 67:37-45, 2020.
[3] Yanming Ren, Zongyao Huang, Lingling Zhou, Peng Xiao, Junwei Song, Ping He, Chuanxing Xie. Ran Zhou, Menghan Li, Xiangqun Dong, Qing Mao, Chao You. Jianguo Xu, Yanhui Liu, Zhigang Lan, Tiejun Zhang, Qi Gan, Yuan Yang, Tengyun Chen, Bowen Huang, Xiang Yang, Anqi Xiao, Yun Ou, Zhengzheng Su, Lu Chen, Yan Zhang, Yan Ju, Yuekang Zhang, and Yuan Wang. Spatial transcriptonies reveals niche-specific enrichment and vulnerabilities of radial glial tike cells in malignant gliomas. Nature Communications, 14(1), February 2023.
10

---

[4] Anna Christina Schaar, Alejandro Tejuda-Lapuerta, Giovanni Palla, Robert Gutgesell, Lennard Halle, Mariia Minaeva, Larsen Vornholz, Leander Dony, Francesca Drummer, Mojtaba Bahrami, and Fabian J Theis. Nicheformer: foundation model for single-cell and spatial omics. bioRxiv, pages 2024-04, April 2024.
[5] Song Liu, Xiaoyan Li, Zhiyue Gu, Jiayu Wu, Shuangzheng Jia, Jinghua Shi, Yi Dai, Yushi Wu, Hailan Yan, Jing Zhang, Yan You, Xiaowei Xue, Lulu Liu, Jinghe Lang, Xiaoyue Wang, and Jinhua Leng. Single-cell and spatial transcriptomic profiling revealed niche interactions sustaining growth of endometriotic lesions. Cell Genomics, 5(1):100737, January 2025.
[6] Reuben Moncada, Dalia Barkley, Florian Wagner, Marta Chiodin, Joseph C. Devlin, Maayan Baron, Cristina H. Hajdu, Diane M. Simeone, and Itai Yanai. Integrating microarray- -aase spatial transeriptomics and single-cell RNA-seq reveals tissue architecture pancreatic ductal adenocarcinomas, Nature Biotechnology, 38(3):333-342, January 2020.
[7] Ludvig Larsson, Jonas Frisén, and Joakim Lundeberg. Spatially resolved transcriptomics adds new dimension to genomics. Nature Methods, 18(1):15-18, January 2021.
[8] Ao Chen, Sha Liao, Mengnan Cheng, Kailong Ma, Liang Wu, Yiwei Lai, Xiaojie Qiu, Jin Yang. Jiangshan Xu, Shijie Hao, Xin Wang. Huifang Lu, Xi Chen, Xing Liu, Xin Huang. Zhao Li, Yan Hong, Yujia Jiang, Jian Peng, Shuai Liu, Mengzhe Shen, Chuanyu Liu, Quanshui Li, Yue Yuan, Xiaoyu Wei, Huiwen Zheng, Weimin Feng. Zhifeng Wang. Yang Liu, Zhaohui Wang, Yunzhi Yang. Haitao Xiang, Lei Han, Baoming Qin, Pengcheng Guo, Guangyao Lai, Pura Muñioz-Cánoves, Patrick H. Maxwell, Jean Paul Thiery, Qing-Feng Wu, Fuxiang Zhao, Bichao Chen, Mei Li, Xi Dai, Shuai Wang, Haoyan Kuang, Junhou Hui, Liqun Wang, Ji- Feng Fei, Ou Wang. Xiaofeng Wei, Haorong Lu, Bo Wang, Shiping Liu, Ying Gu, Ming Ni, Wenwei Zhang, Feng Mu, Ye Yin, Huanming Yang. Michael Lisby, Richard J. Coraall, Jan Mulder, Mathias Uhlén, Miguel A. Esteban, Yuxiang Li, Longqi Liu, Xun Xu, and Jian Wang. Spatiotemporal transcriptomic atlas of mouse organogenesis using DNA pantball-patterned arrays. Cell, 18(110)1177-1192222) May 2022.
[9] Patrik L. Ståhl, Fredrik Salmén, Sanja Vickovic, Anna Lundmark, José Fernández Navarro, Jens Magnusson, Stefania Giacomello, Michaela Asp, Jakub Westholm, Mikael Huss, Annelie Mollbrink, Sten Linnarsson, Simone Codeluppi, Åke Borg, Fredrik Pontén, Paul Igor Costea, Pelin Sahlén, Jan Mulder, Olaf Bergmann, Joakim Lundeberg, and Jonas Frisén. Visualization and analysis of gene expression in tissue sections by spatial transcriptomics. Science, 353(4)))): 78 82, July 2016.
[10] Samuel G. Rodriques, Robert R. Stickels, Aleksandrina Goeva, Carly A. Martin, Evan Murray, Charles R. Vanderburg, Joshua Welch, Linlin M. Chen, Fei Chen, and Evan Z. Macosko. Slide- seq: A scalable technology for measuring genome-wide expression at high spatial resolution. Science, 363(6434): March 2019.
(11) Chee-Huat Linus Eng, Michael Lawson, Qian Zhu, Ruben Dries, Noushin Koulena, Yodai Takei, Jina Yun, Christopher Cronin, Christoph Karp, Guo-Cheng Yuan, and Long Cai. Transcriptome scale super-resolved imaging in tissues by RNA seqFISH+. Nature, 568(7751):235-239, March 2019.
[12] Lambda Moses and Lior Pachter. Museum of spatial transcriptomics. Nature Methods, 19(5): 534-546, March 2022.
[13] James A. Briggs, Caleb Weinreb, Daniel F. Wagner, Sean Megason, Leonid Peshkin, Mare W. Kirschner, and Allon M. Klein. The dynamics of gene expression in vertebrate embryogenesis at single-cell resolution. Science, 360(6392), June 2018
[14] Daniel F. Wagner, Caleb Weinreb, Zach M. Collins, James A. Briggs, Sean G. Megason, and Allon M. Klein. Single-cell mapping of gene expression landscapes and lineage in the zebrafish embryo. Science, 360(6392), June 2018.
[15] Jingyi Ren, Haowen Zhou, Hu Zeng, Connie Kangni Wang, Jiahao Huang, Xiaojie Qiu, Xin Sui, Qiang Li, Xunwei Wu, Zuwan Lin, Jennifer A. Lo, Kamal Maher, Yichun He, Xin Tang, Judson Lam, Hongyu Chen, Brian Li, David E. Fisher, Jia Liu, and Xiao Wang. Spatiotemporally
11

---

resolved transeriptomics reveals the subcellular RNA kinetic landscape. Nature Methods, 20(5): 695-705, April 2023.
[16] Michael Papanicolaou. Amelia L Parker, Michelle Yam, Elysse C Filipe, Sunny Z Wu, Jessica L Chitty, Kaitlin Wyllie, Emmi Tran, Ellie Mok, Audrey Nadalini, et al. Temporal profiling of the breast tumour microenvironment reveals collagen xii driver of metastasis. Nature communications, 13)477887 2022.
[17] Tamim Abdelaal, Laurens M Grossouw, Jeroen Pasterkamp. Boudewijn PF Lelieveldt, Marcel JT Reinders, and Ahmed Mahfouz. Sirv: Spatial inference of rna velocity at the single-cell resolution. NAR genomics and bioinformatics, 6(3):lqae100, 2024.
[18] Wenxin Long, Tianyu Liu, Lingzhou Xue, and Hongyu Zhao. spvelo: Rna velocity inference for multi-batch spatial transcriptomics data. bioRxiv, pages 2025-03, 2025.
[19] Xunan Shen, Lulu Zuo, Zhongfei Ye, Zhongyang Yuan, Ke Huang, Zeyu Li, Qichao Yu, Xuanxuan Zou, Xiaoyu Wei, Ping Xu, et al. Inferring cell trajectories of spatial transcriptomics via optimal transport analysis. Cell Systems, 16(2). 2025.
[20] Dominik Klein, Giovanni Palla, Marius Lange, Michal Klein, Zoe Piran, Manuel Gander, Laetitia Peng-Papaxanthos, Michael Sterr, Lama Saber, Changying Jing, et al. Mapping cells through time and space with moscot. Nuture, pages 1-11, 2025.
[21] John Peterson Bryan, Samouil L Farhi, and Brian Cleary. Accurate trajectory inference in time-series spatial transeriptomics with sructurally-constrrineed optimal transport. bioRxiv, pages 2025-03, 03, 2025.
[22] Xingchao Liu, Chengyue Gong, and Qiang Liu. Flow Straight and Fast: Learning Generate and Transfer Data with Rectified Flow. In The Eleventh International Conference on Learning Representations, 2023.
[23] Yaron Lipman, Ricky T. Q. Chen, Heli Ben-Hamu, Maximilian Nickel, and Matthew Le. Flow matching for generative modeling. In The Eleventh International Conference on Learning Representations, 2023.
[24] Michael Samuel Albergo and Eric Vanden-Eijnden. Building normalizing flows with stochastic interpolants. In The Eleventh International Conference on Learning Representations, 2023.
[25] Floor Eijkelboom, Grigory Bartosh, Christian A. Naesseth, Max Welling, and Jan-Willem van de Meent. Variational flow matching for graph generation. In The Thirty-eighth Annual Conference Neural Information Processing Systems, 2024.
[26] Alexander Tong, Kilian FATRAS, Nikolay Malkin, Guillaume Huguet, Yanlei Zhang, Jarrid Rector-Brooks, Guy Wolf, and Yoshua Bengio. Improving and generalizing flow-based genera- tive models with minibatch optimal transport. Transactions on Machine Learning Research, 2024.
[27] Alexander Y. Tong, Nikolay Malkin, Kilian Fatras, Lazar Atanackovic, Yanlei Zhang, Guillaume Huguet, Guy Wolf, and Yoshua Bengio. Simulation-free Schrödinger bridges via score and flow matching. In Proceedings of The 27th International Conference Artificial Intelligence and Statistics, volume 238 of Proceedings of Machine Learning Research, pages 1279-1287. PMLR, 02-04 May 2024.
[28] Lazar Atanackovic, Xi Zhang, Brandon Amos, Mathieu Blanchette, Leo J Lee, Yoshua Bengio, Alexander Tong, and Kirill Neklyudov. Meta flow matching: Integrating vector fields on the wasserstein manifold. In The Thirteenth International Conference on Learning Representations, 2025.
[29] Dominik Klein, Théo Uscidda, Fabian Theis, and Marco Cuturi. Genot: Entropic (gromov) wasserstein flow matching with applications to single-cell genomics. Advances in Neural Information Processing Systems, 37: 37:103977-103944, 2024.
12

---

[30] Doron Haviv, Aram-Alex Pooladian, Dana Pe'er, and Brandon Amos. Wasserstein flow matching: Generative modeling over families of distributions. arXiv preprint arXiv:2411.0069,, 2024.
[31] Sichen Zhu, Yuchen Zhu, Molei Tao, and Peng Qiu. Diffusion generative modeling for spatially resolved gene expression inference from histology images. In The Thirteenth International Conference Learning Representations, 2025.
[32] Xiaomeng Wan, Jiashun Xiao, Sindy Sing Ting Tam, Mingxuan Cai, Ryohichi Sugimura, Yang Wang, Xiang Wan, Zhixiang Lin, Angela Ruohao Wu, and Can Yang. Integrating spatial and single-cell transcriptomies data using deep generative models with spatialscope. Nature Communications, 14(1):7848, 2023.
[33] Doron Haviv, Ján Remšfk, Mohamed Gatie, Catherine Snopkowski, Meril Takizawa, Nathan Pereira, John Bashkin, Stevan Jovanovich, Tal Nawy, Ronan Chaligne, et The covariance environment defines cellular niches for spatial inference. Nature Biotechnology, 43(2):269-288) 2025.
[34] Kongming Li, Jiahao Li, Yuhao Tao, and Fei Wang. stdiff: a diffusion model for imputing spatial transcriptomics through single-cell transcriptomics. Briefings in Bioinformatics, 25(3): bbae171, 2024.
[35] Amir Akbarnejad, Lloyd Steele, Daniyal Jafree, Sebastian Birk, Marta Rosa Sallese, Koen Rademaker, Adam Boxall, Benjamin Rumney, Catherine Tudor, Minal Patel, et al. Mapping and reprogramming human tissue microenvironments with mintflow. bioRxiv, pages 2025-06, 2025.
[36] Tingyang Yu, Chanakya Ekbote, Nikita Morozov, Jiashuo Fan, Pascal Frossard, Stéphane d' Ascooli and Maria Brbić. Tissue reassembly with generative ai. bioRxiv, pages 2025-02. 02, 2025.
[37] Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. Advances in neural information processing systems, 33:6840-6851,, 2020.
[38] Duy Pham, Xiao Tan, Brad Balderson, Jun Xu, Laura F Grice, Sohye Yoon, Emily F Willis, Minh Tran, Pui Yeng Lam, Arti Raghubar, et al. Robust mapping of spatiotemporal trajectories and cell-cell interactions in healthy and diseased tissues. Nature communications, 14(1):7739, 2023.
[39] Peter Halmos, Xinhao Liu, Julian Gold, Feng Chen, Li Ding, and Benjamin Raphael. Dest-ot: Alignment of spatiotemporal transcriptomies data. Cell Systems, 16(2), 2025.
[40] Titouan Vayer, Lactitia Chapel, Rémi Flamary, Romain Tavenard, and Nicolas Courty. Fused gromov-wasserstein distance for structured objects. Algorithms, 13(9):212, 2020.
[41] Aram-Alexandre Pooladian, Heli Ben-Hamu, Carles Domingo-Enrich, Brandon Amos, Yaron Lipman, and Ricky Q. Chen. Multisample flow matching: Straightening flows with minibatch couplings. In Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 28100-28122 PMLR, 23-29 Jul 2023.
[42] Xiaoyu Wei, Sulei Fu, Hanbo Li, Yang Liu, Shuai Wang, Weimin Feng, Yunzhi Yang, Xiawei Liu, Yan-Yur Yun Zeng, Mengnan Cheng, Yiwei Lai, Xiaojie Qiu, Liang Wu, Nannan Zhang, Yujia Jiang, Jiangshan Xu, Xiaoshan Su, Cheng Peng. Lei Han, Wilson Pak-Kin Lou, Chuanyu Liu, Yue Yuan, Kailong Ma, Tao Yang. Xiangyu Pan, Shang Gao, Ao Chen, Miguel A. Esteban, Huanming Yang, Jian Wang, Guangyi Fan. Longqi Liu. Liang Chen, Xun Xu, Ji-Feng Fei, and Ying Gu. Single-cell stereo-seq reveals induced progenitor cells involved in axolotl brain regeneration. Science. :acb6440 2022.
[43] Eric D. Sun, Olivia Zhou, Max Hauptschein, Nimrod Rappoport, Lucy Xu, Paloma Navarro Negredo, Ling Liu, Thomas A. Rando, James Zou, and Anne Brunet. Spatial transcrip- tomic clocks reveal cell proximity effects in brain ageing. Nature, 638(8049):160-171 Feb 2025.
13

---

[44] Kok Hao Chen, Alistair N Boettiger, Jeffrey R Moffitt, Siyuan Wang, and Xiaowei Zhuang. Spatially resolved, highly multiplexed ma profiling in single cells. Science, 348(6233)aaux6090 2015.
[45] Titouan Vayer, Nicolas Courty, Romain Tavenard, and Rémi Flamary. Optimal Transport for Structured Data with Application on Graphs. In Proceedings of the 36th International Conference Machine Learning, volume 97 of Proceedings of Machine Learning Research, pages 2284. PMLR, 2019.
[46] Haoqiang Fan, Hao Su, and Leonidas Guibas, A point set generation network for 3d object reconstruction from single image. In Computer Vision and Puttern Recognition, 2016.
[47] Facundo Mémoli. Grosevewasserstein distances and the metric approach object matching. Foundations of Computational Mathematics, 1 1(4))417-487, Aug 2011.
[48] Alexander Tong, Jessie Huang, Guy Wolf, David Van Dijk, and Smita Krishnaswamy. Trajecto- ryNet: A dynamic optimal transport network for modeling cellular dynamies. In Hal Daumé IIII and Aarti Singh, editors, Proceedings of the 37th International Conference Machine Learning, volume 119 of Proceedings of Machine Learning Research, pages 9526-9536. PMLR, 2020.
[49] Guillaume Huguet, Daniel Sumner Magruder, Alexander Tong, Oluwadamilola Fasina, Manik Kuchroo, Guy Wolf, and Smita Krishnaswamy. Manifold interpolating optimal-transpon flows for trajectory inference. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho, editors, Advances in Neural Information Processing Systems, 2022.
[50] Zhenyi Zhang, Tiejun Li, and Peijie Zhou Learning stochastic dynamics from snapshots through regularized unbalanced optimal transport. In The Thirteenth International Conference on Learning Representations, 2025.
14

---

A Broader impacts
This work addresses fundamental challenges in spatial transcriptomics by modeling complex spatial and compositional changes in developing tissues. We demonstrate how efficient representations of nigh-dimensional spatial cellular data can advance the understanding of developmental trajectories and microenvironment dynamics. We anticipate releasing NicheFlow as an open-source, user-friendl, tool to enable broad application in spatial biology studies. Given its use with biological data, NicheFlow may also be applied in sensitive contexts involving clinical or patient information.
B Limitations
Our approach relies on fixed OT feature weighting by a parameter during training (see Sec. 4.2), limiting flexibility at inference and potentially constraining certain biological analyses. Moreover, radius-based niche definitions may also be sub-optimal for small or irregularly shaped microenviron- ments, where the radius captures excessive spatial context and does not allow fine-grained modeling of the functional region's 's evolution.
The model assumes that spatial slides can be aligned with respect to each other in time and requires normalization-bait
pre-processing. Future work will be directed towards rotational and translational invariant spatial constraints and the incorporation of cell-to-cell communication priors in the neigh- borhood definition. While the learned flow models cell population dynamics, does not explicitly capture biological events such division or death.
Finally, this study focuses in-distribution testing and does not consider generalization unseen slides or full anatomical regions excluded from the training process. To achieve prediction in unseen settings, we foresee the need for technical replicates of the same slide across time points. as the model cannot extrapolate spatial arrangements without prior exposure to the associated region. We leave these analyses future work when spatio-temporal measurements across multiple replicates become increasingly available.
C Mixed-zactorized Variational Flow Matching
C.1 Theoretical aspects of Variational Flow Matching
Variational Flow Matching (VFM) [25] relies on the observation that one can write the time-resolved marginal vector field (x)) in FM as the expected conditional field uz(a: æ1) under the posterior (r x) as:
ue (x(x) = Ep.(#: [u:( xx)]]
(15)
Since ue(2 wi) has a closed form and u, (x) is all that we need to generate the probability path Pt from noise to data, this opens the door to new interpretation of the objective as variational inference problem, where we approximate Pr (x x) with a variational posterior of (x1 x). In
other words, one can optimize the following objective:
(V) - 1,,, (mi) (=)))))) [log qi (æ) [)]
where (a)) is the data distribution and (= mi) straight probability path. When (x ax) is linear in Ø1, this model formulation acquires convenient properties listed below.
Mean parameterization. The expected conditional field under the posterior only depends on the posterior mean:
E ppp(ss (a)) (x x1)] = u (ms [01]).
suggesting that it is sufficient to parameterize the posterior mean simulate data under the marginal flow. The posterior mean can be regressed against real samples æ during training.
Equivalence between posterior and approximate posterior formulation. From the previous point, follows that the expectation of the conditional field the same under the true and approximate posterior, long their first moments match.
Efficient simulation. Given a parameterized posterior mean n. simulating the generative field in Eq. (15) is efficient under the linearity condition. For example, in the standard FM setting with
15

---

straight paths [23], the marginal generative field becomes:
uz ( of (z) (x) (ur(er ax))
(16)
Ut (or(l)))))) (x1)
(17)
/ (x) a
(18)
1
which can be easily simulated in the range E [0, 1].
C.2 Factorized posterior
Similar to Eijkelboom et al. [25], in our work, we use a fully factorized posterior, where individual dimensions can follow different families of distributions with finite moments (see Sec. 4.3), Notably, a factorized approximate posterior over Φ is allowed as a choice for df since the only requirement to simulate u, (æ) is for of (=1 x)) to match the expectation pt (z1 x)) over æ1, irrespectively of higher moments correlations between factors.
In this regard, it is useful to consider the following proposition.
Proposition 1. Let 21 E RDD be D-dimensional larget data point, pe( x)) the posterior probubility path conditioned on a noisy point x ~ pe(x), and ua(a xx) the conditional velocity field. Assume that ux(x) x1) is linear in Φ1. Then, for any dimension d {I...., D}, the following holds:
E (. [xr]] = (e
(19)
Uy ux(x ue (xd (e
(20)
where I refers to the t scalar dimension of the vector Proof. We begin by proving Eq. (19) using marginalization: pe (x) [x-]]
f si pe(ar ax)) day
= f x4 (Dalter a) ma)"" dz
f xt pe(aft x)) dat.
(21)
Next, we prove Eq. (20). Under the assumption that the conditional velocity field ur(x xx)) is linear in æ1, we have:
Uy (x
Eyg (m])
Uy (x.
x
En p. (ang (mn)
Uc (
x
xd xd
(a (mn)
[xr]] xdd
(21 (e) (rg xd
1-t
We (1da (2)
(22)
Here, step (1) follows from the linearity assumption, which ensures that the conditional velocity at rd depends only on x.f-
In other words, the expected value under the posterior at an individual feature d does not depend on the other features and has an influence only on the d-th dimension of the conditional vector field. This flexibility allows each dimension's 's approximate posterior to be chosen from a potentially different distributional family, long the first moment exists and is correctly parameterized.
16

---

C.3 Marginal field derivation in source-conditioned VFM
When applying source conditioning to VFM, the marginal conditional vector field given souree zp is:
Pa(
To)
(x xo)
ue( x1)
dæ
(23)
Pr
where the pa x1) is probability path interpolating observations æı with noise. Note that omit To from the probability path and conditional velocity as they are fully determined by æ1 under linear conditional probability paths. Furthermore, can rewrite the marginal as expectation:
Pa( x1) (( æo) dæ (e, [ur at1)
We((r ai))
(24)
(a m0))
where we used that pr( (x) (x | æo, 21).
C.4
Gaussian and Laplacian Hybrid VFM
We define hybrid Variational Flow Matching (VFM) model using fully factorized variational distribution over individual points in the target microenvironment M Following the mean-field assumption, the variational distribution factorizes over spatial and feature dimensions:
(M' |,
II (oL.))EM.
(ct, M,,
(25)
( M.M') 1 ( d=1
II (CL.
(26)
Akml
In line with Eijkelboom et al. [25], using FM with straight paths enables us to efficiently simulate the marginal generating field using the first moment of the posterior distribution. In other words, for fully factorized posterior, we only need to parameterize mean predictor. In our setting, the mean prediction is neural network ie time-condition function of noisy microenvironment M and source M with outputs:
= ( (M,
where j and P!dd are the expected values for the hthh coordinate and dhh cell feature. Then, choose a parameterization for the variational factors at time & E [0, 1] as follows:
N(r 1),
(27)
Laplace (Je 1),
(28)
Substituting into the negative log- likelihood yields:
log(g((M' M,
(29)
log
II
II k=
(f M,M") II r I
(30)
Σ
log 2 lct F2.t1) +
Σ
log(2n)
1 2 (2 - =Q. )²
(EM
1=1
d=1
Σ
[|c ||a- fi'lt l1 +
1 - rfll3]
+ const w.r.t. 0
(31)
EEM
This results in loss consisting of an é error spatial coordinates and a mean squared error on gene expression features, consistent with the hybrid variational design.
17

---

D
Additional results
D.1 Additional comparisons with moscot on embryonic development
We propose a similar analysis as presented in Sec. 5.2.
Mapping anterior neural crest cells E9.5 E10.5 E11.5
Source
Predicted intermediate
Expected target
E9.5
E10.5
E11.5
Toural creat
Neural crest
Joe und tooch
Cavily
push
push
NicheFlow mapping probability E10.5
moscot mapping probability E10.5
À 0.5
0.5
0.0
Brain
Guily
Serfice ecosder
0.0
Cavelyy
Surface extudion
NicheFlow Flow mapping probability E11.5
moscot mapping probability E11.5
0.5
0.5
0.
Seclace ectodemm
Coots
0.0
Cavity
and tosth
Surtace ectaderm
Figure 4: Comparison of NicheFlow and moscot on the prediction of the anterior neural crest cells' fate. For both models, we take source facial neural crest cells at E9.5, push them to time point E10.5, and show the compositional and density predictions in the middle panel. Then, the predictions 10.5 are used as a source for second trajectory prediction operation from 10.5 to 11.5, for which we inspect again the cell density over the target slide and the cell type probabilities.
In Fig. 4, compare the ability of NicheFlow to predict an entire spatial structure trajectory by pushing an initial source cloud through all the developmental stages. To this end, the trajectory of an initial point cloud is first pushed the next time point, and the model's prediction is used as source for predicting the subsequent time point. An accurate niche trajectory reconstruction signifies that our model can be used to sequentially predict microenvironment evolution by treating its intermediate predictions as inputs, corroborating their accurate reflection of real point clouds.
In Fig. 4, we show that pushing anterior neural crest cells twice from E9.5 to 1.5. through the flow generates realistic target point clouds with a cell composition reflecting the expected cranial structure, mostly made of cavity cells, jaw and teeth (arising E1 1.5 for the first time) and surface ectoderm. Doing the same with moscot oversamples regions outside the cranial structure, thereby incorrectly mapping most of the neural crest density to brain cells.
In Fig. 5, we also show that NicheFlow is more accurate than moscot transporting mass from defined organs like the liver across development. More specifically, while density leaks from the liver to the tract in the mapping produced by moscot, the prediction computed by NicheFlow more accurately retrieves the liver structure at the later time point. Together with previous evidence, our results underscore the importance of accounting for spatial correlations between cells during OT-based trajectory inference to buffer out the noise resulting from single-cell-based predictions.
18

---

Mapping liver E10.5-E11.5
Source
NicheFlow mapping probability
= 0.5
0.4
Livr
0.2
Connective Hone
Gl eret
Hewr
Expected target
moscot mapping probability
0.5
fresuency
Livat
Commective Hopoo
uoed
Hebr
Liper
Figure 5: Comparison between moscot and NicheFlow on mapping the liver structure from E10.5 to E11.5. The liver at time E10.5 is used as source for trajectory prediction using the different models. The left column shows the source and expected target regions highlighted on the respective E10.5 and E11.5 embryos. The middle column displays the density of the prediction obtained by transporting niches from the source to the target slide. On the right, the aggregated cell type proportions according to the density in the middle column (see App. F.8).
D.2 parameter sweep in moscot
The comparison with moscot assesses how well each model captures compositional changes in anatomical structures or migratory patterns, depending on the use case. Qualitatively, we found that the spatial component in the OT problem considered moscot and regulated by the hyperparameter a is comparable to our spatial term a across the 0.1-0.9 range the Fig. 3 analysis.
In Tab. 2, we report the proportion of source density mapped to the correct cell type across values, to support our choice for the results in Sec. 5.2. For the spinal cord (Fig. 3, left), values above 0.75 yielded better qualitative and quantitative results. For neural crest cells (Fig. 3, right), a 0.5 performs best.
More specifically, when mapping fixed structures over time, values below 0.75 caused excessive density dispersion outside the anatomical region. For migration, no value led generally accurate transitions, though 0.5 mapped the highest density to the expected cell type.
Table 2: Effect of the parameter balancing spatial and cell state preservation in moscot. The results in the table indicate the percentage of source density mapped to the correct cell type from the source anatomical structure (the higher, the better).
Tissue
0.1
0.25
0.5
0.75
0.9
Spinal Cord
0.146
0.190
0.712
0.729
0.725
Neural Crests
0.159
0.160
0.162
0.155
0.112
D.3 Additional experiments on the axolotl brain development and aging datasets
We provide additional visualizations of the generated samples on the axolotl brain development dataset, presented in Fig. 6. As can seen from the figure, NicheFlow correctly retrieves the spatial and anatomical characteristics of the brain, including hemisphere formation and cavity.
19

---

Source 44)
SPFTow CFM
SPF GVFM
SPPiow GLVFM
Nicheflow CFM
Nicheflow GVFM
Nicheflow GLVFM
Target
Source (St.
SPflow CFM
SPFlow GVFM
SPFlow GLVFM
Ncheflow CFM
Nichellow GVFM
Nicheflow GLVFM
Target 57)
Source SPFiow GFM
SPFlow GLVFM
NcheFlow NicheFiow GVPM
NicheFlow GLVFM
SPFlow CFM
SPFlow
SPFlow
NichFllww
NicheFlow GVFM
NicheFlow GLVFM
Target (Acull)
(Adut)
Target (Meta)
CMFN
Immeture colN
MCG
VLMC
CNBL2
dpEX
reyiN
dnFX
cckNN
cNBL3
Immature CMPN
Immature mpEX
Oligo
CEGC
CNBL4
mplN
raEGC
INBL
Immature DMIN
Immeture npoxEX
Unknown
cNBL1
cNBL5
WHEGG
Immature
Figure 6: Qualitative comparison of generated samples on the axolotl brain development dataset (Stage 44, 54, 57, Juvenile, Adult, Meta). We show source and target samples alongside predictions from SPFlow and NicheFlow with different objectives. NicheFlow captures the spatial structure and cell-type organization more faithfully across developmental stages.
To further support result, propose application on the axolotl brain development dataset. More specifically, in Fig. 7 we demonstrate that model predicts the formation of crucial anatomical structures like the left and right lobes both spatially and compositionally. This vouches for flexibility in NichFlow's performance, which extends non-trivial topology changes and simultaneously accounts for accurate cell state and coordinate generation in time. Similar results can be observed in Fig. 8, where we showcase the correct prediction of the formation of left lobe cavity, predicting trajectories from immature brain region.
in-depth qualitative analysis of the NicheFlow
Moreover, in Fig. pallium in the axolotl brain development. Following Wei et al. [42], we know that early time points populate the dorsal pallium of immature cell types like ependymoglial cells (EGC), neuroblasts (NBL), and immature neurons. In the left dorsal pallium, these disappear at the juvenile stages (the 3rd) and lead to differentiation into mature neurons (nptxEX) and later EGC (WntEGC, sfrpEGC). This fixed structural development was accurately predicted by NicheFlow when pushing left dorsal pallium cells forward across the trajectory.
predict the compositional and structural time evolution of the left dorsal
Finally, similar to Fig. 2 and Fig. 6, Fig. 10, we qualitatively show that NicheFlow with microenvi- ronment sampling strategy and mixed-factoriaa factorized VFM is the best approach for reconstructing mouse brain trajectorics in time.
20

---

Source
EDGC
EM812
Predicted frequencies
Generated niche
Generated region
St.54
Retwerce
Innature CuN
Inmature dpEx
Inmature maEk
VLMCC
DEGC
(NELL
ONILD
St. 57
AMEZG
Juv.
Aclorenc
UNEL WHEDOC
MSN
UnAnoAn
VLMC
QUEX
mpEx
TOYN
Adult
Generated
Eouo
SKO
WHZGG
Meta.
C7 MDC
RECEEX
Gevenste
STADGG
Unrrown
Figure 7: Prediction of hemisphere formation the axolotl brain development. We evaluate the ability of NicheFlow to generate the anatomical splitting of central brain structure upon the formation of the right and left brain regions. The top panel shows the reference source region and cell type composition. For each stage (St.54, St.57, Juvenile, Adult, and Meta), we show: (1) predicted vs. reference cell type frequencies, (2) the generated niche visualized via 2D embedding colored by cell type, and (3) spatial projection of the generated region onto the anatomical reference (right). We set 0.5.
21

---

Source
difa. chBL11 chBL2
Predicted frequencies
Generated niche
Generated region
St.54
Immmatone uAl
cerented
Inmmature mpEx
dNAL1
eetter
(NAL2
dNBLI
strpEGC
St. 57
felourss
mpex
Inmature roraEx
VINC
dEGC
UNBLL
ONBL2
OURI
stuEGC
WHEEGC
Juv.
SIUN
MON
Reference
Untnown
WHRGGC
VUNC
resiN
maEx
noexEx
Adult
sstN
Reforence Generred
MSN
LNBL
moEx
nocrEX
noyiN
rbbGCC
strEOC
Meta.
NCG
NSN
WITECC
Olige
LiVETOW
VLHC
dOEK
MIEE
niDEE
noyIN
Figure 8: Prediction of cavity formation in the left brain lobe during axolotl development. We assess the ability of NicheFlow model the emergence of cavity structure within the left lobe of the axolotl brain. The top panel shows the reference source region and its cell type composition. For each developmental stage (St.54, St.57. Juvenile, Adult, and Meta), we display: (I) predicted versus reference cell type frequencies, (2) the generated niche visualized in 2D embedding space, colored by cell type, and (3) spatial projection of the generated region onto the anatomical brain reference (right). The progression illustrates the model's 's ability to recapitulate the asymmetric cavity formation localized to the left lobe. set 0.5.
22

---

Compositional evolution left dorsal pallium
Composition
Time poiti
Time point =
Time point =
Time point =
Time point -
Time point -
Iminature Mature
0.0
02
Cel type
Source
Generated region
Figure 9: Structural and compositional prediction the left dorsal pallium during the axolotl brain development. We use the left dorsal pallium region at the St.44 developmental stadium (highlighted on the left) and predict its trajectory over time. In the top row, we show the proportion of different cell types in the predicted region, colored immature (light red) and mature (dark red). At the bottom, we show the structural prediction for the left dorsal pallium overlaid on the true slide. We set 0.5.
Source (3.4m.) SPFlow CFM
SPFlow GVFM
SPFlow GLVFM
NichoFlow CRM
NitheFlow GVFM
NicheFlow GLVFM
Tergei (3.8m.)
Source (3.8m.)
SPF CFM
SPFlow GVFM
SPFlww GLVFM
NicheFlow CFM
NicheFlow GVFM
NicheFlow GLVFM
Targel (4 (4.3m.
Source (4.3m.)
SPFiew CFM
SPFlow GVFM
SPFlow GLVFM
NicheFlow CFM
NicheFlow GVFM
NicheFlow GLVFM
Target (5.4m.)
Sourca (5.4m..)
SPFiew CFM
SPFlew GVFM
SPFlow GLVFM
NichaFlow CFM
NicheFlow GVFM
NicheFlow GLVFM
Target (6.6m.)
Source (6. (8.8m.)
SPFiow CFM
SPFIow GVFM
SPFlow GLVPM
NicheFlow CFM
NicheFlow GVFM
NichoFlow GLVFM
Targer (3.8m.)
Astrocyte
Ependymal
Nesrra-Inhiiitrry
Macrophage
Neumbiast
Neuran-MSN
Oligodendraaytee
VLMC
Endothelial
Microglia
Neuron-Exctattry
Neutrophi
Pericyte
VSMC
Figure 10: Qualitative comparison of generated samples on the mouse brain aging dataset (3.4, 3.8, 4.3, 5.4, 6.6, 9.8 months). We show source and target samples alongside predictions from SPFlow and NicheFlow with different objectives. NicheFlow captures the spatial structure and cell-type organization more faithfully across developmental stages.
23

---

D.4
Conditional generation
Source It:)
Target
Predictions (50x)
Sample Likelinood
-06
6.2
Source it,)
Target (t2)
Predictions (50x)
Sample Likelinood
64
0.6
6.4
622
Source (t:)
Target
Predictions (50x)
Semple Likelinood
6.6
4.4
6.2
Source It.)
larget
Predictions (50k)
Likelinood
-04
0.2
Tareet
Predictions
Likelinood
6.8
0.4
0.2
Source
Predictions
Likelinood
0.8
Usee
Figure 11: Qualitative evaluation of conditional generation with NicheFlow on the embryonic development dataset. Each row corresponds a different source microenvironmentas time t (blue), shown alongside the ground truth target microenvironment at time to (green). The third column displays 50 samples (orange) generated conditionally by the model, while the fourth column visualizes kernel density estimate of the sample likelihood over the spatial domain.

---

To assess whether the model accurately conditions on input microenvironments, we visualize the spatial distribution of generated samples given fixed source microenvironment. Figure 11 illus- trates such cases on the embryonic development dataset, showing the input region at time 11, the corresponding target at time 12, and 50 independently generated samples from NicheFlow. Directly computing the likelihood of ground truth cell coordinates under our generative model is intractable, as it would require evaluating the density of an implicitly defined distribution over point clouds. Instead, we approximate the spatial likelihood using kernel density estimation (KDE) over Monte Carlo samples drawn from the model. Given a set of generated coordinates {c, L.V. likelihood at ground truth location as:
we estimate the
p(c) Σ exp
- 20²
(32)
where a is a fixed bandwidth parameter. The resulting KDE heatmap visualizes the spatial concentra- tion of samples, allowing us to qualitatively assess whether the model produces consistent predictions conditioned on the source microenvironment.
D.5 OT Ablation Study
We begin by emphasizing that the optimal choice of the OT feature-weighting parameter in Eq. (12) depends on the downstream application. This parameter determines the relative importance assigned to gene expression versus spatial coordinates during OT. For developmental processes such as organogenesis or regeneration, higher values of X prioritize transcriptional similarity, facilitating the reconstruction of continuous differentiation trajectories and capturing fate-driven transitions, which may span large spatial distances. In contrast, for applications that aim to monitor changes within fixed spatial region, such as shifts in cell-type composition over time, a lower À is more appropriate. In such cases, emphasizing spatial locality helps avoid spurious long-range transport assignments caused by molecular noise.
Table 3: Ablation study on the feature-coordinate trade-off parameter À in Eq. (12) across mouse embryonic development (MED), axolotl brain development (ABD), and mouse brain aging (MBA). We report INN-F1, PSD, and SPD for each training objective (CFM, GVFM, GLVFM) and model variant (RPCFlow, NicheFlow). All results are averaged over five evaluation runs. Bold indicates the best value per dataset and metric.
MED
ABD
MBA
Model
Olg
INN-FI
PSD 1
SPD
INK-FI
PSD
SPD |
INK-FI 1
PSD |
SPD
RPCFFow
CFM
0.10
0.56::amm
0.981 20004
0.564
20005
0.524
2.05|
11301
1,015
11.00
0,271 :cam:
1,543 :c00
0.810:100)
RPCFliww
CFM
0.25
0.545 st
0.95-ccMM
0.570
2401
0.511
-1301
2.070
:1985
004
0.273
1.535 -t.e
0.818
-1803
RPCHlw
CFM
0.50
0.554 2c.00
0.985_c00
-C001
0.582
-2001
0.507
.6201L
2,077
.1385
1,023
a1002
0,273 caaes
1,566 20000
0.810 -130
RPCFlw
CFM
0.75
0,537
icawe
0.981
0,595
_6002
0.517
_6205
2,033
.ean2
1,026 11000
0,275
11945
1,553 110K5
0,801
.3300
RPCHlw
GLVFM
0.10
2ca00
0.979
0.586.6000
0.554
-1200
2,053
L,03 11005
0,265 cca0s
1,723 21
0.779
-1801
RPCFlw
GLVTM
0.25
0.593 1ca00
0.924
-8924_0201
0.575 _0007
0.555
.62013
2,076
11387
1,06 11304
0.267 11a45
1.728 11305
0.777
RPCFlww
GLMPM
0.50
0.585.1.00
0.934
c00
0.509 -6001
0.551
+800
2,038
.03M5
1,122 4000
0.269 c1000
1,715 +1000
0.763
RPCFlow
GLMPM
0.75
0,593
20000
0,945
c0005
0,570
0.561
+8008
2.014
1,055 #10
0.268 ca0s
1,675
+0000
0,754
1
RPCFlow
GVFM
0.10
0.503
.0000
1,155 -COM
0.578
-800)
0,477
4800
26260
1883 4000
0.249 ca0es
1,753 +0000
0.764
RPCFlw
GVFM
0.25
0.20).amm
1,223 -C0M
0,566 cm0m
0
:00044
2.364- :1395
1.030
0.246 :1:45
st
0.787
RPCFkw
GVFM
0.50
0.521 sco0n
11855
0.569-6003
-6004
0.480-6004
1,025 +1
0.245 m:
1,756 at
0.774
RPCFlw
GVTM
0.75
1c800
1,202 .00014
0,573
0.471 280012
2.476
L,037 11.05
0,248
1831 11005
0,750
-1802
NichoFFow
CFM
0.10
0.609 anns
0.402-8000
0.604 -001
2,086-1-999
0.568 +100
0.283 -1:000
1,557 stan
0.556-19m0
NicheFllow
CPM
0.25
0.569
1C801
0.973
0.425 _0082
0.586.12003
2.106 .tanz
0,565 1100
0.251
11.46
1,46 L1
0.612 ..002
NicheFlow
CTM
0.50
0.551 .caxe
1,051
0.471 _02113
0.585 .1301
2.089 .1.30
0.591 11005
0.283
11942
L,88 11905
0.604 1.el
NidieFlo
CFM
0.75
0.519 ccus
1.103 (150)
0.515
0.571
2.126x6a11
0,592 100
0.278 ca96
1,566 9
0.616 1
Nichellow
GVTM
0.10
0.5951cact
0.991 -0317
0
-6025
0.574
11315
2.230 JEar0
0.594 110M
0 11945
L661 11:03
0,531
NidieFio
GVFM
0.25
0.563ccmm mc
1.051 -0117
0.491
40405
0.571
+800
2,343
40313
0.619 1100:
0.265 camss
L,59) 40002
0,590
NidieFlo
GVEM
0.50
0,533 .coe)
1.034 01442
0,000 -154E
0.556
2.283 483135
0.742 211144
0.269 caas
1,607 .190)
0,605
.000
NicheFlow
GVFM
0.75
0.526xtart
1,121 -1917
0,870
0
2.156-150
0.684-15m2
0,263 ma
16613-cm
0,731
NicheFlow
GLNFM
0.10
0,664
0,883 -10M
0.625
2.079 :1:M3
0,576
0,285
1,554 =1
0,530=19e>
NidhoFlow
GLNFM
0.25
0.629xcam)
ccom
0.923-6514
0,394
0668-8800
101 -sam4
0.577 -1am
0.264
1,509 -10m
0.549
Nichefllow
GLVTM
0.50
0.610
0.919
0.417
2.136
0,579
0,283
1,600
0,546
Nichelllow
GLNTM
0.75
0,592
0,879
0,472
0.600
16204
2.106 1133)
0,592
0,281
1,573
0,595
In Tab. 3, we sweep the value of À across multiple models and training strategies on the same task as in Sec. 5. showing consistent results across settings as reported in Tab. 1.
Furthermore, Figure 12 visualizes the impact of varying on the OT plans described in Sec. 4.2, where intermediate values yield more coherent and biologically plausible couplings.
25

---

0.00
A 0.10
A = 0.25
A- 0.50
A-0.7
A
==000
Figure 12: Visualization of OT couplings computed under varying values of the pooling parameter in eq. (12), which balances spatial coordinates and gene expression in microenvironment matching. Lower values prioritize spatial proximity, resulting more dispersed and less structured alignments, while intermediate values yield tighter, biologically consistent mappings. Very high À settings ignore spatial context and may lead to implausible long-range matches.
D.6 K-Means regions ablation study
To ensure diverse and spatially distributed sampling during training, we partition each tissue section into K spatial regions using K-Means clustering over the 2D cell coordinates (see Sec. 5.1.1). At each training step, microenvironments broad spatial coverage and preventing oversampling of densely populated areas. In this ablation study, we investigate how varying the number of spatial regions, K, affects model performance.
sampled uniformly from within these regions, encouraging
Table 4: Ablation study the number of spatial regions K defined over the datasets. We evaluate NicheFlow with the GLVFM objective across three datasets: mouse embryonic development (MED), axolotl brain development (ABD), and mouse brain aging (MBA). Results are reported as mean standard deviation over five evaluation runs.
MED
ABD
MBA
INN-FL
PSD L (10")
SPD (10")
INN-FL 1
PSD (10")
SPD (10")
INN-FI 1
PSD J (10")
SPD (10")
0.640
400043
0.617 400017
2000s
0.283
16
0.661
000333
0.881
0.393
à 4n00w
0.633 moses
1.936 -mm
0.628 aome
0.283
nocon
1.564 .00m
0.538 omee
32
0.659
0.899
0.391
0.622
1.968
0.640 1000z
0.279
1.600
0.537
64
0.664
28/004
0.883 20,004
0.398 2 000025
0.628
000013
-009220003
0.57624000
0.285 00000
1.554 280021
0.532 : 2,000
As shown in Figure 13, increasing leads to increasingly fnn-ggraned spatial partitions. While moderate values of K help improve spatial resolution, excessively high values (e.g., K = 128 or 256) result in overly small and fragmented regions. This can cause significant overlap between sampled microenvironments and reduce sampling diversity. Moreover, in sparsely populated tissue sections, high values may yield regions with insufficient cells, degrading both representativeness and stability.
26

---

293 810.5 T-8115
393
Figure 13: Visualization of spatial partitions obtained via KMeans clustering with K (4,8, 8, 16, 32, 64, 128, 256} on the embryonic development dataset. For low K, each region cov- ers large heterogeneous areas; for high K. regions become small, dense, and highly overlapping, potentially degrading the diversity and utility of sampled microenvironments.
Conversely, using too few regions (e.g., K = 4) results in broad spatial partitions that may encompass multiple heterogeneous tissue compartments. This undermines the locality assumptions of our model and increases intra-region variability, which can impair the model's ability to learn.
Given the trade-offs outlined above, we focus our evaluation on E (8, 16, 32, 64), which spans range of granularities that preserve both spatial interpretability and sampling robustness. Tab. 4 sum-
27

---

marizes the quantitative results for this ablation across three datasets: mouse embryonic development
(MED), axolotl brain development (ABD), and mouse brain aging (MBA). We observe that using
K 64 consistently yields strong performance, achieving the highest INN-F1 scores on both the
MED and MBA datasets, while also performing competitively on ABD. These findings indicate that
K 64 offers an effective balance between spatial resolution and stability, and we adopt it as the
default configuration in our main experiments.
D.7 Justifying the choice of NicheFlow over RPCFlow
Although RPCFlow achieves competitive performance on spatial metrics such as PSD and SPD, lacks the essential capability of meaningful conditional generation. In RPCFlow, conditioning is performed using randomly sampled point clouds from the spatial regions without explicit microenvi- ronment structure. As shown in Fig. 14, when conditioned such random sources, RPCFlow tends to reconstruct the entire target tissue, rather than capturing local dynamics driven by the source input. This undermines its ability to model spatiotemporal evolution in biologically grounded manner.
Source
Target
Predictions (50x)
Semple Likelhood
0.045
8.000
0,22
0.20
8.15
8010
Source (y)
Target ((g)
Predictions (50k)
Sample Likelihood
0.04
Figure 14: RPCFlow fails to perform meaningful conditional generation: the model generates diffuse reconstruction resembling the entire target tissue rather than conditioning on the provided source points. Each row shows different source-target pair.
In contrast, NicheFlow explicitly designed to push localized microenvironments through time. By conditioning on fixed-radius neighborhoods centered around specific spatial positions, the model learns how cellular contexts evolve, preserving both spatial coherence and transcriptional identity. Figure 15 visualizes this distinction: while NicheFlow consistently generates well-localized pre- dictions aligned with the input microenvironment, RPCFlow fails retain spatial specificity, often diffusing the prediction across broader regions.
From biological perspective, predicting the fate of local tissue region over time is far more relevant than mapping random point sets. Microenvironments encode structured cellular contexts, such as signaling interactions or niche-specific cell states, that are crucial for downstream analysis (e.g.. lineage fate prediction, microenvironment baasedd- intervention simulation). Because RPCFlow lacks this interpretability and fails to enable such downstream tasks, it cannot serve practical generative model in biological settings.
In summary, while RPCFlow may appear performant under some aggregate spatial metrics, only NicheFlow enables localized, conditionally consistent generation of evolving tissue regions. This makes not only superior for evaluation but also for practical use in biological modeling.
D.8 Structure-aware evaluation
While we refer to our primary evaluation metries as PSD and SPD, they correspond to the two asymmetric directions of the Chamfer Distance (CD) widely used metric in point cloud generation [46]. PSD measures how closely the predicted points adhere to the ground truth structure (fidelity), while SPD captures how well the prediction covers the full extent of the target (coverage). We
28

---

[t:)
Target
Semple
0.6
0.4
(a) NicheFlow Sample
Source It.)
Target (t))
Predictions (50x)
Sample Likalihood
0.8
(b) NicheFlow Sample 2
Source
Target ((g)
Predictions (50k)
Sample Likelihood
0.15
0.10
(c) RPCFlow Sample
Source it.)
Target (t))
Predictions (50x)
Sample Likelhoodd
0325
0.00
0.055
0.050
0.025
(d) RPCFlow Sample 2
Figure Comparison between NicheFlow and RPCFlow fixed-source microenvironment setting. Each row shows imput (source), ground truth target, model prediction, and KDE-estimated likelihood. While NicheFlow (top) produces well-localized samples consistent with the input, RPCFlow (bottom) generates less structured and spatially inaccurate predictions.
compute both metrics over collections of microenvironments tiled across the tissue, providing an indirect assessment of global morphological accuracy.
To directly assess the preservation of internal structure within generated microenvironments, we additionally evaluate [47, 40]. These metrics compare the intra-poing pairwise distances in predicted and ground truth point clouds, capturing latent structural relationships. FGW further incorporates transcriptomic similarity, offering a holistic measure of structural and feature- feature-based alignment.
-a (GW) and Fused -aatoov-ssserteennn (FGW) distances
Given the computational cost of GW and FGW, we restrict their evaluation to micrenniionmmen-leee generations. For each source niche, we sample predictions from the trained model and compute GW/FGW against the corresponding real microennironment, averaging results across multiple runs.
29

---

The source-taie target pairs are obtained by computing OT between all microenvironments at time and t ++ 1, yielding a one-to-one coupling used for conditioning and evaluation.
As SPFlow generates individual points and RPCFlow produces randomly sampled, unstructured clouds, these models lack coherent internal structure and not amenable to structural comparison via GW or FGW. We therefore report these metrics only for models that generate full microenvironments.
The results in Tab. 5 show that NicheFlow, particularly when trained with the GLVFM objective, consistently achieves lower GW and FGW distances compared to baseline models. These findings underscore the structural fidelity of our predictions and further validate the modeling advantages of microenvironment-uuee aware generative training.
Table 5: Gromov-Wasserste (GW) and Fused Cromov-Wassersteee (FGW) distances between generated and real microenvironments on mouse embryonic development (MED), axolotl brain development (ABD), and mouse brain aging (MBA). All models are trained with fixed À = 0.1. Results are averaged over five generation runs. Bold indicates the best (lowest) value per column.
MED
ABD
MBA
Model Objective
GW (102)
FGW
GW (10")
FGW
GW (10?)
FGW
NicheFlow
CFM
0.315 2000
3,273 a 0,003
0.783
0.002
3,543 0.00
0.315 0.001
3,531 20000
NicheFlow
GVFM
0,463 10001
3.167 0.00
0.797 1.006
3.414 0.003
0.339
1.000
3.387 20000
NicheHlow
GLVFM
0.224 0.00
3.147 0.007
0.720 0.004
3.420 0.00
0.262 0.000
3.367 0.000
D.9 Wasserstein metries
To complement the local structural assessments with GW and FGW, we additionally include global evaluation metries that are widely adopted in spatial transcriptomics [48-50]. Specifically, we report the 1-Wasserstein (Wi) and 2- 2-Wasserstein (W2) distances between predicted and real cell distributions, computed separately over spatial coordinates and gene expression features. These metrics are calculated at the level of individual cell types and averaged across all timepoints and generated samples. This provides a more global perspective on how well the model reconstructs tissue-wide spatial and transcriptional distributions.
Table 6: Comparison of Wasserstein distances between generated and real tissue structures on the mouse embryonic development (MED) dataset. We compute W1 and W2 distances separately for spatial coordinates (Pos.) and gene expression features (Genes), averaged across all generated cell types and timepoints. All models are trained with fixed value of À = 0.1, and results are averaged over five generation runs. Bold indicates the best (lowest) value per column.
Model
Obj.
Wi Pos.
Wi Genes
W/ Pos.
W2 Genes
SPFlow
CFM
0.634 a 0.004
6.545 10004
0.773 1000
5.585 1000
SPFlow
GVPM
0.612 +0.00
6.517
+0002
0.743 +000s
5.609 +0.113
SPFlow
GLVFM
0.613 a 0.001
6.457
1.0000
0.738 10001
5,546 10.008
RPCFlow
CFM
0.290 0.906
6.303 +0.04
0.460 o007
5.386
RPCFlow
GVFM
0,258 a 2..00
6.134
1 0.004
0.398
0.006
5.292 a 101114
RPCFlow
GLVFM
0,221
a 0.003
6.087
= 0.002
0.365 0.00
5,244 : 0.1113
NicheFlow
CFM
0.253 0.002
6.177 0.000
0.420 0.000
5.314 101111
NicheFlow
GVPM
0.222 0.004
5.985 0.005
0.352 0.006
5.173 0.112
NicheFlow
GLVFM
0.212 a 0.004
5.930 = 0,003
0.342 0.007
5.031 0.003
We evaluate all model variants and training objectives across the three datasets (MED, ABD, MBA) and report the results in Tabs. 6 to 8. NicheFlow with the GLVFM objective consistently achieves the lowest Wasserstein distances across nearly all dimensions in the embryonic development dataset (MED), outperforming both baseline models and alternative objectives.
In the ABD and MBA datasets, which exhibit more gradual spatial evolution or less distinct cell-type boundaries, we observe that RPCFlow sometimes matches or marginally outperforms NicheFlow in isolated metries. However, as we discussed in App. D.7, RPCFlow not a conditionally consistent model: it does not produce meaningful trajectories given a source microenvironment and lacks
30

---

Table Comparison of Wasserstein distances between generated and real tissue structures on the axolotl brain development (ABD) dataset. We report Wi and W2 distances separately for spatial positions (Pos.) and gene expression features (Genes), averaged across all generated cell types and timepoints. All models use À = 0.1, and results are averaged over five generation runs. Bold indicates the best (lowest) value per column.
Model
Obj.
Wi Pos.
Wi Genes
W22 Pos.
W2 Genes
SPFlow
CFM
0.474 0008
6.711 =000
0.573 0.00
6.096 0.00
SPFlow
GVFM
0,538 0000
6,681
2000
0.639 0.00
6.209 40009
SPFlow
GLVEM
0.489 40002
6.533
20,002
0,586 0.003
6.090 0.166
RPCFlow
CFM
0.234 10,003
6.363
0.374 0.005
6.031 0.089
RPCFlow
GVFM
0.247 0.002
6.093 0:006
0.369 0.00s
5.921 0.051
RPCFlow
GLVFM
0,254 0000
6.070 0.084
0.406 cums
5.884 0,04
NicheFlow
CFM
0.267 +0.02
6.347
0.003
0.428
6.107 0.003
NicheFlow
GVFM
0.284 0.004
6.081
0:004
0,448 0.006
5.917
NicheFlow
GLVEM
0.279 0.006
6.085 0:000
0.444 0011
5.923 0.004
Table 8: Comparison of Wasserstein distances between generated and real tissue structures on the mouse brain aging (MBA) dataset. We compute Wi and W2 distances for spatial coordinates (Pos.) and gene expression features (Genes), averaging per cell type and timepoint. All models are trained with À = 0.1, and results are averaged across five generation runs. Bold marks the lowest value in each column.
Model
Obj.
Wi Pos.
Wi Genes
We Pos.
W2 Genes
SPFlow
CFM
0,515 #0001
7.544 20004
0.600 400m
5.237 #0.66
SPFlow
GVFM
0,564 10001
7.461 20.041
0.646 20001
4.961 0.114
SPFlow
GLVEM
0.552 0.002
7.431 =0.000
0.632 0002
4.117 0.19
RPCFlow
CFM
0.385 00000
7.318
0:006
0.480
6.659 0.09
RPCFlow
GVFM
0.416 0.08
7.202 0000a
0.510 0.00
6.346 0.66
NicheFlow
CFM
0.416 0.088
7.321 00012
0.510 a 0.008
6.487 0.084
NicheFlow
GVFM
0.412 0000s
7.102 0.001
0.503
6.403 0.07
NicheFlow
GLVFM
0.431 00008
7.121 0:000
0.527 a 0.007
6.298 0.029
biological interpretability. Therefore, even when RPCFlow achieves slightly lower Wasserstein distances in aggregate, these gains not actionable in practice, as the model cannot be used to infer biologically meaningful transitions or trace the evolution of specific spatial regions.
Taken together, Wasserstein metrics offer complementary, global validation of the structural and semantic accuracy of NicheFlow. They reinforce the quantitative evidence already provided by Chamfer- Chanfer-based metrics (PSD/SPD) and Gromov-based structural evaluations (GW/FGW), confirming that our model accurately captures the tissue-scale organization of both spatial and gene expression patterns.
E
Algorithms
Here, we present the algorithmic procedures underlying NicheFlow. To streamline the exposition, we first define compact notation for representing noisy microenvironments and their interpolations.
We denote single noisy microenvironmenn (simplifying Eq. (11)) as:
M (O, I0+2) 1*k [(c], 2;) of ~ (O, L)) 2; (O Vi 1....,, k).
(33)
where & denotes the number of spatial points per microenvironment, ot 2D cell coordinates, and z; associated feature vectors in R For collection of N microenvironments, define the corresponding set of noisy microenvironments
M ~ I = (Mi ~ No I [i 1, N}.
(34)
Given noisy sample M and a corresponding ground-truth microenvironment M we define their linear interpolation at time t e [0, as: M (1 t)M² + tMM, where the interpolation is
31

---

performed clement-wiss across the matrix rows. For batched data, this generalizes to:
M (1- LM
(35)
Algorithm 1 SAMPLEANDILTERPOLLSS
Input: Number of samples N, feature dimension D, microenvironment size k, OT plan =
Output: Source (M"), larget (M'), ) interpolated (M'), and noisy (M*) microenvinonments
I: (M". M < Sample from K-Means regions
Initial microenvironment pairs
2: (M" M²) < T., Nxk (M)
Resample with OT plan
3: M (N
Noisy initial states
t t~0
Random interpolation time
5: M + (1 t)M + tM
Linearly interpolated states
E.1 OT Conditional Flow Matching
The OT Conditional Flow Matching (OT-CFM) algorithm consists of training and a generation phase. During training, the model learns time-dependent vector field ug conditioned on a source mi- croenvironment M which maps noisy microenvironment M M from Gaussian noise-to its corresponding target microenviromment M The supervision provided via cource-target pairs (M,,M') obtained through EOT over pooled microenvironment representations. At generation time, the learned vector field is integrated starting from M conditioned on given source M to generate the predicted target microenvironment M
We optimize the following loss:
L(M", M M
Σ || (M' -
(36)
MM M'EM
M'EM'
M'EM'
The full pseudocode for both phases is provided in Algs. 2 and 3.
Algorithm 2 OT CFM Training
Input: Number of samples N, feature dimension D, microenvironment size k. OT plan ,.... conditional velocity field uf
Output: Trained parameters of uf
MMMMM. M .M' M (M) SAMPLEANDINTEEPOLAAEEA D, x2,)) Microenvironments (Alg. I) 2: 0 VoL(M M'
I: (M'
b Compute loss (Eq. (36)) & update parameters
Algorithm 3 OT CFM - Generation
Input: Source microenvironment M learned conditional velocity field uf Output: Generated microenvironment M
I: M (0 ID-2) 1xf
Sample noisy sample b Solve ODE
2: M + M + lo ui (M' MM) dt
E.2 OT Gaussian Variational Flow Matching
The OT Gaussian Variational Flow Matching (OT-GVFM) algorithm adopts variational perspective on Flow Matching, Instead of directly learning a time-dependent conditional velocity field, the model learns a factorized variational posterior of (M² .M'MM) over target microenviroomennss M conditioned on an interpolated microenvironment M and source microenvironment M The predicted velocity field then computed as the difference between the posterior mean g (M* and the current state M
32

---

The training objective minimizes the discrepancy between ground truth targets and the predicted posterior means (ft r!)):
C(M M M (4)
Σ M (OLM)JMM M'EM' (f. (EE"(MMM" M'EM'
Σ
(11= - 1 + 1| - -3,1))
(37)
2
At generation time, trajectories are generated by integrating the vector field μ (M' .M M'. starting from a noise-sampled microenvironment M and conditioned on given source M To ensure stability near 1, the vector field is scaled by a time-dependent denominator, in prior VFM formulations.
The pseudocode for both phases is provided in Algs. and 5.
Algorithm 4 OT-GVVFM Training
Input:
Number of samples N, feature dimension D, microeevirrnmment size k OT plan T source-conditioned posterior mean predictor 44
Output: Trained parameters 0 of u I: (M'. M .M', .M²) SAMPLEEANDINTERPPLNNNEE D, =2,)) Microenvirromments (Alg.
2: 0 < VoL(M", M M
Compute loss (Eq. (37)) & update parameters
Algorithm 5 OT-GVVM Generation
Input: Source microenvironment M learned conditocnd posterior mean predictor μ Output: Generated microenvironment M
I: M ~ (O, ID-2) 1xk
b Sample noisy sample
2: M M So
LM (MM) -M") M
dt
Solve ODE
1-t
E.3 NicheFlow: OT Gaussian-Laplace Variational Flow Matching
NicheFlow extends OT-GVFM by modifying the variational posterior: it assumes Gaussian distribu- tion for gene expression features and Laplace distribution for spatial coordinates. This change leads hybrid loss that combines an 2 loss on gene expression and an L' loss on spatial locations:
C(M MM ,²
Σ (MLM)JEM (f. (M) M'EM!
E M M'EM'
(Ic: [la--ft - fo + 2
1 =11=1 = 1
(38)
At generation time, the model integrates the velocity field defined by the difference between the predicted mean μ (MM,, MM) and the current state M starting from noise and conditioned on the source M identical to the OT-GVFM procedure.
F
Experimental setup
F1
Dataset description
We evaluate model three publicly available, time-resolved spatial transcriptomics datasets spanning development and aging processes. Each dataset provides single-cell resolution profiles with matched spatial coordinates and curated cell type annotations. A detailed summary of the organism, tissue, number time points, cell types, total number cells, and acquisition technology for each dataset is provided in Table 9.
F.2
Dataset and microenvironments preprocessing
Dataset preprocessing. All datasets used in our study underwent a preprocessing procedure appropriate for spatial transcriptomic analysis, involving total count normalization, logarithmic
33

---

Table 9: Overview of the reme-vesolved spatial transcriptomics datasets used our experiments. MED: Mouse Embryonic Development, ABD: Axolotl Brain Development, MBA: Mouse Brain Aging. Each dataset varies in organism, tissue type, number of timepoints, cell types, total number of cells, and spatial transcriptomies technology.
Dataset
Organism
Tissue
Timepoints
Cell Types
Cells
Technology
MED
Mouse
Whole embryo
3
54k
Stereo-seq
ABD
Axolotl
Brain
6
33
36k
Stereo-sey
MBA
Mouse
Brain
20
1.5M
MERFISH
transformation, and principal component analysis (PCA). For the mouse embryonic development [20] and axolotl brain development [42] datasets, total count normalization, and logarithmic transformation had already been applied; we, therefore, performed PCA ourselves on the transformed data. For the mouse brain aging dataset [43], we applied all three steps: we first normalized raw count matrices so that each cell had the same total expression, followed by a natural logarithm transformation of the form log( X + 1) to stabilize variance and mitigate the influence of large values. We then computed PCA on the fog-transformed data.
To reduce computational overhead due to high cell counts in the aging dataset, we subsampled the data by factor of 0.2.
Finally, we standardized PCA components across cells to ensure consistent scaling across time points. Spatial coordinates were standardized independently for each time point by subtracting the per-time-pe time-point mean and dividing by the standard deviation. This standardization preserves the relative spatial configuration while accounting for scale and position differences over developmental time.
Microenvironments preprocessing. To facilitate efficient training and enable consistent microen- vironment construction, we precompute all fixed-radius neighborhoods for each dataset using a radius r chosen based spatial resolution. To reduce variability in the number of neighbors and improve batching efficiency, we fix the number of nodes per microenvironment to the most frequent neighbor count observed within each slide. This standardization ensures structural comparability across microenvironments while significantly reducing computational overhead during training, as costly radius nearest-neighhor queries avoided runtime.
F.3 Quantitative evaluation metrics
We assess spatial fidelity using two complementary asymmetric distance measures. The point-to- shape distance (PSD) captures how much predicted cell positions diverge from the actual tissue layout, computed the average squared distance from each simulated point its nearest neighbor in the ground truth. Conversely, the shape-lo-point distance (SPD) quantifies how comprehensively the predicted distribution spans the target structure by averaging the squared distance from each ground truth point to its closest generated counterpart.
Let G denote the set of generated coordinates and Rt the set of ground truth coordinates time t. Define NN (Ct) as the nearest neighbor of c; E G* in Rt, and NN! gen (c:) as the nearest neighbor of
Cj Re in G Then, the two metrics are given by:
PSD
Σ
||: - NNt (N
(39)
G'ES c:E!
SPD
||. - NN Σ
gon (c, )12
(40)
R
RIER CRER
where G U and R R
E.4 Cell-type classification for evaluation
To evaluate cell-type fidelity of generated microenvironments, we train supervised classifier to assign cell type labels based gene expression profiles. We apply the same preprocessing steps
34

---

used for training our generative models: total count normalization, log sng-trassformaioon and PCA (see App. F.2). The resulting low-dimensional embeddings are used input features for a simple multilayer perceptron (MLP), trained to predict discrete cell type labels.
The classifier consists of two-layer feedforward network with ReLU activations and final linear projection to the number of cell types. It is trained using cross-entropy loss and optimized with the AdamW optimizer. We report performance using the weighted Fl-score.
We obtain strong classification results across all datasets. On the mouse embryonic development dataset, the classifier achieves weighted Fi-score of 0.85; axolotl brain development, 0.80; and on the aging dataset, 0.97. These results correlate with the number of input genes and the variance retained in the PCA-reduced space. The aging dataset contains only 300 genes, and 50 principal components explain sufficient variance to accurately distinguish most cell types. In contrast, the embryonic development dataset contains approximately 2,000 genes, and the axolotl brain development dataset includes over 12,700 genes, making classification more challenging due to higher gene expression variability.
E,5 Discretized microenvironments
To ensure consistent and reproducible evaluation across methods and datasets, we construct fixed set of evaluation microenvironments by discretizing the spatial domain each tissue section. For each time point, define regular 2D grid over the tissue and select the closest cell to each grid point as the centroid of a microenvironment. Each centroid is then used to construct a fixed-radius neighborhood, following the microenvironment definition in Section 4.1. This procedure ensures full spatial coverage by verifying that every cell belongs to at least one microenvironment.
E9.5
E10.5
E11.5
336
Figure 16: Discretized grid of microenvironments for the mouse embryonic development dataset. Fach blue point denotes centroid around which microenvironment is constructed. To ensure consistent coverage across tissue sections, additional centroids are randomly sampled such that each slide contains the same number of microenvironments.
Figure 16 illustrates this discretization for the mouse embryonic development dataset. Each blue dot corresponds to a selected centroid. In cases where the number of grid-based centroids falls below target threshold, additional centroids are randomly sampled to match a fixed total count per slide. This augmentation ensures that all slides contribute equally the evaluation and prevents bias from sparse regions.
We apply the same discretization procedure all three datasets used in experiments: mouse embryogenesis, axolotl brain development, and mouse brain aging. By standardizing the evaluation regions spatially and deterministically, we eliminate the need for stochastic region sampling during evaluation, which would otherwise lead to nondeterministic and irreproducible results.
F.6 Microenvironment Transformer
To model the spatiotemporal evolution of local cellular neighborhoods, we design permutation- equivariant transformer-base architecture tailored to structured point cloud data. Our Microenviron- ment Transformer processes local microenviromments- sets of cells with gene expression features and spatial coordinates and predicts time-dependent outputs such as velocity fields future states.
35

---

The model operates on source M and noisy M microenvironments with per-cell features 2; R and 2D coordinates c; IR, The architecture consists of the following components:
1. Input Embeddings:
(a) Feature Embedding: A linear transformation applied the input features æ; of
each cell.
(b) Coordinate Embedding:
Spatial coordinates e, linearly projected and concate- nated the feature embedding.
(c) Time Embedding: For the noisy microenvironnent only, time E [0, 1] encoded using sinusoidal functions cos(wt) and sin(wt) followed by linear projection and concatenation with the input embedding.
2. Transformer Encoder (Source Microenvironment):
(a) Self-Attention: A stack of transformer encoder blocks with multi-head self-attention processes the embedded source microenvironment.
(b) No Time Embedding: Time information not provided the encoder, as encodes the source M
(c) Residual Feedforward: Each block contains feedforward subnetwork with LeakyReLJU activation and residual connection.
(d) Layer Normalization: Applied after both the attention and feedforward layers.
(e) Masking: Binary masks are used to ignore padding in variable-length microenviron- ments.
3. Transformer Decoder (Noisy Microenvironment):
(a) Time Embedding: Temporal context is injected into the decoder by embedding the time and concatenating it to the target point embedding.
(b) Cross-Attention: Decoder layers apply cross-attention between the noisy microenvi- ronment and the encoded source microenvironment.
(c) Self-Attention and Feedforward: Each decoder block includes standard self-attention and residual feedforward layers.
(d) Layer Normalization and Masking: As with the encoder, normalization, and masking applied throughout.
4. Final Output Projection:
(a) Prediction Head: A linear layer maps the decoder outputs to the desired dimensional- ity.
This architecture allows for flexible and expressive modeling of temporal dynamics in cellular point clouds. By encoding only the source and decoding the temporally conditioned target, the model supports variational and flow-based training objectives with explicit temporal conditioning.
F.7 Hyperparameters and Computational Costs
Model hype hyperparameters.
For all experiments, we use the same configuration for the Microenvi- ronment Transformer architecture. The full set of hyperparameters is as follows:
Input feature dimension: 50 PCA-based gene expression features concatenated with one-hot encoding of the time-point, resulting total dimensionality of 50 [T], where T is the number of slides (timepoints) in the dataset.
Input coordinate dimension: 2
Embedding dimension: 128
MLP hidden dimension: 256
Number of attention heads:
Number of encoder layers:
Number of decoder layers:
Dropout rate: 0.1
Output dimension: 52 (gene expressions features coordinates)
36

---

OT and mini-batching. To ensure spatial diversity and computational tractability, we uniformly sample 256 source-targe microenvironment pairs from the spatial clusters obtained via K -MMeass We then compute the entropic OT plan between these sampled pairs and resample 64 source-target pairs from this plan to defin a single training instance. During training, we process 16 such instances per batch.
Optimization. All models are trained using AdamW optimizer with a learning rate of 10-4 and weight decay of 1 10-5. We train each model until convergence.
Computational cost. All models were trained on a single NVIDIA GeForce GTX 1080 GPU with 11GB of memory. Depending on the dataset and training objective (e.g., CFM or VFM), training takes approximately 12-16 hours per model.
F.8 Comparison with moscot
Here, we describe how we conduct the comparisons with moscot, as shown in Fig. 3, Fig. 4, and Fig. For both NicheFlow and moscot, we select source microenvironment that we want to track over time. In the case of NicheFlow, this corresponds to an aggregate of point clouds. For moscot, it refers group single cells spatially located within the region of interest. The same set of cells is used for both methods.
NicheFlow. To generate contour plots over the spatial slide, push forward the selected region and assign the generated points to their nearest real neighbors based on spatial coordinates. We then compute a probability value for each real position by normalizing the number of assigned generated points. In other words, the more generated points that are close to given real point, the higher the probability assigned that location the contour plot. Cell type proportions are computed as the aggregated frequencies of the generated cell types across the slide. Each plot considers 10 independent generation runs from the same niche.
moscot. This baseline is not a generative model, but rather a standard discrete OT framework using a Fused Ceoomo-Waaseerserr cost. As such, does not generate new features or coordinates. Instead, it outputs a transition matrix that assigns matching probabilities between each source slide cell and each target slide cell. We use these transition probabilities to compute contour plots over the target slide and to aggregate cell type probabilities for the bar plots. For the latter, moscot provides a custom method called sell_transition()
For the Appendix figures Fig. 4.Fig. 7, Fig. 8, and Fig. 9, we propagate the initial source region across multiple time steps. In NicheFlow, this is achieved by using the simulated point cloud from step to predict the next state at t + 1, and then feeding this output as the input for the following trajectory step from t ++1 to t and this process. For moscot, pushing the source across time points can automatically done by setting non-subsequent time points in the cell, transiannin tion () function.
Ground truth points are not used intermediate sources during
37

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: The abstract and introduction state the contributions. These claims substantiated through targeted experiments and ablation studies.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
is fine to include aspirational goals as motivation as long as is clear that these goals not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification:
We discuss the limitations in App. B.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
the paper. The authors should use their best
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
38

---

Justification:
We introduce the Lixed-Factorized VFM formulation alongside NicheFlow in Sec. 4, and provide complete theoretical derivations in App.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
Justification: The experimental setup is described in Sec. 5.1.1, with additional details provided in Apps. E and F. We outline the dataset preprocessing, evaluation discretization procedure, backbone architecture, and all hyperparameter and optimizer settings. The appendices also include training and inference algorithms for OT-CFM and OT-VFM under both GVFM and GLVFM objectives.
Guidelines:
The answer NA means that the paper does not include experiments. If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, other means that appropriate to the research performed.
also be provided via detailed
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (c.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (c.g., to registered users), but it should be possible for other researchers to have some path to reproducing verifying the results.
39

---

reproducible, they 5. Open access to data and code
Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, described in supplemental material?
Question:
Answer: [Yes]
Justification: A link to the code is included in our project page, which is reported in the footnotes of page 1.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nin) cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (c.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https: /nips.
gwwcc//gg//e//eeeeeeeeeceeo.seee /aaeee//oooebbiiiiecc (esiinnlccy)) for more details.
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: We report the values used in our experiments in Sec. 5.1.1 and provide additional justification through ablations in Apps. D.5 and D.6 and App.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: We report the mean and standard deviation across five independent evaluation runs for all experiments and evaluation metries.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emoo bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
40

---

The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not show in tables figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification: We provide details on the hardware that has been used in this research in App. F.7.
Guidelines:
The answer NA means that the paper does include experiments.
The paper should indicate type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (c.g., preliminary or failed experiments that didn't 't make into the paper).
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics https http:///eee//eeeeeee https/cccccccccccccccccieeeeeee
Answer: [Yes]
Justification; The research presented in this work conforms with the NeurIPS Code of Ethics.
Guidelines:
The answer NA means that the authors have reviewed the NeurIPS Code of Ethics.
the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (c.g.. if there special consid- eration due to laws or regulations in their jurisdiction..
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification:
We discuss the broader impacts of our work in App. A.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
41

---

If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), faimess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g.. gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Question:
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Answer: [NA]
Justification: The research discussed this paper does not require safeguards to be put in place.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
Justification: Where relevant, acknowledge and cite the original authors of prior work, including associated codebases and datasets, and specify the licenses under which these resources were released.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
42

---

For scraped data from a particular source (e.g., website). the copyright and terms of service of that souree should be provided.
assets are released, the license, copyright information, and terms use in the package should be provided. For popular datasets, paperswi thcode com datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
If this information is available online, the authors encouraged to reach out the asset's s creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification: We use only publicly available datasets, which are cited appropriately. The complete codebase, including documentation and data preprocessing scripts, is public.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification: This work did not conduct research human subjects crowdsourcing experiments.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included the main paper.
According to the NeurIPS Code of Ethics, workers involved in data collection, curation, other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: This work did not conduct experiments where human subjects were involved and therefore does not require IRB approvals.
Guidelines:
43

---

The answer NA means that the paper does not involve crowdsoureing nor research with
human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent)
may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and expect authors adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Does the paper describe the usage of LLMs if it is important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, formatting purposes and does not impact the core methodology, scientific rigorousness, or originality of the research, declaration is not required.
Question:
Answer: [NA]
Justification: LLMs were not used for the development of our core method.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, or non-standard components.
Please refer to our LLM policy (https://neurt)) wwtccc////cooccecccccccoeooo cotp.c/oeeeeeenee/22/5////////// for what should or should not be described.
44


Reviews:
[{"review_id": "PlXdghHCCj", "reviewer": "Reviewer_1m53", "summary": "This paper proposes a framework, NicheFlow, for predicting both the physical migration of cells within the microenvironment and the associated changes in gene expression. The framework combines variational flow matching with optimal transport-based conditioning, and is evaluated on three time-resolved spatial transcriptomics datasets to assess its representational capacity and fitting accuracy.", "questions": "* Appendix D.3 explores how changing $\\lambda$ affects prediction and metrics, but the paper does not explain how to select $\\lambda$ on data not seen during training. Without a principled way to choose $\\lambda$, the method may have limited practical utility.\n* In Fig. 10, NicheFlow appears to better fit niche migration than RPCFlow, but could this difference be due to hyperparameter choices or distribution assumptions (e.g., Gaussian vs. Laplace), rather than model architecture itself?\n* In Fig. 3, NicheFlow seems to outperform Moscot, but this might simply be due to better parameter tuning. While both $\\alpha$ (in MOSCOT) and $1-\\lambda$ (in NicheFlow) control the spatial weighting, their formulations differ (e.g., GW loss vs. distance), so direct comparison at the same values is not necessarily fair.\n\nMinor Comments\n* What do the names RPCFlow and SPFlow mean? Please provide explicit definitions for these models.\n* Table 1 is not referenced in the main text.\n\nI would be open to raising my score if the authors are able to address the critical concerns raised above.", "limitations": "yes", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 2, "originality": 3, "strengths_and_weaknesses": "**Strengths**\n\n* The problem of jointly estimating niche migration and gene expression changes is well-motivated and biologically relevant.\n* The proposed model integrates recent advances, including the variational formulation of flow matching and OT-based alignment of unpaired source and target spaces, into a cutting-edge framework.\n* The paper conducts component-wise validation of the model, including ablation studies, which strengthens the credibility of the empirical findings.\n* The manuscript is well-structured, and the main idea is clearly explained.\n\n**Weaknesses**\n1. Ambiguous notation: The notation is somewhat unclear and potentially confusing. The paper refers to $x_0$ and $x_1$ as source and target variables, respectively, which implies that they live in different, unaligned domains. However, this creates confusion because $x_{t=0}$ and $x_0$ may refer to different variables.\n2. Lack of hold-out evaluation: The most critical issue is the lack of explanation regarding how the hold-out test set was constructed. Throughout the paper, it appears that the model is evaluated directly on the source and target slides used for training, thereby assessing fitting performance rather than generalization. Since the proposed model is generative and intended to handle unseen states, it is essential to evaluate its performance on held-out samples or slides. Without such validation, the model’s generalizability remains unclear.\n3. Model design: It is not evident that the current formulation of NicheFlow is effective in capturing heterogeneous niche dynamics. As described around eq. (13), the model appears to share the same flow functions $f_t$ and residual terms $r_t$ across all samples in $\\mathcal{M}^1$, which could lead to excessive averaging and loss of niche-specific information. Although Table 1 shows some improvements of NicheFlow over RPCFlow, the performance likely depends on proper tuning $\\lambda$, whose selection process is not discussed. This omission makes it difficult to judge whether the observed improvements are robust or incidental.\n4.  Limited evaluation validity: The evaluation in Fig. 3 appears to measure how well the model fits the training data, not its predictive performance on unseen inputs. This makes it difficult to assess the true modeling capability. Moreover, the comparison of hyperparameters such as $\\alpha$ (in MOSCOT) and $1 - \\lambda$ may not be fair due to differences in underlying loss formulations. A comparison under optimal settings for each method is needed.", "was_revised": true, "final_justification": "While I appreciate the new results provided in the rebuttal, I am not fully convinced by the claim that the value of $\\lambda$ can be selected based on the biological question without a principled procedure. The updated evaluation does provide a clearer picture of the sensitivity to $\\lambda$, but the absence of any holdout evaluation (e.g., a validation set or cross-validation) still limits the method’s generalizability and undermines model selection.\nOverall, I will raise my score slightly.", "rebuttal": "We thank Reviewer 1m53 for the constructive and thorough feedback. We appreciate the positive remarks on our model's structure and core idea, as well as the critical comments. Below, we address the remaining concerns.\n\n> Notation\n\nAs noted, $\\mathbf{x}_0$ and $\\mathbf{x}_1$ explicitly refer to source and target samples. When source and target spaces are incomparable, $\\mathbf{x}_0$ serves as a conditioning input, not as the ODE’s initial state, and our notation reflects it.\n\nWe intentionally avoid using $\\mathbf{x}_t$ for the ODE solution, instead denoting it as $\\phi_t^\\theta(\\cdot)$ throughout. Using $\\mathbf{x}_t$ would suggest $\\mathbf{x}_0 = \\mathbf{x}\\_{t=0}$, which is ambiguous when source and target domains are incomparable. By using $\\phi_t^\\theta$, we clarify that $\\phi_0^\\theta(\\mathbf{z} \\mid \\mathbf{x}_0) = \\mathbf{z}$, with $\\mathbf{x}_0$ as a condition rather than the starting point of the trajectory. This is made explicit in L138 to L140. The same convention applies to interpolations, where we introduce $g_t$ (e.g., Eq. (2), L105).\n\n> Hold-out evaluation\n\nWe agree that hold-out evaluation is generally desirable, and this could be implemented by leaving out either an anatomical region or an entire slide. However, both options are challenging due to **data constraints**.\n\nPredicting the trajectory of a region on a new slide requires prior observation of the same region in another sample at the same time point. Without this, the task becomes fully out-of-distribution, which remains difficult for current generative models [1]. Existing datasets lack replicated slides per time point, and to our knowledge, none provide repeated spatiotemporal measurements. In addition, current datasets span only coarse developmental stages and cover a few time points due to experimental costs. Holding out a slide would therefore introduce large temporal gaps, making it hard for the model to infer transitions.\n\nInstead, we showcase the model as follows:\n- ML perspective: We show that our approach improves over baselines at enforcing spatial constraints for spatiotemporal prediction (Fig. 2 and Tab. 1), evaluated on a grid of **test point clouds** (Sec. 5.1.1, App. F.5).\n- Biological perspective: We highlight its use for hypothesis generation, in contrast to exact OT methods, by analyzing spatial and transcriptional cell fate and the compositional evolution of selected anatomical regions (Sec. 5.2, App. D.6).\n\nOur goal is to introduce a scalable deep generative framework for niche-based spatiotemporal analysis, addressing the limitations of exact OT methods, which do not scale to large slides. With advances in dataset generation, we expect neural OT maps to play a key role in enabling generalization and modeling in spatiotemporal and perturbation contexts.\n\n> Flow function/residual \n\nIn NicheFlow, $r_t$ and $f_t$ are neural networks predicting posterior means over denoised target coordinates and cellular features (L213–222). These means are used to infer the velocity field on noisy coordinates and cell states, conditioned on a source point cloud during generation (Sec. 3.3 and Sec. 4.3). Thus, $r_t$ is not a residual, but a posterior mean predictor that enables computing the generative velocity on the coordinates. \n\n$r_t$ and $f_t$ are outputs of a point cloud transformer that attends to **all nodes and positions** in both the source and noisy point clouds with cross-attention, and predicts a velocity for each sample. Samples from $\\mathcal{M}$ are passed through the same global function, but this is **not** a shared MLP over individual points, but a **permutation-invariant** point cloud transformer. This architecture computes the velocities of samples and coordinates, integrating the intra-cloud context as well as the conditioning source point cloud (Sec. 4.4 and App. F.6 for details).\n\n> Justification for Tab. 1\n\nThe metrics in Tab. 1 provide a technical evaluation of the **source-conditioned point cloud generation task**, supporting our architectural and modeling choices over baselines. Since we assess how well the global point cloud structure is recovered by a flow model and training strategy, this comparison was performed using a fixed $\\lambda$ to enforce spatial preservation and the same **biological semantics** across models. For the rebuttal, we followed the suggestion and evaluated RPCFlow across varying $\\lambda$ values (see below).\n\n> $\\lambda$-dependent improvement \n\nTo address the reviewer's concern, we conducted an extensive ablation study where we varied $\\lambda \\in \\\\{0.1, 0.25, 0.5, 0.75 \\\\}$ for both RPCFlow and NicheFlow across all three datasets and across all three training objectives. This resulted in 72 evaluations, each reported as the mean ± std. over 5 generation runs.\n\nDue to the limited space, we present only the best-performing method for each metric on each dataset. We will add the full table to the Appendix. When NicheFlow outperforms RPCFlow on all metrics, we report the overall best RPCFlow configuration.\n\n#### MED\n\n| Mod. | Obj. | $\\lambda$ | 1NN-F1 | PSD $(10^2)$ | SPD   $(10^2)$ |\n|---|---|---|---|---|---|\n| NicheFlow | GLVFM | 0.1 | **0.663±0.002** | **0.897±0.021** | **0.824±0.013** |\n| RPCFlow | GLVFM | 0.25 | 0.593±0.001 | 0.926±0.003 | 1.137±0.002 |\n\n#### ABD\n\n| Mod. | Obj. | $\\lambda$ | 1NN-F1 | PSD $(10^2)$ | SPD   $(10^2)$ |\n|---|---|---|---|---|---|\n| NicheFlow | GLVFM | 0.1 | **0.629±0.001** | 2.074±0.006 | **1.235±0.006** |\n| RPCFlow | GLVFM | 0.75 | 0.561±0.001 | **2.015±0.005** | 1.546±0.005 |\n\n#### MBA\n\n| Mod. | Obj. | $\\lambda$ | 1NN-F1 | PSD $(10^2)$ | SPD   $(10^2)$ |\n|---|---|---|---|---|---|\n| NicheFlow | GLVFM | 0.1 | **0.285±0.001** | 1.556±0.001 | 0.694 ± 0.002 |\n| RPCFlow | CFM | 0.75 | 0.276±0.000 | **1.554±0.001** | 0.708±0.002 |\n| RPCFlow | GVFM | 0.5 | 0.245±0.000 | 1.755±0.001| **0.684±0.001** |\n\nWe conclude that RPCFlow benefits from tuning the $\\lambda$ parameter. However, even with its optimal configurations, it does not outperform NicheFlow with different patterns than when $\\lambda=0.1$. Biologically, different values of $\\lambda$ have distinct meanings (see below). \n\n> Biology of $\\lambda$\n\nWhile Tab. 1 evaluates point cloud generation, the choice of $\\lambda$ carries strong biological significance, which we explore throughout the paper. Specifically, $\\lambda$ controls the trade-off between spatial coherence and feature similarity in the OT coupling and **should be selected based on the underlying question** (L318–323):\n\n1. High $\\lambda$: suited for modeling cell type migration across anatomical regions (Fig. 3, right).  \n2. Low $\\lambda$: appropriate for studying how fixed structures evolve compositionally (Fig. 3, left).\n\nSince the optimal value depends on the biological context, we emphasize a qualitative analysis of $\\lambda$’s impact, rather than relying solely on quantitative metrics. We show examples where a high or low value is necessary for meaningful biological interpretation.\n\nNotably, RPCFlow is excluded from our qualitative evaluation, as Sec. D.5 (Fig. 10) shows it fails to condition predictions on the source. Despite strong generative performance, it **ignores the source microenvironment** and yields biologically implausible mappings.\n\n> Comparison with moscot\n\nMoscot solves an **exact OT problem**, i.e. it computes transition probabilities between all source and target point pairs. Since its transition matrix is fixed, the model cannot generalize to test samples. We therefore compare our model to moscot on its native task of mapping embryonic structures across developmental time. Our results show that the proposed spatially constrained generative model offers a promising alternative to this established method. Additionally, App. D.6.2 presents tasks that exact OT models like moscot cannot handle, such as generating niche-specific anatomical changes and compositional remodeling.\n\n> $\\alpha$ in moscot\n\nThe comparison with moscot assesses how well each model captures compositional changes in anatomical structures or migratory patterns. Qualitatively, we found that the spatial component in the FGW problem is comparable to our spatial term across the 0.1–0.9 $\\alpha$ range in the Fig. 3 analysis. When mapping fixed structures over time, values below 0.75 caused excessive density dispersion outside the anatomical region. For migration, no $\\alpha$ value led to generally accurate transitions, though $\\alpha = 0.5$ mapped the highest density to the expected cell type.\n\nFollowing 1m53's suggestion, we report the proportion of source density mapped to the **correct cell type** across $\\alpha$ values, to support our choice. For the spinal cord (Fig. 3 left), values above 0.75 yielded better qualitative and quantitative results. For neural crest cells (Fig. 3 right), $\\alpha = 0.5$ performs best.\n\n| Tissue | α=0.1 | α=0.25 | α=0.5 | α=0.75 | α=0.9 |\n|---|---|---|---|---|---|\n| Spinal Cord | 0.146 | 0.190 | 0.712 | **0.729** | 0.725 |\n| Neural Crests | 0.159 | 0.160 | **0.162** | 0.155 | 0.112 |\n\n> RPCFlow vs. NicheFlow\n\nThe difference between RPCFlow and NicheFlow is that the former is trained on point clouds sampled across the slide, while the latter uses strictly defined spatial neighborhoods. Both methods are evaluated on the same test niches, and across the same model backbones, OT strategies and posterior losses (Gassian vs. Laplacian). Moreover, we showed above that NicheFlow overcomes RPCFlow for different values of $\\lambda$, making our quantitative results robust to it. \n\n> RPCFlow and SPFlow?\n\nRPCFlow refers to a Random Point Cloud Flow model, where point clouds are sampled randomly across the tissue. SPFlow stands for Single Point Flow, which models trajectories at the single-cell level rather than on microenvironments. While described in Sec. 5.1.3, we will clarify them further in the revision.\n\n> Tab. 1 ref.\n\nWe reference Tab. 1 in L299.\n\n[1] Nalisnick et al. Do deep generative models know what they don't know?, ICLR (2019)."}, {"review_id": "lPYNjAEweD", "reviewer": "Reviewer_yXEk", "summary": "This paper introduces NicheFlow, a generative model designed to capture the spatiotemporal evolution of cellular microenvironments in spatial transcriptomics (ST) data. Unlike existing cell-centric approachStes, NicheFlow models local neighborhoods (“niches”) as point clouds of cells, each comprising gene expression and spatial coordinates, and learns their dynamics over time using a combination of: Optimal Transport (define soft couplings between niches across timepoints), Variational Flow Matching (simulate niche trajectories via learned velocity fields), factorized probabilistic model over expression and spatial modalities, and a Transformer-based point cloud architecture. Empirical evaluations on three datasets (mouse embryonic development, axolotl brain development, and mouse brain aging) demonstrate that NicheFlow outperforms both traditional FM-based models and the widely used MOSCOT framework in terms of spatial and semantic accuracy. Qualitative results highlight NicheFlow’s ability to capture known biological differentiation trajectories.", "questions": "- Can $\\lambda$ be learned or adapted across time or tissue regions? This seems essential given its major role in balancing spatial vs transcriptomic cues.\n- How sensitive is the model to the choice of neighborhood radius $r$? Could learned neighborhoods or graph-based neighborhood definitions be more effective? \n- Could this model be extended to continuous-time predictions or interpolation across multiple timepoints?", "limitations": "While NicheFlow models gene expression trajectories over time, it does so directly in the observed (or PCA-projected) expression space. However, gene expression data is known to lie on a lower-dimensional manifold shaped by biological processes like differentiation or zonation. Many successful models in scRNA-seq (e.g., scVI, DPT, diffusion models) explicitly exploit this structure for denoising, interpolation, and trajectory learning. By ignoring the data manifold, NicheFlow may be susceptible to noise and struggle to distinguish technical variability from biological transitions. Incorporating latent space modeling or manifold-aware trajectories could improve robustness and interpretability.", "rating": 6, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 4, "significance": 4, "originality": 4, "strengths_and_weaknesses": "Strengths:\n- Shifts from modeling single-cell trajectories to niche-level dynamics, which is more aligned with tissue-scale processes like differentiation, migration, and patterning.\n- Uses point cloud representations of microenvironments, which are permutation-invariant, variable in size, and expressive of spatial structure.\n- Integrates EOT and source-conditioned VFM, extending prior FM formulations to handle unpaired, noisy, and variable-sized spatiotemporal point clouds.\n\nWeaknesses:\n- NicheFlow is trained on local microenvironments only — there’s no mechanism to enforce global tissue-wide coordination, which is often important in developmental processes governed by morphogen gradients or topological constraints.\n- The trade-off parameter $\\lambda$ (controlling spatial vs gene similarity in the OT cost) is fixed at 0.1. There’s no exploration of $\\lambda$ sensitivity or adaptation, which is important for different biological regimes (e.g., morphogenesis vs migration).\n- PSD and SPD capture spatial proximity but not structural or topological similarity (e.g., niche boundaries, internal configuration, or inter-cell distances). All metrics are pointwise — no graph or structure-aware evaluation is used.\n- While using the microenvironments concept is good, the key part of utilizing this concept is missing -- the cells actually interact with their neighbors. Without considering cell-cell interactions, grouping cells solely based on their spatial location has limited biological significance.", "was_revised": true, "final_justification": "I thank the authors for their complete response. Their responses address my concerns, therefore I would like to increase my ratings from 5 to 6.", "rebuttal": "We thank Reviewer yXEk for the thoughtful review. We appreciate the positive assessment and the insightful suggestions. Below, we address the remaining questions and points of criticism.\n\n> Global tissue coordination. \n\nWhile we recognize the value of global coordination, our work focuses on individual functional regions in a tissue, enabling users to study the migratory or compositional development of custom biological areas over time. This aligns with the locality paradigm in graph and point-cloud methods for cell communication, which we extend to the spatiotemporal setting.\n\nIn future work, we could incorporate global context by conditioning velocity fields on global niche embeddings from long-range graph or point cloud representation models. This allows biasing the model by the influence of the global cellular context. Still, scaling point-cloud methods like NicheFlow to whole-slide data remains a computational challenge, especially at higher resolutions, and is left for future exploration.\n\n> Tradeoff parameter. \n\nWe note that the role of $\\lambda$ is discussed and explored throughout the main text and Appendix (L314–323):\n\n- **Sec. 5.2.** In our comparison with moscot, we test two scenarios:\n  * For studying compositional changes in fixed spatial structures, we use a low value ($\\lambda=0.1$) to preserve spatial location over time.\n  * For modeling migratory processes (e.g., neural crest to head mesenchyme), we use $\\lambda=0.5$, showing that adding cell-state similarity aids fate mapping across space.\n\n- **App. D.3.** The value $\\lambda=0.1$ is used for the **point cloud generation** task in Tab. 1. We also provide:\n  * An ablation across $\\lambda$ values (Tab. 2),\n  * A visualization of its effect on coupling distance (Fig. 7, ref. in L199).\n\n- **App. D.6.1 and D.6.2.** We use $\\lambda=0.5$ for additional biological experiments, such as mapping contained organs over time and inferring anatomical changes during axolotl brain development.\n\n> Pointwise metrics\n\nPSD and SPD are standard metrics in point cloud generation and correspond to the two directions of the Chamfer distance, quantifying fidelity and coverage [1]. Importantly, we compute them over many aggregated microenvironments to reconstruct the entire tissue, so they indirectly reflect global morphology.\n\nTo directly address this concern, we include Gromov–Wasserstein (GW) and Fused GW (FGW) distances [2, 3] between generated and real point clouds to compare different training paradigms for NicheFlow. These metrics compare intra-point-cloud pairwise distance matrices between source and target, capturing internal structure (GW) and gene expression similarity (FGW). Values are averaged over 5 generation runs.\n\nAs computing GW/FGW on entire slides is computationally infeasible due to the number of points, we evaluate them at the generated niche level, where structural comparisons are tractable. A direct comparison with SPFlow or RPCFlow is not appropriate: SPFlow predicts single points and lacks structure, while RPCFlow produces random point clouds that vary across runs.\n\n#### MED\n\n| Mod. | Obj. | GW $(10^2)$ | FGW |\n|---|---|---|---|\n| NicheFlow | CFM | 0.315 ± 0.003 | 3.273 ± 0.003 |\n| NicheFlow | GVFM | 0.463 ± 0.004 | 3.167 ± 0.007 |\n| NicheFlow | GLVFM | **0.224 ± 0.001** | **3.147 ± 0.007** |\n\n#### ABD\n\n| Mod. | Obj.| GW $(10^2)$ | FGW |\n|---|---|---|---|\n| NicheFlow | CFM | 0.783 ± 0.002 | 3.543 ± 0.004 |\n| NicheFlow | GVFM | 0.797 ± 0.006 | **3.414 ± 0.003** |\n| NicheFlow | GLVFM | **0.720 ± 0.004** | 3.420 ± 0.003 |\n\n#### MBA\n\n| Mod. | Obj. | GW $(10^2)$ | FGW |\n|---|---|---|---|\n| NicheFlow | CFM | 0.315 ± 0.001 | 3.531 ± 0.000 |\n| NicheFlow | GVFM | 0.339 ± 0.000 | 3.387 ± 0.000 |\n| NicheFlow | GLVFM | **0.262 ± 0.000** | **3.367 ± 0.000** |\n\nNicheFlow with the GLVFM objective consistently achieves **the lowest or competitive GW and FGW distances**, confirming the value of our novel mixed-factorized posterior at preserving geometric constraints.  \n\n> Cellular interactions \n\nWhile modeling cell-cell interactions is important, spatial proximity alone can reveal meaningful biology, e.g., how brain regions evolve with disease and treatment, or how processes like fibrosis and inflammation alter tissue structure beyond direct cell-cell signaling. Our goal is to provide a framework for modeling combined spatial and molecular shifts as a response to time.\n\nAdditionally, our qualitative analyses (Fig. 3, App. Figs. 11–14) show that NicheFlow generalizes to biologically annotated regions in the mouse embryo and axolotl brain, demonstrating its flexibility to **handle more sophisticated definitions of source clouds** based on labeling and prior biological knowledge.\n\n> Adapt $\\lambda$\n\n$\\lambda$ controls how source and target batches are paired during training, influencing the modeled dynamics.\n\nIn principle, $\\lambda$ can vary over time (e.g., higher in early development to capture large-scale migrations, lower later to emphasize local changes) or space (e.g., focusing long-range cell migration only in specific regions). However, spatial variation requires biological knowledge and annotation to identify areas where such dynamics are expected, making this approach less easily automatable.\n\nLearning $\\lambda$ is technically feasible, but it requires additional supervision. Simply including it in the gradients of Flow Matching lacks biological interpretability. Biologically motivated proxies (e.g., stemness) could guide learning, but at the cost of added complexity. Instead, we treat $\\lambda$ as a biologically informed prior that shapes learning in an interpretable way.\n\n> Sensitivity to $r$\n\n$r$ defines the spatial resolution at which NicheFlow learns tissue dynamics and **constraints its generated outputs to it**. According to our experiments, the choice of $r$ has an influence on the resulting microenvironments, which highlights the importance of selecting a biologically meaningful radius that matches the desired level of granularity (i.e., resolution of the anatomical regions to study). In practice, $r$ can be chosen to emphasize either coarse or fine structures, and it could even **vary across timepoints** to reflect biological contexts.\n\n> Alternative neighborhoods\n\nWe refer the reviewer to our response above on \"Microenvironments and cellular interactions\", where we discuss the flexibility of NicheFlow in handling non-trivial neighborhoods. As noted there, one could define point clouds using alternative logics, such as pruning or weighting edges based on cell–cell interactions, rather than using a fixed spatial radius. While we opted for radius-based neighborhoods due to their interpretability and prevalence in spatial transcriptomics [4], our framework is not limited to this choice and could integrate graph-based or learned groupings in future work.\n\n> Continuous-time/Multi-timepoint  \n\nNicheFlow enables multi-timepoint generation by sequentially predicting future states: we predict $t{+}1$ from $t$, then use the output as a source to generate $t{+}2$, and so on, without relying on ground truth between steps (see App. D.6). This is demonstrated in App. Figs. 11, 13–15 across various developmental processes in mouse (e.g., neural crest cell propagation) and axolotl (e.g., tracking hemisphere and cavity formation).\n\nA true continuous formulation would require aligning the model’s internal flow with biological time, for instance via a Neural ODE–style vector field. This extension is challenging due to point cloud density changes and the issue of incomparable source and target spaces, but it represents an interesting direction for future work. More crucially, such a formulation would require a more densely sampled temporal dataset to approximate a continuum. This is a technological limitation of nowadays' spatiotemporal data acquisition techniques, which focus on a coarse view of a dynamical process. \n\n> Limitations\n\nOur representation of attributed point clouds enables **flexible modeling of continuous cellular states** and can be combined with predefined or learned embeddings. In this work, we focused on datasets with a single sample per time point, where time is the dominant source of variation in both gene expression and spatial organization. We thus avoided scVI or similar approaches that remove sample-specific batch effects, as they risk discarding meaningful temporal signals. When replicates are available, VAE-based models can help remove technical variation, but in our setting, training a VAE would likely capture similar biological structure as PCA, albeit via a non-linear encoder.\n\nWe also highlight the generative nature of NicheFlow: it uses Conditional Flow Matching to transform noise into spatially organized cell neighborhoods, conditioned on earlier time points. The model integrates spatial geometry via a point cloud transformer, capturing both expression and spatial context for generating spatial neighborhoods.\n\nLastly, unlike pseudotime methods (e.g., DPT) that infer latent trajectories from expression data only, our approach models how spatial cell assemblies evolve over physical time. Our method does not attempt to infer latent cell orderings from expression data, but instead focuses on modeling how local tissue organization changes over time.\n\nWhile manifold learning is relevant, spatially constrained versions thereof are an interesting future direction, but beyond the scope of this work.\n\n**Ref**\n\n[1] Fan, Haoqiang, et al.. \"A point set generation network for 3D object reconstruction from a single image.\" CVPR (2017).\n\n[2] Mémoli, Facundo. \"Gromov–Wasserstein distances and the metric approach to object matching.\" Foundations of Computational Mathematics 11.4 (2011): 417–487.\n\n[3] Vayer, Titouan, et al. \"Fused Gromov-Wasserstein Distance for Structured Objects.\" Algorithms (2020)\n\n[4] Haviv, Doron, et al. \"The covariance environment defines cellular niches for spatial inference.\" Nature Biotechnology 43 (2025): 269–280."}, {"review_id": "7ilbbSNINI", "reviewer": "Reviewer_fydQ", "summary": "This paper introduces NichFLOW, a flow-based generative model designed for modeling time-series spatial transcriptomics data. The authors propose incorporating microenvironments, variational flow matching, and conditional generative modeling to learn correspondences between spatially adjacent cell slices in terms of spatial position and gene expression. The method is validated on several datasets. Overall, the idea is interesting, as while static optimal transport (OT) methods have been widely used in this context, flow matching-based methods are less considered.", "questions": "1.  Can the authors clarify whether NichFLOW can perform interpolation (e.g., interpolate time point 1 from time points 0 and 2) or similar tasks? If so, how were these capabilities evaluated, and could additional experiments demonstrate these advantages?  \n \n\n2.  Could the authors provide insights into the key differences or advantages of NichFLOW over other methods to better position its contribution?  \n\n3. The evaluation metrics focus on cell type performance but lack gene-level assessments, such as W1 or W2 distances for gene expression vectors, spatial positions, or combined metrics.  Would adding such metrics provide a clearer understanding of the method’s generative capabilities? Additionally, why were neural ODE-based methods like TrajectoryNet, MIOFlow, DeepRUOT not included as baselines? \n\n4. The method does not account for cell growth, death, which may impact the accuracy of dynamics, especially in developmental datasets. Why was an unbalanced OT (UOT) or unbalanced Fused Gromov-Wasserstein OT framework, such as those in MOSCOT, not considered to model these effects, and could its inclusion improve the method’s applicability to such datasets?\n\nOverall, I view this work positively, and if the authors can address the weaknesses outlined, incorporate additional experiments, baselines (e.g., TrajectoryNet, MIOFlow, DeepRUOT), gene-level metrics (e.g., W1, W2 distances), and comprehensive discussions, I will consider raising my score.", "limitations": "yes", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n\n1. The paper is generally well-written, includes extensive experiments, and presents a comprehensive study.  \n\n2. The work builds on some concepts from prior studies, such as GENOT and VFM, and introduces the concept of microenvironments to handle time-series spatial transcriptomics data. The integration of these ideas in this context is new.\n\n3.   Reconstructing dynamics from time-series spatial transcriptomics data is an important problem and has gained increasing attention.\n\nWeakness:\n\n1. The paper employs flow matching as a conditional generative model, which appears limited in capturing the continuous dynamics of cells across time points; however, this can be achieved by many dynamical-based optimal transport (OT) methods. Specifically, the method seems less suited for interpolation tasks, such as predicting data at time point 1 given data at times 0 and 2, unless relying on simple linear interpolation. A major weakness to me is that NichFLOW appears to achieve functionality similar to static OT methods (such as WaddingOT, Moscot, etc). The paper does not clearly show whether NichFLOW can accomplish tasks that static OT methods cannot, limiting its broader applications in modeling complex spatiotemporal cellular dynamics.\n\n\n2. The paper lacks a comprehensive discussion comparing the proposed method with other works working on the spatio-temporal ST data. Including such a comparison, e.g., with [1,2], could better clarify the similarities and differences between NichFLOW and other methods. \n\n\n3. While the evaluation includes metrics for cell type prediction, it lacks detailed gene-level assessments, such as Wasserstein-1 (W1) and Wasserstein-2 (W2) distances between generated and true gene expression vectors, spatial positions, or their combined distributions. These metrics could better elucidate the method’s generative performance. Additionally, the paper does not include neural ODE-based methods (e.g., TrajectoryNet [3], DeepRUOT[4]) as baselines for comparison. \n\n\n4. The method does not account for cell growth or death, as acknowledged in the limitations. Given that some tested datasets involve developmental processes, neglecting cell proliferation effects may lead to inaccurate dynamics. The paper could benefit from exploring unbalanced OT (UOT) frameworks, such as those used in WaddingOT/MOSCOT, which could be a straightforward way to address this issue.\n\nReference:\n\n[1] DeST-OT: Alignment of Spatiotemporal Transcriptomics Data, Peter Halmos, Xinhao Liu, Julian Gold, Feng Chen, Li Ding, Benjamin J. Raphael\n\n[2] Inferring cell trajectories of spatial transcriptomics via optimal transport analysis, Xunan Shen, Lulu Zuo, Zhongfei Ye, et. al.\n\n[3] TrajectoryNet: A Dynamic Optimal Transport Network for Modeling Cellular Dynamics, Alexander Tong, Jessie Huang, Guy Wolf, David van Dijk, Smita Krishnaswamy\n\n[4] Learning stochastic dynamics from snapshots through regularized unbalanced optimal transport, Zhenyi Zhang, Tiejun Li, Peijie Zhou", "was_revised": true, "final_justification": "I appreciate the authors' thorough and thoughtful response. The further discussions on existing literatures would enhance readers’ understanding of the proposed approaches and their broader implications. I believe the authors have effectively addressed my concerns, and I think this work is an impactful contribution. So I am very happy to raise my score to 5.", "rebuttal": "We thank Reviewer fydQ for their constructive feedback and positive assessment of our method and manuscript. Below, we address the reviewer’s suggestions and remaining concerns point by point.\n\n> Dynamic with conditional model\n\nOur approach uses source-conditioned CFM to approximate a continuous Entropic Optimal Transport (EOT) coupling between source and target distributions [1]. This makes our method a principled predictor for the next time point, **provided the OT assumption holds.**\n\nUnlike previous scRNA-seq methods that evolve individual cells, we model the joint evolution of spatially defined point clouds, treating each microenvironment as a **rigid object**. Since spatial slides can differ in cell density, matching source and target point clouds becomes an **OT problem across incomparable spaces** with different cell counts.\n\nA standard ODE model like the cited ones would evolve each point individually and would be inadequate in the point cloud setting. Instead, we use the continuous EOT coupling of a source-conditioned CFM model to approximate how a source point cloud evolves over time, capturing both compositional and spatial changes in cellular neighborhoods. As shown in **Fig. 13-15** (Appendix), our model generates spatially resolved point clouds from an initial state at $t=0$.\n\n\n> Functional similarity to static OT approaches. \n\nOur method is intrinsically more flexible than moscot and WaddingtonOT, as it provides a generative model approximating a neural Monge map. NicheFlow predicts how a cellular point cloud evolves from $t$ to $t+1$ in both space and gene expression. In contrast, static OT methods only compute the probability that each earlier cell evolves into each later cell, and do not generate later time points.\n\nFor example, see Fig. 13-15 in the Appendix, where we predict the evolution of a source spatial region to a target point cloud across time, using the generated neighbourhood at $t+1$ as a source to predict the dynamics until $t+2$. All the colored points are purely generated. In other words, our model can be used for discrete time point interpolation and cell state generation, which cannot be achieved by static OT.\n\nThe reason for our comparison with moscot is that it is the most established method performing spatiotemporal trajectory predictions across time, and we believe it is important to show the community how our model represents a promising alternative in an applied setting. With an extra page of results, we will add some results on generative interpolations from the Appendix to the main text. \n\n> Comparison with references. \n\nWe thank the reviewer for highlighting these relevant works and will cite them in our revised manuscript. Below, we briefly summarize the key differences between our approach and the suggested models, DeST-OT and SpaTrack. \n\nDeST-OT aligns spatial slides using semi-relaxed optimal transport, incorporating a loss that preserves transcriptomic and spatial proximity between ancestors and descendants. It also introduces a growth vector to estimate the expected number of descendants per source cell. \n\nSpaTrack uses Fused Gromov-Wasserstein OT, balancing transcriptomic and spatial differences based on spatial autocorrelation of features.\n\nOur work is comparable to SpaTrack and DeST-OT in terms of biological settings: they all aim to capture temporal information in time-resolved spatial transcriptomics. However, there are crucial differences:\n\n* Both suggested methods use discrete OT to probabilistically assign each source cell to an existing target cell, similar to moscot. In contrast, our approach parameterizes a continuouos map in the Monge sense (as in [1-4]), predicting how cell distributions (or point clouds, in our case) evolve over time via a flow.\n* Our model trains an approximate coupling in batches, making NicheFlow scalable to large spatial slides and datasets (e.g., the 4.2 million-cell brain ageing dataset we used). Moreover, parameterized deep models also enable future generalization to large treatment studies or cross-species tasks. \n* The main difference is that NicheFlow models the evolution of spatial microenvironments (neighborhoods), not individual cells. It generates how a spatial region changes across slides, while prior methods focus on discovering single-cell descendants from ancestors. By adjusting neighborhood size and shape, NicheFlow enables both localized and large-scale tissue analysis.\n\nIn summary, with NicheFlow we propose two new perspectives on spatio-temporal modeling: the use of a continuous flow-based generative model to learn spatial trajectories and the formulation of the trajectory problem at the level of point clouds.\n\n> W1 and 2 metrics \n\nWe agree that incorporating gene-level and spatial assessments using Wasserstein distances can offer additional insight into the generative capabilities of the models. In response, we have now computed $\\mathcal{W}_1$ and $\\mathcal{W}_2$ distances for both spatial coordinates (“Pos.”) and gene expression vectors (“Genes”) for all models presented in Table 1.\n\nTo ensure scalability, we compute Wasserstein distances per cell type and report the average across types. Shown below are the best-performing configurations for each model on the mouse embryonic development dataset; additional results will be provided in the appendix due to space constraints.\n\n| Model | Objective | $\\mathcal{W}_1$ Pos. | $\\mathcal{W}_1$ Genes | $\\mathcal{W}_2$ Pos. | $\\mathcal{W}_2$ Genes |\n|---|---|---|---|---|---|\n| SPFlow | GLVFM | 0.613 ± 0.001 | 6.457 ± 0.003 | 0.738 ± 0.001 | 5.546 ± 0.108 |\n| RPCFlow | GLVFM | 0.221 ± 0.003 | 6.087 ± 0.002 | 0.365 ± 0.007 | 5.244 ± 0.113 |\n| NicheFlow | GLVFM | **0.212 ± 0.004** | **5.930 ± 0.003** | **0.342 ± 0.007** | **5.031 ± 0.003** |\n\nThese results reinforce our findings: NicheFlow with the GLVFM objective outperforms other models across all metrics, including the new Wasserstein evaluations, confirming its ability to generate realistic spatial and transcriptional patterns.\n\n> Proposed baselines: MIOFlow, TrajectoryNet and DeepRUOT?\n\nThese methods are designed for single-cell transcriptomic trajectory inference and do not model spatial structure or microenvironment-level dynamics. TrajectoryNet, MIOFlow, and DeepRUOT learn continuous trajectories in gene expression space (or PCA reductions thereof) rather than generating spatially resolved point clouds conditioned on tissue context. \n\nOur work addresses a different problem: modeling how spatial neighborhoods of cells evolve jointly in both space and gene expression composition. Adapting these models to **predict coordinates** cell state simultaneously, as well as generating an increased number of cells in time, would alter their original design and put them at a disadvantage. Intuitively, the lack of spatial information in the prediction prevents the computation of the SPD and PSD spatial metrics on the outputs of the proposed models. \n\nWe will clarify these aspects in the revised manuscript and cite the papers.\n\n> The use of unbalanced OT. \n\nWhile unbalanced OT is a promising direction, integrating it into NicheFlow is non-trivial. Our model tracks the evolution of cell groups as point clouds, not individual cells. Unlike single-cell frameworks, where growth reflects mass conservation, growth at the point cloud level is harder to define. Since our OT formulation models how cell neighborhoods evolve jointly in space, it does not rely on explicit growth modeling or enforce mass conservation at the individual-cell level.\n\nInstead, our model captures changes in cell density by using a conditional generative approach to approximate the EOT coupling, allowing mappings between source and target clouds of differing sizes. The number of sampled cells at each time point defines the noisy cloud size during training, enabling the model to learn spatially resolved density shifts, as proliferation in development or expansion in disease. \n\n> Interpolation capabilities.\n\nGreat question. In principle, NicheFlow can interpolate unseen time points. Given a densely sampled developmental process, one could condition generation on a continuous target time and simulate evolution from a source. If the samples sufficiently approximate a continuum, this extends to unseen slides.\n\nUnfortunately, this approach faces key data **limitations:**\n\n- Current spatial datasets are sparsely sampled in time and do not approximate a continuum. Intermediate slides are often **anatomically distinct** from both their predecessors and successors, making interpolation a fully out-of-distribution (OOD) task.\n\n- One could test interpolation using replicates (e.g., training on one embryo and predicting held-out slides from another sampled at the same timepoints), but to our knowledge, such replicated spatiotemporal datasets are not available.\n\n- Crucially, some intermediate slides contain **unique cell types or anatomical structures** not observed at other time points. In spatial transcriptomics, where transcriptional profiles are constrained by spatial context, this can lead to biologically infeasible predictions.\n\nIn other words, while interpolation is algorithmically possible, its biological feasibility depends on denser temporal sampling or the availability of replicates, which are currently lacking in existing datasets.\n\n**Ref**\n\n[1] Klein et al. \"GENOT: Entropic (Gromov) Wasserstein flow matching with applications to single-cell genomics.\" Advances in Neural Information Processing Systems 37 (2024): 103897-103944.\n\n[2] Tong et al. \"Conditional flow matching: Simulation-free dynamic optimal transport.\" TMLR (2023).\n\n[3] Atanackovic et al. \"Meta flow matching: Integrating vector fields on the wasserstein manifold.\" ICLR (2024).\n\n[4] Wang et al. \"Joint Velocity-Growth Flow Matching for Single-Cell Dynamics Modeling.\" arXiv preprint arXiv:2505.13413 (2025)."}, {"review_id": "LBR4lYXzq5", "reviewer": "Reviewer_yBc1", "summary": "In this manuscript the authors propose a flow matching framework to model the temporal evolution of spatially resolved microenvironments from time-resolved spatial transcriptomics data. \n\nThey develop a method denoted as NicheFlow and apply it to several spatiotemporal datasets involving embryonic and brain development.", "questions": "P2 Lines 39-44 Are the methods mentioned here applied to spatio-temporal datasets, or just to spatial datasets ? \n\nP2 Line 61: independently optimizes.... why is it important to perform this optimization independently ? \n\nP2 Line 75: in modeling cellular data [25-28]. How are references [25-27] relevant to this task ? \n\nP3 L90 identify likely descendants.. sentence is unclear. \n\nP4 L150 NicheFlow predicts the corresponding tissue structure... Why is this interesting from a biological / scientific perspective?  \n\nP5 top: it is a bit unclear what dynamics can be learned by looking at only two consecutive time points. It is also unclear if\nthe two time points one measures exactly the same cells (did the cells not move / where they tracked / no cell divisions occurred) ? \n\nP5 Eq. 12 - it seems a bit simplistic to simply average the gene expression of neighborhood cells? What are typical processes\nin the dynamical evolution of niches and are these easily captured by such averaging ? \n\nP8 Which reference is SPFlow ? \n\nP12 are refs 34 and 38 the same?", "limitations": "The authors do not discuss the limitations of their approach in the main text, but they do mention several of them in the supplementary, appendix B.", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 1, "significance": 2, "originality": 3, "strengths_and_weaknesses": "The main strength of the manuscript is the development of NicheFlow, a method to model the temporal evolution not of individual cells,\nbut of niches.  \n\nThat said, the manuscript has several major weaknesses, as I outline below. \n\nClarity of Exposition and Fit to NeurIPS: The manuscript is not very accessible to a non-expert, and not sufficiently self contained. \nTo me it seems that the target audience of this manuscript are researchers in computational biology. \nAs such, in my opinion, the manuscript would be more suitable for publication in such venues rather than in a ML conference. \nAlso, the organization of the paper could be improved. For example, section 3 with a background on flow matching is a bit out of context\ngiven that the precise formulation of what is measured and what are the tasks of interest appear only later. \n\nQuality and Significance: While I am not an expert in this field, to the best of my understanding, the reason for collecting\nspatial-temporal biological datasets is to allow to gain novel insights into temporal biological processes, such as development\nand disease progression. \nIn this manuscript, the authors present their method of modeling the evolution of Niches.  \nHowever, there is no sufficient biological background, no concrete biological questions considered. \nThe manuscript would benefit if it included an in-depth example where new insight is gained using the proposed model \non a concrete biological problem of major importance. \n\nOriginality and Novelty: The model proposed is to the best of my knowledge novel.  \n\nThe evaluation section would benefit from a specific example where the authors would illustrate how NicheFlow allows\ngaining novel insight into some dynamical biological process. \n\nAdditional comments: \nP2 L56 prediction of spatial coordinates and gene expression... It may be good to provide some concrete biological examples\nwhere this would be useful / important. \n\nP2 L66 NicheFlow enables principled learning of dynamics... Again, it may be good present a concrete biological example\nwith specific scientific findings enabled by NicheFlow.", "was_revised": true, "final_justification": "Based on all referee reports, and the author replies, in my view, in its current form (or with small revisions) this manuscript is much more suitable for publication in a computational biology venue. As I am not an expert in this field, I increased my rating to borderline reject.", "rebuttal": "We would like to thank reviewer yBc1 for reviewing our paper. We value the criticism and feedback as an opportunity to improve the clarity of our submission, and hope that our answers will contribute to a more positive evaluation of our work. \n\n> Fit to NeurIPS\n\nWe respectfully disagree with the assessment that our manuscript does not fit NeurIPS due to its application domain. Applications of ML to computational biology have become increasingly representative at major venues, including papers on spatial transcriptomics accepted at NeurIPS, ICML and ICLR (some examples: [1-3]). Moreover, NeurIPS explicitly lists life sciences as one of the key areas under its Call for Papers, reinforcing its commitment to applied domains. \n\nOur work does not simply apply an existing algorithm but introduces a new methodological perspective on generative Flow Matching (FM), motivated by an applied setting. Specifically, we extend Variational Flow Matching (VFM) to approximate a source-conditioned Entropic Optimal Transport (EOT) model and introduce a factorized posterior to model structured data such as attributed point clouds.\n\nWe describe these novel contributions in the introduction (L51-68).\n\n> Background organization\n\nOur work builds on different aspects of FM that have not previously been explored together: EOT couplings across incomparable spaces, and a variational perspective on the generated space.\n\nWe believe it is important to formally introduce these concepts in a general mathematical setting. This structure allows us to tie the ideas together in a principled way later in the paper, and moves from general background to model desiderata and finally to our point-cloud-based formulation.\n\nWe introduce the mathematical formulation of spatial niches after the general background, as presenting it earlier would require framing the entire background in point cloud terms, which is specific to NicheFlow and better suited to the contribution sections.\n\nWe propose adding a new section before subsection 3.1 to introduce how biological trajectories are typically modeled using OT.\n\n> Biological introduction and scope\n\nIn the introduction, we highlight the importance of **considering cells in their environments**, describe static and dynamic spatial transcriptomics, identify key gaps in existing models, and state our research question explicitly (L45–46). This establishes both the biological context and our main goal. We frame the problem as a **general trajectory modeling task** over microenvironments, aiming to develop an algorithm applicable across **diverse biological systems with spatiotemporal axes**. \n\nThis broad scope motivates our submission to NeurIPS and is reflected by our experiments with generation and biological trajectory inference on development and ageing using widely available data from benchmarks [4].\n\nWe will clarify the relevance of our task in the final version with examples from tumour progression and brain development.\n\n> Biological analysis and findings\n\nSec. 5.2 and the Appendix present predicted dynamics and developmental hypotheses enabled by our model:\n\n* Fig. 3, right panel: we predict that immature cell groups emigrate in the upper cranial region **and** that they differentiate into head mesenchymal cells, simultaneously. This demonstrates our model's combined spatial and transcriptomic insight, unlike standard single-cell RNA-seq methods.\n* App. Fig. 11: We draw similar knowledge studying how facial neural crest cells evolve into jaw and teeth structures, this time predicting across multiple timepoints. \n* App. Fig. 13-14: We show dynamic formation of anatomical regions during the axolotl brain development, capturing spatial and compositional changes.\n* App. Fig. 15: We demonstrate our model predicts a decrease in stemness in the left dorsal pallium.\n\n> L56 - gene expression and coordinates\n\nWe refer to our previous answer for the biological use cases presented in our paper. Predicting both the location and gene expression of cells is important, as it provides an understanding of both **tissue composition** and **organization**. Different from standard scRNA-seq that only specifies how gene expression evolves across cell types, here we study how earlier cell states evolve and where they are positioned anatomically, adding a new layer of information to standard trajectory inference (e.g., prediction of both fate and migration of immature cells in Fig. 3, right). Predicting spatial and state transitions is also clinically relevant, where future model applications could, for example, study metastatic progression in cancer over time.\n\n> P2 L39-44: referenced methods\n\nThe models cited in these lines are applied to spatial transcriptomics and infer trajectories at the single-cell level. In this paragraph, we highlight how most current approaches do not consider the trajectories of entire microenvironments of cells, even though the dynamics of tissue restructuring during temporal processes are strongly determined by the coordinated evolution of cells. \n\n> P2 Line 61: optimization\n\nWe thank the reviewer for pointing this out. Our model **jointly** trains on both spatial coordinates and gene expression features, but the loss function is factorized into two components, one for spatial reconstruction and another for biological fidelity, each with its likelihood and distribution (Laplace for coordinates, Gaussian for gene expression).\n\nWe will revise the sentence as follows:\n\n> \"...that jointly trains on spatial and cell state dynamics using a factorized loss, modeling spatial reconstruction and biological fidelity with tailored distributional assumptions.\"\n\n> P2 Line 75\n\nAll these references point to Flow Matching papers that were applied to learning cellular trajectory from earlier to later time points, which is a common benchmark in OT-based generative models. All these methods predict how cell states evolve into one another, but not how they organize in the anatomical space (L75-76).\n\n> P3 Line 90: descendants\n \nExact-solution OT, as implemented in moscot estimates a transition matrix containing the probability of each target cell to be a descendant of each ancestor cell. We propose to make the sentence clearer as follows:\n\n> ... use discrete optimal transport to link cells across time and infer the evolution of cell states from spatially-resolved gene expression.\n\n> P4 L150 tissue structure\n\nPredicting tissue structure from a source microenvironment is biologically significant as tissue organization reflects biological function, with spatial changes indicating migration, proliferation, or reorganization during development or disease.\n\nAdditionally, spatial reorganization has been crucial for critical phenomena such as immune cell infiltration and tumour microenvironment remodeling in cancer [5], fibrosis progression and matrix remodeling [6].\n\nIn this context, our model provides a generative tool to investigate dynamics at the neighborhood level, offering insights that go beyond single-cell trajectories by modeling how spatially localized groups of cells evolve as structured units.\n\n> P5 top: type of dynamics\n\nIn scRNA-seq and its spatial variant, each cell is destroyed during sequencing. Thus, each slide or sample is an independent measurement. This is a standard setting in single-cell transcriptomics and a key reason why OT has emerged as a framework for trajectory inference across timepoints. Once an appropriate optimality criterion is defined, OT can be used to infer spatio-temporal dynamics on independent snapshots. In other words, computational methods are used to reconstruct a trajectory from successive measurements, even though longitudinal tracking of individual cells is not experimentally possible. In our work, each whole slide, as defined in Eq. (8), contains a different number of cells. We will add a sentence to make this aspect clearer in Section 4.1.\n\n> P5 eq 12: Averaging \n\nWe act on the transcriptional and spatial means to match point clouds across batches. We highlight that this defines a criterion to pair up point clouds across batches, but the dynamical model still learns to flow to all points of a target point cloud, attending all the points and coordinates of a source microenvironment. We agree that using the mean in combination with very large point clouds may mask the signal, but in a more localized setting (see App. Fig. 7), other established studies have demonstrated that cellular microenvironments are well-approximated by the statistics of a Gaussian around a central cell [7]. This supports our matching criterion. \n\n> SPFlow \n\nSPFlow is our internal baseline, introduced in Sec. 5.1.3, which models trajectories at the single-cell level rather than the microenvironment level. While not from a published paper, it was designed to highlight the benefits of modeling point cloud dynamics over individual-cell dynamics. We will add a clarifying sentence in Sec. 5.1.3 stating that SPFlow is introduced in this work.\n\n> refs 34 and 38?\n\nThank you for pointing this out. We will correct the typo in our revised version. \n\n> Limitations \n\nDiscussed in Appendix B and referenced at L227. \n\n**Ref.**\n\n[1] Xie et al., \"Spatially resolved gene expression prediction from histology images via bi-modal contrastive learning\", NeurIPS 2023\n\n[2] Oh et al., \"Global Context-aware Representation Learning for Spatially Resolved Transcriptomics\", ICML 2025\n\n[3] Zhu et al., \"Diffusion Generative Modeling for Spatially Resolved Gene Expression Inference from Histology Images\", ICML 2023\n\n[4] Halmos et al., \"DeST-OT: Alignment of Spatiotemporal Transcriptomics Data\", Cell Sys. 2025\n\n[5] Zhao et al., \"Spatial transcriptomics at subspot resolution with BayesSpace\", Nat. Biotech. 2021\n\n[6] Raredon et al., \"Single-cell connectomic analysis of adult mammalian lungs\", Science Advances 2021\n\n[7] Haviv et al., \"The covariance environment defines cellular niches for spatial inference\", Nat. Biotechnology 2025"}]

Author rebuttal:

Meta-review: The reviewers and I have a near-unanimous positive assessment of the paper: modeling the spatiotemporal evolution of ST data is an important setting and the approach performs favorably in comparison to MOSCOT and other FM approaches. It's a relevant ML problem (arising in computational biology) with appropriate baselines. The original paper lacked some quantitative results (and held-out performance), but the results in the rebuttal address this and show favorable results.
Decision: Accept (poster)