Curly Flow Matching for Learning Non-gradient Field Dynamics
Authors: Katarina Petrović
Lazar Atanackovic
Viggo Moro
Kacper Kapuśniak
Ismail Ilkan Ceylan
Michael M. Bronstein
Joey Bose
Alexander Tong
Venue: NeurIPS 2025 poster ()
Keywords: Single-cell
Trajectory Inference
Schrödinger Bridge
Optimal Transport
Flow Matching
Primary area: machine_learning_for_sciences
TL;DR: We learn non-gradient field dynamics by solving Schrödinger Bridge problem with non-zero reference process drift

Abstract: Modeling the transport dynamics of natural processes from population-level observations is a ubiquitous problem in the natural sciences. Such models rely on key assumptions about the underlying process in order to enable faithful learning of governing dynamics that mimic the actual system behavior. The de facto assumption in current approaches relies on the principle of least action that results in gradient field dynamics and leads to trajectories minimizing an energy functional between two probability measures. However, many real-world systems, such as cell cycles in single-cell RNA, are known to exhibit non-gradient, periodic behavior, which fundamentally cannot be captured by current state-of-the-art methods such as flow and bridge matching. In this paper, we introduce Curly Flow Matching (Curly-FM), a novel approach that is capable of learning non-gradient field dynamics by designing and solving a Schrödinger bridge problem with a non-zero drift reference process---in stark contrast to typical zero-drift reference processes---which is constructed using inferred velocities in addition to population snapshot data. We showcase Curly-FM by solving the trajectory inference problems for single cells, computational fluid dynamics, and ocean currents with approximate velocities. We demonstrate that Curly-FM can learn trajectories that better match both the reference process and population marginals. Curly-FM expands flow matching models beyond the modeling of populations and towards the modeling of known periodic behavior in physical systems. Our code repository is accessible
at: https://github.com/kpetrovicc/curly-flow-matching.git

Curly Flow Matching for Learning Non-gradient Field Dynamics
Katarina Petrović' Lazar Atanackovic2:344 Viggo Moro', Kacper Kapušniak', ismail Ilkan S66 Michael Bronstein" 166 Avishek Joey Bose Alexander 6.7.8 University of Oxford, 2Broad Institute of MIT and Harvard, 3 University of Toronto, 4 Vector Institute, STU Wien, AITHYRA, Mila Quebec Al Institute, Université de Montréal
Abstract
Modeling the transport dynamics of natural processes from population-level observations is ubiquitous problem in the natural sciences. Such models rely on key assumptions about the underlying process order to enable faithful learning of governing dynamics that mimic the actual system behavior. The de facto assumption in current approaches relies on the principle least action that results in gradient field dynamics and leads trajectories minimizing energy functional between two probability measures. However, many real-world systems, such as cell cycles in single-cell RNA, known exhibit non-gradient, periodic behavior, which findamentally cannot be captured by current state-of-the-er methods such flow and bridge matching. In this paper, we introduce CURLY FLOW MATCH- ING (CURLY-FM), a novel approach that is capable of learning non- conggradient field dynamics by designing and solving Schrödinger bridge problem with a non-zero drift reference process-in in stark contrast to typical zero-drift reference processes which is constructed using inferred velocities in addition to population snapshot data. We showcase CURLY-FM by solving the trajectory inference problems for single cells, computational fluid dynamics, and ocean currents with approximate velocities. We demonstrate that CURLY-FM can learn trajectories that better match both the reference process and population marginals. CURLY-FM expands flow matching models beyond the modeling of populations and towards the modeling known periodic behavior in physical systems. Our code repository is accessible at: https::////// /ighhttt.b.. com http://ptttctcoc/ceeeeeeeeeeeeet currc/urll-eeeeell ffoo wow-matchinn.. git.
XM
Mix
(a) Asymmetric circles
(b) CURLY-FM
(c) CFM
Figure 1: Particle trajectories generated between samples drawn from asymmetric circles distri- bution at = and t = 1 and respective to underlying reference velocity field fe(x)). Traditional Row-based models such OT-CFM and CFM cannot capture cyclical patterns in physical systems. CURLY-FN FM is capable of learning non con-gradient field dynamics behavior in the underlying data.
"Correspnndenee to: [katarina, petpeiieaeetttorreessss... ox.ac.
'Equal advising
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

Introduction
Understanding the temporal evolution of multi-bod bodyi systems remains central challenge across many applications in life sciences [Lange et al., 2024, Pan and Zhang, 2024], as such systems are often characterized by complex dynamics governing their evolutionary behavior. For example, biological tissues are complex systems that evolve through tissue differentiation characterized by cell divisions and deaths as well morphological changes. Learning evolutionary behavior in such systems can be formulated as a trajectory inference problem [Hashimoto et al., 2016, Lavenant et al., 2024], where the goal is to recover full trajectories of particles given partial and noisy population- snapshots.
The dominant paradigm in solving the trajectory inference problems in scientific applications involves leveraging tools from computational optimal transport (OT) [Peyré and Cuturi, 2019] to learn neural dynamical systems, e.g. NeuralODE [Chen et al., 201 8], such that sampled trajectories under the model optimize a notion of likeliness of being observed [Bunne et al., 2024b]. For instance, in single-cell trajectory inference, such methods broadly follow a pipeline that first infers "optimal" cell trajectories that follow the gradient of some potential function (often termed a Waddington Landscape [Waddington, 1942]), before searching for important regulators of biological process in development [Schiebinger et al., 2019, Shahan et al., 2022] or disease [Tong et al., 2023, Klein et al., 2025]. Despite the ability to produce approximately optimal trajectories w.r.t. the energy landscape, current methods are limited in their ability to model only gradient-fiele dynamics. Consequently, trajectories inferred under the model are realistic and fail to model crucial system dynamics such as periodic behavior that arises in natural systems, e.g., cell cycling [Riba et al., 2022] where periodic behavior is observed at the scale of months, days, hours, and minutes. These behaviors cannot be addressed by current OT-based methods, as gradient field dynamics cannot capture periodic behavior.
Present work. In this paper, we tackle the modeling of systems governed by non- won-gradient field dynamics. We introduce CURLY FLOW MATCHING (CURLY-FM), novel approach for learning non-gradient field dynamics by solving a Schrödinger bridge problem with a designed reference process that induces the learning of periodic behavior. Specifically, we consider reference processes with non-zero drift-in stark contrast to zero drift processes in approaches such as Diffusion Schrödinger Bridges [De Bortoli et al., 2021b, Shi al., 2024]. Such modification elevates the established (entropic optimal) mass transport problem to new class of problems that require matching the reference drift while also transporting mass between time marginals associated with observations. As a result, solutions to this Schrödinger bridge problem are capable of learning non non-gradient field dynamics and exhibit behaviors such as periodicity as found in cell cycling.
In addition to conventional population snapshot data, we design CURLY- by leveraging approx- imate velocity information that is used to construct the drift of a reference process. Consequently, to model periodic dynamics CURLY-FM solves the Schrödinger bridge problem by decomposing it into a two-stage algorithm. The first stage learns a neural path interpolant by regressing against the drift of our constructed reference process. Unlike straight paths in optimal transform conditional flow matching trajectories, the neural path interpolant exhibits cyclic behavior due to the optimization objective of matching a constructed reference drift. In stage two, we learn, in a simulation-free manner, to construct the generative process that solves the mass transport problem as a mixture of conditional bridges built using optimal transport-baase based couplings that minimize the length of the velocity field of the neural path interpolant. The combination of our two-stage approach enables CURLY-FM to learn dynamics that do not affect the population data, but do affect individual particles, which include observed periodic non -raddidient field dynamics, unlike other methods (see table 1).
We instantiate CURLY-FM on modeling a suite of problems in natural sciences that exhibit known non-gradient field dynamics, such as cell cycle systems found in scRNA-seq data with RNA- velocity [L.a Manno et al., 2018], ocean currents, and computational fluid dynamics for PDEs simulated using the Lagrangian particle discretizations [Toshev et al., 2023]. In each case, the system under consideration is under the influence of known reference process, c.g., RNA-velocity gives an approximation of the instantaneous velocity of a cell, which must be adhered when solving the trajectory inference problem. More precisely, using the reference process drift and population snapshots CURLY- FM solves the Schrödinger bridge problem by searching for a bridge that matches the reference process as closely as possible, but also matches the marginal distributions at both timepoints of the system dynamics. Consequently, we show that previous flow matching approaches fail to model this eyclic behavior as they are not able to take advantage of the additional reference process information. Furthermore, while there exist simulation-ban able to learn the correct dynamics [TTng et al., 2020], we show that in practice CURLY-FM performs
methods that are in principle
2

---

significantly better both in terms of accuracy to match the reference crift-enablingg modeling of cyclic behavior (c.f. fig. whhil being computationally cheaper due its simulation-free training nature.
2 Background and preliminaries
Given two distributions po and pt. the distributional matching problem seeks to find push-forward map 16 : IR_ + Rd that transports the initial distribution to the desired endpoint, - [6] (Po). Such a problem setup is pervasive in many areas of machine learning and notably encompasses the standard generative modeling and optimal transport settings [De Bortoli et al., 2021a, Peyré al., 2019]. In this paper, we consider the setting where each distribution Pu and f1 is an empirical distribution that is accessible through a dataset of observations {~6) po(2o) and {ai)
pI (I1). Thus the modeling task is to learn the (approximate) transport map dh.
2.1 Continuous normalizing flows and flow matching
One common choice for modeling & is as a deterministic dynamic system with a timeddepeneent generator vhe : [0, 1] x R - RR The solution to this dynamical system is ordinary differential equation (ODE) and the learned transport map is known as a continuous normalizing flow (CNF). A CNF is time-indexed neural transport map 1/,, for all time [0, 1]. that is trained to push forward samples from prior μo to a desired target μ1. Specifically, a CNF models the ODE d di (x) = ft (vb(()) with initial conditions t/a (20) = To and ft [0, 1] x Rd + R2 being the dimendependent vector field associated with the ODE and transports samples from μio + f1.
The most scalable way to train CNFs is to utilize a simulation- free training objective which regresses a learned neural vector field ( [0, 1] 12d > IR to the desired target vector field (xx,) for all time. This technique is commonly known as flow-matching [Liu, 2022, Albergo and Vanden-Eijndee 2023, Lipman et al., 2023, Tong et al., 2024a] and has the neural transport map thr, which is obtained through a neural differential equation [Chen et al., 2018] de vr, (z) = v1, (Vr.e Specifically, flow- fow-matching regresses 14 (xx) to the target conditional vector field fx(xxa z) associated to the target flow (x We say that this conditional vector field fe(xx(2), generates the target density μn (I1) by interpolating along the probability path fe (2xl2) in time. We often do not have closed-form 11). we
access to the generating marginal vector field fe(z)). [. Still, with conditioning, e.g., can obtain a simple analytic expression of conditional vector field that achieves the same goals. The conditional fow-matching (CFM) objective can then be stated as a simple simulation-free regression,
(CF (r,,(()),, (0,,, - (z)
(I)
The conditioning distribution ((z) can be chosen from any valid coupling, for instance, the independent coupling q(z) = (z(z)) (21). To generate samples and their corresponding log density according to the CNF, we may solve the following flow ODE numerically with initial conditions 2o = (zu) and log (o(z)), which is the log density under the prior:
d (1)) di log (xr))
(t, [- - (t, -V.00,((t
(2)
In the next section, we outline a different methodology to build a transport map leveraging stochastic dynamics. This allows us frame the mass transport problem as Schrödinger bridge, which is well suited to modeling noisy measurements in applications such as single-cell evolution.
3 Schrödinger Bridge with non-zero reference field
The complex nature of particle dynamics can be captured as a mass transport problem under prescribed reference process. Specifically, stochastic differential equation (SDE), with drift vr,0 : [0, 1] x Rdd 4 Rd, diffusion coefficient yi > 0:
model particle evolution using parametrized
dX,, U2 V,,0 (X.) dt + gtd Xo no, X1
11,
where B, is a standard Brownian motion and by convention time t E [0.1] flows from 0 to
I 1= 1 such that marginal distribution at the endpoints are /o and /). These endpoints are provided as empirical distributions and represent endpoint observations along a transport trajectory. The SDE
in eq. (3) induces path measure in the space of Markov path measures (Pr,0) E 10,
3

---

such that the marginal density Pa evolves according to the following Fokker-Plant equation:
dp - (X,,)(X)), (x.)) + d
gi App (Xt), Po = po: P1 = P1-
2
In addition, our modeling of particle dynamics is informed by reference process which is defined by the following SDE with corresponding drift ft : IR + Rd and diffusion coefficient gt >
dX. (X.) dt + g.dB,
(5)
We denote the induced path measure of eq. (5) as (Q0)
Note further that we assume the diffusion coefficient, 9t. to be the same for both processes Pr and Qu to simplify the setting and facilitate easier exposition of the setting considered in this paper.
Schrödinger bridge with zero-drift. We now state the Schrödinger bridge problem, which finds an optimal path measure F that is the solution to the following KL-divergence minimization problem:
P
argmin [KL (PollQ) Po = No. P1
=11
In settings where eq. (5) is zero-drift and with constant diffusion coefficient-i.ce i.. dX, gr B. the Schrödinger bridge problem [Schrödinger, 1932] devolves into the Diffussion Schrödinger Bridge problem [De Bortoli et al., 2021a, Bunne et al., 2023]. In this special case, the Schrödinger bridge problem admits unique solution and is linked to the entropic optimal transport plan through the seminal result of Föllmer [1 [1988]. Specifically, Tp" is a mixture of conditional Brownian bridges Qu (-(oo, 21) weighted by the entropic OT '' E [pp @ p1) which is valid coupling in the product measure Po @ PL. in other words f [ -) = [o(o)) (1,21) = p(æ1).
f Q+ (-120, (1) dr (xo,
pp= =
[cee, x1 (10, + 2L @ pt). =ell)aam)))
argmin
(8)
This holds when gr. = o for some constant a (i.e. 9t. is "fixed"), which assume for the remainder of this work. In this case, operationally, the conditional Brownian bridges take the form of Normal distribu- tion Qu (-|zzo æ1) = (2;; + +1-t)ao- e(1-t)a²) with the mean given as an interpolation between two endpoints. Furthermore, when + O. the entropic OT problem reduces to the regular OT problem. We note that this Schrödinger bridge problem can be reinterpreted as a stochastic optimal control prob- lem where the control cost is the drift vr.0- That is the stochastic optimal control perspective minimizes average kinetic energy' of the learned process which leads to the following optimization problem:
0/ = min Ex.
Ur X, 12 dt : dX, = 0 (X,) dt + ge,B,, = po: P1 = pI)
(9)
We approximate [P* using mini-batch OT [Fatras et al., 2020, 2021] and simulation-free matching algo- rithms [Tong et al., 2024a,b, Pooladian et al., 2023), iterative proportional and Markov fitting [De Bor- toli et al., 2021a, Shi et al., 2024], and generalized Schrödinger bridge matching [Liu et al., 2023a].
3.1 Schrüdinger Bridges with non-zero drift
We now consider the more general case where the drift of the reference process Q is non-zero. In this case, existing computa- tional approaches, which rely gradient field dynamics, no longer apply. As de- tailed in table 1, we consider the case of a constant ge > and tailor our approach to small de previous work has shown this to be preferred empirically [Tong et al., 2024b]. While it is possible also learn ge in GSBM [Liu al., 2023a], this is computationally expensive; in this work, we focus on the setting of fixed or small and vectorfields with Curl, as motivated by applications in our experiments section 4. In this case, we can approximate the
Table 1: Overview of the properties of different approaches.
Method
Pa
Q
Models Curl
9t
DSBM
V
X
fixed
X
OT-CFM
V
limg.-+0
GSBM
<
learned
CURLY- FM (ours)
fixed
Schrödinggr bridges minimize the relative entropy w.r.t. to Q and kinetic energy in the deterministic case.

---

marginal process using a mixture of conditional bridges, albeit not necessarily Brownian. Specifically, (1).o.)))((0... 21), where we use
we consider the case modeling the marginal process P, = Q (- #1) to denote the stochastic bridge pinned at zo, a1 at times 0 and 1 respectively and π is valid coupling. With this decomposition, remains specify the parameterization of Q (- X1) and x (ro. 21). As modeling choice, we parameterize (2 21) using mixture of Brownian bridges with learnable parameters yj:
21,4 (T, To, x1 dt gdd
(10)
with the idea that,
f
=
min
Ep, (()))n(()))) (x2)) : dx, St,n(X,)) g. Qo = Po: Q1 pi
We approximate the mean of Qr by designing neural path interpolant (Pt. with parameters Ty:
Pe, := 1. (a.)(H))ds = tæ + t)=o + t(1 (To
(11)
We optimize Pt.y by minimizing the following simulation-free objective of the relative kinetic energy:
|| @ 0t
ss(uee)) |
[(n)= = [e
[u()),)) ((Men) fo)
Here Ne (He.n σo) ) is any smooth function, e.g. nearest neighbor based distance kernel Ky (H4,,) Too) || =/ll2/ Σi || - 2bll2. Then we use kernel definition of ft based on data uncom- mon to have access to fr for the practical applications we consider. We further discuss our assumptions for Ki in the appendix $G.4. We also note 0ue,,// nn/O can be computed using automatic differentiation:
Opun ((n( = x1 20 + t(1 - t) @pr. (xo, x1) + (1 2t) æ1).
(12)
Øt
8t
The pseudocode for learning the neural path presented algorithm 1. To approximate F, we next learn to approximate the optimal mixture of conditional bridges P7- = (or [Zt,n]] However, this necessitates the feasibility of computing the OT-plan, which is defined below:
n" (Hos 11) argmin / c( x1 (dT()) 21), c(zo: [ (oT)p0)))))
f. 0pun
11 Ot
feluen)
dt
s.t. / T(zo. =oo(ro), / T(,,x1) = p(x1).
Indeed, this optimal transport cost c(2o, 11) can be computed through simulating the entire trajectory. However, we opt for using a stochastic estimator of the cost with K samples:
c(zo, x1) Be
Oærn 06
1
Σ
K
Ozt.n
Se (
(13)
K
at
This cost ensures that we choose coupling which minimizes the total cost in eq. (10). We highlight that the optimal plan x" is intractable and we instead use biased minibatch approximation of the plan (see $G.3). We use the cost in eq. (13) to estimate a transport plan x (Ir, x1) to construct the approximated mixture of conditional bridges P, which is needed to learn the drift (z,) of eq. (3),
(ape.. 01
LRow (0) -E (00=
U2 (Z).)
.ddaaccchh ()
2
The procedure for this marginal (flow) matching objective is presented in algorithm 2. In the that St := a > 0, for constant Ø, we also need to learn the marginal score $1,8 R V log pt (zz). This can be learned using a conditional score matching objective where x 22/1(1 - t)/a and
5

---

Zt.n = Me.n + de t)
with e <~ (O 1) [Tong et al., 2024b],
(eeSL + c12
Lucrr (0) = (
In totality, the combined loss given by thus ((0) = Cnow (0) (
Algorithm Training algorithm for neural path Algorithm 2 Marginal Score and Flow Matching
interpolant network
Require: Marginals (a( and P (T1), network fp
Require: Marginals po (20) and pi(x1), network yen reference drift ft.n
vector field network 24.0-
1:
while Training do
1: while Training do
2:
Sample (Io, 2i, (23) ~ po (o( p1 (O.
2: Sample (Jo, (1, t) (2 (1 1)
3:
Cy (aá, x/) = R.
Ope.n
3: Pe,n (1 ()w0 + Lxi + 1(1 D)
(100
4: 0pi -1- x0-
0 +(1-
4:
2o, I (2o, x1 +- OT x1, On)
(20.
5:
~ ((0, 1), NO,,
5: ((m)= 0pt. fenn(Me.n)
6:
((() = Cn (0) + Cu)
7:
f- Update (0, VeL(O))
6: + Update(n, v,,((n))
return 17g
return ypr.
We highlight that while L(0) seeks to match v1,0 to velocity of the neural 8ps. and the optimal velocity vi + ft. since the reference process Qn does 88t not necessarily transport Po to pt. More precisely, Qn does not have constraints the endpoints, that Qo = Po and Qn = p1. which are required from our learned process Po and its drift vt,0-
Remark
Experiments
We investigate the application of CURLY-FM on multiple applications which exhibit non-gradien field dynamics including a simple toy example, an occan currents modeling application, computational fluid mechanics dataset, and an application to single-cell trajectory inference. We benchmark CURLY-FM against both simulation- free flow matching approaches: Conditional flow matching (CFM) [Liu et al., 2023b, Peluchetti, 2023, Lipman et al., 2023, Albergo et al., 2023], optimal transport conditional flow matching (OT-CFM) [Tong et al., 2024a] and when possible metric flow matching [Kapuániak et al., 2024] which cannot model non-zero drift dynamics, as well as iimulation-based methods in TrajectoryNet and SBIRR [Shen al., 2025] which can model non-zero drift dynamics but are much slower [Tong et al., 2020) and numerically unstable.
We evaluate CURLY-FM using metries both on held out samples (2- (2-Wasserttein (W2)) as well as metrics which directly measure how well the learned drift fo field matches the reference drift (Cosine distance and Le cost). We note that in many cases, it is not possible to match the reference drift exactly as the model is forced to match the marginals.
4.1 Synthetic Experiments
We start our experimental study of learning cyclical patterns from population-level observed popula- tions by considering a synthetic example. We construct source and target distributions on asymmetri- cally arranged circles ( fig. la), each with higher particle population density on one side, Given circular reference velocity field ()) with constant rotational speed, the goal is to learn the velocity-fiee (x,) and trajectories vr. (xx) for E [0, 1]. We find that previous flow matching methods with zero-reference field fi result in straight paths between source and target distributions, thereby failing to capture cycling patterns in the underlying data (see 1b and fig. 1c).
4.2 Modeling Ocean Currents
We model ocean currents in the Gulf of Mexico using a resolution of km of bathymetry data from HYbrid Coordinate Ocean Model (HYCOM), which allows us to obtain a reference field. We follow the data processing pipeline of Shen et al. [2025] and observe 111 particles per time-point (see $D for exact dataset details). We report our quantitative results in table 2 and observe that across the left-out time points, CURLY-FM obtains the best results for the majority of the reported metrics,

---

Data
OT-CFM
Curly
25.20
-95.00
-00.00
26,00
-00.75
20.3
M0.0
10.5
30.5
90.0
90.5
80.5
90.0
90.5
Longitude
Longitude
Longitude
Figure 2: Visualization of ground truth data and vectorfield (left), OT-CFM predicted trajectories (center), and CURLY-FM predictions (right). Curly fits the vortex much better than OT-CFM.
Table 2: Quantitative metrics on left out test timepoints for oceans.
numbers taken from Shen al. [2025]
Metric
Method
t2
44
f6
fs
OT-CFM
0.148 + 0.004
0.227
0.008
0.191 + 0.012
0.250 + 0.018
MFM
0,107 + 0.014
0.056
0.014
0.052 + 0.011
0.070 + 0.021
EMD
Vanilla-SB"
0.270 + 0.058
0.300 + 0.056
0.420 + 0.056
0,410 + 0.048
SBIRR Shen et al. [2025]*
0.073 + 0.020
0.072 + 0.012
0.120 +: 0.029
0.094 + 0.023
CURLY-FM
0.019 + 0.003
0.045 + 0.005
0.027 + 0.001
0.030 + 0.006
OT-CFM
0.229 + 0.004
0.121 0.008
0.034 + 0.002 + 0.001
0.067 + 0.007
Cos. Dist.
MFM
0.179 + 0.010
0.011
0.001
0.004 + 0.002
CURLY-FM
0.231 + 0.004
0.017 + 0.001
0.002 + 0.000
0.002 + 0.000
OT-CFM
0.167 0.014
0.144 0.014
0.095 0.005
0.250 0.023
L2
MFM
0.203 0.011
0.067 0.011
0.101 0.015
0.141 0.018
CURLY-FM
0.151
0.001
0.135 1 0.010
0.178 + 0.017
and also outperforms the previous state-of-the-art SBIRR [Shen et al., 2025] on the EMD metric.
Moreover, we note that CURLY-FM is computationally fast and achieves these results in minutes compared 4hrs for the simulation- based SBIRR. These findings also substantiated in fig. 2, where we see trajectories that look more natural at modeling periodic behavior than OT-CFM.
4.3 Experiments on Single-Cell Data
To show that CURLY-FM is effective in learning dynamic behavior in single-cell data, we leverage two biologically rich datasets consisting of cell cycles in human cell fibroblasts [Riba et al., 2022] and erythroblast development mouse [Pijuan-Sala et al., 2019]. We aim to learn cell state trajectories and development paths considering the respective RNA-velocity fields, providing information about cell cycling, lineage bifurcation, and transcriptional dynamics.
(a) RNA-Velocity Field
(b) Cell Cycles
Figure 3: Ground truth data.
Cell cycle dynamics in human fibroblasts. We study the nature of cell cycling human fibroblasts and reconstruct cyclical patterns in spliced-unspliced RNA space for single genes. We leverage RNA- velocities in figure 3a construct cell state transition paths figure by estimating RNA velocity field between marginals using k-nn algorithm. Further dataset details are included in $E Figure shows learned velocity fields Ve,e (xx)) and trajectories Vr. (x;) between cell cycle distributions at = 0 and f 1. In table 3, show results on the trajectory inference task comparing CURLY-FM to CFM, OT-CFM, and TrajectoryNet. Given the underlying cell cycle process, the aim is to learn circular trajectories resulting from divergence free velocity field. While traditional methods successful in generating end points near ground truth, they fail at learning cyclic patterns, as shown in figures 4f.
Our results show that considering non-zero reference field and velocity inference captures non-gradient dynamics in data. Figure 4a and fig. 4d show learned behavior using CURLY-FM We observe that the trajectory (r,)(x;) inferred with CURLY-FM closely matches expected eycling patterns in the fibroblast dataset, in contrast trajectories inferred using CFM and OT-CFM. This

---

is quantified in table 3, where we can see the cosine distance to the reference field is significantly lower for CURLY-FM.
(a) CURLY- FM (2+)
(b) OT-CFM (2,)
(c) CFM vr,((1)))
Xixt
K1o
Rxt
ĀO
(d) CURLY-FM FM (1zz)
OT-CFM (
(I) CFM (2,)
Figure 4: Vectorfields (top) and trajectory traces (bottom) learned using CURLY-FM (left) OT-CFM (center) CFM (right). CURLY-FM is the only method able learn the cell cycle.
Table 3: Quantitative results for cell cycle trajectory inference task. We report the mean result for a metric with standard deviation over three seeds. CURLY-FN FM performs the best across matching inferred velocity field to the reference process (cosine distance) while maintaining comparable predictive quality.
Datascts
2
d
Algorithm
W2 L
Cos. Dist
W2
Cos. Dist
W22
Cos. Dist J
CFM
0.294 + 0.030
1.065 0.080
0.606 + 0.059
1.001
0.037
1.227
0.013
1.007 + 0.010
OT-CFM
0.248 + 0.030
0.800 + 0.309
0.586 + 0.041
1.008 + 0.039
1.183
0.015
0.978 + 0.125
Trajectory Net
0.531 + 0.021
1.077 + 0.031
0.853 + 0.059
0.979 + 0.064
CURLY-FM (Ours)
1.199 1: 0.177
0.295 0.040 0.930 + 0.024
0.300 1: 0.058
1.261 +: 0.077
0.249 + 0.024
Reconstructing cell differentiation in mouse Erythroid development. Mouse erythroid cells develop in curved trajectory time. We show that FM can adhere this developmental path purely based on velocity data for the first time. Earlier works have used manifold-based penalties follow curved structures. We show that this is no longer necessary with clever usage of velocity information. We observe 9,815 erythroid cells undergoing differentiation and partition the data into three temporal snapshots, withholding the central marginal assess trajectory inference (see $E.2).
We visualize trajectories in figure 5 showing that CURLY-FM clearly follows the devel- opmental pathway of mouse erythroid cells, whereas OT-CFM fails to capture dynamics between marginals. To assess CURLY-FM performance, we measure cosine-distance and L2 norm between learnt and ground truth ve- locities as well as W2 distance between points on left-out marginal. Our quantitative results show that CURLY-FM outperforms OT-CFM reconstructing the underlying RNA-velocity field and cell trajectories in majority of se- lected dimensions. MFM continues to achieve lower W2, indicating stronger adherence to the underlying manifold. Conversely, CURLY-FM attains
Table 4: Erythroid dataset results dimension.
Metrie
OT-CFM
MPM
CURLY-FM (Ours)
Dimension
Cos. Dist
0,146 + 0.001
0.014 = 0.001
0.009 0.000
L2
2.704 + 0.019
1.999 0.014
1663 0.293
W/
0.646 + 0.006
0.269 = 0.004
0.369 0.090
Dimension
Cos. Dist
0.489 + 0.001
0.495 0.001
0.488 0.001
Lg (x10*)
I 1,885 + 0.020
1627 = 0.040
1.721 0.035
6.103 + 0.074
4,855
0.052
6.124 0.027
Dimension = 50
+ 0.0 0.000
0.494
0.000
Lg (x109)
2.215 1 0.022
1.971
0.023
2.045 0.073
WA
7.969 + 0.029
6.727
0.022
7.729 0.046
superior cosine similarity the ground-trutt velocity field, consistent with its objective emphasizing faithful velocity alignment which contributes challenge exactly matching end-point marginals.

---

UMAP2
Thain
UMAP1
UMAP1
UMAP1
Figure 5: Visualization of ground truth data and vectorfield (left), OT-CFM predicted trajectories (center) and CURLY-FM predictions (right). Curly fits the ground truth much better than OT-CFM.
Table 5: Quantitative results for the CFD trajectory inference task. Metrics reported on held-out particles
from the test set for all marginals. Error bars show standard deviation.
Method
Cos. Dist. 4
MSE
Prec. @5 1
Prec. @ 10 1
Prec. @25 f
CFM
0.254 + 0.003
0.085 + 0.002
0.079 + 0.004
0.164 + 0.006
0.337 + 0.016
OT-CFM
0.248 + 0.011
0.095 + 0.001
0.303 + 0.002
0.388 + 0.004
0.496 + 0.001
CURLY-FM
0.189 = 0.027
0.095 + 0.003
0.489 0.010
0,522 + 0.009
0.628 + 0.010
4.4
Experiments on Computational Fluid Mechanies Data
We evaluate CURLY-FM on a particle-based PDE dataset generated by a Lagrangian solver. Unlike grid-based (Eulerian) methods, Lagrangian approaches discretize the fluid as set of particles that move with the flow. These particles evolve under the dynamics of the PDE which provides the particles' positions over time; other quantities of interest, such as velocity or energy, are then computed from these positions. We use data from LagrangeBench [Toshev et al., 2023], specifically the two-dimensional decaying Taylor-Gree Green vortex (2DTGV) dataset(see $F).
In table 5, report quantitative results on left-out particles for each marginal from test set. We evaluate performance using (i) the cosine distance between the learned velocity field and the reference field; (ii) mean squared error (MSE) between the predicted particle positions marginal I and the ground truth positions, using the known coupling (ordering) between particles across marginals; and (iii) precision @ k. measuring how often the predicted position is among the & nearest neighbors of the corresponding ground truth particle. CURLY-FM outperforms baselines in terms of cosine distance and precision @ k while matching or outperforming the baseline methods on MSE. The smaller cosine distance for CURLY-FM shows that CURLY-FM produces velocity fields that better align with the reference field. Equal or lower MSE paired with higher precision @k k shows that CURLY- FM more accurately recovers the true particle coupling and generates more faithful trajectories (c.f. fig. 10).
4.5 Further analysis of CURLY-FM performance
On higher stochasticity. We extend our discussion in section 3.1 on stochasticity levels to consider an ablation where > 0. Despite our work being motivated in the little no stochasticity regime, we demonstrate that considering o > 0 does not impact CURLY-FM efficiency. We find, similar to previous work [Tong et al., 2024b), low values of a perform the best on all metrics. As result, we recommend setting to zero unless some reference o value is known. Therefore, all of our experiments are performed under = St = 0 assumption. We consider ocean currents dataset and find that larger stochasticity monotonically decreases performance on our tasks, thus justifying our choice of for our empirical work.
Table 6: Ablation stochasticity o.
a
Cos. Dist. J
1
W2 1
0.01
0.061 + 0.003
0.141 + 0.009
0.028 + 0,066
0.10
0.062 1 0.002
0,145 1 0.011
0.066 + 0.008
1.00
0.145 1 0.009
0.474 1 0.058
0.871 1: 0.048
On computational efficiency. We further provide the computational cost in wall clock time for TrajectoryNet, SBIRR and CURLY-FM in table 7 (see $H.2 for further baseline comparison). We observe that CURLY-FM completes the Ocean currents problem in minutes with higher accuracy in trajectory and velocity field inference task, while SBIRR and TrajectoryNet are the order of multiple hours and unlike CURLY-FN FM are iimllttionbased.

---

5 Related Work
Flow matching, Flow matching [Lipman al., 2023], also known as rectified flows [Liu, 2022, Liu et al., 2023b] or stochastic interpolants [Albergo and Vanden- Vanden-Fijnden,, 2023, Albergo et al., 2023], has emerged as the default method for training continuous normalizing flow (CNF) models [Chen et al., 2018, Grathwohl et al., 2019]. However, FM can lead to unnatural dynamics less, and therefore many works attempt to derive methods for using minimum energy [Tong et al., 2024a, Pooladian et al., 2023] and more flexible conditional paths [Neklyudov et al., 2024, Kapušniak et al., 2024].
Table 7: Compute cost.
Method
Hours
TrajectoryNet
7.44
SBIRR CURLY-FM (Ours)
4.67
0.06
Schrödinger bridges with deep learning. To tackle the Schrödinger bridge problem in high dimensions many methods propose simulation-base
[De Bortoli et al., 2021b, Chen al., 2022, Koshizuka and Sato, 2023, Liu et al., 2022] and simulation-free [Shi et al., 2024, Tong et al., 2024b, Pooladian and Niles- Weed, 2023, Liu et al., 2023a] set-ups with various additional components incor- porating variable growth rates [Zhang et al., 2025, Pariset et al., 2023, Sha et al., 2024), stochasticity, and manifold structure [Huguet et al., 2022] proposed based on neural ODE and neural SDE [Li et al., 2020, Kidger et al., 2021] frameworks. However, very few methods are able to incorporate approximate velocity data, and either match marginals using simulation [Tong et al., 2020], or do not attempt to match marginals [Qiu et al., 2022]. Finally, Schrödinger bridges with non-zero reference field have also been considered by Bartosh et al. [2024] and concurrently by Bartosh et al. [2025] and Shen et al. [2025]. however, they do not employ a two-stage simulation-fru free approximation as CURLY-FM. We include further details on related work comparison in appendix SB.
RNA-velocity methods on discrete manifolds. common strategy to regularize and interpret RNA-velocity [La Manno et al., 2018, Bergen et al., 2020] is to restrict it to Markov process on a graph of cells representing a discrete manifold or compute higher-level statistics on it [Qiu et al., 2022]. However, these approaches are not equipped to match the marginal cell distribution over time. CURLY-FM can be seen as a method that unites these approaches with m matching approaches.
6 Conclusion
In this work, we introduced CURLY-FM, a method capable of learning non-gradient field dynamics by solving Schrödinger bridge problem with a non-zero reference process drift. In contrast to prior work, CURLY-FM is simulation-free greatly improving numerical stability and efficiency. We showed the utility of this method in learning more accurate dynamics in cell cycle system with known periodic behavior, computational fluid dynamics under Lagrangian solvers, and ocean currents. CURLY-FM opens up the possibility of moving beyond modeling population dynamics with simulation-free training methods and towards reconstructing the underlying governing dynamics [Xing, 2022]. Nevertheless, CURLY-FM is currently limited in its ability discover the underlying dynamics by accurate inference of the reference field, which is an inherently difficult problem, especially over longer timescales. Exciting directions for future work involve additional verification of trajectories through lineage tracing [McKenna and Gagnon, 2019, Wagner and Klein, 2020], and improved modeling across non-stationary populations with the additional incorporation of unbalanced transport or multiomics datatypes [Baysoy et al., 2023].
Acknowledgments
The authors acknowledge funding from UNIQUE, CIFAR, NSERC, Intel, and Samsung. The research was enabled in part by computational resources provided by the Digital Research Alliance of Canada (https://alliancee)
ca), the Province of Ontario. companies sponsoring the Vector Institute ai/partners/> ), Mila (https: /ttss///mla.... quebec), and NVIDIA. KK is supported by the EPSRC CDT in Health Data Science (EP//02428X//I). AJB and LA are partially supported by NSERC Post-doc fellowships. LA is supported by the Eric and Wendy Schmidt Center at the Broad Institute of MIT and Harvard. This research is partially supported by EP- SRC Turing AI World-Leading Research Fellowship No. EP/X040062/1 and EPSRC AI Hub on Mathematical Foun- dations of Intelligence: An "Erlangen Programme" for AI No. EP/Y028872/1. We thank Renato Berlinghieri, author of SBIRR, for valuable discussions and for generously sharing his ocean currents code and data, which made our ocean current experiments possible.
10

---

References
M. S. Albergo and E. Vanden-Fijnden. Building normalizing flows with stochastic interpolants. International Conference Learning Representations (ICLR), 2023.
M. S. Albergo, N. M. Boffi, and E. Canden-Eijnden. Stochastic interpolants: A unifying framework for flows and diffusions. arXiv preprint 2300308797,, 2023.
L. Atanackovic, X. Zhang, B. Amos, M. Blanchette, L. J. Lec, Y. Bengio, A. Tong, and K. Neklyu- dov. Meta flow matching: Integrating vector fields on the wasserstein manifold. International Conference on Learning Representations (ICLR), 2025.
V. Bansal, S. Roy, P. Sarkar, and A. Rinaldo. On the wasserstein convergence and straightness of rectified flow. arXiv preprint, 2025.
G. Bartosh, D. Vetrov, and C. A. Naesseth. Neural flow diffusion models: Learnable forward process for improved diffusion modelling. Neural Information Processing Systems (NeurIPS 2024). 2024.
G. Bartosh, D. Vetrov, and C. A. Nacsseth. Sde matching: Scalable and simulation-free training of latent stochastic differential equations. In International Conference on Machine Learning (ICML), 2025.
A. Baysoy, Z. Bai, R. Satija, and R. Fan. The technological landscape and applications of single-cell multi-omics, Nature Reviews Molecular Cell Biology, 24(10))695-713, 2023.
V. Bergen, M. Lange, S. Peidli, F. A. Wolf, and F. J. Theis. Generalizing ma velocity to transient cell states through dynamical modeling. Nature biotechnology, 38( 1408-11444,, 2020.
C. Bunne, Y.P.. Hsieh, M. Cuturi, and A. Krause. The schrödinger bridge between gaussian measures has a closed form. In International Conference on Artificial Intelligence and Statistics, pages 5802-5833, PMLR, 2023.
C. Bunne, Roohani, Y. Rosen, A. Gupta, X. Zhang, M. Roed, T. Alexandrov, AlQuraishi, P. Brennan, D. B. Burkhardt, A. Califano. J. Cool, A. F. Dernburg, K. Ewing, E. B. Fox, M. Haury, A. E. Herr, E. Horvitz, P. D. Hsu, V. Jain, G. R. Johnson, Kalil, D. R. Kelley, S. O. Kelley, A. Kreshuk, T. Mitchison, S. Otte, J. Shendure, N. J. Sofroniew, F. Theis, C. V. Theodoris,
S. Upadhyayula, M. Valer, B. Wang. E. Xing, S. Yeung-Levy, M. Zitnik, T. Karaletsos, A. Regev,
E. Lundberg, J. Leskovee, and S. R. Quake. How build the virtual cell with artificial intelligence: Priorities and opportunities. 187(25):7045-700 7063, 2024a.
C. Bunne, G. Schiebinger, A. Krause, A. Regev, and M. Cuturi. Optimal transport for single-cell and spatial omics. Nature Reviews Methods Primers, 4(1 4(1):58, 2024b.
R. T. Q. Chen, Y. Rubanova, J. Bettencourt, and D. Duvenaud. Neural ordinary differential equations. Neural Information Processing Systems (NeurIPS), 2018
T. Chen, G.-H. Liu, and E. A. Theodorou. Likelihood training of Schrödinger bridge using forward- backward SDEs theory. In International Conference on Learning Representations (ICLR), 2022.
V. De Bortoli, J. Thornton, J. Heng, and A. Doucet. Diffusson schrödinger bridge with applications to bcore-based generative modeling. Advances in Neural Information Processing Systems, 34: 17695-17709, 2021a.
V. De Bortoli, J. Thorton,, Heng. and A. Doucet. Difffusion Schrödinger bridge with applications to score-based generative modeling. Neural Information Processing Systems (NeurIPS), 2021b.
K. Fatras, Y. Zine, R. Flamary, R. Gribonval, and N. Courty. Learning with minibatch Wasserstein: Asymptotic and gradient properties. Artificial Intelligence and Statistics (AISTATS), 2020.
K. Fatras, Y. Zine. S. Majewski, R. Flamary, R. Gribonval, and N. Courty. Minibatch optimal transport distances; analysis and applications. arXiv preprint 2101.01792 2021.
H. Föllmer. Random fields and diffusion processes. Lect. Notes Math, 1362:101-204 1988.
11

---

M. Gao, C. Qiao, and Y. Huang. Unitvelo: temporally unified ma velocity reinforces single-cell trajectory inference. Nature Communications, 13(1):6586, 2022.
W. Grathwohl, R. T. Q. Chen, J. Bettencourt, Sutskever, and D. Duvenaud. Ffjord: Free-form continuous dynamics for scalable reversible generative models. International Conference Learning Representations (ICLR), 2019.
T. Hashimoto, D. Gifford, and T. Jaakkola. Learning population-level diffusions with generative rnns. In International Conference on Machine Learning (ICML), 2016.
D. Haviv, A.-A. Pooladian, D. Pe'er, and B. Amos. Wasserstein flow matching: Generative modeling over families of distributions. In International Conference on Machine Learning (ICML), 2025.
G. Huguet, D. S. Magruder, A. Tong, O. Fasina, M. Kuchroo, G. Wolf, and S. Krishnaswamy. Manifold interpolating optimal-transport flows for trajectory inference. Neural Information Processing Systems (NeurIPS), 2022.
K. Kapušniak, P. Potaptchik, T. Reu, L. Zhang, A. Tong, M. Bronstein, A. J. Bose, and F. D. Giovanni. Metric flow matching for smooth interpolations the data manifold. Neural Information Processing Systems (NeurIPS), 2024.
P. Kidger, J. Foster, X. Li, H. Oberhauser, and T. Lyons. Neural SDEs as Infinit--Dimensionall GANs. International Conference on Machine Learning, 2021.
D. Klein, G. Palla, M. Lange. M. Klein, Piran, M. Gander, L. Meng-Papaxandie M. Sterr, L. Saber, C. Jing, et al. Mapping cells through time and space with moscot. Nature, pages 1-11, 2025.
T. Koshizuka and Sato. Neural lagrangian schrödinger bridge: Diffusion modeling for population dynamics. Internationcal Conference on Learning Representations (ICI.R), 2023.
G. La Manno, R. Soldatov, A. Zeisel, E. Braun, H. Hochgerner, V. Petukhov, K. Lidschreiber, M. E. Kastriti, P. Lönnerberg, A. Furlan, et al. Rna velocity of single cells. Nature, 560(7719):4944499, 2018.
M. Lange. Z. Piran, M. Klein, B. Spanjaard, D. Klein, J. P. Junker, F. J. Theis, and M. Nitzan. Mapping lineage-traced cells across time points with moslin. Genome Biology, 25(1):277, 2024.
H. Lavenant, S. Zhang, Y.-H. Kim, and G. Schiebinger. Towards a mathematical theory of trajectory inference. The Annals of Applied Probability, 34(1A):428-500 2024.
X. Li, T.-K. L. Wong, R. T. Q. Chen, and D. Duvenaud. Scalable gradients for stochastic differential equations. International Conference on Artificial Intelligence and Statistics, 2020.
Y. Lipman, R. T. Q. Chen, H. Ben-Hamu, M. Nickel, and M. Le. Flow matching for generative modeling. International Conference on Learning Representations (ICLR), 2023.
G.-H. Liu, T. Chen, O. So, and E. A. Theodorou. Deep generalized Schrödinger bridge. In Neural Information Processing Systems (NeurIPS), 2022.
G.-H. Liu, Y. Lipman, M. Nickel, B. Karrer, F. A. Theodorou, and R. T. Chen. Generalized schr odinger bridge matching. arXiv preprint arXiv:2310.02233. 2023a.
G.-H. Liu, Y. Lipman, M. Nickel, B. Karrer, E. A. Theodorou, and Q. Chen. Generalized schröodinger bridge matching. International Conference on Learning Representations (ICLR), 2024.
Q. Liu. Rectified flow: marginal preserving approach optimal transport. arXiv preprint arXiv:2009.14577, 2022.
X. Liu, C. Gong, and Q. Liu. Flow straight and fast: Learning to generate and transfer data with rectified flow. International Conference on Learning Representations (ICLR), 2023b.
A. McKenna and J. A. Gagnon. Recording development with single cell dynamic lineage tracing. Development, 146( (2)d112)eee660 169730, 2019.
12

---

K. R. Moon, J. S. Stanley, D. Burkhardt, D. van Dijk, G. Wolf, and S. Krishnaswamy. Manifold learning-base based methods for analyzing single-cell ma-sequencing data. Current Opinion in Systems Biology, 7:36-46, 2018.
K. Neklyudov, R. Brekelmans, A. Tong, L. Atanackovie, Q. Liu, and A. Makhzani. A computa- tional framework for solving wasserstein lagrangian flows. International Conference on Machine Learning (ICML), 2024.
X. Pan and X. Zhang. Studying temporal dynamics of single cells: expression, lineage and regulatory networks. Biophysical Reviews, 16(1):57-67, 2024.
M. Pariset, Y.-P. Hsieh, C. Bunne, A. Krause, and V. D. Bortoli. Unbalanced diffusion schrödinger bridge. Uncertainty in Artificial Intelligence (UAI), 2023.
S. Peluchetti. Non-denoising forward-time diffusions. arXiv preprint arXiv:2312. 14589, 2023.
G. Peyré and M. Cuturi. Computational optimal transport. Foundations and Trends in Machine Learning, 1 1(6-6))355-6077, 2019.
G. Peyré, M. Cuturi, et al. Computational optimal transport: With applications to data science. Foundations and Trends@ in Machine Learning, 11(5-6):355-601) 2019.
B. Pijuan-Sala, A. Griffiths, C. Guibentif, T. W. Hiscock, W. Jawaid, F. J. Calero-Nieto, C. Mulas, X. Ibarra-Soria, R. C. Tyser, D. L. L. Ho, et al. A single-cell molecular map of mouse gastrulation and early organogenesis. Nature, 566(7745):490-495 2019.
A.AA Pooladian and J. Niles-Weed Plug-in estimation of schrödinger bridges. Internationcal Conference on Learning Representations (ICLR), 2023.
A.-A Pooladian, H. Ben-Hamu, C. Domingo-Fnrich, B. Amos, Y. Lipman, and R. T. Chen. Multi- sample flow matching: Straightening flows with minibatch couplings. International Conference on Learning Representations (ICLR), 2023.
X. Qiu, Y. Zhang, J. D. Martin-Rufinn C. Weng. S. Hosseinzadeh, D. Yang. A. N. Pogson, M. Y. Hein, K. H. J. Min, L. Wang, et al. Mapping transcriptomic vector fields of single cells. Cell, 185 (4):690-711, 2022.
A. Riba, A. Oravecz, M. Durik, S. Jiménez, V. Alunni, M. Cerciat, M. Jung. C. Keime, W. M. Keyes, and N. Molina. Cell cycle gene regulation dynamics revealed by ma velocity and deep-learning Nature communications, 13(1):2865, 2022.
G. Schiebinger. Reconstructing developmental landscapes and trajectories from single-cell data. Current Opinion in Systems Biology, 27:100351, 2021.
G. Schiebinger, J. Shu, M. Tabaka, B. Cleary, V. Subramanian, A. Solomon, J. Gould, S. Liu, S. Lin, P. Berube, L. Lee, J. Chen, J. Brumbaugh, P. Rigollet, K. Hochedlinger, R. Jaenisch, A. Regev, and E. Lander. Optimal-transport analysis single-cell gene expression identifies developmental trajectories in reprogramming. Cell, 176(4):928-94) (018-993.222..... 2019.
E. Schrödinger. Sur la théorie relativiste de l'électron I' l'interprétation de la mécanique quantique. Annales de l'Institut Henri Poincaré, 2(4))269-310, 1932.
Y. Sha, Y. Qiu, P. Zhou, and Q. Nie. Reconstructing growth and dynamic trajectories from single-cell transcriptomics data. Nature Machine Intelligence, 6(1):25 2024.
R. Shahan, C.-W. Hsu, T. M. Nolan, B. J. Cole, I. W. Taylor, L. Greenstrect, S. Zhang, A. Afanassiev, A. H. C. Vlot, G. Schiebinger, et al. A single-cell arabidopsis root atlas reveals developmental trajectories in wild-type and cell identity mutants. Developmental cell, 57(4):543-55) 2022.
Y. Shen, R. Berlinghieri, and T. Broderick. Multi-marginal schrödinger bridges with iterative reference refinement, 2025.
Y. Shi, V. De Bortoli, A. Camphell, and A. Doucet. Diffusion schrödinger bridge matching. Advances in Neural Information Processing Systems (NeurlPS), 2024.
13

---

A. Tong, J. Huang, G. Wolf, D. van Dijk, and S. Krishnaswamy. Trajectorynet: A dynamic optimal transport network for modeling cellular dynamics. International Conference on Machine Learning (ICML), 2020.
A. Tong, M. Kuchroo, S. Gupta, A. Venkat, B. Perez San Juan, L. Rangel, B. Zhu, J. G. Lock, C. Chaffer, and S. Krishnaswamy. Learning transcriptional and regulatory dynamics driving cancer cell plasticity using neural ode-based optimal transport. bioRxiv preprint 2023.
A. Tong, K. Fatras, N. Malkin, G. Huguet, Y. Zhang, J. Rector-Brooks, G. Wolf, and Y. Ben- gio. Improving and generalizing flow-based generative models with minibatch optimal transport. Transactions on Machine Learning Research (TMLR), 2024a.
A. Tong, N. Malkin, K. Fatras, L. Atanackovic, Y. Zhang, G. Huguet, G. Wolf, and Y. Bengio. Simulation- schrödinger bridges via score and flow matching. AISTATS, 2024b.
A. Toshev, G. Galletti, F. Fritz, S. Adami, and N. A. Adams. Lagrangebench: A lagrangian fluid mechanies benchmarking suite. In Neural Information Processing Systems Datasets and Benchmarks Track, 2023.
C. H. Waddington. The epigenotype. Endeavour, 1:18-20, 1942.
D. E. Wagner and A. M. Klein. Lineage tracing meets single-cell omics: opportunities and challenges. Nature Reviews Genetics, 21(7):410- 427, 2020.
Xing. Reconstructing data-driven governing equations for cell phenotypic transitions: integration of data science and systems biology. Physical Biology, 19(6):61001, 2022.
Z. Zhang, T. Li, and P. Zhou. Learning stochastic dynamics from snapshots through regularized unbalanced optimal transport. In International Conference on Representation Learning (ICLR), 2025.
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 main claim is to introduce a novel generative modelling framework for learning non-gradee oragrdddiet field dynamics in natural sciences by designing and solving Schrödinger bridge problem with non-zero reference drift process. Abstract and intro- duction clearly state our contributions, confirming and empirically evaluating these in the remainder of the main text.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in 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 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.
It is fine to include aspirational goals as motivation, long it is clear that these goals are not attainable by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification:
We describe limitations of our work in Section 6.
14

---

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 in 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 on 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 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 or images are taken in low lighting. Or speech-to-text system might not be used reliably provide closed captions for online lectures because it fails 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 faireess
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full of assumptions and a complete (and correct) proof?
Answer: [NA]
Justification: This paper does not provide and include any concrete theoretical results and/or proofs.
Guidelines:
The answer NA means that the paper does include theoretical results. All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or 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 it affects the main claims and/or conclusions the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification: We have disclosed full details of our experimental set-up in Section 4, supple- mentary materials, and throughout the paper.
Guidelines:
The answer NA means that the paper does not include experiments.
15

---

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 a 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 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 are appropriate to the research performed.
also be provided via detailed
While NeurlPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibilityy 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) If the contribution is new model (e.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 instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome deseribe 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 (e.g., registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in the supplemental material?
Answer: [Yes]
Justification: We have provided open access the data and code in the supplementary materials as well as references to the datasets used in the paper which are all publicly available for full reproducibility of our results.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https:/nnips.. cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g.. for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https: / //nips. cc/ppbblc/gguieee httpc//guies/eeeeeeeee Sttps//oooSSummmmeee iisisinnn insiisssooooo onPolicy) for more details.
The authors should provide instructions on data access and preparation, including how to 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 are reproducible, they 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).
16

---

Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental setting/details
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?
Question:
Answer: [Yes]
Justification:
We have provided full training and test details in supplementary materials.
Guidelines:
The answer NA means that the paper does include experiments.
The experimental setting should be presented in the core the paper to a level detail that is necessary 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 I-sigma standard deviation across three chosen seeds across all experimental results in the main text and supplementary materials.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results accompanied by enror bars, confi- dence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper.
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 library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
should be clear whether the error bar is the standard deviation or the standard error of the mean.
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 errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (e.g. negative еттог rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: We have provided full details on used compute resources and model configura- tions in supplementary materials.
Guidelines:
The answer NA means that the paper does not include experiments.
17

---

The paper should indicate the 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 (e.g., preliminary or failed experiments that didn't It make it into the paper).
9. Code ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https: http:s:///euriig nnwrrps....ppp/ cc htc/ppllic/thhiscuuueeeeeeeeeeee
Answer: [Yes]
Justification: We confirm that our research adheres NeurlPS code of ethics. We rely exclusively on publicly accessible datasets with permissible licenses, acknowledge all third- party resources, and uphold fairness, transparency, and reproducibility throughout our model development and evaluation.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics.
The authors should make to preserve anonymity (e.g., if there is 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 include a broader impacts section thoroughly assessing positive and negative societal impacts of our work in supplementary materials.
Guidelines:
The answer NA means that there no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness 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 not 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 improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that 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 negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility ML).
11. Safeguards
Question: Does the paper describe safeguards that have been put in place for responsible release of data models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
18

---

Answer: [NA]
Justification:
Our work does not pose high risks of misuse.
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 to 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: We have provided citations and credits to all authors of publicly used datasets and baselines that are used to empirically validate our work in main text and supplementary materials.
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 (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license dataset.
For existing datasets that are re-packaged, both the original license and the license of derived asset (if has changed) should be provided. this information is not available online, the authors are encouraged to reach out to the asset's creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
No new assets are released in this research.
Guidelines:
The answer NA means that the paper does not release assets.
Researchers should communicate the details the dataset/code/moodll as part 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 anonymized URL or include an anonymized zip file.
19

---

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 participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification; Our work does not include any form of crowdsourcing or research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in 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 data collection, curation, or 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
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 an equivalent approval/review based on the requirements of your country or institution) were obtained?
Question:
Answer: [NA]
Justification: Our work does not include any form of crowdsourcing or research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing 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 we expect authors to 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.
Declaration LLM usage
Question: Does the paper describe the usage of LLMs if it is important, original, 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.
Answer: [NA]
Justification: Core method of our research has not been developed with use of LLMs.
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: //t////eee...//ee/ wee/oonferences.. /LLM) for what should or should be described.
20

---

Appendix
A
Complementary Orthogonal Work on Single Cell
There are many recent works tackling the single-cell trajectory inference problem on single-cell RNA-sequencing data. In this work we focus on incorporating a reference velocity field and learning non-gradient field dynamics. There are many other related works that may be used in conjunction with ideas in this paper. Here we detail some of these complementary works and how they could be combined with ideas in CURLY-FM. Specifically number of areas have been identified as improving how cells are modeled by flow-based networks.
Optimal transport / minimal energy: Cells are modeled by optimal transport over short enough time scales and is therefore desirable almost all applications of single-cell trajec- tory inference [Hashimoto et al., 2016, Schiebinger et al., 2019, Tong et al., 2020, 2024a].
Density or manifold assumptions: Cells are also known to lie on low dimensional manifold within gene space [Moon et al., 2018]. This knowledge has been exploited in a number of works that are both require simulation during training [Tong al., 2020, Huguet et al., 2022] and more recent methods which do not [Neklyudov et al., 2024, Kapušniak et al., 2024]
Unbalanced transport (modeling cell birth and death): By default, flow models assume conservation of mass over time through the continuity equation. Over long time scales this is not good model of population of cells as [Zhang al., 2025]
Stochasticity: Cells move stochastically based unobserved factors. This has led to number of methods that attempt to model cells with stochastic dynamics [Koshizuka and Sato, 2023, Tong al., 2024b, Schiebinger, 2021], and in particular with Schrödinger bridges, as the stochastic extension of dynamic optimal transport. While there have been many works on learning Schrödinger bridges simulation-fri there is little work on efficient learning of general reference drift functions which we tackle here.
Velocity estimates: RNA-velocity [La Manno et al., 2018] exploits particularities about the RNA collection data to measure both older and newer RNA transcripts a single timepoint. This allows the approximation of RNA-velocity of the RNA expression of each cell. First used in Tong et al. [2020], velocity estimates are relatively underutilized in the trajectory inference problem as data is more scarce and more difficult to process. To our knowledge, CURLY-FM is the first method to provide a simulation-free method for incorporating velocity information into learned flow.
the approximate instantaneous change
Distribution conditioning Where most flow-based frameworks model cells a non- interacting point cloud, recent work has also considered flows which include terms for the interaction between cells [Atanackovic et al., 2025, Haviv et al., 2025]. This allows the modeling of more complex interactions between cells which is an extremely prevelent dynamic in real cell systems.
While CURLY-FM focuses the incorporation of optimal transport, stochasticity, and in particular velocity, it does not address problems of unbalanced transport or manifold structure of the single-cell datatype. Future work will incorporate ideas from CURLY-FM in combination with existing ideas on how to model density and unbalanced transport for more expressive, accurate, and useful models of cell dynamics towards virtual cells [Bunne al., 2024a].
B Further Details on Related Work
Iterative Algorithms. Shen et al. [2025] propose an iterative bi-level algorithm, first solving a Schrödinger bridge and estimating the forward-backwar drift given a current guess of the reference drift. Given the solution to step 1, trajectories are simulated to obtain an updated estimate of the reference drift. This bi-level approach is considered since the reference drift belongs to a family of possible reference drifts rather than a single prescribed one. In contrast, CURLY-FM employs two- stage algorithm that is not iterative and also does not search over family of reference drifts. More critically, CURLY-FM is also simulation-free which achieved by making the modeling assumption that the conditional mixture of bridges is modelled Brownian bridges. In practice, this means that
21

---

CURLY-FM is sigmificantly more efficient to train shown in table 15. Finally, note, both Shen et al. [2025] and CURLY-FM can model non gradient field dynamies, but only the latter is simulation-free
Latent SDEs. Bartosh et al. [2024] and Bartosh et al. [2025] consider Latent SDEs for learning stochastic dynamics with simulation-free training without explicitly solving Schrödinger Bridge problem considered by CURLY-FM framework. This means that unlike CURLY-FM, work Latent SDEs does not consider minimizing KL divergence to find an optimal path measure with respect the reference process. Further, SDEs CURLY-FM are not Latent SDEs. Whilst building computational approach to solving the Schrödinger Bridge problem with latent SDEs may be possible, this is an orthogonal research direction our work.
C Algorithmic Details
In this section we detail python code for training and inference for reproducibility.
for range (n_iter))
sample_x0 (batch_size)
x1 sample_x (batch_size)
x0, x1 coupling (x0, x1)
rorch.raud(batchhsttze
type_as (x0)
eps torch. rand (x0)
lambda_t (2* torch.sqrt (t (1-t))) / si sigma
xt, xt_ =
x0, x1, geodesic_modei] model)
= epa sigma drift_model (zt,
11
st score_model (xt, t)
12
flov_loss = torch... mean ((vt xt_dot) ** 2) score _oorelloss__ = torch mean ((lambda_t st + eps) ** 2) loss flow_loss
13
acore_loss
Listing 1: Python implementation of CurlyFM Score and Flow Training algorithm.
torch randn (batch_size size, dim)
for t in torch. linspace (0, 1, 100) [: 11:
drift drift_model (x, t) score model (x, score_model
drift * dt + sigma torch sqrt (dt) torch fandn_like (x)
Listing 2: Python implementation of CurlyFM inference algorithm.
D Ocean Currents Dataset
We assess performance of CURLY-FM on real ocean currents dataset acquired from HYbrid Coordinate Ocean Model (HYCOM) reanalysis released by US Department of Defense.
Dataset. This data consists of real ocean currents measurements in Gulf of Mexico acquired at 1km bathymetry, providing hourly ocean currents velocity fields for the geographic region between 98E and 77B in longitude and I8N 32N in latitude each day since January 1st, 2001.
Experimental Set-up. We leverage experimental set-up presented in Shen et al. [2025], and focus on a specific time point at 17:00 UTC on June 1st, 2024, extracting ocean surface velocity field that contains a vortex. From point near its center, radius 0.05 and evolve them across nine time steps such that At = 0.9, computing velocities by the nearest grid node for 111 observations per time step.
uniformly draw 1,000 initial positions whithin
Ablations. Tables 8 and 9 show ablations using coupling cost in algorithm 2 for trajectory inference between drifter observations ocean currents experiments. We observe that there limited effect in integrating coupling cost over increased numbers of time steps between marginals.
22

---

Table 8: Comparison of CURLY-FM with and without the coupling in algorithm. 2. Without coupling refers to an independent coupling. The results are averaged over three seeds.
Metric
Method
48
CURLY-FM (without coupling)
0.230 0.003
0.017
0.001
0.002
0.000
0.002 + 0.000
Cos. Dist.
CURLY-FM
0.231 0.004
0.017
0.001
0.002
0.000
0.002 = 0.000
CURLY-FM (without coupling)
0.152 + 0.003
0.099 +: 0.001
0.132 + 0.007
0.187 =: 0.012
L2 cost
CURLY-FM
0,151 + 0.004
0.098 4: 0.001
0.135 = 0.010
0.178 = 0.017
Table Comparison of different numbers of times (n) used to compute the coupling cost in algo- rithm 2. The times are equispaced except for n = 1, where they are drawn uniformly at random. The results are averaged over three seeds.
Metric
Method
f2
44
f6
fg
CURLY-FM (n=1)
0.231 + 0.004
0.017 1 0.001
0.002 + 0.000
0.002 + 0.000
Cos. Dist.
CURLY-FM (n=3)
0.230 + 0.002
0.017 + 0.001
0.002 + 0.000
0.002 + 0.000
CURLY-FM (n=5)
0.229 + 0.005
0.018 + 0.001
0.002 + 0.000
0.002 + 0.000
CURLY-FM (0=10)
0.227 + 0.005
0.017 = 0.001
0.002 + 0.000
0.002 + 0.000
CURLY-FM (n=1)
0.151 + 0.004
0.098 1 0.001
0.135 1 0.010
0.178 1 0.017
L2 cost
CURLY-FM (n=3)
0.153 + 0.002
0.101 + 0.002
0.132 1 0.004
0.184 1 0.011
CURLY-FM (n=5)
0.153 + 0.003
0.104 + 0.002
0.131 + 0.001
0.193 + 0.011
CURLY-FM (n=10)
0.150 = 0.004
0.101 = 0.017
= =1100
0.191 = 0.014
E
Single Cell Datasets
F.1
Human Fibroblasts Dataset
We consider the human fibroblasts dataset [Riba et al., 2022] that contains genomic information about 5,367 cells observed across a fibroblast cell cycle. Cell data further contains information about cycling genes, and more specifically, their RNA velocities, which are used to estimate the reference RNA velocity field. Figure 7 shows a distribution cell rotations and phases during a cell cycle process.
Pre-processing. Data is pre-processed by selecting top d variable genes from the data. Fur- ther, we use scvelo package [Bergen et al., 2020] to compute imputed unspliced (Mu) and spliced (Ms) expressions as well as velocity graph. We construct a cell-cell k-nn graph us- ing sCv. pp..pp.... neihhbrbss (adata) in the joint spliced and unspliced expression space using default scvelo hyperparameters. For each gene and cell, we compute averaged fits moment and the centered second moment of spliced and unspliced genes as well as relared attributes. RNA velocities are computed using scv. t1. velocity
fitting the stochastic tran- scriptional dynamics model [La Manno et al., 2018J. Finally, we compute low-dimensional embeddings for velocities using UMAP with scv. tl. velocity_embedeee our visualizations in figure 4d. For selecting top d highly variable genes, we use sc. sc.pp. wwwp.hihhly.vriiabe..g... (adataa, n_top_genes-d)
(adata) for
E.2 Mouse Erythroid
We consider dataset showing mouse gastrulation subset to erythroid lineage [Pijuan-Sala et al., 2019], representing the developmental pathway during which embryonic cells diversify into lineage-specific precursors, evolving into adult organisms. The data consists of 9,815 cells evolving through five lineage stages as shown in Figure 8, and it is available through scvelo package API.
Pre-processing. Data is pre-processed using scvelo and uni tvelo [Gao et al., 2022] package to compute imputed unspliced (Mu) and spliced (Ms) expressions as well as the velocity graph. With unitvelo, construct the cell latent time used to approximate the differentiated time experienced by cells. Mouse erythroid velocity field is computed using unitvelo. run_ sun_model () following instructions from unitvelo documentation.
23

---

(a) Cell Rotations 7
(b) Cell Phases
(a) Cell Stages
(b) RNA Velocity Field
Figure 7: Human Fibroblasts Dataset.
Figure Mouse Erythroid Dataset.
CFM
OT-CFN
Figure 9: Trajectories of particles for CFD for Figure 10: Learned CFD velocity field for CFM different methods.
(left), OT-CFM (middle), and CURLY-FM (right).
Filtering RNA Velocities. To address tail-effects of noisy single-cell data, filter RNA velocities by dowghwiig distant neighbors in the k-NN estimate at o; and injecting small Gaussian noise. We construct new estimate fł by interpolating between k-nn velocity estimate ft and noise such that f7 = why (x+)) ft + (x,) * (0, 0.1). The weight w. (xx) ) E [0, 1] is given by a sigmoid of the distance between the k-NN distance and threshold hyperparameter y, so that larger distances yield larger wi (x:) and thus stronger penalization of distant neighbors.
Computational Fluid Dynamics Dataset
Experimental details for CFD. We conducted our CFD experiments using the two-dimensional decaying Taylor-Green vortex (2DTGV) dataset provided by LagrangeBench [Toshev et al., 2023]. We subsampled 2000 particles and considered five equispaced marginals (snapshot distributions over particle positions). The goal was to perform trajectory inference from unordered population snapshots. Since instantaneous velocity is not directly observed, we constructed the reference drift field using finite differences. This aligns with how derived physical quantities such as velocity or energy-are computed from particle positions in Lagrangian PDE solvers.
We considered a dataset split of (80,, 20%] across the train and test sets, respectively, For the 2000 particles, this resulted in 1600 being used for training and 400 for testing. All other hyperparameters of FURLY-FM were the same as for the other experiments.
F1 Ablations on CFD
Figure 9 illustrates the trajectories of 25 particles under CFM, OT-CFM, and CURLY-FM. Both CFM and OT-CFM tend produce straight, relatively short paths, most notably in the case of OT-CFM, indicating a preference for minimal transport effort. In contrast, CURLY-FM learns longer, more intricate trajectories that better resemble the expected fluid dynamics.
Figure 10 visualizes the velocity fields inferred by each method. While CFM and OT-CFM yield smoother and simpler velocity patterns, CURLY-FM captures richer and more structured field. This complexity reflects closer alignment with the reference field and suggests improved physical fidelity, The resulting velocity field from CURLY-FM more accurately models the underlying dynamics, adhering to both data-driven transport and the governing reference flow.
Tables 10 and 11 show CURLY-FN without coupling and for different number of times used to evaluate the coupling cost, respectively. In particular, table shows that the coupling based on minimizing the kinetic energy only marginally improves performance for CFD experiments. Furthermore, table shows that there little no benefit in using additional times to approximate the coupling cost in algorithm 2.
24

---

Table 10: Comparison of CURLY-FM with and without the coupling algorithm. 2. Without coupling refers to an independent coupling. The results are averaged over five seeds.
Method
Cos. Dist. +
MSE
CURLY-FM without coupling CURLY-FM
0.214 = 0.007
0.092 = 0.007
0.189 = 0.027
0.095 = 0.003
Table 1I: Comparison of different numbers of times (n) used to compute the coupling cost in algorithm 2. The times are equispaced except for n where they are drawn uniformly random. The results are averaged over five seeds.
Method
Cos. Dist. 4
MSE +
CURLY-FM (n=1)
0.189 + 0.027
0.095 +: 0.003
CURLY-FM (n=3)
0.185 = 0.027
0.092 0.003
CURLY-FM (n=5)
0.179 + 0.020
0.091 +: 0.005
CURLY-FM (n=10)
0.182 + 0.029
0.095 + 0.006
G Further Experimental Details
G.1 mingle-marginal Set-up
Pt.n (zo, 21, .1,() .) and (r,0(are; n)) design. Both (t,,, (), and vr,o(1z, 17) are designed as MLP models with layers. We select MI.P dimensions based on the number chosen highly variable genes d. For Pt.. (r.n .1.() f)
we choose din 2 x d and dout = d. For ( 17,, we choose dimensions of din = d. The dataset is split as [80%, 10%, 10%] across training, validation, and test.
Training. All CURLY-FM and baseline experiments are run using Ir 10-4 learning rate and Adam optimizer with default 31, 32, and e values across three seeds and with 1,000 epochs split into 500 epochs train Pr., followed by 500 epochs to train 21.0-
Baselines. TrajectoryNet was with 250 epochs with the Euler integrator with 20 timesteps per timepoint. We use 250 epochs to limit the experimental time and the number function evaluations to roughly 5x that of simulation-free methods. We use a batch size of 256 samples. We Fermand-Prince 4-5 (dopri5) adaptive step size ODE solver to sample trajectories with absolute and relative tolerances of 10-4
Compute. All experiments conducted using a mixture of CPUs and A10 GPUs.
G.2 Multi-marginal set-up
In the multi-marginal setting, we train by randomly sampling interpolation times & within intervals corresponding to each adjacent pair of marginals-foor example these intervals will be [0,1] and [1,2] if our marginals lie at E {0, 1, 2}. For each interval [ta, t4+1), we sample random time t, then compute both the neural interpolant 2; and its derivative Σ.. We also compute global time e 4, + t, if say t 0.5, the global time for the second marginal pair will be f = + 0.5 = 1.5. This effectively ensures that neural interpolant and derivative match the global axis of time. We then connect all the marginals by concatenating neural interpolant, its derivative and global times into single tensors and use these as inputs to train our neural path interpolant in algorithm and the drift in algorithm 2.
PPn, (20,011 1) and (2,, n) design. We design (r,, Ii, () and Ur.o (xt, 11) MLPs as shown in table 12 across all multi-marginal experiments. CFD experiments additionally use four residual connection blocks with no dropout. We select MLP dimensions based on the number of chosen highly variable genes d. For ( æ1,t) we choose d 2 x d and dout = d. For Uz (Zt, m). we choose dimensions of de. = d.
Training. We show training details for multimarginal set-up in table 12 and compute test metrics on left-out marginals. We select = 20 neighbours for ocean currents and CFD experiments and k = 30 neighbors for mouse erythroid experiments compute ground-truth velocities.
Compute. All experiments were conducted using a mixture of CPUs and A10 NVIDIA GPUs.
25

---

Table 12: Overview of model design and training hyperparameters across multimarginal experiments.
Mouse Erythroid
Ocean Currents
CFD
Pt.,, (To. f1,
(tt,, (No, k1,t) t)
(2,)
(t,n((uu 21, 21,t)
U2 (2t)
Channels
256
256
64
64
64
64
Batch size
256
256
64
64
256
256
Epochs
2k
3k
5k
3k
1.5k
1.5k
Learning rate
10-
10-4
10
10
10
10
G.3
Cost of OT plan
The cost of the optimal plan x² is intractable as it would in general require computing costs over stochastic paths. Consequently, we make several approximations to these couplings that enable faster throughput offered in simulation-free training. In particular, we make three approximations the cost c( zo. 21) between two points, and one on the coupling given this cost:
Algorithm has learned path or set of paths for the means of Gaussians in proximity to the optimal (x (r, #1)
We consider low stochasticity a setting, and therefore the cost is close to the distance of the means travel
The integral of the squared length of the curve is approximated empirically using Monte Carlo estimate
After approximating c( 11), ), we use we use mini-batch OT to approximate the entropic OT problem following set-up in Tong et al. [2024b]. Interestingly, in the low stochasticity setting, Tong et al. [2024b] found that mini-batch OT with no stochasticity was empirically a better approximator of " in the Gaussian case. This is because mini-batch transport adds some amount of "entropy" to the plan due to its approximation. We follow this approximation in our setting.
G.4 Further details on ground-truth velocity estimate
On use of kernels. We highlight that in practical settings we investigate, continuous ground truth reference field does not exist. In our applications, we only have access to the reference field at discrete points in space and time that correspond to ground truth data at the different time marginals. Consequently, we use a kernel to build a continuous reference field f((2t.))
Reference drift velocity estimate ft- In our considered setting, we assume access to the ground truth reference drift at samples To Ho and xı ~ μ41 part of the problem setup. Note that in the high-impact application domains we consider, such as trajectory inference single-cell data, we have access to the RNA velocity [Riba et al., 2022, Bergen et al., 2020], which is assumed to be a reasonable estimate (up to a scaling factor) the SDE velocity [Tong et al., 2020J. As we need construct smoothed
to estimate CURLY-FM everywhere in space and time during training, version of the reference drift by using a kernel Kg. This allows us to construct the reference drift f((xx)-using knowledge from the ground truth reference drift in the existing datasct-in places where there are no ground truth samples. Consequently, we take fi in these intermediate points as ground truth reference drift. We further highlight that kernel estimate is not used in cases where ground-truth velocities are given on a continuous domain.
Ablation on kernel estimate accuracy. To show that estimate accuracy does not effect our findings, we conduct an ablation study by constructing noisy reference drift fruise for various noise levels 3 E [0, 1]. We show results in table 13 for the Ocean Currents dataset. The noisy reference drift is obtained as a linear combination of the ground truth reference drift f. and noise from standard gaussian distribution (fmonse = (1- 3) * f ). We find that the performance between / = 0 (no noise, ie. regular CURLY-FM) and 8 0.25 are similar whilst the performance for 8 0.5 and higher gradually becomes worse, as the noise dominates over the ground truth reference drift in fmoise This shows that CurlyFM is robust to moderate amounts of noise added to the ground truth reference drift.
Reconstructed field. We provide comparison to ground truth velocity field and k-nn estimate for the human fibroblasts data. From figure 12 clear that our approach faithfully reconstructs ground-truth RNA velocities.
26

---

Table 13: Noisy reference drift ablation.
Cos. Dist J
L/
W2 J
0.00
0.062 0003
0.143 + 0,010
0.034
0.006
0.25
0.057 0.03
0.021
+ 0.036
0.051 0,030
0.50
0.087 0.047
0.301
+ 0.085
0.091 0.046
0.75
0:261 + 0.123
0.381 + 0.120
0.145 + 0.062
1
(a) Ground truth velocity
(b) k-nn velocity estimate
Figure 12: Ablation and estimate of ground truth velocity.
G.5 On the tradeoff between directly computing OT plan iterative refinement
Two main strategies exist for approximating the static (possibly regularized) optimal transport coupling. In this work we primarily use the mini-batch approximation. However, there are also iterative-refinemeent type approaches that unbiased in the infinite limit like those presented in De Bortoli et al. [2021b], Shi et al. [2024]. These methods work by creating an "outer loop" where the bridge simulated in direction, then matched in the other create iterative refinement to match marginals. While this is possible to incorporate into our framework, it is not as suitable for our application domain, where we assume small stochasticity level. Iterative approaches are not suitable for small stochasticity levels because the number of iterations of iterative proportional fitting or iterative Markov fitting (IPF IMF) for suitable convergence depends inversely proportional to the stochasticity [Tong et al., 2024b, Shi et al., 2024]. Indeed, at the zero noise limit, the IMF approach collapses to a rectified flow, which only solve the OT problem in limited domains e.g. ID [Liu et al., 2023b], and a few select Gaussian settings [Bansal et al., 2025]. For low noise levels the convergence rate significantly slower.
H Supplementary Results
H.1 Additional Synthetic Experiments
CURLY-FM robustness. To further assess CURLY- robustness, we designed a toy experiment with analytical Schrödinger Bridge solution and non-gradient dynamics based on [Tong et al., 2024b, De Bortoli et al., 2021b, Shi et al., 2024]. We bridge two Gaussians in the presence of spiral reference field across various dimensions and stochasticity levels Specifically, we initialize two Gaussians in 20 dimensions centered at μuu [-0.1, 0, O.... 0] and J4 = [0.1, 0, .,....0, with standard deviations do σi , 1]. We also define a ground truth transport field, which unlike the standard OT field, has an additional rotational component. We note that this has equivalent marginal probability distributions, but also rotates around the origin in the second and third dimensions. This allows us to test how well our method works in simple toy setting.
9t 0.
(a) Source target distribution
(b) Velocity field
Figure 13: Example of chosen synthetic setting in case of 3-dimensions.
27

---

Table 14: CURLY-FM ablation in synthetic setting and comparison to OT-CFM, SF2M and SB-CFM.
(a) CURLY-FM
(b) OT-CFM
Dim
K1-(p1 9x)
Mean KI.
Wa(Ps: 9t.
Dim
K1-(p) 9x)
Dist.
Wa(Ps 9.)
0.043
0.044
0.001
0.192
0.033
0,044
0.909
0.369
0.040
0.017
0.038
0.60
1.0
0.049
0.049
0.519
0.633
1.0
1.000
0.879
0.1
0.041
0.018
0.018
0.880
0.1
1.0
0.065
0.038
0,522
0.887
1.0
20
0.00
0.258
0.178
0.022
3,887
20
0.0
0.201
0.165
0.997
3.921
20
0.1
0.250
0.169
0.029
3
20
0.1
(c) SF2M
(d) SB-CFM
KL[P1 91)
KL.
Wa(Ps 9..
KL 91)
Dist.
Wy(Ps 91.
0-480)
0.999
0.024
0.038
0.509
0.003
1.0
0.028
0.015
0.998
0.418
D.1
0.038
0.040
0.516
1.118
0.1
0.026
0.012
1.000
0,002
1.0
0.108
0,058
0,518
1,127
1.0
0.024
0.017
0.997
0,878
20
0.1
20
0.0
20
0.1
0.537
0.481
0.515
4.194
30
0.1
0.231
0.156
0.995
3,850
20
4.347
0.165
Table 15: Computational efficiency comparison
Method
DM-SB
Vanilla-SB
TrajectoryNet
SBIRR
CURLY-FM (ours)
Hours
15.44
0.43
7.44
4.67
0.06
def pei_xt(xt):
rot_speed np
velocity torch. stack
(mu_1 mu_0) torch. oneo ones_like (xt [..., 0]),
rot_speed xt []].
rot_speed -xt[....1]
1, dim=-1)
xt. shape [-1] 3:
extra_zeros = torch, zeros_like [..., 3:] velocity torch. cat ([velocity, extra_zeros], dim=-1) return velocity
Listing 3: Python implementation of ground truth vector field.
We compare against the closed-form solution and OT-CFM, measuring KI. divergence, Wasserstein distance, and the cosine distance to the ideal rotational angle. Results table confirm CURLY-FM performs best in low-dimensional, low-stochasticity settings, achieving high cosine similarity.
Comparison to baselines in synthetic setting. We further compare CURLY-FM robustness in syn- thetic setting to SF2M and SB-CFM baselines. We find that CURLY-FM outperforms both SF2M and 0.1 on cosine distance. In case of a 0.1, CURLY-FM outperforms SB-CFM, and achieves similar cosine distance with SF2M while providing better W/ metrics.
SB-CFM stochasticity
H.2
Computational Efficiency
We provide comparison across DM-SB, Vanilla-SB, TrajectoryNet and SBIRR baselines in terms of computational efficiency evaluated on the ocean currents experiments (shown in table 15). We further report results for the 2D and 10D settings for TrajectoryNet and CURLY-FM FM in table 16. We find that TrajectoryNet considerably slower than CURLY-FM and is not scalable approach. This is evident in that training TrajectoryNet takes, average, 11 times longer than CURLY-FM in the 2D case, and 17 times longer in the 10D case. Moreover, as dimensionality increases from 2D to 10D, CURLY-FN FM incurs increase in computational cost by 1.7x, whereas TrajectoryNet incurs
28

---

Table 16: Computational efficiency comparison
Method
d d= (seconds)
d = 10 (seconds)
TrajectoryNet CURLY-FM
17005 + 110
*43080 1:65
1429 = 31
2471 #10
Table 17: GSBM Control Points and Efficiency
Method
Cos. Dist. 4
L2 |
W2 4
Train (s)
GSBM 2
0,279 0.006
0.395 + 0.008
0.337 + 0.009
37
GSBM 10
0.158 = 0.166
0.130 + 0.060
0.247 + 0.144
68
GSBM 15
0.075 + 0.017
0.134 + 0.030
0.248 + 0.045
68
GSBM 30
0.088 + 0.034
0.077 + 0.023
0.225 + 0.085
68
GSBM 60
0.269 0.126
0.200 + 0.124
0.400 + 0.065
68
CURLY-FM
0.062
0.143
0.034
176
at least a 2.5x increase in computational cost. We use (*) to denote runs that did not finish within allotted resource allocation time.
H.3
Further Baseline Comparisons
H.3.1 Generalized Schrödinger Bridge Matching
In this section we provide a more extensive comparison of CURLY-FM to Generalized Schrödinger Bridge Matching (GSBM) [Liu et al., 2024].
Experimental Set-up. We make two modifications to improve GSBM in our setting for the fairest comparison CURLY-FM: We use our loss instead of GSBM loss in equation 6a control splines and to incorporate the signal from reference drift directly. We fix to the target a as we assume a constant o and this is the optimum for each bridge. We keep the iterative algorithm as found in the original GSBM paper. We follow the formulation in Tong et al. [2024b] which separately learns the deterministic flow and stochastic score functions, which then can be added together to calculate the stochastic drift or used separately to integrate deterministically. We note that we do this so that we can experiment with different stochasticity and simulate forward or backwards in time.
Number of control points. We initially set the number of control points to We provide an additional ablation over the number of control points in GSBM, noting that the default number in GSBM is 15. We also find 15 provides the best tradeoff. We find that more or fewer control points do not improve performance. Otherwise we keep the same hyperparameters as CURLY-FM in terms of iterations, kernels, batch size, learning rate, models, etc. for fair comparison. Table 17 compares CURLY- FM to GSBM on the Ocean currents dataset with varying number of control points
On loss used in baseline comparison. We clarify that the loss (6a) in GSBM cannot take into account non-zero reference drift directly. Loss 6a only considers potential functions (.x)) and the kinetic energy of ue, and therefore is theoretically problematic in our setting. We therefore use CURLY-FM loss function, which we believe is more fair comparison. For completeness we further include GSBM baseline with loss 6a found in the GSBM paper. We find that CURLY-FN FMM is able to outperform both variations of GSBM on our oceans dataset across the three main quantitative metries of interest. We further find our modification of GSBM's loss allows it to better match the non-zero reference drift as evidenced by a lower cosine distance, while also performing slightly better in W2 distance.
Table 18: GSBM loss comparison
Method
Cos. Dist. 1
L2 +
W2 1
GSBM (our loss)
0.075 + 0,017
0.134 + 0.030
0.248 + 0.045
GSBM (6a)
0.083 = 0.037
0.134 0.029
0.201 - 0.058
CURLY-FM
0.070 = 0.001
0.107 + 0.003
0.052 - 0.004
29

---

H.3.2 Metric Flow Matching
For completeness of our analysis, we additionally report results for Metric Flow Matching (MFM) [Ka- pušniak et al., 2024]. MFM introduces geometric bias by enforcing interpolations that remain close to the underlying data manifold, effectively learning smooth geodesic paths that reflect intrinsic geometric structure. The two-stage learning strategy in CURLY- FM is directly adapted from MFM replacing the manifold-constrained interpolations with regression against the reference non-gradient dynamics. In other words, CURLY-FM can be viewed introducing an alternative inductive bias on trajectories: whereas MFM constrains paths to lie on the manifold defined by data, CURLY-FM enforces a bias on velocities, aiming to learn aeference-consisten vector fields.
Since the core algorithmic structure of CURLY- FM (Algorithm 1) mirrors that of MFM, the two for- mulations are fact compatible and can be combined-leveraging manifold-constrained interpolants together with velocity-bal basee biases. We leave this promising direction for future work.
30


Reviews:
[{"review_id": "x9hmbF7L34", "reviewer": "Reviewer_Tr4n", "summary": "The paper introduces Curly Flow Matching (Curly-FM), a simulation-free method for learning non-gradient, possibly periodic dynamics from snapshot data and approximate velocity measurements.  The authors cast trajectory inference as a Schrödinger-bridge problem whose reference process has non-zero drift built from the given velocities.  Because (most of the) existing diffusion/flow-matching techniques assume zero drift, the paper proposes a new method to account for non-zero drift via a two-stage approximation: (1) find a neural path interpolant such that its instantaneous velocity matches the reference drift, and (2) marginal flow matching to learn a drift that marginally modifies the neural path interpolant to make sure that the induced stochastic differential equation satisfies the marginal constraints at the initial and final times. The authors show promising results on empirical and real data experiments with rotational flows, improving over state-of-the-art baselines.", "questions": "In this section, I will get in more details of the comments provided in the strengths and weaknesses section above, as well as provide some additional minor comments and questions. \n\n**Major Comments and Questions**\n\n- When I first encounter $\\eta$, it is not clear what this represents and why we need it. I think it would be great if the notation is setup and clearly described why we need this before using it. Similarly (and possibly relatedly) I am not sure I follow how $Q_{\\eta}$ (line 158) is actually constructed. I get it that is parametrized as a mixture of Brownian bridges, but above you are saying that these bridges need not be Brownian (Line 154). Is this making a modeling/approximation choice? I think it would be much clearer if the mixture were explicitly written out, with each component described both formally and intuitively.\n- After line 159, $f*_{\\eta}$ seems to correspond to the minimum value of the integral. Is that correct? Or is it instead meant to refer to one of the $f_{t,\\eta}$ terms? I am also not sure that the star notation is used anywhere later in the paper, so I don't follow why we need this.\n- What is the role of $\\kappa$, and how is it chosen in practice? Is it the same for each experiment? Why is it introduced below line 162? Similarly, what is the precise function of the neural path interpolant? How does it relate to model drift and the approximated velocity field? My understanding is that with the neural path interpolant I am trying to make the path $x_{t,\\eta}$ curly, and I am doing it by finding the best interpolant that matches a reference velocity field. And then once I found the interpolant I can use it in the derivative of the path to approximate the drift in the original SDE (equation (3)) with the loss after Line 172.\n- Overall, I find Section 3.1 quite difficult to follow, as it lacks both formal clarity and intuitive explanation. A clear breakdown of the two main steps, along with an explanation of what each quantity represents (both formally and intuitively), would greatly improve readability. I understand the high-level two-step procedure from the introduction, but I find it challenging to connect that to the formalism presented in this section.\n- More broadly, I think I am confused by how the reference velocity field is obtained. Is it computed by deriving the velocity along the stochastic path interpolant and then fitting a neural network to that data, using the loss after equation 162? Or are the velocities approximated using finite differences outside the optimization loop and then used in a regression setup? Is approximating the velocity field necessary at all? Could the method be adapted to scenarios where the true velocity field is known? And how does the accuracy of this approximation affect the performance or theoretical guarantees of the method? I know these are a lot of questions, but I am asking because I think this is the crucial innovation of the method and I am wondering whether it can be spelled out better. I am also asking because for most of these data there is no notion of observed trajectory for a particle, since the particles are unpaired, so I am not sure if this approach is impacted by this or not.\n- Lines 265–267: What is meant here? Could you provide a citation or further explanation? This relates to my broader confusion about the role of the velocity field (as noted above). It’s also not clear whether the method assumes access to RNA velocity fields for this experiment, as seems implied in Line 230.\n- How does the method generalize to more than two time points? In the ocean current experiment, there are more than two time points.  But in the method part the discussion is focused on just two end points at time 0 and 1. But from the experiment results it doesn't look like you're learning the drift piecewise between each pair of points (since the connections look smooth enough), so what is exactly going on? Does the need to approximate multiple velocity fields introduce cumulative bias?\n\n---\n\n**Minor Comments and Questions**\n\n- I really like the inclusion of Figure 1 on the first page! Could you increase the font size for the legend labels? Also, it would be helpful to spell out “OT-CFM” and “CFM” in the figure caption, since the acronyms are not introduced in the main text at that point. Many readers will look at the figure together with the abstract and will not yet know what those terms mean.\n- Line 43: I suggest adding a “Most of” qualifier here, since, as noted above and in the extended related work section, there are existing methods that can model non-gradient dynamics.\n- Line 68: The meaning of this sentence is unclear — could you rephrase or clarify?\n- Lines 140–141: A citation is needed.\n- Equation (9): Unless this particular characterization of the Schrödinger Bridge problem is used later in the paper, it may not be necessary to include it here.\n- Line 146: The phrase “We approximate...” makes it sound like these are the algorithms used in the current paper, which is not true. A more general phrasing would be preferable.\n- Lines 152–153: Is there a double negation issue here? It may be worth rephrasing for clarity.\n- Line 170: Is Equation (10) the correct reference? If so, what exactly do you mean by “total cost” in this context?\n- It would be helpful to explain how boldface is used in the tables. Currently, it seems that only the proposed method is bolded when it achieves the best performance, while other methods are not highlighted even when they perform best. This could be misleading unless explicitly stated.\n- Not a big deal, but note that Figure 5 currently appears before Figures 3 and 4.\n\nMy score is mostly driven by the points raised in the main comments and weaknesses sections. I will be very happy to revise my score if the authors can clarify my confusion and doubts.", "limitations": "Yes, the authors discuss some limitations of their work in the discussion section. In particular, they briefly note that the accuracy of the inference process depends on how well the reference family can be reconstructed. As mentioned above, I believe it would be very valuable to explore this more formally — for example, by investigating whether certain assumptions on the data or on the quality of the derivative approximation could guarantee a desired level of inference accuracy.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 2, "originality": 2, "strengths_and_weaknesses": "**Quality.** The submission is generally technically sound. It does not contain any formal propositions or main theoretical results, but most of the claims are well supported by the experimental findings. The authors do a good job discussing important background ideas, although in some parts it would be helpful to include either more intuitive explanations or more precise references (see the *Questions* section below for details).\n\n**Clarity.** I find most of the paper to be clearly written. I appreciate that the authors clearly motivate their work and set up the problem before presenting potential solutions. However, the paper becomes harder to follow in Section 3.1, where additional details and intuition would be very helpful for the reader to fully understand the method and be able to reproduce it. I will provide more specific feedback in the *Questions* section, but in my view, this is the first major weakness of the paper.\n\n**Significance.** The results have the potential to be impactful for the community. I believe that some researchers will find the experimental contributions appealing and may be motivated to pursue theoretical developments based on them. Additionally, practitioners working with similar data may be interested in applying this method. That said, the paper would benefit from an analysis of how crucial it is to accurately approximate — or have access to — the velocity field in order to obtain good results, since this data may not be readily available in all practical settings. More broadly, I think the paper should include a more thorough discussion of the role of the velocity field in the proposed approach.\n\n**Originality.** The paper presents a new method that clearly differs from most of the prior works cited. However, the authors repeatedly claim (e.g., Lines 43 and 151) that existing methods can only model gradient field dynamics, which is not accurate in general. For example, the work by Shen et al. (2024), which the authors use in one of their experiments, supports more general dynamics, The two works by Bartosh et al. (2024 and 2025) are also relevant in this context—they are briefly mentioned in Section 5 (line 284), but the authors do not provide any intuition about how these methods differ from their own (beyond a short comment about the “simulation-free” aspect) and do not compare against them experimentally. I view this as the second major weakness of the paper. If the proposed method is very similar to Bartosh’s work and its main advantage is being simulation-free, then I would expect some discussion or analysis of computational speed or efficiency. Without such analysis, it is difficult to justify a new method that does not clearly differentiate itself from recent work on a very similar task. Overall, I believe the work is original and constitutes an interesting contribution, but the discussion of related work should be expanded and the motivation for the proposed approach strengthened.", "was_revised": true, "final_justification": "As explained in my reviews, I found the paper an interesting submission tackling an important problem. I had some initial concerns (raised in my review) and I appreciated the clarifications and the authors’ efforts to address them. I also appreciated the follow-up comment. For these reasons, I increased my score.", "rebuttal": "We thank the reviewer for their time and effort in reviewing our paper. We appreciate that the reviewer found our paper “technically sound” and that the majority of our claims “are well supported by the experimental findings,” which have the potential to be impactful to the community. We are also heartened to hear that the reviewer finds our paper to be “clearly written” and that the work is clearly motivated. Finally, we thank the reviewer for acknowledging that CurlyFM “presents a new method that clearly differs from most of the prior works”. We now address the key points of clarification raised in the review, grouped by theme.\n\n## Questions regarding the clarity of presentation\n\nWe value the reviewers' feedback on improving the presentation clarity of our paper. We now answer the specific pain points raised in the questions.\n\n> When I first encounter $\\eta$, it is not clear what this represents …\n\n$\\eta$ represents a general mixture of conditional bridges $\\mathbb{Q}_{\\eta, t}$. As a modelling choice, we then set  $\\eta$ as the parametrized neural path interpolant, which we use to parametrize **a mixture of Brownian bridges** as stated in line 154. We understand that this point may not have been clear, and we will update the notation to improve clarity.\n\n> After line 159, $f^\\star_\\eta$ seems to correspond to the minimum value of the integral. Is that correct? …\n\nIndeed, this is the correct interpretation. We use the star notation to denote the solution to an optimization problem, such as $f^\\star$ or the optimal coupling $\\pi^*$, which leads to the Schrödinger Bridge.\n\n> What is the role of $\\kappa$\n\nThe kernel allows us to construct the reference process everywhere in space and time by smoothening out ground truth instances of the reference drift found in the data samples from each time marginal. This is part of the problem setup and is considered given.\n\n> what is the precise function of the neural path interpolant?\n\nThe neural path interpolant is used to match the reference drift $f_t$, which in general is no longer a straight path between $x_0$ and $x_1$. Intuitively, this makes our paths bend and can be thought of as the mean of the conditional mixture of Brownian bridges. To select which $x_0$ is paired with $x_1$ we use the neural path interpolant as part of Algorithm 2 to compute the OT cost and consequently solve the SB problem.\n\n## Differences with related work (e.g. Shen et. al 2024 and Bartosh et. al (2024 and 2025))\n\n### Shen et. al 2024\n\nShen et. al 2024 use an iterative bi-level algorithm wherein step 1) they solve a Schrödinger bridge and estimate the forward-backward drift given a current guess of the reference drift, 2) Given the solution to step 1, they simulate trajectories to estimate a new best reference guess. This bi-level approach is considered since the reference drift belongs to a family of possible reference drifts rather than a single prescribed one. In contrast, CurlyFM employs a 2-stage algorithm that is not iterative and also does not search over a family of reference drifts. More critically, CurlyFM is also simulation-free, which is achieved by making the modeling assumption that the conditional mixture of bridges is modelled as Brownian bridges. In practice, this means that CurlyFM is significantly more efficient to train. We report in Table R1 the computational cost in wall clock time for each approach and notice that CurlyFM completes the Ocean currents problem in a matter of minutes with better performance (see Table 2 in main paper), while the method proposed in Shen et. al 2024 is in the order of multiple hours— thus clearly demonstrating the computational benefit and empirical caliber of CurlyFM. Finally, note, both Shen et. al 2024 and CurlyFM **can model non-gradient field dynamics**, but only the latter is simulation-free. \n\nTable R1: Compute cost in hours\n| Method                   | Hours |\n|--------------------------|--------|\n| DM-SB                   | 15.44  |\n| Vanilla-SB              |  0.43  |\n| TrajectoryNet           |  7.44  |\n| SBIRR (Shen et al. 2024)|  4.67  |\n| CurlyFM (Ours)          |  **0.06**  |\n\n\n\n### Bartosh et. al (2024 and 2025)\n\nBoth papers from Bartosh et. al consider Latent SDEs with simulation-free training. The key difference between these papers and CurlyFM is that they do not solve a Schrödinger Bridge problem as considered in our work. Furthermore, SDEs in CurlyFM are not Latent SDEs. More precisely, Latent SDEs learned in Bartosh et. al 2025, have no theoretical basis to converge to a Schrödinger Bridge—i.e. the learning process does not minimize the KL between a reference process, and they do not learn a bridge in their current presentation of the paper. Consequently, we argue that the approaches and goals of these papers are different from CurlyFM and cannot be meaningfully compared experimentally, despite both approaches learning an SDE in a simulation-free manner. Whilst it may be possible to build a computational approach to solving the SB problem with latent SDEs, this is an orthogonal research direction to this paper. We will revise our related work section to include a more detailed commentary.\n\n## Handling more than two time points?\n\nWe acknowledge the reviewer's question. Please note we handle the multi-marginal case as done in prior work (e.g., MFM [1]), where by we construct CurlyFM between two time marginals by skipping an intermediate time marginal, for example, $t=1$ is skipped when $t=0$ and $t=2$ data is used. After which we proceed sequentially by constructing CurlyFM (instantiated from the previous iteration) on the next two time marginals and skipping the next time marginal, i.e., $t=2$ is skipped and $t=1$ and $t=3$ are used.\n\n## Access and role of the reference velocity field \n\nIn our considered setting, we assume access to the ground truth reference drift at samples $x_0 \\sim \\mu_0$ and $x_1 \\sim x_1$ as part of the problem setup. Note that in the high-impact application domains we consider, such as trajectory inference in single-cell data, we have access to the RNA velocity [2,3]. As we need to estimate CurlyFM everywhere in space and time during training, we construct a smoothed version of the reference drift by using a kernel $\\kappa_t$. This allows us to construct the reference drift $f_t (\\mathbf{X}t)$---using knowledge from the ground truth reference drift in the existing dataset— in places where there are no ground truth samples. Consequently, we take $f_t$ in these intermediate points as our ground truth reference drift.\n\n> Could the method be adapted to scenarios where the true velocity field is known?\n\nYes! If the reference is known everywhere, we would simply plug it in for $f_t$ rather than use our smoothed version via $\\kappa_t$.\n\n> How does the accuracy of this approximation affect the performance\n\nThis is a great question. To ablate this, we now construct a noisy reference drift for various noise levels in the tables below (for the Ocean Currents dataset). The noisy reference drift $f_{t}^{\\text{noise}}$ is obtained as a linear combination of the ground truth reference drift $f_{t}$ and noise from a standard gaussian distribution ($ f_{t}^{\\text{noise}} = (1 - \\beta) * f_{t} + \\beta * \\text{noise}$). We find that the performance between $\\beta = 0$ (no noise, i.e. regular CurlyFM) and $\\beta = 0.25$ (25\\% noise) are similar whilst the performance for $\\beta = 0.5$ and higher gradually becomes worse, as the noise dominates over the ground truth reference drift in $f_{t}^{\\text{noise}}$. This shows that CurlyFM is robust to moderate amounts of noise added to the ground truth reference drift.\n\nWe hope that our answers and ablation in Table R2 help alleviate the reviewers' concerns. We will also include this ablation in the updated draft of the manuscript.\n\nTable R2: Ablation with noisy reference drifts\n| $\\beta$ | Cos. Dist.     | L2 cost       | W1            |\n|--------:|----------------|---------------|---------------|\n| 0.00    | 0.062 ± 0.003  | 0.143 ± 0.010 | 0.034 ± 0.006 |\n| 0.25    | 0.057 ± 0.033  | 0.021 ± 0.036 | 0.051 ± 0.030 |\n| 0.50    | 0.087 ± 0.047  | 0.301 ± 0.085 | 0.091 ± 0.046 |\n| 0.75    | 0.261 ± 0.123  | 0.381 ± 0.120 | 0.145 ± 0.062 |\n| 1.00    | 0.428 ± 0.157  | 0.445 ± 0.121 | 0.237 ± 0.079 |\n\n\n## Concluding remarks\n\nWe thank the reviewer for their time and effort in reviewing our work, and we hope our efforts to clarify the main points, along with the additional ablations, allow the reviewer to consider improving their score, as they initially said they would be open to. We are also more than happy to answer any other questions that arise.\n\n## References\n\n[1] Kapusniak, Kacper, et al. \"Metric flow matching for smooth interpolations on the data manifold.\" Advances in Neural Information Processing Systems 37 (2024): 135011-135042.\n\n[2] A. Riba, A. Oravecz, M. Durik, S. Jiménez, V. Alunni, M. Cerciat, M. Jung, C. Keime, W. M. Keyes, and N. Molina. Cell cycle gene regulation dynamics revealed by rna velocity and deep-learning. Nature communications, 13(1):2865, 2022\n\n[3] V. Bergen, M. Lange, S. Peidli, F. A. Wolf, and F. J. Theis. Generalizing rna velocity to transient cell states through dynamical modeling. Nature biotechnology, 38(12):1408–1414, 2020."}, {"review_id": "xcljF0S7Ly", "reviewer": "Reviewer_BnED", "summary": "An approach is presented that allows learning flows that include drifts rather than only gradient-field dynamics. The authors first learn an interpolant (neural path) and then fit a velocity field based on the coupling provided by the interpolant.", "questions": "- It is unclear how the trade off between matching the drift and marginals is handled. What knob needs to be turned in Algorithm 1 and 2 to trade off these two learning objectives? A more detailed discussion should be added because this question comes up several times in the paper (e.g., Intro on page 3 top and line 186-187)\n\n- There are other approaches than simply pairwise connecting marginals with OT such as action matching by [Neklyudov et al., 2024]. These approaches incur often lower training costs than solving many OT problems. Can this be used as a first stage in the presented approach as well to avoid some of the limitations of pairwise interpolants?\n\n- Minor questions: The authors state that TrajectoryNet is much slower and unstable (line 183). Can the authors quantify this for one of the examples?", "limitations": "Yes, the authors have addressed in the Conclusion the limitation of first learning the reference field, which is the main limitation of this approach.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n+ Learning non-gradient field dynamics is indeed an important problem and little work in this direction has been done. So this work is timely and useful.\n\n+ The proposed two-stage approach of first learning a neural path interpolant and then fitting a velocity field is straightforward but reasonable and sound. The approach is well described.\n\n+ There are extensive numerical experiments that demonstrate that the approach captures drift dynamics whereas standard conditional flow matching fails\n\nWeaknesses:\n- There are extra training costs incurred by the additional steps (e.g., the OT step) that can be high. No discussion about it is provided.  \n\n- There is a limitation by first fitting the neural path interpolant, which is not addressed well. For example, if the initial neural path is poorly fitted, how will it affect the second step? The authors do mention briefly in the Conclusion that the approach hinges on inferring the reference field well, but this requires more discussion.", "was_revised": true, "final_justification": "There were no major issues before the rebuttal. The answers by the reviewers clarified some minor comments and overall confirmed by Accept rating.", "rebuttal": "We are grateful to the reviewer for their time and effort reviewing our paper, as well as their positive appraisal of our work. We appreciate that the reviewer agrees with us that learning non-gradient field dynamics “is indeed an important problem and little work in this direction has been done”. We are also heartened to hear that the reviewer finds our two-stage approach to be “reasonable and sound” and that our approach is “well described.” Finally, we thank the reviewer for finding our set of numerical experiments “extensive” and that they demonstrate settings where CurlyFM captures drift dynamics while “standard conditional flow matching fails”. We next answer the main questions raised in the review.\n\n## Tradeoff between matching the drift and marginals\n\nWe acknowledge the reviewer's question regarding the possibility of an inherent trade-off in matching the reference drift versus matching the marginals. We believe there may be a potential cause for confusion, as there is no single hyperparameter that negotiates a trade-off between Algo 1 and Algo 2. Indeed, CurlyFM is a two-stage method which first (1) matches the reference drift as closely as possible, and then (2) uses the learned reference drift to learn the optimal marginals that solve the mass transport problem. The hyperparameter and simulation-free nature of CurlyFM is a significant advantage over other methods, e.g., TrajectoryNet [1].\n\nOur comments on the top of page 3 and lines 186-187 were simply to illustrate the fact that learning the reference drift perfectly **does not necessarily solve the mass transport problem as a Schrödinger Bridge**, and this is exactly the reason we need Algo 2. As a result, the perfect final model should be a valid Schrödinger Bridge, which is close to the reference as much as is required to actually solve the mass transport problem and no more. We understand that this point may not have been sufficiently clear in the main paper, and we will update the draft to include more commentary on this point.\n\n## Other mechanisms for connecting marginals\n\nWe thank the reviewer for their interesting question regarding alternative approaches to connecting the marginals without solving an OT problem, such as Action Matching. Unfortunately, to the best of our knowledge, methods such as Action Matching are more computationally expensive than what is presented in Algo 2. This is because Algo 2 remains simulation-free due to an unbiased stochastic estimate of the OT cost—rather than computing the entire cost using simulation. In practice, this is quite fast and is not a significant computational overhead in all the experimental settings we consider in the paper. At a broader level, approaches that do not connect marginals via OT run the risk of not being valid solutions to the Schrödinger Bridge problem, which is the thesis of our work. We hope that this answers the very interesting question raised by the reviewer.\n\n## Computational speed of TrajectoryNet\n\nThank you for your great suggestion. We have now included the computational cost of TrajectoryNet (and Curly-FM) in Table R1. We report these results for the 2D and 10D settings for the experiments reported in Table 4.\n\nTable R1: Computational efficiency comparison between CurlyFM and Trajectory Net\n| Method           |        d = 2 (seconds)    |        d = 10 (seconds)   |\n|--------------------|-------------------------------|-------------------------------|\n| TrajectoryNet  |    $17005 \\pm 110$     | *$43080 \\pm 65$         |\n| CurlyFM         |     $\\mathbf{1429 \\pm 31}$        |  $\\mathbf{2471 \\pm 10}$           |\n\n\nAs observed in the above table, we find that TrajectoryNet is considerably slower than CurlyFM and is not a scalable approach. This is evident in that training TrajectoryNet takes, on average, 11 times longer than Curly-FM in the 2D case, and 17 times longer in the 10D case. Moreover, as dimensionality increases from 2D to 10D, Curly-FM incurs an increase in computational cost by ~1.7x, whereas TrajectoryNet incurs at least a 2.5x increase in computational cost. We use (*) to denote runs that did not finish within allotted resource allocation time. We will include these results in the updated draft of the paper.\n\n## Concluding remarks\n\nWe thank the reviewer again for their time and effort during this review process. We hope that our rebuttal here was successful in answering all the great questions raised by the reviewer and allows the reviewer to continue positively endorsing our paper. We are also happy to answer any further questions the reviewer may have during the remainder of the rebuttal period.\n\n\n## References\n[1] Tong, Alexander, et al. \"Trajectorynet: A dynamic optimal transport network for modeling cellular dynamics.\" International conference on machine learning. PMLR, 2020."}, {"review_id": "g1OyeICw0F", "reviewer": "Reviewer_9hEa", "summary": "In this paper, the authors propose a numerical methodology to solve the standard dynamic Schrödinger Bridge (SB) problem ($\\min KL(P|Q):P_0=\\pi_0, P_1=\\pi_1$), in the particular case where the reference process $Q$ has non-zero (or at least non-linear) drift. This is motivated by the fact that various real-world dynamics exhibit non gradient, periodic behavior, which cannot be modeled by simple reference dynamics. In theory, the solution to this problem writes as the mixture of the induced entropic OT plan and the reference \"bridge\" process $Q_{|0,1}$ (i.e., the reference process pinned at terminal times 0 and 1). The goal of the method is to approximate the drift of this solution by a neural network to be able to generate dynamics afterwards. A popular and simulation-free approach to compute the SB problem is bridge matching, as proposed in [1], which relies on two requirements wrt $Q$ : (i) being able to sample from the conditional distribution $Q_{t|0,1}$ at any time $t$ and (ii) having a tractable expression of the drift of $Q_{|0,1}$ (which is Markovian). In the current setting, these two requirements are not satisfied.\n\nTo circumvent this issue, the authors propose to approximate the reference bridge by a neural interpolant (i.e., an interpolant that relies on a neural network) before applying the standard bridge matching technique, in a similar way to [2]. This ends up in a two stage procedure coined \"Curly Flow Matching\", which requires two neural networks (one to approximate the reference bridge, one to approximate the drift of the SB solution). Then, they present the results of numerical experiments on various types of dynamics (single-cell trajectory, fluid dynamics...) and compare their method to classic flow matching and mini-batch OT-based flow matching. \n\n[1] Diffusion Schrödinger Bridge Matching, Shi et al., Neurips 2024. \n\n[2] Generalized Schrödinger Bridge Matching, Liu et al. ICLR 2024.", "questions": "**Questions**\n- from what I understand, the first stage of Curly-FM consists in approximating the marginal distribution of the reference bridge (i.e. $\\mathbb{Q}_{t|0,1}$) by a deterministic neural interpolant. Since $\\sigma_t$ (the volatility of the reference process) is not zero, I have the feeling that there is a mismatch between what needs to be approximated (a probability distribution) and the approximation itself (a Dirac mass). In this case, the resulting SB solution is not correct, since it is not stochastic. Do I miss something ? In my opinion, a more adapted parameterization of the marginal distribution of the reference bridge would be a Gaussian distribution whose mean is given by the proposed deterministic interpolant with parameterized variance. Of course, this would change the matching losses (since the drift of the corresponding bridge process would not simply given by the time derivate of the Gaussian mean).\n- I do not understand how the loss of Alg. 1 is built : putting aside my previous point (ie considering the deterministic case only), the objective function $\\mathbb{E}[\\|\\partial_t X_t,\\eta-f_t (X_{t, \\eta})\\|^2]$, where the expectation is taken wrt $(X_0, X_1, X_{t, \\eta})$ and $f_t$ is exactly the drift of the reference process, is perfectly legitimate in my opinion to approximate the reference bridge. In the stochastic case, this would correspond to minimize wrt $\\eta$ the objective $KL(\\mathbb{Q^\\eta}|\\mathbb{Q})$, where $\\mathbb{Q^\\eta}=\\int \\mathbb{Q^\\eta}(\\cdot|x_0, x_1)d\\pi(x_0,x_1)$,  as also proposed in [1] for the general SB setting. Why justifies the use of kernels here ?\n- I do not understand how the cost of the OT plan can be deduced. Regarding the formulation of the SB problem, the optimal plan $\\pi^\\star$ is given by $\\min KL(\\pi|Q_{0,1})$ (with the right marginals), where $Q_{0,1}$ is the joint distribution of the reference process at times 0 and 1 (which is not tractable). Even when using the approximated bridge (defined via the neural interpolant), this joint distribution remains untractable and an additional entropy term is missing in the loss. How do you obtain this OT cost ? How do you explain that you obtain a non-regularized OT cost, while the SB problem naturally leads to an entropy-regularized OT problem when reducing to the couple (0,1) ? \n- do you consider different values for $\\sigma_t$ and $g_t$ ? How these values are set in the experiments ? How does the algorithm behave when changing these values ?\n- the tradeoff between directly computing the OT plan (for example with Sinkhorn algorithm or mini-batch OT) + applying bridge matching VS doing an iterative refinement of bridge matching loss, by taking the current outputs of the model (see [1, 2]) has not been obviously solved in the past years: have you tried to turn Curly-FM into an iterative algorithm ?\n\n\n**Suggestions**\n- I would suggest the authors to write \"non-linear drift\" rather than \"non-zero drift\", which is more accurate since the DSBM approach also provides expressions for the first case, which is more general (see Appendix B.2)\n- I would recommend the authors to put more highlight on the mathematical formulation they want to solve, in particular on the choice of the volatilities $\\sigma_t$ and $g_t$ : is it any positive couple ? is the same value for both volatilities ? do you always consider the zero-limit ?\n- I would also suggest the authors to rewrite Section 3.2 in a clearer manner (for example, explicitly define two parts : (i) learning the bridge and (ii) learning the drift of the SB solution) to include the computation steps that justify the way the algorithm is built\n\n[1] Generalized Schrödinger Bridge Matching, Liu et al. ICLR 2024.\n\n[2] Schrödinger Bridge Flow for Unpaired Data Translation. De Bortoli et al, 2024\n\n**Update**: these questions have been addressed by the authors in the discussion and imply changes in the manuscript.", "limitations": "I think that the paper lacks detailed ablation studies on the proposed algorithm to evaluate its robustness (toy setting with increasing dimension, various values of volatility...). For example, it would be interesting to have an experiment, where the reference bridge can be analytically computed (i.e linear reference SDE) and compared with the output of Algorithm 1. This would enable to have some empirical guarantee on the first stage of Curly-FM.\n\n**Update**: these questions have been addressed by the authors in the discussion and imply changes in the manuscript.", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 1, "significance": 2, "originality": 2, "strengths_and_weaknesses": "**Strengths**:\n- the motivation of tackling the SB problem with non-zero drift is well explained\n- the range of the numerical experiments is interesting and based on real-world systems\n\n**Weaknesses**:\n- in my opinion, the mathematical explanation of the Curly-FM algorithm (Section 3.2) lacks crucial details on the manner the losses are defined (I could not find any additional details in the appendix), which makes the reading very hard. While the authors claim to solve the SB problem (which is stochastic by nature, due to the consideration of non-zero volatilities), the proposed method seems to only solve the deterministic-limit case (i.e. when the volatilities go to zero)... In particular, I have serious doubts about (i) the loss to approximate the reference bridge and (ii) the way the OT plan is approximated (see section \"Questions\") to solve an arbitrary SB problem with non-linear drift.\n-  in [1], the authors propose a very close approach to Curly-FM in a general SB setting, where the reference process does not have a trivial drift too (i.e., (1) approximation of the reference bridge and (2) application of bridge matching). The key advantage in comparison with the current method is that the reference bridge is approximated by a Gaussian Stochastic interpolant (which has since been extended to interpolants with Gaussian mixtures) which does not rely on any neural parameterization (rather a spline parameterization). This method should be compared with in the numerical experiments, especially via the trade-off expressivity/computational budget.\n\n[1] Generalized Schrödinger Bridge Matching, Liu et al. ICLR 2024.", "was_revised": true, "final_justification": "The authors have addressed most of my concerns (described in my final answer in the discussion). I increase my score based on the integration of the elements given by the authors in the manuscript.", "rebuttal": "We thank the reviewer for their time and effort in reviewing our work. We are glad the reviewer found our work “well explained” and found the range of experiments interesting. We now address the main concerns raised in the review, grouped by theme.\n\n## Motivation for CurlyFM \n\nOur work is motivated by practical applications, which are similar to Schrödinger bridge problems and motivated by them, but have different criteria for success. Namely, we aim to model dynamical systems well empirically in a variety of diverse settings, where the reference drift may only be observed at measurements rather than everywhere in space and in time. This necessitates the need to devise a learning framework that is efficient to scale to these problems and to iterate quickly, two points that we have sought to demonstrate empirically through our selection of benchmarks and results.\n\nIn more detail, our work attempts to find a computationally efficient solution to model our high-impact application settings, such as trajectory inference in single-cell RNA with RNA velocity, rather than computing the best numerical approximation to the Schrödinger bridge problem in a general setting. Given the vast and storied literature on SBs with many numerical approximations that require simulations to solve, we position CurlyFM as a simulation-free method that takes key ideas from the SB problem and allows modeling of complex non-gradient field dynamics found in real-world systems, which may have relatively low noise relative to the scale of the drift. \n\nAlso of note, in these applications, previous work [1] has found low or no stochasticity to work the best empirically. Consequently, most of our approximations are optimized assuming that we will likely be in the low or no-stochasticity settings, as motivated by our choice of experimental benchmarks. We understand that this motivation may not have been clear in our original draft, and we will update the paper accordingly. \n\n## Approximation of SB problem in the stochastic case\n\nWe thank the reviewer for their insightful question and opportunity to clarify our work. Despite our work being motivated in the little to no stochasticity regime, we now consider an ablation below where $\\sigma > 0$ to demonstrate that this is a simple modification to CurlyFM that does not lose any of its efficiency. \n\n### Behavior with respect to $\\sigma_t$ and $g_t$\n\nWe find, similar to previous work [1], that low values of $\\sigma$ perform the best on all metrics. As a result, we recommend setting $\\sigma_t$ to zero unless some reference $\\sigma$ value is known. Therefore, all of our experiments in the current manuscript are with $\\sigma_t=g_t=0$. To show this in our setting, we perform a new ablation of $\\sigma_t$ and $g_t$ on the oceans application. We find that larger sigma monotonically decreases performance on our tasks, justifying our choice of $\\sigma_t=0$ for our empirical work.\n\nTable R1: Ablation on stochasticity for CurlyFM\n|$\\sigma$|Cos. Dist.|L2 cost|W1|\n|-|-|-|-|\n|0.01|0.061 ± 0.003|0.141 ± 0.009|0.028 ± 0.066|\n|0.10|0.062 ± 0.002|0.145 ± 0.011|0.066 ± 0.008|\n|1.00|0.145 ± 0.009|0.474 ± 0.058|0.871 ± 0.048|\n\nWe further note that bridge matching approaches based on iterative refinement do not converge well with small values of $\\sigma$ and therefore are only really practical in the large $\\sigma$ setting [1]. In settings where a large $\\sigma_t$ is required, we believe that these methods may be more well-suited to the task, but present a different set of tradeoffs with better theoretical approximation of the bridge, but in contrast, may incur larger computational complexity due to simulation.\n\n## Comparison with GSBM\nWe appreciate the reviewer’s suggestion on the inclusion of GSBM as a baseline. We now include this baseline adapted to the Ocean currents dataset below using spline neural interpolant $x_t$ (as in GSBM) as well as a comparison of the training times between GSBM and CurlyFM. \n\nTable R2: Comparison with GSBM\n| | Cos. Dist. | L2 cost | W1 | Train time (s) |\n|-|-|-|-|-|\n| **GSBM** | 0.279 ± 0.006 | 0.395 ± 0.008 | 0.337 ± 0.009 | **37** |\n| **CurlyFM**| **0.062 ± 0.003** | **0.143 ± 0.010** | **0.034 ± 0.006** | 176 |\n\nWe see that CurlyFM performs significantly better than GSBM, albeit at a higher cost which represents the expected tradeoff between expressivity and performance. We will include these baselines in our update to the paper.\n\n## Use of kernels\nThe reviewer is correct that if we have direct access to the reference $f_t(x)$ which is continuous in time and space, then kernels are unnecessary. However, in the practical settings that we investigate, a continuous ground truth reference field does not exist $f_t(x)$. In our applications, we only have access to the reference field at discrete points in space and time that correspond to ground truth data at the different time marginals. Consequently, we use a kernel to build a continuous reference field $f_t(x)$. We will clarify this further in the main text and thank the reviewer for their comments on how to improve the readability of our work, and we will add motivation for each step in the methods section.\n\n## Cost of OT plan\nThe reviewer makes an insightful observation that the optimal plan $\\pi^* = \\min KL( \\pi | Q_{0,1})$ is indeed intractable. Consequently, we make several approximations to these couplings that enable faster throughput as offered in simulation-free training. In particular, we make three approximations to the cost $c(x_0, x_1)$ between two points, and one on the coupling given this cost. These are the following:\n\n- 1) Algorithm 1 has learned a path or set of paths for the means of Gaussians not too far from the optimal $Q_t(x_t | x_0, x_1)$. \n- 2) The diffusion term $\\sigma_t$ is small, and therefore the cost is close to the distance the means travel.\n- 3) The integral of the squared length of the curve is approximated empirically using a Monte Carlo estimate\n \nAfter approximating $c(x_0, x_1)$, we use mini-batch OT to approximate the entropic OT problem following [1]\n\nInterestingly, in the low stochasticity setting [1] found that mini-batch OT with no stochasticity was empirically a better approximator of $\\pi^*$ in the Gaussian case. This is because mini-batch transport adds some amount of “entropy” to the plan due to its approximation. We follow this approximation in our setting. We hope that this clarifies some of the approximations made in the OT computation for our method, and will add a detailed explanation of these choices to the appendix, as well as a brief explanation in the main text.\n\n## Tradeoff between directly computing OT plan vs. iterative refinement\nThis is an interesting question. We have not tried turning CurlyFM into an iterative refinement-based algorithm. However, we argue that this is a fairly easy adaptation, but not suitable for our application domains, where we assume small stochasticity levels. We also note that iterative refinement-type algorithms all require extensive simulation for each refinement step, and hence are not really simulation-free methods, which is a key goal here. However, we do acknowledge that iterative refinement-based approaches can still be significantly faster than methods that require backpropagation through simulation (e.g. [2,3]). We will update this within our related work and also discuss the tradeoffs depending on stochasticity levels in the appendix. \n\n## Other suggestions around clarifications\nWe appreciate the reviewer's thoughtful suggestions, and we will update the manuscript accordingly, taking into account all three suggestions. \n\n## Ablation to evaluate its robustness in a toy setting\nWe thank the reviewer for their insightful suggestion. To assess robustness, we designed a toy experiment with an analytical SB solution and non-gradient dynamics based on [1,4,5]. We bridge two Gaussians in the presence of a spiral reference field across various dimensions $d$ and noise levels $\\sigma_t​=g_t$​. We compare against the closed-form solution and OT-CFM, measuring KL divergence, Wasserstein distance, and the cosine distance to the ideal rotational angle. Results confirm CurlyFM performs best in low-dimensional, low-stochasticity settings, consistently achieving high cosine similarity.\n\nTable R3: CurlyFM Ablation\n|Dim|$\\sigma$|KL$(p_1, q_1)$|Mean KL|Cos. Dist.|$W_2(p_1, q_1)$|\n|-:|-:|-:|-:|-:|-:|\n|3|0|0.043|0.024|0.001|0.192|\n|3|0.1|0.04|0.017|0.038|0.61|\n|3|1|0.049|0.049|0.519|0.613|\n|5|0|0.021|0.017|0.01|0.871|\n|5|0.1|0.041|0.018|0.018|0.88|\n|5|1|0.065|0.038|0.522|0.877|\n|20|0|0.258|0.178|0.022|3.887|\n|20|0.1|0.25|0.169|0.029|3.891|\n|20|1|0.358|0.234|0.541|3.749|\n|OT-CFM Baseline|\n|3|0|0.033|0.024|0.999|0.369|\n|5|0|0.031|0.021|1|0.879|\n|20|0|0.201|0.165|0.997|3.921|\n\n## Concluding remarks\nWe thank the reviewer again for their time, effort, and thoughtful suggestions, which have enabled us to clarify technical details and also include ablations and baselines that strengthen the empirical results. We further hope that this rebuttal was successful in answering all the reviewer's great points and alleviating the concerns shared by the reviewer. If the reviewer is satisfied, we invite the reviewer to consider a fresher evaluation of our work with this rebuttal in context. We are also more than happy to answer any further questions the reviewer might have. Please do let us know.\n\n## References\n[1] Tong, Alexander, et al. \"Simulation-free schrödinger bridges via score and flow matching.\" AISTATS. 2024.\n\n[2] Tong, Alexander, et al. \"Trajectorynet: A dynamic optimal transport network for modeling cellular dynamics.\" ICML. 2020.\n\n[3] Koshizuka, Takeshi, and Issei Sato. \"Neural lagrangian schrödinger bridge: Diffusion modeling for population dynamics.\" ICLR. 2023.\n\n[4] De Bortoli, Valentin, et al. Diffusion Schrödinger bridge with applications to score-based generative modeling. NeurIPS. 2021\n\n[5] Shi, Yuyang et al. Diffusion Schrödinger bridge matching. NeurIPS. 2023"}, {"review_id": "PnoQkIlN5q", "reviewer": "Reviewer_qrFN", "summary": "The authors proposed a generalization of flow matching method for trajectory inference without assuming a gradient field reference drift. The main idea is by decomposing Schrodinger bridge's solution into mixture of endpoint couplings and conditional bridges, and construct them in a simulation free manner separately. The authors tested their method in several different setups and seen their proposal out perform existing baselines.", "questions": "1) In eq (9), with general diffusion term $g_t$, should the KL divergence being $E_{P}\\int_t g_t^{-1}||v||_2^2$? Similar question holds for the unnumbered equation below eq. (10).\n\n2) I might miss understood it, as I understood the authors modeled $\\mathbb{Q}$ using a mixture of Brownian bridges, with potentially different diffusion terms and $\\mathbb{P}$ to have $g_t$ as diffusion term. Will the KL divergence between $\\mathbb{P}$ and $\\mathbb{Q}$ constructed this way being finite? KL divergence between SDEs not sharing diffusion terms might be infinite since they are not absolute continuous to each other. \n\n3) The authors applied their method in a RNA velocity dataset treating RNA velocity as similar velocity in the SDE model. I am not sure to what extend they are connected. RNA velocity take advantage of splicing happens after transcription and use the ratio to infer whether genes are being transcribed more or less. These velocity would depends on each gene's splicing rate and mRNA's degradation rate. I wonder how they author think about these? In the ocean example one might indeed have some sparse velocity measurements but need interpolation. I encourage the authors to make some comments on how good the velocity measurements needs to be.", "limitations": "Yes mostly. But I encourage the authors to provide some comments on how velocity measurements were taken in each experiments and assumptions on them.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "## Strength\n1) The algorithm is simulation free and performs well.\n2) Can handle none potential dynamics and leverage velocity information. \n\n## Weakness \n1) I am not fully convinced that some of the model choice is compatible with Schrodinger bridge formulation esp. when the diffusion term is different. See questions for details. \n2) The notations is a bit hard to follow -- I understand this is partly due to Schrodinger bridge problems involving many different SDEs.", "was_revised": true, "final_justification": "The author was able to address all my concerns and I keep my positive evaluation of this paper.", "rebuttal": "We thank the reviewer for their time and effort in reviewing our work. We are pleased that the reviewer found our method to “perform well” and that the reviewer agrees that CurlyFM can handle non-gradient field dynamics by leveraging velocity information. We now address the main questions posed in the review. \n\n\n## Compatibility of the model with Schrödinger Bridge\n\nWe acknowledge the reviewer's questions regarding the compatibility of our diffusion terms $g_t$ in various places of our Schrödinger Bridge formulation. We now answer each question in turn.\n\n> In eq (9), with general diffusion term , should the KL divergence being $E_P \\int_t g_t^{-1} \\|v\\|^2_2$? Similar question holds for the unnumbered equation below eq. (10).\n\nThe reviewer is indeed correct in their observation. However, for our problem, we can drop the scaling factor as it doesn’t affect the optimization problem, so we drop this term for simplicity. We will add this as a note to clarify this point.. \n\n> Compatibility with different diffusion terms\n\nThe reviewer is again correct in their statement that SDEs in the SB problem that do not share diffusion terms could lead to problems. However, in this work, we consider the KL divergence between SDEs with matching diffusion terms, and always construct $\\mathbb{Q}$ as a mixture of Brownian bridges with the same diffusion term as $\\mathbb{P}$. We hope that this clarifies this technical point.\n\n## Velocity data in transcriptomics and oceans datasets\n\nWe thank the reviewer for their helpful suggestions here. We are happy to provide further discussion on velocity measurements and the assumptions on them in each application. The reviewer is correct to point out that the assumptions on velocity measurements may be quite different between applications, and the quality of the velocities may vary greatly as well. \n\n### Relating RNA velocity to velocity in the SDE model \n\nWe agree with the reviewer that RNA velocity attempts to compute the rate of transcription of each gene, which is not directly analogous to the velocity used in the SDE model. We follow [1] in assuming that the RNA velocity estimate is a reasonable estimate (up to a scaling factor) of the SDE velocity. There are a number of discrepancies between these two types of velocities, but they nevertheless seem correlated enough in practice to provide improved estimates of the SDE velocity field as demonstrated in a number of subsequent works [2, 3, 4]. \n\nIntuitively, this interpretation makes sense, as assuming we have a perfect estimate of transcription and degradation rates and at an infinite number of locations, then following the continuous vector field in RNA state space defined by these rates should be a reasonable model of cell dynamics (at least at the RNA level). These, of course, are unrealistic assumptions, but they give intuition on how RNA velocity estimates may be useful in SDE models of cell dynamics.\n\n### Sparse velocity interpolation\n\nYes, we thank the reviewer for their insightful comment on velocity interpolation. This is indeed a necessary step in both the oceans and single-cell experimental setups—the most realistic setting—we consider in CurlyFM, as most velocities will be measured at sparse locations both in time and space. \n\nIn this work, we primarily explore how to use a continuous velocity prior in flow matching models and leave detailed exploration of this interpolation to future work. More precisely, we rely on using simple kernel smoothing following [5] for all applications with sparse velocities. \n\n### How good do the velocity measurements need to be?\n\nWe thank the reviewer for their insightful remark on how the quality of velocity measurements might affect the performance of our model. To investigate this, we ablate the velocity quality by adding various levels of noise. In particular, we construct a noisy drift  $f_{t}^{\\text{noise}}$ as a linear combination of the ground truth reference drift $f_{t}$ and noise from a standard gaussian distribution ($ f_{t}^{\\text{noise}} = (1 - \\beta) * f_{t} + \\beta * \\text{noise}$). Note that $\\beta = 0$ (no noise) corresponds to the setting of CurlyFM used in the main paper. The results of this ablation are presented in the Table R1.\n\nWe find that the performance between $\\beta = 0$ (no noise, i.e. regular CurlyFM) and $\\beta = 0.25$ (25\\% noise) are similar whilst the performance for $\\beta = 0.5$ and higher gradually becomes worse, as the noise dominates over the ground truth reference drift in $f_{t}^{\\text{noise}}$. This shows that CurlyFM is robust to moderate amounts of noise added to the ground truth reference drift.\n\nTable R1: Ablation with noisy reference drifts\n| $\\beta$ | Cos. Dist.     | L2 cost       | W1            |\n|--------:|----------------|---------------|---------------|\n| 0.00    | 0.062 ± 0.003  | 0.143 ± 0.010 | 0.034 ± 0.006 |\n| 0.25    | 0.057 ± 0.033  | 0.021 ± 0.036 | 0.051 ± 0.030 |\n| 0.50    | 0.087 ± 0.047  | 0.301 ± 0.085 | 0.091 ± 0.046 |\n| 0.75    | 0.261 ± 0.123  | 0.381 ± 0.120 | 0.145 ± 0.062 |\n| 1.00    | 0.428 ± 0.157  | 0.445 ± 0.121 | 0.237 ± 0.079 |\n\n\n## Concluding remarks\n\nWe would like to thank the reviewer again for their time reviewing our paper. We hope that our rebuttal here was successful in answering the reviewer's great points and allows the reviewer to continue positively endorsing our paper. We are also more than happy to answer any further questions that arise in the discussion phase.\n\n## References \n\n[1] Tong, Alexander, et al. \"Trajectorynet: A dynamic optimal transport network for modeling cellular dynamics.\" International conference on machine learning. PMLR, 2020.\n\n[2] Koshizuka, Takeshi, and Issei Sato. \"Neural lagrangian schr\\\" odinger bridge: Diffusion modeling for population dynamics.\" arXiv preprint arXiv:2204.04853 (2022).\n\n[3] Neklyudov, Kirill, et al. \"A computational framework for solving wasserstein lagrangian flows.\" arXiv preprint arXiv:2310.10649 (2023).\n\n[4] Qiu, Xiaojie, et al. \"Mapping transcriptomic vector fields of single cells.\" Cell 185.4 (2022): 690-711.\n\n[5] Shen, Yunyi, Renato Berlinghieri, and Tamara Broderick. \"Multi-marginal schr\\\" odinger bridges with iterative reference refinement.\" AISTATS (2025)."}]

Author rebuttal:

Meta-review: This paper proposes a two-stage method for learning non-gradient field dynamics. Authors appreciated the novelty of the and significance of the method. The main remaining issue were concerns regarding clarity. In particular, reviewers Tr4n and 9hEa bring up some significant issues of clarity / opportunities for in the manuscript. The authors have committed to make these changes. These are quite substantial, to the degree that I consider it borderline if it might be appropriate for the paper to undergo another round of review. But I am trusting the authors to make the promised changes.
Decision: Accept (poster)