Sequence Modeling with Spectral Mean Flows
Authors: Jinwoo Kim
Max Beier
Petar Bevanda
Nayun Kim
Seunghoon Hong
Venue: NeurIPS 2025 poster ()
Keywords: sequence modeling
time series
hidden Markov models
mean embeddings
linear operators
maximum mean discrepancy
gradient flows
Primary area: deep_learning
TL;DR: We propose a new algorithm for deep generative modeling of sequence data in continuous spaces based on a novel adaptation of operator theory for probabilistic dynamical systems.

Abstract: A key question in sequence modeling with neural networks is how to represent and learn highly nonlinear and probabilistic state dynamics. Operator theory views such dynamics as linear maps on Hilbert spaces containing mean embedding vectors of distributions, offering an appealing but currently overlooked perspective. We propose a new approach to sequence modeling based on an operator-theoretic view of a hidden Markov model (HMM). Instead of materializing stochastic recurrence, we embed the full sequence distribution as a tensor in the product Hilbert space. A generative process is then defined as maximum mean discrepancy (MMD) gradient flow in the space of sequences. To overcome challenges with large tensors and slow sampling convergence, we introduce spectral mean flows, a novel tractable algorithm integrating two core concepts. First, we propose a new neural architecture by leveraging spectral decomposition of linear operators to derive a scalable tensor network decomposition of sequence mean embeddings. Second, we extend MMD gradient flows to time-dependent Hilbert spaces and connect them to flow matching via the continuity equation, enabling simulation-free learning and faster sampling. We demonstrate competitive results on a range of time-series modeling datasets.

Sequence Modeling with Spectral Mean Flows
Jinwoo Kim KAIST
Max Beier TU Munich
Petar Bevanda TU Munich
Nayun Kim KAIST
Seunghoon Hong KAIST
Abstract
A key question in sequence modeling with neural networks is how represent and learn highly nonlinear and probabilistic state dynamics. Operator theory views such dynamics as linear maps on Hilbert spaces containing mean embedding vectors of distributions, offering an appealing but currently overlooked perspective. We propose new approach sequence modeling based an operator-theorets view of a hidden Markov model (HMM). Instead of materializing stochastic recurrence, we embed the full sequence distribution tensor in the product Hilbert space. generative process is then defined as maximum mean discrepancy (MMD) gradient flow in the space of sequences. slow sampling convergence, we introduce spectral mean flows, novel tractable algorithm integrating two core concepts. First, we propose new neural architecture leveraging spectral decomposition of linear operators to derive scalable tensor network decomposition of sequence mean embeddings. Second, we extend MMD gradient flows to time-dependent Hilbert spaces and connect them to flow matching via the continuity equation, enabling simulation-free learning and faster sampling. We demonstrate competitive results range of time-series modeling datasets
overcome challenges with large tensors and
Introduction
A fundamental question in sequence modeling with neural networks is how to represent and learn highly nonlinear and probabilistic state dynamics throughout a sequence. A popular method is to model each step with stochastic nonlinear network [3] 23
47 62, 98). Since this requires serializa- tion over steps at test time, there has been recent interest in linear recurrence [36] 37 77] that offers more parallelizability. Yet, their underlying theory often only covers deterministic linear dynamics, requiring probabilistic dynamics to be handled in a post hoc manner via autoregression 36 37 or
stepwise variational inference 82 104 72 This leaves room for new method that computationally
uses linear recurrence while natively handling the nonlinear and probabilistic nature of state dynamics.
Operator theory suggests appealing approach in this context: embedding probability distributions Hilbert spaces in which nonlinear dynamics become linear [68]. Specifically, consider Markovian state space Z. The theory states that, with proper choice of Hilbert space G, the distribution of each state Z" can be uniquely embedded as a vector μ2z in G. Furthermore, the state dynamics r [
dynamics (2"),.
can be represented linear operator U G G such that:
pz Upz
(1.1)
This enables the use of standard linear algebraic tools, such as spectral decomposition of linear operators, to implicitly but correctly represent and learn nonlinear probabilistic dynamics. This approach has been extensively explored in kernel algorithms for learning dynamical systems [5 However, its application neural sequence modeling has been limited.
In this work, we propose new approach to neural sequence modeling that leverages this perspective. We begin with the classical notion of hidden Markov model (HMM), which is probabilistic model
Code available at baattps://giie httpa//llshlpp/ggggtgta.oo/////e -eatsttrrrllmmm--lfooooo
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

for sequences X XI:N_ = (X) .'....,,", )) based on hidden states Z = Z under Markovian dynamics PI and local observation model 2"*) While HMMs offer a structured and universal way to model sequences, directly modeling their nonlinear stochastic recurrence can be challenging. This motivates the aforementioned operator- theoretic perspective.
Specificallyy, our first step is to embed the distribution of the entire sequence X1:N as an element in a tensor product Hilbert space. This provides a basic setup for sequence modeling in the operator- theoretic context, allowing us to exploit linearity (T.I). The setup also allows us to define a generative process as a maximum mean discrepancy (MMD) gradient flow in the space of sequences. Then, to overcome practical challenges associated with nigh-dimensional tensors and slow convergence of standard MMD flows, we introduce spectral mean flows, a novel method with two core contributions:
We propose new neural architecture that utilizes spectral decomposition of the linear operator underlying an HMM. This leads to scalable tensor network decomposition of the sequence distribution embedding, making computations tractable.
We propose an extension of MMD gradient flows to time-dependent Hilbert spaces, which may of independent interest. This connects to flow matching (60] via the continuity equation, enabling end-to-end, simulation-free learning and faster sampling convergence.
We test spectral mean flows range of time-series datasets, demonstrating competitive results.
The remainder of this paper is organized as follows. In Section 2 we introduce necessary background. Section introduces spectral mean flows. Section [4 discusses related work, and Section presentss experiments on synthetic data and time-series modeling datasets. Section draws conclusions. All theoretical arguments are formally stated and proven in the appendix.
2 Background
In this section, we introduce key mathematical tools that enable us work with probability distribu- tions using linear algebraic methods. More details can be found in Appendix
Reproducing kernel Hilbert spaces (RKHSs) offer framework for treating functions as vectors in Hilbert spaces [83] Consider a domain X' C Rd, and let HH be an RKHS induced by kernel k: : X' x X + R. The space IHI contains functions ʃ : X' + R and is equipped with an inner product (., >/ and a corresponding norm || ||. Its canonical feature map Φ : X -> IHT is defined by Φ(x) == x),, and satisfies the reproducing property: f fx) (ø(x), f for any f H..
To handle dependencies between variables, we also consider a second domain C R" with its own RKHS G induced by a kernel 1: 2 2 + IR with canonical feature map yp Z - G.
Mean embeddings provide a way to represent probability distributions vectors in RKHS 770 For a random variable X with distribution π on X, the mean embedding is defined as:
|1 EX)( H.
(2.1)
We write μi = My when the distribution is clear from context. For a characteristic RKHS 191] the map + 24 is injective, meaning each distribution has unique embedding. This property underlies the maximum mean discrepancy (MMD) distance of distributions, MMD(2, T) = ||,, - |4x #x||a l| (34).
Conditional mean embeddings extend the above to conditional distributions 70 Consider random variable Z on 2, defining PP X Z). Under standard assumptions (68, (7), the mean embedding of X
given Z= = can be expressed with the conditional mean embedding (CME) operator Uxiz G - H:
#x|z== == [x)(() Z = z] = Uxiz ppz) E Hl.
(2.2)
Notably, Uxiz is a linear operator. The linearity allows the of spectral methods, and also yields useful properties such as ,²x = Uxiz 1zz (1.T) via the law of total expectation.
MMD gradient flows enable sampling from distributions specified by mean embeddings (4) Given target distribution स with embedding Mr. the flow is driven by a time-dependent vector field (v:)
(x() = -V (Mp. Mxx (x) = - Vx (d(x), Hpp Mx/H:
(2.3)
where (p: ):2o is a probability path induced by the continuity equation a,p, + div(p.vt) = 0. Starting at any initial distribution Po. the path converges towards the target distribution: pt - as - xo. In practice, one can generate samples x ~ æ by starting at an initial sample xo ~ Po and numerically integrating the ordinary differential equation (ODE) dx, vr 0x(x,) dt to obtain X Pt-
2

---

3 Spectral Mean Flows
Problem setup We formulate sequence modeling in an operator-theoretic context using hidden Markov models (HMM). We consider sequence distribution X IINN p. where X E X, and model it as an HMM with hidden states Z ZI:N where Z" E Z. The hidden states evolve under iiveiinvara Markovian dynamics P[Z Iz),, and determine X I:N
X via local observation model P[X" The HMM can express any sequence distribution the hidden Markov chain sufficiently expressive.
Our approach begins by embedding the entire sequence distribution X1:N p in an appropriate Hilbert space, such that MMD flow be defined the space of sequences X.A.. We leverage the fact that tensor product HHN = Hi @ @ H of a characteristic RKHS HH remains characteristic under mild assumptions (Lemma B. l.. The canonical feature map of this tensor product space is I::N + o(x') @ ... @ o(x*), ), which naturally leads to the sequence mean embedding μo = #xxiN
Up El(((+)) co. ............ @o ((x*) E TH
(3.1)
With the above concepts in hand, we can extend the MMD flow (2.3) to operate on entire sequences.
Specifically, we define a time-dependent vector field (v+)e20 that induces a probability path (pr (pr) of sequences XP:VV Pt. converging towards the target distribution X1:N ~ p over time
Up (x = -V (o(x²) X
@ o(x"), Hp. - Mp)
(3.2)
Samples x!::N P. can then be generated by solving the ODE dx! = (x!!:N) dt from xb po-
While theoretically sound, the idea in its naïve form faces practical challenges. The mean embedding Mp (3.1) is a higher-order tensor whose size grows exponentially with the sequence length N. making training and sampling intractable. In addition, sampling with MMD flow typically slow, because its convergence pt + P is only guaranteed asymptotically as + 0o (].
To overcome the challenges, easy to train, and converges fast at sampling. In Section 3.1 we derive a scalable tensor network decomposition of mean embeddings (3.1) using spectral decomposition of linear operators underlying the HMM. In Section 3.2 we present an extension of MMD flow (3.2) to time-dependent RKHS that connects flow matching and converges within finite time I. Section [3. integrates these ideas into a neural sequence architecture and learning algorithm, and discusses implementation.
propose spectral mean flows, novel approach that is tractable,
3.1 Tractability with Spectral Decomposition
To make the MMD flows 3.2 with sequence mean embeddings 3..)) computationally tractable, we exploit the linear operator structure underlying the HMM by applying spectral decomposition.
Linear operators underlying HMMs underlying our HMM for sequences X XI:N, we can define the respective transition and observation CME operators U G + G and O G + Hi as follows (see Section 2 and Appendix AA.3)
For the conditional distributions P[Z*+1|Z" and P[X"|
Uę(z)) = El(p(Z Z" = z E G,
(3.3)
Os(Zz) E(o(X") 2" z] E H.
(3.4)
For the MMD flow generative process (XI:N) Je20 (..2). ) we make a modeling choice that the hidden
dynamics P| remain fixed over time and only the observation model P P[X?1 [ evolves
over t. Then we obtain time-dependent observation CME operator (St) where S, G + H:
Sip(2) = E(x(X)) | Z" = 2] H..
(3.5)
Mean embedding decomposition To make MMD flow (.2) tractable, our key idea is to decompose mean embeddings Up and Uy, into tractable forms by leveraging lincarity of CME operators U, O, S..
As first step, we decompose the mean embeddings Me and Up, into their respective observation operators O and S,, along with the mean embedding pzuv of the hidden Markov chain defined as:
U2IN == ) 8 ... & (( E ONN
(3.6)
We have two different notions of time: I...N over sequence length, and / 2 over generative process.
3

---

Ob(x²)
Ow(x* x²)
13
Ow(xx)
14
Ow(x*)
45
Ou(x²)
Figure I: A tensor diagram of the inner product (Ø(x*) @
6(xx))
Mp)
3.10 for N
showing how its evaluation involving a tensor Mp decomposes into matrix and vector multiplications.
Defining the N-fold tensor product operator by OSN = =0 10.00... :0: G +
H
we obtain:
He E @ 8 ((x*)]
(mean embedding 3.1
E [E(((X)) 0........ Blo(X*) |
(HMM structure)
=E[(O0((')) 8...8 (0+(Z)))
(CME operator
3.4
=O02) @ ..... x((
(linearity)
O UZI.N.
(3.7)
Similarly, we have Up. S Mzin with S G&N + H Formal proof is in Proposition B.2
The remaining challenge is decomposing the hidden chain embedding pzuxx We use the fact that the linearity of the transition operator U G + G enables a decomposition the following form:
U
Σ iEN
Aeh @ gh.
(3.8)
In practice, we rank-r decomposition by restricting to i [r]. where controls expressiveness.
possible choice of 3.8 is singular value decomposition (SVD) with singular values A and left/right singular functions Js. h4 Z + R. Or, assuming U is normal with a discrete spectrum, we can use eigenvalue decomposition (EVD) with eigenvalues X and eigenfunctions gi, 9i, While both are valid, we use EVD that efficiently encodes long-range interactions via oscillations 12
C.
Then, to decompose fzzı.& GON draw inspiration from result in Koopman operator theory that eigenfunctions are closed under pointwise multiplication, i.c., they form a multiplicative algebra, with products yielding new eigenfunctions whose eigenvalues multiply [66] [2] [0]. As we show, encoding this closure in finite-rank representations makes instantiating G SNN redundant, as nonlinear or higher-order interactions are already encoded in the algebra generated by the span of eigenfunctions.
Specifically. we make a technical assumption that the RKHS G is closed under pointwise multi- plication. f.g E G = f g G. This holds for Sobolev spaces HRR") of order > d/2 (.) which are characteristic and thus suitable for our framework. Under this assumption and leveraging a decomposition of U (3.8) with fi, h4 G, derive the following novel decomposition of pznx:
MZIN Σ big & Ar Wi
Xix 1ix & Aix hix:
(3.9)
$2.....
where b;, Wi.j are fixed elements in G determined by gts hy and the initial state distribution pzz-. Complete details and proof are given in Proposition B.6
Tractability With the decompositions in 3.7 and 3.9 we can convert the intractable compu-
tations in MMD flow (3.2) into tractable form. Specifically, in (3.2). the inner product between feature of sequence data x 1:N and mean embedding "p dramatically simplify (Corollary B.7
(6(x') 8 ... ((xN), (((x)) ... xx 0 0 pzzi)))))
Σ
Obi (xt) dis Owi (x >ix Ohiw (x*), (3.10)
$2...
where we use the shorthand (..:) (r)) For Mp: it is only required to substitute O with S.

---

Importantly, the reformulation in (3.10) eliminates the requirement to materialize exponentially large tensors in H when computing (3.2). Instead, the computation reduces a standard tensor network contraction (illustrated in Figure achieving linear complexity in sequence length N and polynomial complexity in rank r. This makes the approach tractable for realistic sequence lengths.
3.2 Faster Convergence with Time-dependent RKHS and Flow Matching
Having established tractable mean embeddings, we now address the remaining challenge of slow sampling convergence in standard MMD flows.
Recall MMD gradient flow (pt):20 for a target distribution н on X' within an RKHS HI (Section 2) The flow follows the steepest direction of the MMD between the current distribution Pt and target ı measured in H. [ Since the RKHS HI determines the geometric structure underlying the gradient flow via MMD, convergence behavior is essentially tied to the fixed H..
Time-dependent RKHS
To improve the convergence of MMD flow, natural way is to use time-dependent RKHS HH, [30]. This creates time-evolving geometry where the steepest direction of the MMD changes allowing more flexibility of the flow. Let (Hl))22 be a time-dependent family of characteristic RKHS, with canonical feature map (ó1)120 where Φi X + H,. The mean embedding of distribution æ is:
time,
Ux,t [x): HH,.
(3.11)
The MMD distance Hi, is accordingly given as MMD, T) 11 Px. Pe.
Then, MMD flow with time-dependent RKHS THI, can be obtained from standard construction 27 The flow defines probability path 120 driven by vector field (uz) via the continuity equation O.P. + div(p: = 0 that, at each time f, follows the steepest direction of the MMD measured in Hy,
(x) = Vx (xr Pp... ().
(3.12)
and samples x, ~ Pr. can be obtained via the ODE dx, = (x(x)) dt. More details are in Appendix B.3
Connection to flow matching With time-dependent RKHS, 330 used discriminator learning of de empirically improve the convergence of MMD flow. Instead, by connecting flow matching 66) we propose way to achieve convergence at finite time = rather than asymptotically as co.
We hinge the fact that we often find an analytic probability path qi and vector field u, defined in E [0, 1] that satisfy the continuity equation 0,q, + d = O and also qi = π This is the idea of flow matching (FM) (60), which allows us borrow established pairs of (gt, ut).
That is, for a choice of (gt, ue), we aim to identify the components of the MMD flow (.12) such that its vector field U2 regresses ut- If we have 8 U, and Po Go. it follows from the continuity equation that pt = qr and the MMD flow converges closely to the target p1 = q1 & N at
Here, an important caveat is that our vector field U2 is always a gradient field with zero curl. For the target equality , ue to be well-specified, the target 2e must also be gradient field. While this is not always true, we fortunately find that many commonly used vector fields in FM are gradient fields. This includes the optimal transport (OT) field and the diffusion-typee as we prove in Appendix B.4
Sequence instantiation We now apply the aforementioned ideas HMM for sequences X1:N p based on hidden states 21:N We consider time-dependent RKHS (Hl))e00,1I) By extending (3.12). the vector field in of the corresponding MMD flow X Pt can be written as follows:
Ny (x I:N) V (xx(x)) s @ th Mpi.t - Ppp)ER.
(3.13)
Then, assuming that the transition operator U: G >> G for TP - [ does not depend on time f, we extend the observation operators O,, S, C IHI for P| ~ and P[XP|Z"], respectively, as:
Ogp(z) E[o h =z] H4,
(3.14)
Szp(z) = Ex)((X)) Z" = z] E HI.
(3.15)
For example, gi is smoothing of with Gaussian filter with sufficiently small standard deviation.
5

---

With this, by extending the decomposition in 3.7 we obtain:
Ppp O M Mpr.t S Pzz.N.
(3.16)
And by extending the decomposition in (310 we obtain:
(xx(x') @ .. (x((x ( Mpp)
Σ (x')
Aig OeWigids (x²) ... Xix Oshin (x*) (3.17)
where we use the shorthand (.,-) (.- mus.. For Mpr.ts it is only required to substitute O₂ by St. This retains the tensor network structure, and hence the tractability of the MMD flow as in Section 3.1
3.3 Neural Parameterization and Learning
Having established spectral mean flows in RKHSs, we now develop a neural sequence architecture and learning algorithm. In (3.13) and (3.17). we can see that the operators O2, : + H,, together with the (Sobolev) RKHS elements b,, Way h2 G and X for i, [r]. determine the MMD flow vt-- We describe how to parameterize and learn these components as neural networks to pursue 8t- = u4-
Neural parameterization To develop a neural parameterization, we adopt neural tangent kernel (NTK) theory [42] as a theoretical footing that connects neural networks and Sobolev spaces, e.g., G. Specifically, the function space defined by multi-layer perceptron (MLP) in the infinite-width limit norm-equivalent to the RKHS a Sobolev kernel [8) In addition, any function in a Sobolev space can be approximated to arbitrary precision with a width MLP, where the theoretical scaling law predicts decreasing approximation error as the network size increases [86] Theorems I and 3].
These results justify. c.g., parameterizing each b, Wyy. h4 E G as scalar-valued MLP. Yet. instead of separating O,, St G -> IH,, we use an end-to-end parameterization (39 of each O,b,, Seb, ... E THI, as f-conditioned scalar- scalar-valued MLP on X. While this f-conditioning is inspired by [30], a difference that we directly parameterize elements of H. without using linear kernel approximation.
Learning and sampling Under the aforementioned neural parameterization, the MMD flow (313 becomes neural gradient field (=(:;)) 0) with learnable components 0 = (O,b,, S,b,, OpWay- Ax)-
Then, using a known pair (gt, un [e) of probability path qe and vector field u, satisfying qı 23 p.
we can use flow matching (FM) objective (60) to achieve the equality oz(:; = (.(-) (Section 3.2)
8 arg min (1 [xle 1:N () ue 1:N (x: [[]]
(3.18)
0
While 3.18) is intractable, thanks to compatibility with FM (Section 3.2 we can use an equivalent
tractable objective named conditional FM 600 using conditional counterparts ge I:N) (xx(-xx)))))
0" arg min px )). (x (xfo))))))
[xxx,(" 1:N, O) u, 1:N
I:N x
1|
(3.19)
This leaves us with a tractable learning algorithm for the neural MMD flow (-; 0). The algorithm is simulation- free, i.e., no ODE solving with (= () is necessary during training, and uses double back- propagation, i.e., ve(:; 9) internally invokes backpropagation during forward pass. In our experiments, we use the OT flow path and field (qt, ue) for training (Proposition [B.I] due to their simplicity.
After training, we can sample xt:N
②1 Rs p via the ODE dx!:N = (x..NN 1:N 00),, x1: Qo-
We discuss implementation of spectral mean flow used in our experiments. More details are in Appendix I.I First of all, we use rank-r EVD as discussed in Section 3.1 requiring complex eigenvalues E C and complex-valued MLPs that parameterize O,b,, S,b, ... C.
Implementation
To parameterize Obbe, O.wsy. Oyh, for j E [r]. we employ a shared feature extractor MLPo( X + Cd with feature dimension dy and combine it with readout heads B, L, Cry
Oxb,(x) R [B]; O,wiy (x) Rs [Ldiag(h))]j O,hx(x) RJ [hH]; = MLPo(x,)) t). (3.20)
We note that similar approaches are found in operator-theoretic learning dynamical systems (57. To parameterize Syb,, Stway, S,h,,
only change the feature extractor, denoted MLP
X + Cd2, We design the two MLPs to share most of their parameters through switching scheme:
MLPo(:,t) t) = MLP(., t) MLP s( t) == MLP(:, s, s,t).
(3.21)

---

where MLP(., : R + C is shared and o, R are trainable switching parameters.
The design of readout heads in (3.20) offers computational benefit. With the linearity of the tensor network (3.17). we can rearrange matrix multiplications to avoid materializing r x r matrix states Ldiag(h R. With this, we achieve O(Nr r²) space complexity of (3.17, avoiding the Nr² term.
To parameterize and initialize the complex-valle from the universal spectrum argument of 7 Proposition 3], which proposes to take eigenvalues from the uniform distribution on the complex unit disk. We employ the exponential parameterization of 775 Lemma 3.2] for this, which has an additional benefit of stabilizing the training.
parameters (A, B, L, R, H), we take inspiration
We lastly discuss the design of MLP( t). Since MMD flow (3.13) uses its gradient with respect to input, non-differentiable components such as ReLU can cause discontinuities [26). We design (,,((:,;t) to be fully differentiable with squared ReLU activation (87 and root mean square (RMS)
normalization [03]. For i-conditioning, we use sinusoidal embedding [96]. We find these design choices to work robustly while being simpler over alternatives such as SwiGL.U [85 and adaL.N 77]
4 Related Work
Operator theory Operator theory lifts the idea of vector spaces from points to functions [80]. The framework was initially developed for solving equations of infinitely many variables. This point of view is particularly advantageous in high-dimensional spaces, maps of functions describe the map of every points simultaneously instead of each individual point.
used the theory of linear operators to describe physical evolution equations [24] of probability measures and their dual observable functions. a modern context, linear operator theory for evolution equations often used to find practical linear algebra techniques to minimize, typically unsupervised, learning objectives. Following [40]. we can categorize them into four goals: revealing underlying structure 20 learning representations [64 57 81), modeling data distributions for prediction [43
[52 54
55), and lastly, generating points from a distribution over sequential data. To the best of our knowledge, the fourth task of generating new points from distribution over sequences has not been approached yet.
28
68
Mean embeddings Hilbert space embeddings of statistical quantities are long-term staple in statistical learning (6. In particular, RKHS mean embeddings have proven to be powerful tools [70]. Regressing the embedding of conditional expectations is studied with conditional mean embed- dings 88 35 94). This framework generalizes linear regression in RKHS from the space of vectors to the space of distributions, naturally giving rise so-called conditional expectation operators that map distributions to their expectations of function. They appear in the kernel Bayes rule (29] and are used in filtering [65], where a recursive application of conditional expectation operators over sequences first appeared in (88). Recently, conditional expectation operators were connected to the MMD [68], the transfer operators of Markov models [3 88 (5), and have been integrated with the linear operator theory for evolution equations [55] Again, none of the methods studies generation
with transfer operators, especially generating an entire coherent sequence, which is of our interest.
Generating samples from mean embeddings A key underlying idea of our work is generating samples from a distribution specified as a mean embedding. A classical approach for this is matching a mean embedding against those of a parametric family of distributions and recovering the minimizing distribution [90]. Another approach to estimate the solution to an inverse problem from data to obtain operators that directly recover densities from mean embeddings [84). Our work is more closely related approaches that evolve an empirical distribution to match the mean embedding, which
includes kernel herding [14 (58), and particularly MMD gradient flows 4 that form the basis of our
work. A downside is that fast convergence only guaranteed for a small class of kernels 44 30 which we improve by connecting to flow matching [60] via the continuity equation.
5 Experiments
We demonstrate spectral mean flows on two synthetic setups and generative modeling on a range of time-series datasets. Details of the experiments and supplementary results are in Appendix C
7

---

wis met antwel orgrotiee #f lowar setwork dicorpattion
16
MMD Gradient Flow RLB kernel
10P
42-
MMD Gradleen Flon w/ time-dependent RBF lorrel
11
Spectel
109
Sequence sength
Spectral Mron (Oum)
Figure 2: Peak GPU memory of
Figure 3: 2D checkerboard experiment. Left: Intermediate distribu-
inner product (x' -xx pi)
tions over sampling timesteps (zoom in for better view). Right: MMD between the intermediate and target distributions over sam- pling timesteps, measured with an RBF kernel of bandwidth
depending the use of tensor network decomposition (3117.
5.1 Synthetic Experiments
We first verify our claims on tractability (Section 3.1 focusing on how the tensor network decompo- sition in equation 13.17) and Figure contributss tractability and scalability under the hood. We run a numerical experiment of evaluating the inner product (x* 8 ... 6 x .1)) between rank-one tensor product x' @ x of vectors x" R and a higher-order tensor Rdd for d = 32. We measure the peak GPU memory usage across sequence lengths N, depending on the availability of tensor network decomposition in the form of (3.17). The results in Figure showing that without the tensor network decomposition, it is almost impossible to evaluate the inner product for N > 4.
Then, we verify our claims on faster convergence than MMD gradient flows (Section 3.2). For ease of visualization, we use 2D checkerboard datasct with scale -4.5, +4.5). We compare against MMD flows based on radial basis function (RBF) kernel k(x, x') exp(-0.5)
x' /=²) with bandwidth o = 1, and time-dependent RBF kernel (kg je [30] with (t)) = (1 t) 0.1t. For MMD flows, use 10,000 particles and 100 sampling steps with step size 1. For spectral mean flow, we treat each data as a sequence of length = 2 and use 100 sampling steps with a midpoint ODE solver. The results are in Figure 3 MMD flows explicitly minimize MMD measured by RBF kernels, and their samples move rapidly towards the target near / = 0. However, after some steps they face stagnation, suggesting much longer sampling required. In contrast, spectral
flow accurately converges to the target at f = 1. This supports our claims in Section 3.2 that spectral mean flow can converge arbitrarily close to the target at fixed time = 1, faster than naive MMD flow that requires t + Co. We note that, in the context of generative modeling, convergence at = 1 is sufficient and we need not worry about near i where spectral mean flow slightly increases MMD initially.
5.2 Time-Series Modeling
We demonstrate spectral mean flow on unconditional generation of time-series data. The baselines
include autoregressive models B3 15), generative adversarial networks (GANs) 67
25
100
variational autoencoders (VAEs) 21 104 72), and diffusion models
53
16
P01
71.
Regular time series For the first experiment, we follow [0O 101 use four real-world datasets Stocks, ETTh, Energy, fMRI and two simulated datasets Sines, MuJoCo of length-24 time series. We use four existing metrics to measure the quality of generated sequences. Context-Fréched distance (context-FID) score [44] measures the discrepancy between the distributions of features of the real and generated data encoded by a pre-trained TSVec model [02]. Correlational score [59] measures the absolute error between corsecoorel matrices of the real and generated data. Discriminative score [011] trains GRU classifier to distinguish between the real and generated data, and measures how close its accuracy is the chance level 50%. Predictive score [101] trains GRU next-step predictor on the generated data, and measures the mean absolute error (MAE) on the real data.
The main results are in Table 1 Spectral mean flow achieves the best metric in 18 out of 24 cases, showing that it is competitive with the state-of-the-art Diffusion-TS [OT]. The gain the quality metric is often significant, c.g., approach improves previous best context-FID, correlational score,

---

Table Time-series generative modeling.
Metric
Methods
Sines
Stocks
ETTh
MuloCo
Energy
fMRI
Ours
0.004+0001
0,008+.003
0.058=.007
0.018+0002
0,051+.009
0.116+.004
Ditfusion- TS
0.013+.001
0.169#.021
0.126: 007
0.015
0.113+0011
0.1184.000 007
miffT
0.006+.001
0,236+.074
0.299+0044
0.188+.028
0.279+.045
0.340+:015
Context-FID Score |
wallwave
0.014+.002
0.232+.032
0.873+0061
0.3934.041
1.0
0.244+.018
TinctiAN
0.101 0.103+.013
0.300#0013
0.563++052
0.767#. 103
1.292
VAF
0.307+0060
0.215+.035
0.805 186
0,251+.015
1.631 142
14.449+.969
Col-GAN
1.337=.068
0.408
0.980=.071
1.094+.0799
1.039+.028
7.813+.550
Ours
0.027+0012
0.010# 007
0.040=0015
0.173+.016
0.732:.107
0.7374.021
Ditfusion- TS
0.016=.005
0.010#0009
0.049=0013
0.188+.035
0.788+.075
1.252+.070
millT
0.017+0004
0.006+.002
0.067+0005
0.2184.031
1.158 095
1.501+.04
Conrelational
Score
Dillwave
0.022+0005
0.030+0020
0.175=0006
0.579+0018
5.001
..9224.049
TimeCiA
0.045+0010
0.063+.005
0.210=0006
0.886*.039
4.010#. 104
23.502+
VAT
0.131+0010
0,095+.008
0.111=020
0.388+0041
1.688+2226
17.296
Col-GAN
0.049#.010
0.087+.004
0.249=0009
1.042+0007
3.1642.061
26.82A2.449
Ours
0.0062.006
0.022+.013
0.027=.010
0,00S=.004
0.161# 021
0,136+.207
Diffusion- Ts
0.030#.006
0.085+.026
0.075=0007
0.012++006
0.1542.012
0,158+.020
mil'I
0.013+0006
0.097+.016
0.100=0007
0.154
0,445 004
0,245+ 051
Discriminative
Diflwave
0.017+0008
0.232+.061
0.190=2008
0.2033+.09
0.493++004
0.402#.029
Score 4
TimeGAN
0,011+0008 0.102+.021
0.114+0055
0.238
0.236+.012
0.484+0042
Time VAL
0.041#.044
120
0,209 058
0.230 102
0.499#.000
0,476 044
CoI-GAN
0.254+.137 0,230+.016
0.325=.099
0,426#.022
0.498#
0.492+.018
RNN-AR
0.495+.001
0.226++035
0.483+.004 004
Ours
0.093+.000
0.037+0000
0.125+0005
0.008+0001
0.251=0000
0.100++000
Diffusion-
0.095+.000
0.037+.000
0.121#0002
0.0072.001
0.251
0.100:.000
Tiee
0.093+.000 000
0.038 001
0.121=.004
0.010++001
0.252+ 000
0.100=0000
Diffwave
0.093=.000
0.047+0000
0.130=t001
0.013+0000
0.251#.000
0.101+0000
Predictive Score J.
IimeGAN
0.093+.019
0.038+.00 001
0.124#0001
0.025+.003
0.273+.004
0.1 0.12
Time VAL
0,093+.000
0.039+.000 000
0.126:2004
0.012+0002
.029220000
0.113+0003
Col-GAN
0. 1,0000000
0.047 001
0.1 0.00
0.068+t009
0.2599.000
081
RNNN-K
0.150+.022
0,038+ 001
0.315+.000 005
Original
0.094+.001
0.036*.001
0.121
0.007+0001
0,250+.003
0.090+0001
Table 2: Time-series modeling in larger model regime.
Table 3: Long time-series modeling.
Metric
Methods
Sines
Stocks
MuloCo
Metric
Methods
FRED-MD
NNS Daily
Ours
0.002+.000
0.004s.001
0.01 001
Ours
0.019
0.006imma..
Context-FID Score
SDFormer-Ae AR
0.008&0001
0.006: 001
0.0082.000
Marginal Score J
Imagen Time
0.022tn.a.
0.009
0.010+0022
0.034e.008
0.030+.003
0.022ana.
0.9074m.a.
Imagei
0.009
0.011a.002
0.017..002
SaSIIM-AR
0.048en.a
0.020tn.a..
Ours
0.0072.008
0.012s.013
0.019:1009
Ours
1.33&s753
0,950ts.257 257
Discriminative
SDFomee--AN
0.016+.010
0.0062:006
0.009=0006
Classification Score 1
TmagentTime
0.7551.343
0.560a.
Scare
S
SDFormee-hh
0.008+.004
0.020:.011
0.025:.007
S4
0.544mm.a.
0.636tm.a.
0.016+:010
0.010=.007
0.011+0088
SaStIME-AR
0.00
0.045tm.a.
Qurs
0.093=.000
0.037s.000
0.008 001
Ours
0.030+.006
0.539+. 196
Predictive
SDFomer-AR
0.0932.000
0.037s:000
0.008::002 002
Predictive Score
Imagen Time
0.034+.020
0.5841. 188
Soore
SDFW
0.093+.000
0.037=0000
0.0017+0011
0.037am.a.
0.241am.0.
Inagon Hane
0.095:.000
0.0372.000
0.033:.002
SANDM-AAR
0,232en.a.
0.8494m.a.
and discriminative score in ETTh by around 50% on average, while predictive score does not lag much behind Diffusion- TS. The results are notable considering that previous methods often make use of carefully designed auxiliary losses specialized for time-series, e.g., Diffusion-TS uses Fourier contrast, our models only use standard flow matching loss (3.19) for training, which illustrates the advantages of principled architecture design. Specifically, spectral expansion induces the use of complex-valued parameters, and tensor network decomposition induces a multiplicative parameter sharing in the length direction (Section 3.I). offering a natural inductive bias to handle periodicity without specialized loss or components.
domain loss to improve the handling of periodicity.
To provide further comparisons in a larger parameter regime, we conduct an additional experiment following [15] [ comparing against larger models than ones in Table We use Sines, Stocks, and MuloCo considering resource constraints. The results are in Table Spectral mean flow achieves the best in 6 out of cases, showing that it is competitive with the state-of-the-art SDFormer-AR [[5] and generally outperforms SDFormer-M 15 and imagenTime (71. Similarly as in Table I we note that SDFormer-AR uses a sophisticated two-stage training of vector-quantized codebook and then an autoregressive model on top of it, while our models are end-to-end trained with flow matching.

---

Table 4: Irregular time-series modeling based on Stocks dataset, evaluated with discriminative score 4.
Table 5: thysics-informen model-
Tisk
Methods
0% Drop
30% Drop
50% Drop
70% Drop
ing of a nonlinear pendulum.
Ours
0.001920008 00S
0.020+.011
0.019a.018
0.015+.0077
Koopman VAE
0,021 022
0.109
0.0674.038
0.0494.052
Methods
Corr. Score
-+ Regular
Irregular
GLAGN
RUUAN
ImeGAX
0.102s.021
0.0774.081
0.1962.027 0.43620064
0.4
0.251+0097
0.477a.021
0.478=009
0.265=0073
0.485a.022
0.2811.086
0.2304.053
Ours w/ stability loss
KoVAE w/ stability loss
0.0005=.0000
0.030000004
CANN-GAN
0,399
0.500sa000
0,5002.000
0.500040000
Ours wlo stability loss
0.0029t.0008
RNN-AR
0,2262.035 0.305
0.30820010
0.3172..19
KoVAE w/o stability loss
0.0040L0000
Irregular
Ours
0,0092.0088
0.049
0.044:.017
0.1
Imegular
Koopman VAF
0.021+:022
0.227+0096
0.211=0078
0.187++075
Long time series For additional demonstrations of modeling longer time series than in Table 1 we
use FRED-MD and NNS Daily from the Monash repository (32), containing time series of lengths 728
and 791, respectively. Our models are designed as in Section 3 equipped with time-delay observables
used in operator-theoretic methods for dynamical systems 49 We use three metrics following [7T.
Marginal score [59] measures the absolute error between the empirical probability densities of the
real and generated data. Classification score [04 trains an S4 (37] classifier to distinguish between
the real and generated data, and measures its loss on the generated data. Predictive score [T04] trains
an S4 10-step predictor on the generated data, and measures the prediction error on the real data. The results are in Table Spectral mean flow shows strong result, achieving the best metric in 5 out of cases and generally outperforming the previous best methods Imagen'imme 77 and LS4 [04
Irregular time series For further demonstrations of generality, we run experiments on irregularly sampled time series. We obtain 3 irregularly sampled datasets from Stocks by randomly dropping 30%, 50%, and 70% of the observations, following 45 72.. We consider two tasks, irregular > regular; generating the full time series including the missing timesteps, and irregular - irregular: generating only the irregularly sampled time series. The former is standard in literature 45] 72]. and the latter provides proxy task of modeling informatively sampled time series where sampling can be determined by the system state. Baselines include Koopman VAE (72). state-of-ther
operator- based method for irregular time series. The results are in Table 4 Spectral mean flow achieves the best discriminative metrie, supporting its generality to handle irregular time series in both settings.
Physics-informer modeling Lastly, we demonstrate a benefit of being an operator method: incor- poration of physics-based prior knowledge. While interpreting eigenfunctions is not direct due to our end-to-end design, this can be done via explicit spectral regularization thanks linearity of the tensor network. We test this using problem from [72], where we consider a nonlinear pendulum governed by an ODE of angular displacement f from an equilibrium, Ø + 9.8 sin 0 0, á())- 0. As the pen- dulum stable and conservative, the physical knowledge is that eigenvalues of underlying operator have non-pon non-positive real part and some have zero real. In Koopman VAE (72). this stability constraint built-in using a loss on matrix A for hidden states Z Az,, specifically [s| 112 + |s2 - 1/2 where (s1, $2) are the largest eigenvalues of A. As our model is governed by products of matrix hidden states O,w(x"), S,w((") 10. (3.20). we incorporate this in an end-to-end manner via the same spectral loss on the matrix states. The results are in Table 5] With stability loss, we see clear improvement over Koopman VAE that uses the same constraint. This implies incorporating physical knowledge into our model is possible, leading to better match with the true data- settreenng process.
Conclusion
We proposed a new algorithm for sequence modeling based on operator- theoreetc interpretation of a hidden Markov model. Instead of implementing its stochastic recurrence directly, we considered embedded distributions in Hilbert spaces, which enabled us to use powerful linear-algebraic tools including spectral decomposition to derive a tensor-network based neural architecture, as well as a sampling procedure based on time-dependent MMD gradient flow paired with flow matching. On synthetic setups and time-series modeling datasets, we verified our theoretical claims and observed performances competitive with the ttate-of-the-art.
10

---

Acknowledgments This work was in part supported by the National Research Foundation of Korea RKS-2024-00351221 and RS R8-2024-004361155) and the Institute of Information & Commu- nications Technology Planning & Evaluation (IITP) (RS-2024-005099279 RS-2022-11200926, and RS-2022-11220959) funded by the Korean government (MSIT). Max Beier is supported by the DAAD programme Konrad Zuse Schools of Excellence in Artificial Intelligence, sponsored by the German Federal Ministry of Education.
References
[1] R. A. Adams and J. Fournier. Sobolev spaces, volume 140. Elsevier, 2003. (pages 22
[2] T. D. Ahle. The tensor cookbook, 2024. Version:
September. 2024. (page
[3] Ail-Sahalia, L. Hansen, and J. A. Scheinkman. Operator methods for continuous- markov processes. Handbook of financial econometricx: tools and techniques, pages 1-66, 2010. (page 37)
[4] M. Arbel, A. Korba, A. Salim, and A. Gretton. Maximum mean discrepancy gradient flow. In NeurIPS,
2019. (pages 2
27
[5] Y. Bengio, R. Ducharme, P. Vincent, and C. Jauvin. A neural probabilistic language model. Journal of machine learning research, 137 2003. (page T
[6] A. Berlinet and C. Thomas- Agnan. Reproducing kernel Hilbert spaces in probability and statistics. Springer Science & Business Modia, 2011. (page 7
[7] P. Bevanda, M. Beier, A. Lederer, S. Sosnowski, E. Hüllermeier, and S. Hirche. Koopman kernel regression. In NeurlPS, 2023. (pages
[8] A. Bietti and R. Bach. Deep equals shallow for relu networks kernel regimes. In ICLR, 2021. (pages
(2
[9] B. Boots, G. J. Gordon, and A. Gretton. Hilbert space embeddings of predictive state representations. In UAI, 2013, (page 7
[10] N. Boullé and M. J. Colbrook. Multiplicative dynamic mode decomposition. SIAM Journal on Applied Dynamical Systems, 24(2): (122):19:: 194--1688,..19999911111191119998 2025. (page 4
[11] J. C. Bridgeman and C. I. Chubb. Wavd-w and interpretive dance: an introductory course on tensor networks. Journal of physics A: Mathematical and theoreticul, 50/22)222001. 2017. (page 4
[12] L.. Brunton, M. BudiSić, E, Kaiser, and N. Kutz. Modern koopman theory for dynamical systems.
arXiv, 2021. (pages
27
[13] C. Carmeli, F. De Vito, A. Toigo, and V. Umanitá. Vector valued reproducing kernel hilbert spaces and universality. Analysis and Applications. 8(01):1961,, 61, 2010. (page 23)
[14] Y. Chen, M. Welling, and A. Smola. Super-samples from kernel herding. In UAI, 2010. (page
[15] Z. Chen, F. SHIBO, Z. Zhang, X. Xiao, X. Gao, and P. Zhao. Sdformer: Similarity discrete transformer for time series generation. In NeurlPS, 2024. (pages R 9 32)
[16] A. Coletta, Gopalakrishnan, D. Borrajo, and Vyetrenko. On the constrained time-series generation problem. In NeurlPS, 2023. (page(8)
[17] G. Daniel, J. Gray, al. Opt'_einsuuma python package for optimizing contraction order for einsum like expressions. Journal of Open Source Software, 3(26):753, 2018. (page 30)
[18] S. Das and D. Giannakis. On harmonic Hilbert spaces compact abelian groups. Journal of Fourier Analysis and Applications, 29(1):12, 2023. (page 22
[19] S. Das, D. Giannakis, and M. Montgomery. Correction to: On harmonic Hilbert spaces on compact abelian groups Journal of Fourier Analysis and Applications, 29(6)677, 2023. (page[2)
[20] M. Dellnitz O. Junge. On the approximation of complicated dynamical behavior. SIAM Journal on Numerical Analysis, 36(2):491-515, 1999. (page
[21] A. Desai, C. Freeman, Z. Wang, and I. Beaver. Timevae: A variational auto-encoder for multivariate time
series generation. arXiv, 2021. (page B
[22] K. Diederik and B. Jimmy. Adam: A method for stochastic optimization. arXiv, 2014. (pages 31 32
[23] J. L. Elman. Finding structure in time. Cognitive science, 14(2):179-211, 1990. (page D
[24] K.-J. Engel and R. Nagel. One-Parameter Semigroups for Linear Evolution Equations, volume 194 of Graduate Texts in Mathematics. Pnriner-Verrlg. New York, 2000. (page
[25] C. Esteban, S. L. Hyland, and G. Rätsch. Real-valued (medical) time series generation with recurrent conditional gans. arXiv, 2017. (page
11

---

[26] C. Etmann. A closer look at double backpropagation. arXiv,, 2019. (page
[27] L. C. Ferreira and J. C. Valvaria-Guevara. and applications to pdes. Monatshefte für Mathematik, 185:231-268, 2018. (pages
Guevara. Gradient flows of time-dependent functionals in metric spaces
[28] Froyland, Junge, and P. Koltai, Estimating long-term behavior flows without trajectory integration: The infinitesimal generator approach. SIAM Journal on Numerical Analysis, (1(1): 2013. (page
[29] K. Fukumizu, L. Song, and A. Gretton. Kemel bayes" rule. arXiv, 2011. (pages
[30] A. Galashov, V. D. Bortoli, and A. Gretton. Deep MMD gradient flow without adversarial training. arXiv,
2024. (pages S 6
[31] A. H. Galib, P.N. Tan, and I., Luo. Fide: Frequency-inflated conditional diffusion model for extreme- aware time series generation. In NeurlPS, 2024.
[32] R. Godahewa, C. Bergmeir, G. Webb, R. J. Hyndman, and P. Montcro-Manso. Monash time series
forecasting archive. arXiv, 2021. (pages TO
32)
[33] K. Goel, A. Gu, Donahue, and C. Ré. It's raw! audio generation with state-space models. In International conference on machine learning, pages 7616-7633, PMLR, 2022. (page
[34] A. Gretton, K. M. Borgwardt, M. Rasch, B. Schölkopf. and A. J. Smola. A kerne two-sample test. J. Mach. Learn. Res., 13:723 773, 2012. (pages 23
[35] S. Grünewälder, G. Lever, L. Baldassarre. S. Patterson, A. Gretton, and M. Pontil. Conditional mean embeddings regressors. In ICML, 2012. (page
[36] A. Gu and T. Dao. Mamba: Lincar-time sequence modeling with selective state spaces. arXiv, 2023. (page D
[37] A. Gu, K. Goel, and C. Ré. Efficiently modeling long sequences with structured state spaces. In ICLR, 2022. (pages 1 10)
[38] J. Hertrich, Wald, F. Altekrüger, and Hagemann. Generative sliced MMD flows with riesz kernels. ICLR, 2024, (page())
[39] D. Z. Huang, N. H. Nelsen, and M. Trautner. An operator learning perspective on parameter-to-observvee maps. Foundations of Data Science, 7(1):163-225, Mar. 2025. Publisher: Foundations of Data Science. (pagelo) (6)
[40] A. Hyvärinen, I. Khemakhem, and H. Morioka. Nonlinear independent component analysis for principled disentanglement in unsupervised deep learning. Putterns, 4(10), 2023. (page
[41] P. Inzerilli, V. Kostic, K. Lounici, P. Novelli, and M. Pontil. Consistent Long- Term Forccasting of Ergodic Dynamical Systems. arXiv, 2023. (page
[42] A. Jacot, F. Gabriel, and C. Hongler, Neural tangent kernel: Convergence and generalization in neural networks. In NeurlPS, 2018. (page(6)
[43] H. Jaeger. Observable operator models for discrete stochastic time series. Neural Computation, 12:1371- 1398, 2000. (page 7
[44] Jeha, M. Schneie,, P. Mercado, R.-. Nirwan, S. Kapoor, Flunkert, Gasthaus, and
T. Januschowski. Psa-gan: Progressive self attention gans for synthetic time series. arXiv, 2021. (page S
[45] J. Jeon, J. Kim, H. Song, S. Cho, and N. Park. Gt-gan: General purpose time series synthesis with generative adversarial networks. In NeurlPS, 2022. (pages I0
[46] K. Jordan, Y. Jin, V. Boza, Y. Jiacheng, F. Cesista, I.. Newhouse, and J. Bernstein. Muon: An optimizer for hidden layers in neural networks, 2024. (page
[47] M. I. Jordan. Attractor dynamics and parallclism in connectionist sequential machine. In Proceedings of the Annual Meeting of the Cognitive Science Sociery, volume 8, 1986. (page
[48] O. Kallenberg. Foundations of Modern Probability, volume 99 of Probability Theory and Stochastic Modelling. Springer International Publishing, Cham, 2021. (page 23)
[49] M. Kamb, F. Kaiser, L.. Brunton, and N. Kutz. Time-delay observables for koopman: Theory and
applications. SIAM Journal on Applied Dynamical Systems, 19(2)886 -177 2020. (pages) 10 33
[50] T. Karras, M. Aittala, J. Lehtinen, J. Hellsten, T. Aila, and Laine. Analyzing and improving the training
dynamics diffusion models. CVPR, 2024. (pages 31 32)
[51] S. Klus, I. Schuster, and K. Muandet. Eigendecompooitioons of transfer operators in reproducing kernel hilbert spaces. J. Nonlinear Sci., 30 3:23333333 2020. (pages [
[52] A. Kolmogoroff. Über die analytischen methoden in der wahrscheinlichkeisrrchnnng Mathematische Annalen, 104:415-455 1931. (page
12

---

[53] Z. Kong, W. Ping. J. Huang, K. Zhao, and B. Catanzaro. Diffwave: A versatile diffusson model for audio synthesis. arXiv, 2020. (page
[54] B. O. Koopman. Hamiltonian systems and transformation in hilbert space. Proceedings of the National Academy of Sciences, 17(5):315 -18. 1931. (page 7
[55] Kostic, P. Novelli, A. Maurer, C. Ciliberto, L. Rosasco, M. Pontil. Learning dynamical systems via Koopman operator regression in reproducing kernel hilbert spaces. In NeurIPS, 2022. (pages |
[56] V. R. Kostic, K. Lounici, H. Halconruy. T. Devergne, P. Novelli, and M. Pontil. Laplace transform based
low-complexity learning of continuous markov semigroups. arXiv, 2024. (page 37
[57] V. R. Kostic, P. Novelli, R. Grazzi, K. Lounici, and M. Pontil. Learning invariant representations of time- homehhoogeeneos stochastic dynamical systems. In ICLR, 2023. (pages
[58] S. Lacoste-Julien, F. Lindsten, and F. R. Bach. Sequential kernel herding: wrllk-wolfe optimization for particle filtering. In AISTATS, 2015. (page 7
[59] Liao, Ni, L. Szpruch, Wiese, M. Sabate- Vidalee, and B. Xiao. Conditional sig- wasserstein gans for time series generation. arXiv, 2023. (pages [0]
[60] Y. Lipman, R. T. Q. Chen, H. Ben-Hamu, M. Nickel, and M. Le. Flow matching for generative modeling. In ICLR, 2023. (pages 2
29
[61] Y. Lipman, M. Havasi, P. Holderrieth, N. Shaul, M. Lc, B. Karrer, R. T. Chen, D. Lopez-Paz, H. Ben-Hamu, and I. Gat. Flow matching guide and code. arXiv, 2024, (page 30)
[62] J. Liu, M. Saleh, E. Pot, B. Goodrich, R. Sepassi. L.. Kaiser, and N. Shazeer. Generating wikipedia by summarizing long sequences. arXiv, 2018. (page [
[63] I. Loshchilov and Hutter. Decoupled weight decay regularization. arXiv, 2017. (page 32
[64] A. Mardt, L. Pasquali, H. Wu, and F. Noć. Vampnets for deep lcarning of molecular kinetics. Nature communications, 91) 2018. (page
[65]
O' Callaghan, and F. Ramos. Multi-modal estimation with kernel embeddings for learning motion models. In IEEE International Conference on Robotics and Automation, 2013. (page 7
McCalman,
[66] I. Mezić. Applied koopmanism. Nonlinear Dynamics, 41(1):533-546, 2013. (page
[67] O. Mogren. C-run-gan: Continuous recurrent neural networks with adversarial training. arXiv, 2016.
(page 8
[68] M. Mollenhauer and Koltai. Nonparametric approximation of conditional expectation operators. arXiv, 2020. (pages ! 2
27
[69] M. Montgomery and D. Giannakis. An algebra structure for reproducing kernel Hilbert spaces. Banach Journal of Mathematicul Analysis. 19(1):11, 2025. (page 22)
[70] K. Muandet, K. Fukumizu, B. Stiperumbudur, B. Schölkopf, al. Kernel mean embedding of distributions:
A review and beyond. Foundations and Trends in Machine Learning, 10(1-2):1 141, 2017. (pages
23
[71] I. Naiman, N. Berman, 1. Pemper, I. Arbiv, G. Fadlon, and O. Azencot. Utilizing image transforms and
diffusion models for generative modeling of short and long time series. In NewrlPS, 2024. (pages 8 9 10
32
33
[72] I. Naiman, N. B. Erichson, P. Ren, M. W. Mahoney, and O. Azencot. Generative modeling of regular and irregular time series data via koopman vacs. arXiv,, 2024. (pages] I
10 33)
[73] Nishiyama, A. Afsharinejad, S. Naruse, B. Boots, and L. Song. The nonparametric kernel bayes
smoother. In Artificial Intelligence and Statistics, pages 547-555. PMIR, 2016. (page 37
[74] Novelli, M. Pratticò, M. Pontil, and C. Ciliberto. Operator world models for reinforcement learning. arxiv, 2024. (pages) 23
[75] A. Orvieto, S. L. Smith, A. Gu, A. Fernando, Ç. Gülççhhre, R. Pascanu, and S. De. Resurrecting recurrent
neural networks for long sequences. In ICML, 2023. (pages 7 27 3T)
[76] J. Park and K. Muandet. A measure-theoretic approach to kernel conditional mean embeddings. In NeuriPS, 2020. (page 23)
[77] N. Parnichkun, S. Massaroli, A. Moro, T. H. Smith, R. M. Hasani, M. Lechner, Q. An, Ré,
H. Asama, S. Ermon, T. Suzuki, M. Poli, and A. Yamashita. State-free inference of state-space models:
The transfer function approach. In JCML, 2024. (page
[78] A. Paszke, S. Gross, F. Massa, Lerer, Bradbury, G. Chanan, T. Killeen. Lin, N. Gimelshein.
L. Antiga, A. Desmaison, A. Köpf, E. Yang, Z. DeVito, M. Raison, A. Tejani, S. Chilamkurthy, B. Steiner, L. Fang, J. Bai, and S. Chintala. Pytorch: An imperative style, nigh-performance deep learning library. arXiv, 2019. (page 30)
13

---

[79] W. Peebles and S. Xie. Scalable diffusion models with transformers. In ICCV, 2023, (pagel
[80] W. Rudin. Functional analysis. International series in pure and applied mathematics. McGraw-Hilll New York, 2nd ed edition, 1991. (pages 7 27)
[81] J. Ryu, X. Xu, H. S. M. Erol, Y. Bu, L. Zheng, and G. W. Wornell. Operator svd with neural networks via nested low-rank approximation. ICML, June 2024, (page 7
[82] F. Schmidt and T. Hofmann. Deep state space models for unconditional word generation. In NearIPS,
2018. (page D
[83] B. Schölkopf and A. J. Smola. Learning with kernels: support vector machines, regularization, optimiza- tion, and beyond. MIT press, 2002. (page 2
[84] I. Schuster, Mollenhauer, Klus, Muandet. Kernel conditional density operators. In AISTATS. 2020. (page
[85] N. Shazeer. Glu variants improve transformer. arXiv, 2020. (page
[86] J. W. Siegel. Optimal approximation rates for deep relu neural networks on sobolev and besov spaces. Journal of Machine Learning Research, 24(357):1-52, 2023. (pages 6 22
[87] D. R. So, W. Mańke, H. Liu, Z. Dai, N. Shazeer, and Q Le. Primer: Searching for efficient transformers for language modeling. arxiv, 2022. (pages
30)
[88] L. Song, B. Boots, S. M. Siddiqi, G. J. Gordon, and A. J. Smola. Hilbert space embeddings of hidden markov models. In JCML, 2010. (page 7
[89] I.. Song, J. Huang, A. Smola, and K. Fukumizu. Hilbert space embeddings of conditional distribations with applications to dynamical systems. In ICMI, 2009. (pages
[90] Song. Zhang, A. Smola, A. Gretton, and B. Schölkopf. Tailoring density estimation via reproducing kernel moment matching. In JCML, 2008. (page
[91] B. K. Sriperumbudur, K. Fukumizu, and G. R. G. Lanckriet. Universality, characteristic kernels and RKHS embedding measures. J. Mach. Learn. Rex., 12:2389-2410, 2011. (pages
[92] I. Sutskever, J. Martens, and G. Hinton. Generating text with recurrent neural networks. In JCML, 2011. (page 33)
[93] Z. Szabó and B. K. Sriperumbudur. Characteristic and universal tensor product kernels. J. Mach. Learn. Res., 18:233:1 11:33::-233:22 233:299, 2017. (page 24
[94] Szabó, K. Sriperumbudur, Póczos, and Gretton. Learning theory distribution regression. Journal of Machine Learning Research, 17(152):1-44 40, 2016. (page
[95] F. Takens. Detecting strange attractors in turbulence. In Dynamical Systems and Turbulence, Warwick
/ 1980: proceedings of a symposium held at the University of Warwick /979/80, pages 366-381. Springer, 2006. (page 33)
[96] A. Vaswani, N. Shazeer, N. Parmar, J. Uszkorcit, L. Jones, A. N. Gomez, L. Kaiser, and I. Polosukhin. Attention is all you need. In NeurlPS, 2017. (pages 30)
[97] C. Villani et al. Optimal transport: old and new, volume 338. Springer, 2008. (page
[98] R. Williams and D. Zipser. A learning algorithm for continually running fully recurrent neural networks. Neural computation, 1()):270-280, 1989. (page
[99] T. Xu, L. K. Wenliang, M. Munn, and B. Acciaio. Cot-gan: Generating sequential data via causal optimal transport. In NeurlPS, 2020. (page
[100] J. Yoon, D. Jarrett, and M. van der Schaar. Time-series generative adversarial networks. In NeurlPS,
2019. (pages 8
32
[101] Yuan
Qiao. Diffusion-Is: Interpretable diffusion for general time series generation. arXiv, 2024.
(pages 8
[102] Z. Yue, Y. Wang. J. Duan, T. Yang, C. Huang, Y. Tong, and B. Xu. Ts2vec: Towards universal representation of time series. In AAAI, 2021. (page X
[103] B. Zhang and R. Sennrich. Root mean square layer normalization. arXiv, 2019. (pages
30
[104] Zhou, M. Poli, Xu, S. Massaroli, and S. Ermon. Deep latent state space models for time-series generation. In /CML, 2023. (pages] T0p
[105] W. Zucchini and I. MacDonald. Eruptions of the old faithful geyser, 2009. (page 33
14

---

NeurIPS Paper Checklist
1. Claims
Question:
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: The abstract and introduction summarizes the contributions and the scope of the work, including both theoretical and empirical aspects.
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 contri- butions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings,
is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification:
The limitations of the work are discussed in Appendix D
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 a separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g.. independence assumptions, noiseless settings, model well-specificattin, asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. 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, a facial recognition algorithm may perform poorly when image resolution is low images taken in low lighting. Or specch-to-text system might not be used reliably to provide closed captions for online lectures because fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimess.
While the authors might fear that complete honesty about limitations might be used by reviewers grounds for rejection, worse outcome might be that reviewers discover limitations that aren't it acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play important role developing norms that preserve the integrity of the community. Reviewers will be specifically instructed not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?
Answer: [Yes]
15

---

While we state the theoretical results informally in the main text for presentation, the full set of assumptions, formalization of the results, and complete proofs are provided in the appendix.
Justification:
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced 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 short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix 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 experimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
Justification: We included the implementation details and descriptions of the experiments in the main text. More details are provided in Appendix
Guidelines:
The answer NA means that the paper does not include experiments.
the paper includes experiments, 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.
the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible 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, if the contribution necessary to either make it possible for others to replicate the model with the same dataset, or provide access the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to hosted model (e.g., the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
specific model and empirical evaluation,
may be
While NeurIPS does not require releasing code, the conference does require all submissions provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily new algorithm, the paper should make it clear how to 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 a new model (e.g., large language model), then there should either be way access this model for reproducing the results way reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reprodueibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but should be possible for other researchers to have some path to reproducing verifying the results.
5. Open access to data and code
16

---

Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, described in supplemental material?
Answer: [Yes]
Justification:
We provide open access to data and code.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https::/nips.cc///
guiiees/CodeeubniissiooPPP....
for more details.
While we encourage the release of code and data, we understand that this might not be possible, "No" acceptable answer. Papers cannot be rejected simply for including code, unless this is central to the contribution (e.g., for a 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 hatess/oooeSmmbisiinoooorro
/
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 proposed method and baselines. If only a subset of experiments are reproducible, they should state which ones are omitted from the script and why.
experimental results for the new
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (c.g., data splits, hyperparame- ters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: The details of the experiments are provided in the main text and Appendix
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or 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 error bars from repeated tests in Section
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confidence 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 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 a library function, bootstrap, etc.)
The assumptions made should be given (c.g., Normally distributed errors).
17

---

It should be clear whether the error bar is the standard deviation or the standard error of the mean.
It is OK to report I- I-sigma error bars, but one should state it. The authors should preferably report 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not show in tables figures symmetric error bars that would yield results that are out of range (e.g. negative error 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 in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification:
The resources are described in Appendix
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type 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 well 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 make it into the paper).
9. Code of ethics
Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics
Question:
Answer:
[Yes]
Justification:
We conform with the NeurIPS code of ethics.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there is special consideration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [No]
Justification: The work mainly contains theoretical results along with time-series generation experiments, and immediate societal impacts are not expected.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
the authors answer NA or No, they should explain why their work has no socictal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess considerations (e.g.. deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
18

---

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 any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, 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 as 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 or unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (c.g.. gated release models, providing defenses in addition to attacks, mecha- nisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Question:
Does the paper describe safeguards that have been put place for responsible release of data models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Answer: [NA]
Justification:
The paper does not release data or models that have high risk for misuse.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recogmize that providing effective safeguards is 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
Question: 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?
Answer: [Yes]
Justification:
We cited the original sources of all datasets, evaluation setups, and models.
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 is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from particular source (c.g., website), the copyright and terms of service of that source should be provided.
assets are released, the license, copyright information, and terms use in the package should be provided. For popular datasets, paperswi thcode. com/ com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset.
For existing datasets that are re-packaged, both the original license and the license of the 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
19

---

Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification: Details will be provided within the code.
Guidelines:
The answer NA means that the paper does not release new assets. Researchers should communicate the details of the dataset/code/model part of their sub- missions 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 is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text instructions given to participants and sereenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor 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 contribution the paper involves human subjects, then as much detail as possible should be included in the main paper.
According 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
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements your country institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor 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. 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 to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Does the paper describe the usage of LLMs if it an important, original, or non- standard component of the methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Question:
Answer: [NA]
20

---

Justification:
The core method development does not involve LLMs as any important, original, or non-standard components.
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. http://urripssccc/ooo/eeeeeeeeee 2025 /LLM) for
Please refer to our LLM policy what should or should not be described.
21

---

The appendix is organized as follows. In Appendix A
we introduce the mathematical background. In Appendix B, we formally state and prove theoretical results in the main text. In Appendix C we we discuss limitations and future work.
provide details of experiments. In Appendix D
A
Mathematical Background
We provide overview of the necessary mathematical background, and refer the reader to Mollen- hauer& & Koltai (2020) [68] and Arbel et al. (2019) 4 for more details.
A.1 Problem Setup
For a Polish space S, we define P(S) as the set of probability distributions on S. For any т E P(S),
we define I,2 (xr) the set of real-valued Lebesgue square integrable functions on S with respect to т. Any closed subset X' C IRd, and its product x for finite N E N. are Polish spaces.
For closed X' C R4 and Z Rm we consider sequence of random variables X1:N (X1 X on X' and model it as hidden Markov model (HMM) with hidden states Z = (Z) ..... Zd) on Z. The joint distribution of hidden states Z1:N is determined by a time-invariant state transition model P[Z"+1|Z") and determines X1:N through time-invariant local observation model |[X" Z").
Formally, the conditional distributions [ and [[X]|Z"] are specified by the respective Markov kernels u and o that define u(z,
E ((2) and o(2,) :) E P(()) for every z E 2, such that:
PP E EA|
nn = z) = u(z, dz') u(z,
(A.1)
P[X" E B
Z = z] = J
o(z, dx) = olz, B):
(A.2)
for all measurable sets A C 2 and B C X.
A.2 Reproducing Kernel Hilbert Space (RKHS)
Let IH be an RKHS induced by a positive semi-definite kernel function k: X' x X + R. The space HI is vector space consisting of functions f : X + R, equipped with an inner product (-, -) and corresponding norm || lla given by ||A = (f, f) The canonical feature map o X + H of the RKHS is defined by Φ(x)= = k(,x0 x). The reproducing property of H states that f(x) = (Ø(x), f)
for any f E H. We denote the space of Bilbert-Schminn operators from HI to H as S2 (H).
To handle conditional dependencies within HMM, also consider a second RKHS G induced by a kernel / : 2 x 2 - R with the corresponding canonical feature map ş Z - G. We denote the space of -cchittt- operators from G to H as S2 S2, IHI).
We use some assumptions on the RKHS H (or G) and the corresponding kernel (or l). The first five assumptions allow universal approximation of conditional distributions via linear operators on H [68 Section 4.5). The sixth one allows decomposing an HMM embedded in into a tensor network.
(A) (Separability). The RKHS H is separable. This holds if X' is Polish and is continuous. (B) (Measurability)) The canonical feature map Φ X - HI is measurable. This holds if k(x,)) :) is measurable for all x E X.
(C) (Existence of second moments). It holds that Exx<[llxxx) < oo for all N E P(X). This is satisfied if supxe. xxx,x) < 0G. (D) (Co-kernel). H E (X'), where Cx())) 3 is the space of continuous real functions on X' vanishing at infinity. This holds x + k(x, x) is bounded on X' and k(x, -) E Co(X') for all x e X
(E) (L unniveiella kernel). HH is dense in L (T)) for all probability measures E P(X).
(Forming an algebra). H is closed under pointwise multiplication. f.g E H + f f.g. g E H. with the multiplication map (f.g + f bounded.
An example RKHS satisfying all the assumptions Sobolev space H C dd of order 8 > d/2 [1
Theorem 4.39), motivating natural connection to neural tangent kemels g 86 Generally, class of
RKHSs satisfying (M)) is known as reproducing kernel Hilbert algebras (RKHAs) [69 Definition 2.1].
Recent results have shown that (F) may hold under more abstract regularity conditions [1] 19 69
22

---

A.3 Mean Embeddings of (Conditional) Distributions
Under Assumptions (A) to O for any X ~ π E P(X),, the following expectation yields an element Me in H which is called the mean embedding of [70]:
ux = / 6 dx = [X)] E[(()) H.
(A.3)
We write μxx = fo when the distribution clear from context. We call the RKHS H characteristic
the mean embedding map + μ4 is injective, that is, each mean embedding uniquely encodes distribution. The Assumptions (D) and [E) imply that HI is characteristic [13 91.
With a characteristic HHI, we can define a distance metric between distributions in P(2) called the maximum mean discrepancy (MMD), defined follows (4)
/ [(x) du(x) [(x)d((x)
f [xx dr (x) = 11 - M4 IIs
(M, = sup fe
(A.4)
which we can interpret as measuring differences between mean embeddings.
To handle distributions in P(Z), we analogously define mean embeddings and MMD for RKHS G.
We can extend the concept of mean embeddings to the conditional distributions P P 2n) and
[[] of an HMM. This is formalized by conditional mean embeddings [70, 89 76.. defined as:
NZwii/Zm=z / (pxz' u(z, dz') = [l(p(Z=+1) zm =z]: z).
(A.5)
2x= Sxx o(x) ol2, dx) [x)
2 2).
(A.6)
In kernel-based inference, an approximation of the above is often achieved using linear operators
U : G + G and O: : G - Hi called the conditional mean embedding (CME) operators [89 29
68]:
Uz(z)= =
Z 1)
z).
(A.7)
Oc(z) = E
X
Z = z).
(A.8)
For A.7)). if the RKHS H satisfies Assumptions
(A) to E. it is known that the approximation in A.7
can be made arbitrarily precise for a choice of the operator norm [68] Theorem 3.3). Furthermore, for well-specified cases where the map
12" =z is identified an element of S (G), choice
of U (G) exists which achieves zero error Uz((z) E[p(Z*+)) Z
[68 Corollary 5.5].
For (A.8), we need an extension of 68 Corollary 5.5] to different domains X' and 2 (74. This possible, since the result applies as long as there is a Markov kernel connecting the underlying measure spaces, which requires mild regularity conditions
48 Chapter 8, Theorem 8.5 and preliminaries).
Under these conditions, and if IH and G satisfy Assumptions (A) to EE) we consider well-specified
cases where a choice O E (2, HH) exists such that (:p(z) ([o(X")
2d z]..
Remark A.1. Assuming well-specifees CME operators is a standard choice in kernel-basee statistical learning 668 74y. Even if the assumption not exactly met, close approximations can still be achieved.
e.g., representation learning can yield a feature space that is well-specified for the problem at hand.
A.4 MMD Gradient Flow
Let P2(X') be the set of distributions on X' with finite second moment equipped with the 2- Wasserstein metric. Under Assumptions (A) to [E] the MMD (A.4) provides natural process for generating samples π E P22(')) from an embedded distribution Mxx
A.3).
Specifically, the gradient flow of the MMD defines a continuous probability path (p:) (p:) which starts at any initial distribution po and converges towards the target + π as t + co. Let us fix æ, and define F(p)) == MMD² (pe, T), which measures distance between pt and æ. It is known that the time-dependent vector field (t1)120 corresponding the gradient flow of F(pt) is given by 4
(x( = - - Mz (x) = (x())) - px)n:
(A.9)
which generates the path (pt) via the continuity equation dep, div = 0. If & is continuously differentiable on with Lipschitz gradient, for any initial distribution Pn E P2(X'), there exists a unique process (x))22 from Xo Po satisfying d.X, = (X,)d, where the distribution pt of X, is
the unique solution of the continuity equation and F(P.) decreases in time 4 Proposition 1, With
further regularity conditions, it can be shown that Pt + asymptotically t + oo 4 Section 3].
23

---

B Theoretical Results
We formally state and prove all theoretical arguments in the main text.
B.1 Tensor Product Mean Embedding (Section 3
we embedded the joint distribution X1: in the tensor product H@N H .................. H of characteristic RKHS H.. To establish that H is itself a characteristic RKHS, use the following:
In Section 3
Lemma B.1. Let Assumptions (D) and E) be satisfied HL Then H is characteristic.
Proof. From 68
Remark 4.7], Assumption (D) implies: H is T.2-universal
H is Co-universal,
i.e. H is dense in Co(2)) with respect to the supremum norm. By Assumption E)] we have that HI is 2 univverral, and hence Co-universal. Then, from [93 Section 4], we have that HH is Co-universal H@N is characteristic. Therefore, under Assumptions and EE)
H & is characteristic.
B.2 Mean Embedding Decomposition (Section 3.1
In Section 3.1 we derived tensor network decompositions of the sequence mean embeddings μ,, and Up, through (377 and (399.. We provide the respective proofs in Proposition B.2 and Proposition B.6
Proposition B.2. Let Assumptions (A) to (E) be satisfied for G and H. Assume that for P[X"|
and TP| X v 2 O, well-specified CME operators O G + H and S. G + exist, that is:
Ou((z) E(Ø(X") 2" = H
(B.1)
Sip(z) E(o(X?)
2" = H. Hp 08N 1zz:N,
(B.2)
for all z E Z. Then the sollowing holds:
(B.3)
Pp. = S.ON PzzIN.
(B.4)
Proof. The proof for of observations X
B.3 is follows. We use two properties of HMMs: conditional independence
2 Z1:N for * j and locality of observation || = =P(X"
Mp [[x(X*) @ o(xx)]
(by definition)
= [B(o(X*) ... 8 o (xx 1I:N]]
(law of total expectation)
E (B(d(X+))2z4) 11) [|
@ B(
(conditional independence)
E ( 0.. ((((( |2M]]
(locality)
[[0 3 .0 [
(CME operator (B.
= E [OSN (x ... (((
(tensor product of linear operators)
O Blp (Z R ((
(linearity)
N UzI.N.
(by definition)
The proof for B.4 is identical, by substituting O with S..
We now prove tensor network decomposition of pzu. 3.9 We show some useful lemmas:
Lemma B.3. Let be an RKHS on a set Z with norm || ll.. Suppose Assumption holds, so there exists constant C such that
|| gllG S C llalic for all f.g G. Then the bilinear map (f.g) + f extends to a bounded linear operator
T"e G G
G, T" ((f@g) f
with |T S C. Consequently its adjoint
T = (T*)* G GoG
exists and satisfies ||T'||
T C.
24

---

Proof. Define B G:G
G 1 G by B(f.g) = By assumption
[ =If glle s C Ilflla llgll..
so is bounded bilinear map. By universal property of the Hilbert-space tensor product, B extends uniquely to bounded linear operator
T G 2 - G, ( = f.
with | C.
Since 7 is bounded between Hilbert spaces, its adjoint T (7 exists and ||] = |T S C.
n
Corollary B.4. Under the setup of Lemma B.3 we have Ty((z) (zzz) @ ((z) for all E 2.
Proof. This immediately follows from the definition of and the reproducing property:
(Tv(z). @ g)s (p(z), T g))e = (f(z), f g) (f g)) = f())((z) (p(z) @ (zz), f @ g)ccc for all f.g E G and E Z.
(B.5)
Since the collection of elementary tensors {f @ f.g E G} spans the tensor product space G G, it follows that (z) @ pzz) Tz) for all z.
Lemma B.5. Let Assumptions (A) to (F) be satisfied for G. Assume that for P(z"+1Iz") well- specified CME operator : C exists, i.e., it satisfies the following for all z E Z:
Uu((z) Elp(ZZ Z =z] G.
(B.6)
Then, there exists bounded linear operator T G G & G that satisfies the below for all E 2:
= =llo2>)= 10 (2
z z) E G G.
(B.7)
Proof. Let T be the bounded linear operator from Lemma
B.3
Then:
(l 6 (z(=+1)| z" (Zx+1)
z" = .
(Corollary B.4
TE(+(2
2 =
(linearity)
TUp.
(CME operator (B.6)
This completes the proof.
We are now ready to prove the main result (3.9
Proposition B.6. Under the setup of Lemma B.5 and a choice of decomposition of the CME operator
U = Ex Adhs @ gi with hi, E G, the following holds:
HZLN
Σ bis @ Xi Winiig 2 AiN Win- @ Aixhix:
(B.8)
where b, = (Tuz) E G and waj = [Thi]g G.
Proof. We start with the fact that singular value decomposition (SVD) of U (B.6) always exists and written as U = =YEN Xyh; & 9i [68] Section 4.1] where X; E R are singular values and h; 9. E G are left and right singular functions, respectively. We proceed with SVD, which is sufficient for the proof, and discuss eigenvalue decomposition (EVD) at the end of this section.
Define by Tuzz 99: and wi. == [The]g lgg. which are elements of G from the standard identification of tensors as Tilbert-Schmide operators & & (2) [68 Section 4.1]. The following holds:
Eip(Z~)
N Σ
New Jin (Z NN "hixx
(B.9)
iN
Elip(Z")) Minn (Z²) 2n-1] | Σ
Aingi ( (7n-1) Wi
(B.10)
in
Elpp(Z+) )0,,
big
(B.11)
25

---

To avoid clutter, we only write out the derivation of (B.10):
(1p((Z))mmn+) (z*)
-" = =[((((>) @ ((Z)))()))))))))))))))))))))))) 1 12--1]
(reproducing property)
El((Z"))8)) 8 n" 199n41
(linearity)
-TUppZZ"')]0
(Lemma B.5
- - X1. (Z (Z"-1) ()inn 19hn1
(SVD)
In
(z" 19tt1
(linearity)
bn
We now construct recurrent decomposition of Mlz.N using the Markov property of Z!:N The following reduction rule useful. Let K₂ = [Fly(( (² )intt (Z*) zn for n 2 2, then:
Bip(Za-1) @ Kn
zm- = El x
Asn Di. (Zn-1) Win 2m-2
(by B.10
Σ Elsp( zm 19 (z"-1) zm- Ain Windail
(Ineariity)
bn
Σ Kn-1 2 Xin Winidnt1'
(by definition)
The full decomposition is then given as follows:
Pzi. = Bl(((z)) 8 ... @ )
(by definition)
Blp(Z") @ Zz) @ El. Bip(Z~-1) @ Blp( 11 12 2 |21]
B.9
(Markov property)
Bljp(")) 8
din IoN (ZN-1)hin 22.-1....1 Eip(Z=-1) 2 Σ
iN
-lp((")) Σ
Blp
Jii (ZN-1) >| N2N-N 2 @dix ...]
(linearity)
iN
Kx-1
= Blo(Z') 8 ... Bip(ZN-2 N-2) @ Σ KN-1 8 Aix 12 33.31....
iN
= Blp(Z*) 8 ... Elp((N-2) & KN-1 I 8d hin ...]
(linearity)
iN
reduction rule
= Blp(Z') 8 ... Σ
KN-2 2 dea
@ Aix hix ...
Elp(Z*) 8
Σ
@ KN 2 2 @din WiN- @ Dag hig Blp(ZN~3)
1:
(linearity)
Σ
(Z) 6
K₂
② Xig Wig 2 ... 8 Aix hix
(recurrent reduction)
B
El Wig (Z*))
60 60 Wig.is @ dig Wigres @ 8 Aix hix
B.II
biz Aig Wigris & Xig Wighss .... hiy-
This completes the proof.
With Propositions B.2 and B.6 the following is obtained from linearity and reproducing property:
26

---

for Corollary B.7. Let Assumptions (A) to E be satisfied for G and H. and let Assumption (F) be
satisfied for G. With well-specified CME operators in Propositions B.2 and B.6 the inner product
between the feature of data x N E X and the mean embedding Up is given follows:
(o(x') 8 ... 1..00-0..66 0(xx)) ),p = (ő(x') @ ... 00.'0 SNN H2 Jaa
Σ
(xig(x)) Aig OWigiss(e (x²) ... dix Ouiy (x*).
12
(B.12)
For Upp. one only needs to replace O by St.
We conclude the section with discussion on eigenvalue decomposition (EVD), which we adopt in practice. Assuming the CME operator U (B.6) a normal operator with discrete spectrum, we can invoke spectral theorem to obtain an eigendecompossiion U = =, X,h; @ 9: with eigenvalues X; E C
and cigenfunctions 9s, h; : X + C. While SVD already suffices for Proposition B.6 EVD provides practical benefit efficiently captures long-range dependencies via oscillations 75 7 12).
Remark B.8. While U operates on real functions, for EVD we include complex-valued eigenfunctions to ensure sufficient expressiveness. This is standard in the analysis of dynamical systems, as operators on Hilbert space (e.g., Sobolev spaces H* C R) with real spectra self-adjoint and thus can only capture time-reversal invariant dynamics. Normal or bounded operators require complex spectra and eigenfunctions to construct spectral expansions, 180 Definition 12.17 and following).
Remark B.9. As the spectral expansion only used to parameterize the operators, the results in Mollenhauer & Koltai (2020) (68] and Arbel et al. (2019) [ underlying our work are only affected in the sense that, as an additional assumption, a spectral expansion has to exist. complete extension of the results therein to complex-valued RKHS is out of the scope of this paper.
B.3
MMD Flows with Time-dependent RKHS (Section 3.2
In Section 3.2 we introduced an extension of MMD gradient flows using time-dependent RKHS. The
derivation is direct application of existing work on gradient flows of time-dependent functionals 27
We use the fact that, for fixed in (IRd), the time-dependent functional Fi()) = AMMD (M, x) in time interval t E [0, 1] admits the following free-energy expression
4
/ vx(x)dx()) + /
f Wi y) dv(x) dv(y) + C.,
Fi
(B.13)
where Ve, We, C2 are time-dependent confinement potential, interaction potential, and constant:
1 2/2 Vi(x)= - xx(xxx)) dx) xx(x,x)) xx,xx'),,, C. = 2
he kx dx((x) da (x').
(B.14)
Let us assume regularity conditions for V,, We as in 22 Theorem 6.4, Theorem 6.6 and Remark 6.7].
Then, follows from the results that there a continuous path (pe)eejo,1) Pa (Rd) from a given Po
associated time-dependent vector field (we)se[o.1] via the continuity equation +
0.
The vector field is obtained as (x()) - V f + W, (x)(x) 27 Theorem 6.4 and Remark 6.7]
where * is classical convolution, which yields Uy (x) = - VV ((4... p=,,) (x) as in [4] Section 2.1]. This provides minimal derivation of the flow which is sufficient in our problem context, and
leave in-depth investigation of the regularity conditions as future work.
B.4 Gradient Fields in Flow Matching (Section
In Section 3.3
we used vector fields developed in flow matching [60] to define our training objectives, under the premise that they are gradient fields. We provide formal proofs for all vector fields in 66)
Consider converting an initial distribution qu to a target distribution म. by following a continuous path (ge (o in P2(X') with boundary condition 91 π. Such path can be identified by time- dependent vector field (ue)eejo.1 which satisfies the continuity equation dru, + div(guu) = 0 (97].
27

---

In flow matching, the probability path and vector field are defined as marginalization of the conditional counterparts qe (-(xi) and uq (4(-(x)) 660 Eq. (6). (8)], under the assumption that qi (x) >
q: (x(x) = q (x(x(xx)))(XX dx1.
(B.15)
(x))) / ur (x(x)) Te (x (opxx))((( (x11
dx1-
(B.
(x(x)
Known choices include optimal transport (OT), exriance-exploding (VE), and variance-preservinng (VP) conditional vector fields [60] Section 4; Examples and 2J. We show that their marginal vector fields are gradient fields. To ensure the existence of all integrals, we assume that π (x) is decreasing to zero at a sufficient speed as ||x||
and 0t is bounded. From now on, denote Vx by V.
We first show useful lemma:
Lemma B.10. If qt marginal probability path, qe(-(xx) bility path and πt (x1) is the distribution of X1, then the following holds:
(x1) ) the corresponding conditional proba-
f 26 (x(()x))<(x1)) dx1.
1=
(B.17)
(e(xx
/ V log (x(x1)
9o (x(x (e(x(x))(((X)) dx
V log (x()
(B.18)
qt (x)
Proof. (B.17) directly follows from B.15
For B.18 by noting that xxx)) is independent of x in (B15), we have:
Vg((x) = /
/ Vv
dx1
(B.19)
Then, from (x(xx = V Vlog (x|x1) (x(xx), we have:
Vx.(x)
(x(x) = / Vlogg (x/x)) (x(()))((x)) dxi.
V Vlogg((x)
(B.20)
(x(x)
This completes the proof.
now show the main results:
Proposition B.11. Consider the optimal transport (OT) conditional probubility path and vector field from 160 Section 4; Example 2 and Theorem 3]:
q. (x =x(x)uu(xx)) (x (x((xi, (1 (t)
(B.21)
mmin) )x
u. (x(xx))
(B.22)
Omin)
Then the marginal vector field ux(x) is a gradient field.
Proof. By the definition of the Gaussian distribution,
(x(x)) (x (1 (1 (t) exp (- ||x lx1
/Z,
(B.23)
(1 Omin )t}
where 7 a normalizing constant that does not depend on x.
We separate u (x) into two gradient fields: the score function of qe (x(x;) and a linear map.
Let us write a, := (1 mtint. The score function of qi (x|x1) is
txi
X txi
V log de (x(x1) (1 (1 Jmin )0)3
(B.24)
of
28

---

Next, remember that 0z = 1 (1 minu is independent of x, X1 and note that
(1-00)
(1 0 min X
X1
tx oz)
(x(xxx1)
1 (1 Fmin )t
dt
to;
X
X tx1
x
0t
x tx1
t
to
t
of
x
+
Tt
log (x(xx).
(B.25)
t
Or Let t By the definition of marginal vector field, we get
Then a is independent x, x1 and un (x(x +
We (x(xxx Ge (x x, (T1 dxi
ux(x)=
(x(x)
/
x ( + a log (x(xx)) 2i (x(x)x)) dx (x(xx)
/ (((x(x))((XI)) dxi + /Vloogg V ( 2 (x(xx))((X) dxi
2a (x(x)
2e (x(x)
+ a log qt (x);
(B.26)
where the last equality holds by Lemma B.10
Here, x } is linear map on x so it is gradient field, and x ++ a log ge (x) ) is by definition a gradient field. Since un (x) an addition of two gradient fields, it is clearly gradient field. This completes the proof.
The marginal vector fields of the variance-exploding (VE) and variance-preeserving (VP) diffusion probubility paths in 160 Section 4; Example 1] are grudient fields.
Proposition B.12.
Proof. The proof based on a stochastic differential equations (SDE) formulation Consider SDE of the standard form, defined in interval & [0, 1]:
60
Appendix D).
dy fidt + g:dw,
(B.27)
with drift h,, diffusion coefficient Gt. and standard Wiener process dw.
The corresponding marginal vector field is given as follows 660 Appendix DJ:
Üh x = fx(x) 5 v log (x)).
(B.28)
First, it is known that the SDE for the VE path is
dy V d dl
oldw,
(B.29)
where and a, + oo as + 1. This SDE moves from data at / = 0 to noise at / = 1.
Thus the SDE coefficients are fe(x) 0, (x) v o and the following vector field en (x)is
1 d o²V log log (x),
ūx(x) = 2 dt
(B.30)
and the vector field from noise at = 0 to data at t is
1 d
u, (x(x) = 2 dt -
log q (
(B.31)
Since 01-t is independent of x, (x)) is clearly gradient field. Next, it is known that the SDE for the VP path is
(') dy 2 + (t)
(B.32)
29

---

where Tt)) = S0 (ssdd and 3 is noise scale. This SDE moves from data at = 0 to noise at = 1.
Thus the SDE coefficients are Se(x) = T(
x, D6 gx V (1), and the vector field ūx(X) is
2
(t)
1 (1) log (x),
ü (x)
x
(B.33)
2
2
and the vector field from noise at 0) to data at i 1 is
ux(x) = (1
2 x
T" - t) log ()
(B.34)
Here, x + (1-(() 2 x is linear map on x so it is gradient field, and 7 () is independent of x x - - (' 1) log ( (x)xxx) ) is gradient field. Since uz(x)) is an addition of two gradient fields, it is clearly a gradient field. This completes the proof.
C Experiment Details
For all experiments, we implement spectral mean flows as in Section 3.3 PyTorch (78). We use 32 bit- precireciiinn and use opt_elll einsumm package [17] to optimize the tensor network contractions. For flow matching training, we use the OT flow and path (Proposition B.II] with min = 5. Sampling is done with a midpoint ODE solver from flow_mar www_wwthhingggggoooooeeoeeemeene package [61]. Each experiment is done with single NVIDIA RTX A6000 GPU with 48GB and Intel Xeon Gold 6330 CPU 2.00GHz.
C.1 Implementation Details (Sections 3.3 and
We provide the implementation details. Recall that model is parameterized by shared feature extractor MLP(s,,, (): : X' x Rd >> Cdd and parameters o, 8 E IR À E C", B, L, R, H E Cdty xr
We denote trainable linear layer from X' to y by Linx-y; and denote Linxx == Linx We denote squared ReLU activation + ReLU(()) (x)
2 88 by ReLU², RMS normalization [[03] by RMSNorm,
1,2 normalization by L2Norm, and d-dimensional sinusoidal positional encoding 96 by PEga.
For given number of layers L and hidden dimensions dp, d/,, the feature extractor MLP(;, .tt))
X x Rd
+ Cd is designed as follows. We obtain conditioning feature of using sinusoidal PE:
C₁ Lingo ReLU Lingo, ReLU Lingaae_ee PEgas (t).
(C.1)
Then, the full feature extraction process for given input E X' is written as follows:
[0] T (x) + q Ct:
(C.2)
z (1] = z "-1) + Lin 2 -ndn o ReLU o Lin 2 RMMNorm (z a1 (-1)),
(C.3)
( L2Norm ReLU2 Ling 40% (2
(C.4)
MLP((x, q.() = co ..... (41) E (o.s)
(C.5)
In we apply ReLU2 to complex features by separately operating on real and imaginary parts. As shown in C.S). instead of using only the output of the last layer, the MLP concatenates features from all layers to construct the output, i.e., dj = Ldl.
Then, parameterize L. R. H Cdy xr, we ret Ld, for some d, and, to prevent O(L2d, parameter complexity, impose the following block-diagonal structure on each of D E R.
di
D = o
E c(lds)x(Edd) di Cdixd,
(C.6)
d,
In addition to reducing the parameter count by a factor of L, this allows space-efficient, layerwise evaluation tensor network contraction based on (r, (rr))())))) @ .l,) = (r . or,l, In practice, simply set d, dn/2 unless stated otherwise.
30

---

Table 6: Hyperparameters for Table
Hyperparameter
Sines
Stocks
ETTh
MuJoCo
Energy
fMRI
Hidden dimension da, d'
64
64
96
96
Layers L
10
10
10
16
14
16
Batch size
256
128
256
256
128
256
Training steps
12000
10000
18000
28000
25000
15000
Sampling steps
500
500
500
1000
1000
1000
Table Wall-clock training times for Table
Method
Sines
Stocks
ETTh
MuJoCo
Energy
[MRI
Ours
18min
15min
27min
60min
44min
35min
Diffusion-TS
16min
15min
30min
25min
60min
48min
Table 8: Wall-clock sampling times per 1, 1,000 samples for Table
Method
Sines
Stocks
ETTh
MuJoCo
Energy
fMRI
Ours Diffusion-TS
30sec
30sec
29sec
109sec
95sec 202sec
109sec
66sec
66sec
67sec
133sec
262sec
Table 9: Parameter counts for Table
extending
101 Table 7].
Method
Sines
Stocks
Energy
Ours
356,736
356,800
1,086.432
Diffusion-TS TS
232,177
291,318
1,135,144
Diffwave
533,592
599,448
1,337,752
Cot-GAN
40,133
52,675
601,539
TimeGAN
34,026
48,775
1,043,179
TAm
97,525
104,412
677,418
The eigenvalues and nonzero blocks of B, L, R. H. use exponential parameterization of (75), which initializes each entry uniformly between rings of radius mmin and Tiax on the complex plane:
= exp(-v i0)
Log(uz(1") 722 miin min) 0 2 U1 U2 ~ 410 1] (C.7)
When initializing on the complex unit disk, we find necessary to set mmin e and Fmax for small constant = 10-12, , to avoid numerical issues related to logarithm.
B. L. R, H. apply proper scaling prevent the activations from exploding After initializing A, B,
vanishing. We find that scaling À by /2 and the rest by v2d, 1/4 provides the desired stability.
C.2 Synthetic Experiments (Section 5.1
We provide details of the 2D checkerboard experiment. We use MLP(,,,&)) : X² CLds with L = 4 layers and hidden dimensions da = d', = 128. For training, we use Adam optimizer [22] with hyperparameters (51, (32) = (0.9, 0.999). The model is trained for 20k iterations with learning rate le-3 and batch size 10,000, with 1k steps linear learning rate warmup, and gradient norm clipping 1.0. For sampling, keep track of exponential moving average (EMA) of the model parameters during training [50], updating every 10 iterations with decay rate 0.995.
C.3 Time-Series Modeling (Section 5.2)
We provide details of unconditional time-series generation experiments in Section 5.2
Regular time series We first detail the datasets in Table Sines has 10,000 5-channel sequences where each channel is sinusoidal signal x(n) sin sin(2nwn 0) with independently determined
31

---

Table 10: Parameter counts for Table
Table 11; Parameter counts for Table 3
Method
Sines
Stocks
MuJoCo
Method
FRED-MD
NN5 Daily
Ours
32ll
12lk
32,
Ours
100,603k 151,720k
100,569k 151,720k
SDForme--A
44.971k
44,971k
44,97IK
Imagen Time
ImagenTime
57,560k
57,562k
14,430k
Table 12: Parameter counts for Tables 4 and
Method
Stocks, 0% Drop
Stocks, 30-70% Drop
Pendulum
Ours
157,888 42,270
157,888
45,544
Koopman VAE
33,410
37,978
frequency w ~ [4[0, 1] and phase 0 ~ ~=[[[-T, .T. Stocks contains 3,773 6-channel sequences from daily Google stock price data from 2004 to 2019. The channels are volume and high, low, opening, closing, and adjusted closing prices. ETTh contains 17,420 7-channel sequences from electricity transformers recorded every 15 minutes from 2016 to 2018, including load and oil temperature. MuJoCo has 10,000 14-channel sequences from DeepMind dm_ control physics simulation. Energy contains 19,711 28-channel sequences from appliances energy use in a building recorded at a noisy 10-minutes interval, including house temperature and humidity. fMRI contains 10,000 50-channel sequences from a realistic simulation of blood oxygen level dependent signal in functional MRI.
We closely follow the protocol of [01 for hyperparameter selection and training. Detailed hyper- parameters are in Table 6 The ranges considered for MTP((, : () IRd Reds - CLdd. are batch size {128, 256}. layers L E {10, 14, 16), and hidden dimension d,, d', E {64, 96}. Other choices mostly follow [I01]: For training, we use Adam optimizer [22] with (31, B2) = (0.9, 0.96). The model is trained for the same iterations per dataset as in [111] with learning rate Se-4 and gradient clipping at 1.0, with 500 steps of linear warmup and then decaying by 0.5 on plateau with patience 2,000. An exception is MuJoCo, where longer training is necessary for convergence. For sampling, we use EMA of model parameters [50] following [IOI], updating every 10 iterations with decay rate 0.995. We use midpoint ODE solver with the same sampling steps as in [IOI] In T we adopt baseline scores from 100 101 except for Diffusion-TS reproduced with the official implementation.
We provide supplementary results of Table From Figure 4 to 9 we visualize real and generated data from our model in Table [ showing that generations visually realistic. We further perform detailed complexity analysis by measuring wall-clock training and sampling times as well as parameter counts. The results are in Tables 8 and showing that our models for Table | computational complexity, and their generation is faster than the state-of-the-art Diffusion-TS.
have reasonable
We now provide the details of additional comparisons in Table 2 For this setup, we consider two recent state-of-the art methods SDFormer [3] and Imagen Time 77 that operate in larger parameter regime, roughly using 11 1 the parameters of the largest model Diffwave in Table I We use Sines, Stocks, and MuloCo datasets, considering both resource constraints and reproducible datasets in the codebases of 71 We test larger spectral mean flows with comparable number of parameters to the baselines without other architectural changes. Specifically, we use MLP((, ....(,,..).
t) Rdd x IRdh - CLd, with L 12 layers and dimensions dn, = 512, d' 1024. The parameter counts are in Table 10 Other hyperparameters are chosen as follows. Following [1), we train our models using AdamW optimizer
66 with (31, 32) = (0.5, 0.9) and weight decay 1e-6 for 100k iterations and measure the discriminative score with respect to the training dataset every 2,500 iterations for validation. Without hyperparameter tuning, we use batch size 512, learning rate 1c-4, and gradient clipping at 1.0. For sampling, we use EMA of model parameters. updating every 10 iterations with decay rate 0.995. We use midpoint ODE solver with 100 sampling steps for validation and 500 steps for testing. In Table 2 we reproduce SDFormer-AR and Imagen Time using the official implementations.
Long time series We provide details of long time-series modeling in Table 3 We use two datasets from the Monash repository (32]. FRED-MD has 107 time series of length 728 from mccno-cc indicators from the Federal Reserve Bank. NN5 Daily contains 111 time series of length 728 from daily cash withdrawals from ATMS in UK. We choose these datasets as their moderate sizes allow experimentation within resource constraint. Following (7]]. we normalize each trajectory to adhere to zero-centered normal distribution, and use 80% of the data for training and 20% for testing.
32

---

We design our architectures as in Appendix C.I and equip them with time-delay observables used in operator- theoretic methods for dynamical systems 49 995] which we found help stabilizing training. Specifically, for time series xI:N we define a new observable x xI:N' as x == concat (x" ....., ,...,'*""" x for delay a and dimension b. When < b, there is no loss of information and the map x + X is invertible. Hence, we model & using spectral mean flow and convert the result to x using the inverse map. We use 16 and b = 64 for FRED-MD, and a = 20 and b = 20 for NN5 Daily.
We include state oftt--of-the-arr ImagenTime 71 for comparison, and test spectral mean flow with a comparable number of parameters. Specifically, we use MLP(....t) R- x Rd + Cz-da with == 12 layers and hidden dimensions d, 768, d', = 3072. The parameter counts are in Table II Other hyperparameters are chosen as follows. Following (7]. we train our models for 1 epochs with batch size 32 and weight decay 1e-5, measure the marginal score every 100 iterations for validation, and for sampling use EMA of model parameters with 100 warmup steps and decay rate 0.9999. We use Adam W optimizer with ($1. 32) = (0.9, 0.95) and learning rate 3e-4 for all parameters except for (A, B, L, R. H). For these spectral parameters, we use Muon optimizer [46] with learning rate le-3, which slightly improved the performances. While we leave a deeper investigation as future work, we conjecture this is because Muon is second-order optimizer, which are expected to improve the learning of multiplicative matrix states in general [92] Section 3.3]. We use gradient clipping at 1.0. For sampling, we use midpoint ODE solver with 100 sampling steps. In Table 3 baseline scores from [ except for ImagenTime reproduced using the official implementation.
we adopt the
Irregular time series We provide details of irregular time-series modeling in Table We recall that we consider two tasks, (1) irregular - regular: generating the full time series including the missing timesteps, and (2) irregular + irregular: generating only the irregularly sampled time series. We apply spectral mean flow to the tasks as follows. For (1), we handle missing timesteps during training by simply interpolating the observations, which we find to work well in practice. For (2), based on the HMM formulation we use a classical approach [105] of treating the sampling interval itself as part of data. That is, we consider hidden states Z under dynamics P[Z*+I|z") and local observation model ",[ Ar"IZ"
that jointly determines the observation X and sampling interval Ar" between X" X"11 and X", If expressive enough, this HMM can express any joint distribution pi(X*, Ar)), ..... (X*, AT ), i.e., general time series under irregular or informative sampling. The only change in implementation is using sampling intervals A as an additional data dimension.
Detailed hyperparameters are as follows. We use MLP(:, .t) RR x Rdh + CTde with L = 4 layers and hidden dimension dn d' = 64. This leads smaller model than the previous experiments, although larger than Koopman VAE (Table 12p. For other hyperparameters, we closely follow (72) and train all models using Adam optimizer with learning rate 7e-4 for 600 epochs with batch size 64. For sampling, we use EMA of model parameters, updating every 10 iterations with decay rate 0.995, and use a midpoint ODE solver with 100 sampling steps. In Table 4 the baseline scores are from [4] [2) except for Koopman VAE tested using the official implementation. To run Koopman VAE on task (2), we adopt the hyperparameters provided in the official implementation for task (1).
Physics-informed modeling Lastly, we provide details of physics-informed modeling in Table 5
As discussed, we consider a nonlinear pendulum governed by 0+ + 9.8 sin 0 = 0, 0(0) = 0. To obtain a dataset, we simulate 2,000 independent trajectories with initial condition (0) Unif 2.7] over the time interval [0, 10) and sampling interval 0.25, treating the pair (0(t), 0(t) as the observation. Following (72), we add Gaussian noise with standard deviation 0.08 to simulate real-world noise and standardize each trajectory to the range [0. 1). We consider Koopman VAE (72) for comparison, as is a state-of-the- art operator-based method for time series. We use correlational score for evaluation.
To focus on the impact of spectral regularization on the matrix states described in Section 5.2 in the feature extractor we only use the last layer feature (4) (C.). leading to MTP((.,, Rd x Rdd
Cd A Cd, and B. L, R. H E Cd, xd, without block-diagonal structure. We use L 8 layers
and dimension dj = d', = 32, Following (72), we fix the spectral dimension as r = d,- = 4, as well as the latent dimension of Koopman VAE as 4. This leads to similarly sized models (Table 12) with spectral regularization on 4 x 4 matrices (if used). We train all models with Adam optimizer, learning rate 7 for 600 epochs, with batch size 64. For sampling of our model, use EMA of parameters updating every 10 steps with decay 0.995, and use a midpoint ODE solver with 100 sampling steps.
33

---

real
generated
1.0
1.0
0.8
0.8
15
20
15
20
0.8
0.8
0.6
20
20
0.6
0
5
10
15
20
0
10
15
20
1.0
1.0
0.8
0.6
5
15
20
15
20
1.0
0.8
0.8
0.6
20
20
Figure 4: Examples of real and generated data for the Sines dataset.
generated
0.2
0.3
0.1
0.1
0.2
5
10
10
15
20
0.3
0.2
0.2
0
5
10
15
20
10
15
20
0.2
0.2
15
0.2
0.5
0.0
15
20
Figure 5: Examples of real and generated data for the Stocks dataset.
34

---

generated
0.75
0.50
0.25
0.75
0.50
0.25
10
15
20
10
15
20
0.75
0.75
0.30
0,25
0.25
10
15
20
10
15
20
0.75
0.75
0.50
0.25
0.25
Figure 6: Examples of real and generated data for the ETTh dataset.
generated
1.0
0.5
0.25
10
15
20
5
15
0.75
0.5
0.50
0.0
0.75
0.50
5
10
15
20
5
10
15
20
0.75
0.75
0.50
0.50
0.25
0.25
10
15
20
0
10
15
20
1.0
0.5
10
15
Figure 7: Examples of real and generated data for the MuJoCo dataset.
35

---

real
generated
1.0
0
15
20
0
10
15
20
0.5
0.5
0.0
0.0
20
20
0
5
15
20
15
20
1.0
0.0
15
20
15
20
1.0
1.0
0.5
0.0
20
20
Figure 8: Examples of real and generated data for the Energy dataset.
generated
1.0
0.75
0.50
10
15
20
5
10
15
20
0.75
0.0
0.75
0.75
0.50
0.25
10
15
20
5
10
15
20
1.0
0.5
10
15
20
5
10
15
20
0.75
0.75
0.25
0.25
10
15
Figure 9: Examples of real and generated data for the fMRI dataset.
36

---

D Limitations and Future Work
Our current main limitation is the question of how to condition the generative process on partial observations in a principled manner. We believe studying kernel Bayesian inference [29] and Laplace transform of CME operators [56) in our problem context could be fruitful. The reasoning as follows.
Suppose we have partial observation x Im and want to generate x I:N ~ P[X I:N
[ l:n
x [.]
This can be done with MMD flow if we have conditional mean embedding MX X lin. In kernel methods for HMMs, this is precisely estimated via kernel Bayesian inference [2
73 As an
alternative, a type of Laplace transform of Koopman operator has been recently used for kernel-based
forecasting x Iim + x"* 1=N of deterministie nonlinear dynamics [Z]. Given the proximity between the two operators, we conjecture its extension to our setup may enable principled conditioning.
In addition, while in Section 5.2 we showed one possible approach to perform physics-informed modeling using our method via spectral regularization, more direct approach of identifying the operator eigenfunctions from the learned model would strengthen the interpretability of our method.
Also, while we adopted simple strategies handle irregularly or informatively sampled dynamical systems in Section 5.2, another possibly valid approach is using the continuous-time CME operator and its spectral decomposition U = Σ; exp(s;Ar) hi di with continuous-time eigenvalues 8; We leave further investigation of this direction as future work.
U (z) [f 2 (TT)t)]
Lastly, room for future work lies in engineering aspect of the method, including better architectural components and stabilization strategies for long sequences, choices vector fields for flow matching training. well empirical investigations of scaling laws.
While the above directions are exciting, we believe their proper investigation warrants separate work. Thus, we focused our effort on laying the foundations of neural method based on operator theory.
37


Reviews:
[{"review_id": "ZYznlYlhbs", "reviewer": "Reviewer_k6rF", "summary": "This paper introduces a sequence modeling algorithm grounded in an operator-theoretic interpretation of Hidden Markov Models. Rather than directly implementing stochastic recurrence, the authors embed distributions into Hilbert spaces. The model further incorporates a time-dependent MMD gradient flow for sampling and employs a flow-matching learning algorithm trained via double backpropagation.", "questions": "Q1) Do the authors have any preliminary insights on how kernel Bayesian inference or Laplace transforms of conditional mean operators could be practically incorporated into their framework for conditioning on partial data?\n\nQ2) Have the authors explored methods to recover or approximate interpretable components from the end-to-end model?", "limitations": "yes", "rating": 2, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 2, "significance": 2, "originality": 2, "strengths_and_weaknesses": "Strengths: The paper is clearly written and easy to follow. The motivation behind the work is well-articulated, and the related literature is discussed. The proposed method is explained in detail, with clear examples of its integration into state-of-the-art graph learning models.\n\nWeaknesses: \n1) The paper lacks discussion and experimental analysis on the robustness of the proposed method with respect to real word datasets. Also the results in Table 1 is very marginal.\n2) The framework currently lacks a principled approach for conditioning the generative process on partial observations.\n3) The learned model makes it difficult to extract identifiable components like the observation operator or left eigenfunctions, which limits its interpretability.", "was_revised": true, "final_justification": "While some points were clarified, my concerns about the completeness of the evaluation remain. Therefore, I will be keeping my current score.", "rebuttal": "We thank the reviewer for the thoughtful feedback. We address the main concerns below.\n\n>  W1. “The results in Table 1 are very marginal.”\n\nWe would like to remind that our gain in Table 1 is often significant, e.g., Ours improves previous best context-FID, corr., and disc. scores in ETTh by ≈50% (Lines 293-295). We believe this is notable as other methods often use specialized Fourier losses [1] while we only use standard flow matching loss.\n\nIn addition, we did a careful revision of our model (Section 3.4) and were able to improve on high-dimensional Energy and fMRI. The key changes included dropping SwiGLU and adaLN in favor of simpler squared ReLU [2] and sinusoidal embedding [3].\nWe are happy to elaborate on the details during the discussion period.\n\n**Table A.** Updated Table 1. Error bars omitted due to space limits. The best are **bold** and second-best are in [brackets].\n| Metric | Methods | Sines | Stocks | ETTh | MuJoCo | Energy | fMRI |\n|--|--|--|--|--|--|--|--|\n| Cont.-FID ↓ | Ours | **0.006** | **0.020** | **0.051** | [0.024] | **0.058** | [0.115] |\n| | Diffusion-TS | [0.013] | 0.209 | [0.133] | **0.014** | [0.085] | **0.106** |\n| | DiffTime | **0.006** | 0.236 | 0.299 | 0.188 | 0.279 | 0.340 |\n| | Diffwave | 0.014 | 0.232 | 0.873 | 0.393 | 1.031 | 0.244 |\n| | TimeGAN | 0.101 | [0.103] | 0.300 | 0.563 | 0.767 | 1.292 |\n| | TimeVAE | 0.307 | 0.215 | 0.805 | 0.251 | 1.631 | 14.449 |\n| | Cot-GAN | 1.337 | 0.408 | 0.980 | 1.094 | 1.039 | 7.813 |\n| Corr. ↓ | Ours | 0.025 | **0.002** | **0.032** | [0.214] | **0.746** | **0.750** |\n| | Diffusion-TS | 0.032| 0.009| [0.049] | **0.191** | [0.809] | [1.411] |\n| | DiffTime | **0.017** | [0.006] | 0.067 | 0.218 | 1.158| 1.501|\n| | Diffwave | [0.022] | 0.030| 0.175| 0.579| 5.001| 3.927|\n| | TimeGAN | 0.045 | 0.063 | 0.210 | 0.886 | 4.010 | 23.502 |\n| | TimeVAE | 0.131 | 0.095 | 0.111 | 0.388 | 1.688 | 17.296 |\n| | Cot-GAN | 0.049 | 0.087 | 0.249 | 1.042 | 3.164 | 26.824 |\n| Disc. ↓ | Ours | [0.011] | **0.015** | **0.014** | **0.010** | **0.120** | **0.108** |\n| | Diffusion-TS | **0.007** | [0.080] | [0.073] | [0.014] | [0.122] | [0.175] |\n| | DiffTime | 0.013 | 0.097 | 0.100 | 0.154 | 0.445 | 0.245 |\n| | Diffwave | 0.017 | 0.232 | 0.190 | 0.203 | 0.493 | 0.402 |\n| | TimeGAN | [0.011] | 0.102 | 0.114 | 0.238 | 0.236 | 0.484 |\n| | TimeVAE | 0.041 | 0.145 | 0.209 | 0.230 | 0.499 | 0.476 |\n| | Cot-GAN | 0.254 | 0.230 | 0.325 | 0.426 | 0.498 | 0.492 |\n| | AR-RNN | 0.495 | 0.226 | - | - | 0.483 | - |\n| Pred. ↓ | Ours | **0.093** | **0.037** | [0.124] | [0.010] | [0.251] | **0.100** |\n| | Diffusion-TS | **0.093** | **0.037** | [0.124] | **0.008** | **0.250** | **0.100** |\n| | DiffTime | **0.093** | [0.038] | **0.121** | [0.010] | 0.252 | **0.100** |\n| | Diffwave | **0.093** | 0.047 | 0.130 | 0.013 | [0.251] | [0.101] |\n| | TimeGAN | **0.093** | [0.038] | [0.124] | 0.025 | 0.273 | 0.126 |\n| | TimeVAE | **0.093** | 0.039 | 0.126 | 0.012 | 0.292 | 0.113 |\n| | Cot-GAN | [0.100] | 0.047 | 0.129 | 0.068 | 0.259 | 0.185 |\n| | AR-RNN | 0.150 | [0.038] | - | - | 0.315 | - |\n| | Original | 0.094 | 0.036 | 0.121 | 0.007 | 0.250 | 0.090 |\n\nThe updated Table 1 is in **Table A**. Our renovated model achieves the best metric in 16 out of 24 cases.\n\n> W2. “The paper lacks discussion and experiments on real-world datasets.”\n\nWe would like to clarify that 5 out of 6 datasets in Table 1 are real-world or realistic simulations for real-world applications (robotics, medical).\nNevertheless, we ran additional experiments on two challenging real-world setups: (1) long time series and (2) irregularly sampled time series, comparing with recent SoTA ImagenTime [4] and KoVAE [5].\n\nFor long time series, we used FRED-MD and NN5 Daily with 728- and 791-length data [6]. Their modest sizes allowed testing within the rebuttal period. Our models are designed similarly to **Table A**, equipped with time-delay observables used in operator methods [7]. We use marginal, classification, and prediction metrics [4].\n\n**Table B.** Long time series. Baselines are from [4], except for ImagenTime which was reproduced using official code.\n| | FRED-MD | | | NN5 Daily | | |\n|--|--|--|--|--|--|--|\n| | marg. ↓ | class. ↑ | pred. ↓ | marg. ↓ | class. ↑ | pred. ↓ |\n| Ours | **0.0145** | **1.547** | 0.0261 | **0.00453** | **1.07** | 0.427 |\n| ImagenTime | 0.0292 | 1.154 | 0.0279 | 0.00745 | 0.935 | 0.496 |\n| LS4 | 0.022 | 0.544 | 0.037 | 0.007 | 0.636 | **0.241** |\n| SaShiMi | 0.048 | 0.001 | 0.232 | 0.020 \t| 0.045 | 0.849 |\n| SDEGAN | 0.084 | 0.501 | 0.677 | 0.085 | 0.085 | 1.01 |\n| TimeGAN | 0.081 | 0.029 | 0.058 | 0.040 | 0.001 | 1.34 |\n| Latent ODE | 0.042 | 0.327 | **0.013** | 0.107 | 0.000 | 1.04 |\n| ODE2VAE | 0.122 | 0.028 | 0.567 \t| 0.211 | 0.001 | 1.19 |\n| GP-VAE | 0.152 | 0.016 \t| 2.05 | 0.117 | 0.002 | 1.17 |\n| RNN-VAE | 0.132 | 0.036 | 1.47 | 0.137 | 0.000 | 0.967 |\n\nThe results are in **Table B**. Our method shows strong results, achieving the best metric in 4 out of 6 cases and outperforming the previous best (ImagenTime, LS4) in 5 out of 6 cases.\n\nFor irregular time series, we used 3 datasets derived from Google stock data by omitting 30%, 50%, and 70% of observations [5]. The task is generating full signal. Our models are designed similarly to **Table A**. To handle missing values at training, a simple idea of interpolating observed values worked well. Based on the code of [5], we use discriminative metric.\n\n**Table C.** Irregular time series. 30%, 50%, and 70% dropped. Baselines are from [5] except for KoVAE which was reproduced using official code.\n| | Stocks, no drop | 30% | 50% | 70% | | |\n|--|--|--|--|--|--|--|\n| disc. ↓ | | | | |\n| Ours | **0.0150** | **0.1109**| **0.0130**| **0.0628**|\n| KoVAE | 0.0203 | 0.1336 | 0.0899| 0.0740|\n| GT-GAN | 0.077| 0.251 | 0.265 | 0.230 |\n| TimeGAN | 0.102 | 0.411 | 0.477 | 0.485 |\n| RCGAN | 0.196 | 0.436 \t| 0.478 | 0.381 |\n| C-RNN-GAN | 0.399 | 0.500 \t| 0.500 | 0.500 |\n| T-Forcing | 0.226 | 0.305 | 0.308 | 0.365 |\n| P-Forcing | 0.257 | 0.341 | 0.388 | 0.317 |\n| GP-VAE | 0.152 | 0.016 \t| 2.05 | 0.117 |\n| RNN-VAE | 0.132 | 0.036 | 1.47 | 0.137 |\n\nThe results are in **Table C**. Our method achieves the best metric in all 4 cases.\n\n> W3. “It is difficult to extract identifiable components, which limits interpretability.”\n\nWhile interpretation of eigenfunctions is not direct due to our end-to-end design, we believe that a closely related benefit of operator methods is retained to some degree: incorporation of physics-based prior knowledge. Thanks to linearity of tensor network, this can be done via explicit spectral regularization. We test this using a challenging physics-constrained problem from [5].\n\nWe consider a nonlinear pendulum, governed by an ODE of angular displacement $\\theta$ from an equilibrium [5]: $\\frac{d^2}{dt^2}\\theta + 9.8\\sin\\theta = 0, \\dot{\\theta}(0)=0$.\nAs the pendulum is stable and conservative, the physical knowledge is that eigenvalues of underlying operator have non-positive real part and some have zero real.\n\nSince our model is governed by a linear product of matrix hidden states $Ow(x^n)$ (Eq. 3.12 and Figure 2), we can implement this knowledge in an end-to-end manner.\nIn KoVAE [5], this stability constraint is built-in using a loss on matrix $A$ for hidden states $z_{t+1} = Az_t$, specifically $|s_1 - 1|^2 + |s_2 - 1|^2$ where $(s_1, s_2)$ are the two largest eigenvalues of $A$. We compare our model against [5] using the same spectral loss for each matrix hidden state’s two largest eigenvalues.\nWe use correlational score for eval. [8]\n\n**Table D.** Physics-informed nonlinear pendulum modeling.\n|  | w/o stability loss | w/ stability loss |\n|--|--|--|\n| corr. ↓ | | |\n| Ours | 2.298e-3 | 2.320e-4 |\n| KoVAE | 2.164e-3 | 6.954e-4 |\n\nThe results are in **Table D**. Without physics-based loss, ours performs on par with KoVAE. When the physics-based loss is incorporated, we observe a clear improvement, outperforming KoVAE which uses the same stability constraint. This supports incorporating physical knowledge into our model is possible, resulting in a better match with the true data-generating process.\n\n> W4. “The framework lacks a principled approach for conditioning on partial observations.” “Do you have preliminary insights?”\n\nWe do have preliminary insights. Say we have a trajectory $x^{1:n}$ and want to sample $x^{n+1:N} \\sim P(X^{n+1:N}|X^{1:n} = x^{1:n})$. This can be done with MMD flow if we have conditional mean embedding $\\mu\\_{X^{n+1:N}|X^{1:n} = x^{1:n}}$.\n\nIn kernel methods for HMMs, this is precisely what kernel Bayesian inference is used for [9, 10]. As an alternative, in deterministic settings, a type of Laplace transform of Koopman operator (a variant of CME operator in our work) was used for kernel-based forecasting $x^{1:n} \\mapsto x^{n+1:N}$ of nonlinear dynamics [11]. We conjecture its extension to CME operator may give us conditional mean embedding.\n\nWhile the directions are exciting, we believe their proper investigation warrant a separate work. Thus, we focused our effort on laying the foundations of a neural method based on operator theory.\n\n---\n\n[1] Diffusion-TS: Interpretable diffusion for general time series generation, ICLR 2024\n\n[2] Primer: Searching for efficient transformers for language modeling, NeurIPS 2021\n\n[3] Attention is all you need, NeurIPS 2017\n\n[4] Utilizing image transforms and diffusion models for generative modeling of short and long time series, NeurIPS 2024\n\n[5] Generative modeling of regular and irregular time series data via Koopman VAEs, ICLR 2024\n\n[6] Monash time series forecasting archive, NeurIPS 2021\n\n[7] Time-delay observables for Koopman: Theory and applications, SIAM 2020\n\n[8] Conditional Sig-Wasserstein GANs for time series generation, arXiv 2020\n\n[9] The nonparametric kernel Bayes smoother, AISTATS 2016\n\n[10] Kernel Bayes’ rule: Bayesian inference with positive definite kernels, JMLR 2013\n\n[11] Koopman kernel regression, NeurIPS 2023"}, {"review_id": "uOkwKEXAyn", "reviewer": "Reviewer_Mu8H", "summary": "The paper utilizes an operator-theoretic interpretation of hidden Markov models to propose a novel sequence modeling algorithm. The proposed method, spectral mean flows, utilizes a neural architecture that leverages the spectral decomposition of linear operators for scalability and extends the maximum mean discrepancy (MMD) gradient flows to time-dependent Hilbert spaces. The authors validate the method on both simulated and real-world time-series datasets.", "questions": "Could you provide more explanation and intuition on the right panel of Figure 2? More specifically, why does the curve look as it does (initially it goes up)? Moreover, why does it take much longer for the MMD of spectral mean flow to decrease as much as TD-RBF?\n\nCould you provide visualizations for some of the generated sequences in Section 5.2? I think that could help me understand the results from Table 1 better.\n\nThe paper makes claims on faster convergence and scalability of the spectral mean flow. Could you run ablation studies that show how each part of your method contributes to this efficiency?", "limitations": "Yes.", "rating": 5, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 2, "originality": 3, "strengths_and_weaknesses": "Quality: The submission is technically sound. The background and methods sections are especially technically sound, going into the details. However, I am not fully convinced that the spectral mean flow improves much from other baselines. For example, the predictive scores of spectral mean flow are on par or worse than other baselines across all datasets.\n\nClarity: The paper is mostly clearly written and well organized. It has clear methods, background, and related work sections. However, I think the figure and table captions can be improved, providing more explanations or intuition, especially for Figure 1.\n\nSignificance: Although the idea is novel, and the technical details in the paper are sound, the paper lacks experimental results to back up the proposed method. I think a more careful demonstration of the advantage of using spectral mean flow over other baseline methods will improve the significance of the paper.\n\nOriginality: The work provides new insight through an operator-theoretic view of a hidden Markov model.", "was_revised": true, "final_justification": "The authors addressed my concern about the experiments, which was the main justification for my original score. The additional results that the authors shared during the rebuttal are convincing. I thus decide to increase the score to 5.", "rebuttal": "We thank the reviewer for the constructive feedback. We address the main concerns and questions below.\n\n> W1. “I am not fully convinced that the spectral mean flow improves much from other baselines. For example, the predictive scores are on par or worse than other baselines across all datasets.”\n\nWe would like to point out that predictive score has saturated in many datasets in Table 1 almost to best possible value, e.g., in Sines and Energy, evaluating real data gives predictive scores 0.094 and 0.250, close to scores of our method (0.093 and 0.251). We believe it is more reasonable to assess all four metrics holistically.\n\nIn addition, we did a careful revision of our model (Section 3.4) and were able to improve on high-dimensional Energy and fMRI. The key changes included dropping SwiGLU and adaLN in favor of simpler squared ReLU [1] and sinusoidal embedding [2].\nWe are happy to elaborate on the details during the discussion period.\n\nDue to space limits, the updated Table 1 can be found in **Table A** for reviewer k6rf.\nOur renovated model achieves the best metric in 16 out of 24 cases.\n\n> W2. “The paper lacks experimental results to back up the proposed method. A more careful demonstration of the advantage over other baselines will improve the significance.”\n\nTo address the concern, we ran additional experiments on two challenging real-world setups: (1) long time series and (2) irregularly sampled time series, comparing with recent SoTA ImagenTime [3] and KoVAE [4].\n\nFor long time series, we used FRED-MD and NN5 Daily with 728- and 791-length data [5]. Their modest sizes allowed testing within the rebuttal period. Our models are designed similarly to **Table A**, equipped with time-delay observables used in operator methods [6]. We use marginal, classification, and prediction metrics [3].\n\n**Table B.** Long time series. Baselines are from [3], except for ImagenTime which was reproduced using official code.\n| | FRED-MD | | | NN5 Daily | | |\n|--|--|--|--|--|--|--|\n| | marg. ↓ | class. ↑ | pred. ↓ | marg. ↓ | class. ↑ | pred. ↓ |\n| Ours | **0.0145** | **1.547** | 0.0261 | **0.00453** | **1.07** | 0.427 |\n| ImagenTime | 0.0292 | 1.154 | 0.0279 | 0.00745 | 0.935 | 0.496 |\n| LS4 | 0.022 | 0.544 | 0.037 | 0.007 | 0.636 | **0.241** |\n| SaShiMi | 0.048 | 0.001 | 0.232 | 0.020 \t| 0.045 | 0.849 |\n| SDEGAN | 0.084 | 0.501 | 0.677 | 0.085 | 0.085 | 1.01 |\n| TimeGAN | 0.081 | 0.029 | 0.058 | 0.040 | 0.001 | 1.34 |\n| Latent ODE | 0.042 | 0.327 | **0.013** | 0.107 | 0.000 | 1.04 |\n| ODE2VAE | 0.122 | 0.028 | 0.567 \t| 0.211 | 0.001 | 1.19 |\n| GP-VAE | 0.152 | 0.016 \t| 2.05 | 0.117 | 0.002 | 1.17 |\n| RNN-VAE | 0.132 | 0.036 | 1.47 | 0.137 | 0.000 | 0.967 |\n\nThe results are in **Table B**. Our method shows strong results, achieving the best metric in 4 out of 6 cases and outperforming the previous best (ImagenTime, LS4) in 5 out of 6 cases.\n\nFor irregular time series, we used 3 datasets derived from Google stock data by omitting 30%, 50%, and 70% of observations [4]. The task is generating full signal. Our models are designed similarly to **Table A**. To handle missing values at training, a simple idea of interpolating observed values worked well. Based on the code of [4], we use discriminative metric.\n\n**Table C.** Irregular time series. 30%, 50%, and 70% dropped. Baselines are from [4] except for KoVAE which was reproduced using official code.\n| | Stocks, no drop | 30% | 50% | 70% | | |\n|--|--|--|--|--|--|--|\n| disc. ↓ | | | | |\n| Ours | **0.0150** | **0.1109**| **0.0130**| **0.0628**|\n| KoVAE | 0.0203 | 0.1336 | 0.0899| 0.0740|\n| GT-GAN | 0.077| 0.251 | 0.265 | 0.230 |\n| TimeGAN | 0.102 | 0.411 | 0.477 | 0.485 |\n| RCGAN | 0.196 | 0.436 \t| 0.478 | 0.381 |\n| C-RNN-GAN | 0.399 | 0.500 \t| 0.500 | 0.500 |\n| T-Forcing | 0.226 | 0.305 | 0.308 | 0.365 |\n| P-Forcing | 0.257 | 0.341 | 0.388 | 0.317 |\n| GP-VAE | 0.152 | 0.016 \t| 2.05 | 0.117 |\n| RNN-VAE | 0.132 | 0.036 | 1.47 | 0.137 |\n\nThe results are in **Table C**. Our method achieves the best metric in all 4 cases.\n\n> W3. “I think the figure and table captions can be improved, providing more explanations or intuition, especially for Figure 1.”\n\nThank you for the constructive feedback. We will improve the caption for Figure 1 by writing:\n\n**Figure 1.** A tensor network diagram of the inner product $\\langle \\phi(x^1) \\otimes\\cdots\\otimes \\phi(x^N), \\mu\\_\\rho \\rangle$ (Eq. 3.12) for $N = 5$, showing how its computation involving a higher-order tensor $\\mu\\_\\rho$ decomposes into a series of matrix and vector multiplications.\n\nWe will also add a compact guide in the Appendix on how to read tensor diagrams, and leave pointers to the related material [7, 8].\n\n> W4. “Could you provide more explanation and intuition on the right panel of Figure 2? Why does the curve look as it does (initially it goes up)? Why does it take much longer for the MMD of spectral mean flow to decrease as much as TD-RBF?”\n\nThank you for raising this point. There are two reasons for this.\n\nOne reason comes from our notion of “convergence” of the time-dependent distribution $(p_t)_{t\\geq 0}$. We have a target distribution $\\rho$, and, when we say “fast convergence” in our work, it means $p_t$ converges arbitrarily close to $\\rho$ at a fixed time point $t = 1$ (while naïve MMD flow requires $t \\to\\infty$). This notion of convergence is suitable for generative modeling, where we would like to generate samples in a bounded time.\n\nAs long as this boundary condition is satisfied, we have valid convergence at $t = 1$ and hence do not worry about the behavior, e.g., around $t = 0$. This convergence at $t = 1$ is what the optimal transport (OT) vector field (Eq. 3.16) guarantees, and it is precisely reflected in the right panel of Figure 2, where the MMD reaches $\\approx 0$ at $t = 1$ while TD-RBF saturates around $0.1$.\n\nThere might be other choices of vector fields that guarantee convergence at $t = 1$ while offering more rapid convergence at early timesteps, but in this work, we used the OT vector field since it was the simplest and easiest to use. This choice has the possibility for improvement in future work. We will make this point clear in the next revision.\n\nIn addition, the technical reason TD-RBF appears to decrease rapidly in early timesteps is as follows. In the right panel of Figure 2, the y-axis MMD is measured using an RBF kernel with bandwidth $\\sigma = 1.0.$ By definition, TD-RBF at each time-point t tries to minimize the MMD measured with an RBF kernel with bandwidth $\\sigma(t)$. Since we defined $\\sigma(t) = (1 - t) + 0.1t$, it starts with a bandwidth $\\sigma(0) = 1.0$, which gradually decreases over time. This is one reason for its rapid decrease near $t = 0$: it explicitly decreases the y-axis MMD metric of Figure 2.\n\nYet, despite its early behavior, we remind that TD-RBF stagnates around $t = 1$. This is because it does not have the convergence guarantee at $t = 1$ (unlike Ours), which we consider more important in the context of generative modeling.\n\n> W5. “Could you provide visualizations for some of the generated sequences in Section 5.2?”\n\nThank you for the suggestion. While we can confirm that the generated sequences appear visually realistic, this year’s NeurIPS reviewing and discussion policy unfortunately prevents us from including any figures in the rebuttal [9]. We will include a side-by-side comparison of real and generated sequences for all six datasets in Table 1 in the next revision.\n\n> W6. “The paper makes claims on faster convergence and scalability of the spectral mean flow. Could you run ablation studies that show how each part of your method contributes to this efficiency?”\n\nThank you for the constructive comment. We address each aspect below.\n\nOur claim on faster convergence is in comparison to the standard MMD gradient flow (Lines 138-139), and its naïve extension to time-dependent RKHS (Lines 148-149). Thus, the experiment in Section 5.1 and Figure 2 are precisely the ablation study testing improvements compared to these two baselines.\n\nOur claim on the scalability is in comparison to when the tensor network decomposition (Figure 1, Eq. 3.11) is not available. To illustrate how this idea is contributing to scalability under the hood, we ran a simple numerical experiment for a simplified version of Eq. 3.13: evaluating the inner product $\\langle x^1 \\otimes\\cdots\\otimes x^N, T\\rangle$ between a rank-1 tensor $x^1 \\otimes\\cdots\\otimes x^N$ of vectors $x^n\\in\\mathbb{R}^d$ and a higher-order tensor $T\\in\\mathbb{R}^{d^N}$, for $d=32$.\n\n**Table D.** Peak GPU memory usage (GB) of a simple tensor inner product operation $\\langle x^1 \\otimes\\cdots\\otimes x^N, T\\rangle$ depending on the availability of tensor network decomposition of $T$ (Eq. 3.11). Measured on a single NVIDIA RTX A6000 GPU with 48GB.\n| Sequence length $N$ | w/ tensor network | w/o tensor network |\n|--|--|--|\n| 3 \t| 0.00795 \t| 0.00990 \t|\n| 4 \t| 0.00797 \t| 0.13297 \t|\n| 8 \t| 0.00803 \t| OOM \t|\n| 16 \t| 0.00815 \t| OOM \t|\n| 32 \t| 0.00840 \t| OOM \t|\n| 64 \t| 0.00890 \t| OOM \t|\n| 128 \t| 0.00989 \t| OOM \t|\n| 256 \t| 0.01187 \t| OOM \t|\n| 512 \t| 0.01584 \t| OOM \t|\n| 1024 \t| 0.02381 \t| OOM \t|\n\nThe results are in **Table D**, indeed showing that without the tensor network decomposition, it is almost impossible to evaluate the inner product.\n\n---\n\n[1] Primer: Searching for efficient transformers for language modeling, NeurIPS 2021\n\n[2] Attention is all you need, NeurIPS 2017\n\n[3] Utilizing image transforms and diffusion models for generative modeling of short and long time series, NeurIPS 2024\n\n[4] Generative modeling of regular and irregular time series data via Koopman VAEs, ICLR 2024\n\n[5] Monash time series forecasting archive, NeurIPS 2021\n\n[6] Time-delay observables for Koopman: Theory and applications, SIAM 2020\n\n[7] The tensor cookbook, 2024 https://tensorcookbook.com/\n\n[8] An introduction to graphical tensor notation for mechanistic interpretability, arXiv 2024\n\n[9] https://neurips.cc/Conferences/2025/PaperInformation/NeurIPS-FAQ"}, {"review_id": "goTizGdwmU", "reviewer": "Reviewer_EcXj", "summary": "The authors provide a novel method for generating sequences using a technique called spectral mean flows. This uses a neural network parameterization of a spectral decomposition of linear operators to define a generative model over the space of sequences. They demonstrate this on toy data sets and multiple timeseries dat sets showing good distribution match between ground truth and generated data.", "questions": "Following from the potential weakness above, can this method incorporate information from an informatively sampled system?\n\nThe authors comment that “G is closed under point wise multiplication“ and “in the infinite width limit, multi-layer perceptrons (MLPs) are norm equivalent to Sobolev spaces”. How well do finite width and depth MLPs approximate this space? Does this lead to unique scaling performance with MLP size?", "limitations": "The authors recognize a main limitation that methods to generate sequences conditioned on a partial sequence are not provided.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\nThe authors provide adequate background grounding their extensions. This method makes no assumptions on characteristics of the timeseries being modeled (in comparison to others which use Fourier transforms to explicitly parameterize periodicity). The methods are original.\nThe authors recognize a main limitation that methods to generate conditioned on a partial sequence are not provided.\nWeaknesses:\nThe authors do not comment on whether a timeseries produced through irregular or informative sampling would be handled well by this approach. For example, a system with a measurement generator function M(h) which emits measurements of value x at time t according to the value of some hidden state h such that  both the value x and the time difference between x_t and x_{t+1} are dependent on h.", "was_revised": true, "final_justification": "The authors appropriately responded to my comments and I have kept my original rating.", "rebuttal": "We thank the reviewer for the constructive feedback and positive assessment of our work. We address the main concerns and questions below.\n\n> W1. “The authors do not comment on whether a timeseries produced through irregular or informative sampling would be handled well by this approach.”\n\nThank you for raising this point. Thanks to the HMM formulation, our method can be applied to time series under irregular or informative sampling using the classical approach [1] of treating the sampling interval itself as a part of the time series data. We explain the reasoning behind this and perform an experimental demonstration. There is another potential direction based on continuous-time operators, which we discuss at the end of this response.\n\nLet us consider HMM hidden states $Z^{1:N}$ under Markovian dynamics $P(Z^{n+1}|Z^n)$, and a local observation model $P(X^n, \\Delta \\tau^n|Z^n),$ which jointly determines (emits) the measurement of the value of $X^n$ and the sampling time interval $\\Delta \\tau^n$ between $X^{n-1}$ and $X^n$. The setup mentioned by the reviewer is a special case when $\\Delta \\tau^n$ has no randomness given $Z^n$. If the hidden Markov chain is expressive enough, this HMM can express any joint distribution $P((X^1, \\Delta \\tau^1), …, (X^N, \\Delta \\tau^N))$, i.e., general time series under irregular or informative sampling. To incorporate this into our method, the only change necessary is to use sampling time intervals $\\Delta\\tau^n$ as an additional data dimension.\n\nTo experimentally demonstrate this, we consider three irregular datasets commonly used in literature [2] derived from Google stock price data by randomly omitting 30%, 50%, and 70% of the observations. The task is generative modeling of time series values along with their irregular sampling intervals, $P((X^1, \\Delta \\tau^1), …, (X^N, \\Delta \\tau^N))$. Along with our model, we test Koopman VAE [2], a recent SoTA method for operator-theoretic generative modeling of regular and irregular time series, using the recommended hyperparameters for the stock datasets. Based on [2], we use a discriminative metric for evaluation.\n\n**Table A.** Irregular time series modeling (joint modeling of $P((X^1, \\Delta \\tau^1), …, (X^N, \\Delta \\tau^N))$, where 30%, 50%, and 70% of the observations are dropped.\n| | Stocks, no drop | 30% | 50% | 70% |\n|--|--|--|--|--|\n| disc. ↓ | | | | |\n| Ours \t| **0.0150 ± 0.0101** \t| **0.0498 ± 0.0776** \t| **0.1453 ± 0.1390** \t| **0.0350 ± 0.0147** \t|\n| KoVAE \t| 0.0203 ± 0.0199 \t| 0.3070 ± 0.0587 \t| 0.2030 ± 0.0750 \t| 0.1892 ± 0.0631 \t|\n\nThe results are in **Table A**. Our method achieves strong performances, achieving the best metric in all 4 cases. The result demonstrates the potential of scaling our method to irregular time series, and also provides a comparison against an additional baseline. We will add the results and discussion to the main text in the next revision.\n\nIn addition to the above, another potential direction for handling irregular or informative sampling intervals is using the continuous-time CME operator $U^{\\Delta\\tau}\\phi(z) = \\mathbb{E}[\\phi(Z^{\\tau+\\Delta\\tau}) | Z^\\tau = z]$ and its spectral decomposition $U = \\sum_i \\exp(s_i \\Delta\\tau) u_i \\otimes v_i$ with continuous-time eigenvalues $s_i$ [3]. We leave further investigation of this direction as future work.\n\n> W2. “The authors recognize a main limitation that methods to generate sequences conditioned on a partial sequence are not provided.”\n\nThank you for the comment. While we acknowledge the limitation and leave conditioning as future work, following another reviewer k6rF’s request, we discussed preliminary insights on kernel Bayesian inference or Laplace transforms of CME operators (Section 7) and would like to include the response here as well.\n\nLet’s say we have a trajectory $x^{1:n}$ and want to sample $x^{n+1:N} \\sim P(X^{n+1:N}|X^{1:n} = x^{1:n})$. This can be done with MMD flow if we have conditional mean embedding $\\mu\\_{X^{n+1:N}|X^{1:n} = x^{1:n}}$.\n\nIn kernel methods for HMMs, this is precisely what kernel Bayesian inference is used for [4, 5]. As an alternative, in deterministic settings, a type of Laplace transform of Koopman operator (a variant of CME operator in our work) was used for kernel-based forecasting $x^{1:n} \\mapsto x^{n+1:N}$ of nonlinear dynamics [6]. We conjecture its extension to CME operator may give us conditional mean embedding.\n\nWhile the directions are exciting, we believe their proper investigation warrant a separate work. Thus, we focused our effort on laying the foundations of a neural method based on operator theory.\n\n> W3. “The authors comment that “G is closed under point wise multiplication“ and “in the infinite width limit, multi-layer perceptrons (MLPs) are norm equivalent to Sobolev spaces”. How well do finite width and depth MLPs approximate this space? Does this lead to unique scaling performance with MLP size?”\n\nThank you for the comment. The approximation results in [7, 8] provide a theoretical footing for our final method by connecting neural networks and Sobolev spaces. Performance scaling with MLP size is analyzed in [8]. Indeed, [8, Theorem 1,3] shows that approximation of functions in Sobolev space by MLPs scales like $\\inf_{f_{NN}\\in NN (\\theta)}\\|f- f_{NN}\\| \\leq C\\theta^{-2s/d}$, where $s$ is the smoothness if the function to be approximated, $d$ is the input dimension, and $\\theta$ the number of parameters tied to a fixed, large enough layer width. As such, we (at least theoretically) expect performance to scale with MLP size. Whether this leads to unique scaling law in practice is a more complicated question, since in practice, other factors such as numerical precision and optimizers are involved. We leave the investigation of scaling laws of our model as future work.\n\n---\n\n[1] Zucchini & MacDonald, Eruptions of the Old Faithful geyser, 2009\n\n[2] Naiman et al., Generative modeling of regular and irregular time series data via Koopman VAEs, ICLR 2024\n\n[3] Aït-Sahalia et al., Operator methods for continuous-time Markov processes, Elsevier 2010\n\n[4] Nishiyama et al., The nonparametric kernel Bayes smoother, AISTATS 2016\n\n[5] Fukumizu et al., Kernel Bayes’ rule: Bayesian inference with positive definite kernels, JMLR 2013\n\n[6] Bevanda et al., Koopman kernel regression, NeurIPS 2023\n\n[7] Bietti & Bach, Deep equals shallow for ReLU networks in kernel regimes, ICLR 2021\n\n[8] Siegel, Optimal approximation rates for deep ReLU neural networks on Sobolev and Besov spaces, JMLR 2023"}, {"review_id": "JK4Phk0yJa", "reviewer": "Reviewer_j6aF", "summary": "The paper proposes Spectral Mean Flow (SMF), a new sequence modeling framework that combines an operator‑theoretic view of Hidden Markov Models, a tensor‑network decomposition of the sequence mean embedding derived from the spectral decomposition of the underlying operator, and a time‑dependent Maximum Mean Discrepancy (MMD) gradient flow that is trained with the Flow‑Matching objective. The resulting model is fully differentiable, avoids simulation during training, and samples by solving an ODE. Empirically, SMF attains competitive scores on six multivariate time‑series data sets compared to diffusion, VAE, and GAN baselines.", "questions": "See Weaknesses", "limitations": "See Weaknesses", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 2, "significance": 2, "originality": 3, "strengths_and_weaknesses": "**Strengths:**\n\n* The paper unifies three previously disjoint ideas (operator theory, tensor networks, flow matching) into a coherent generative framework. In particular, the extension of MMD flows to time‑dependent RKHSs (Eq. 3.15) and its connection to flow matching appear original.\n* The authors discuss complex‑valued parameterization, low‑rank sharing, and stability tricks (SwiGLU, RMS‑LN, sharing most of the computations), which will be valuable for practitioners. \n\n\n**Weaknesses:**\n\n* The paper does not convincingly articulate the practical advantages of this new paradigm over standard flow-matching or diffusion models. Because the operator is not modeled explicitly, potential benefits such as interpretability through mode analysis remain out of reach, and the added complexity may not translate into clear gains.\n* All empirical studies are restricted to fixed-length sequences of 24 time steps. The authors do not investigate scalability to longer or variable-length sequences, nor do they discuss any practical limitations in this regard.\n* A detailed complexity analysis is missing: the authors do not report wall-clock training or inference times, FLOPs, or parameter counts, nor do they compare these metrics with competing baselines.\n* Important recent baselines, such as SDFormer [1], ImagenTime [2], FiDE [3], KoVAE [4] and other recent works, are omitted, making it hard to judge the proposed method’s relative performance.\n\n[1] Chen, Zhicheng, et al. \"Sdformer: Similarity-driven discrete transformer for time series generation.\" NeurIPS 2024\n\n[2] Naiman, Ilan, et al. \"Utilizing image transforms and diffusion models for generative modeling of short and long time series.\" NeurIPS 2024\n\n[3] Galib, Asadullah Hill, Pang-Ning Tan, and Lifeng Luo. \"FIDE: Frequency-Inflated Conditional Diffusion Model for Extreme-Aware Time Series Generation.\" NeurIPS 2024.\n\n[4] Naiman, Ilan, et al. \"Generative modeling of regular and irregular time series data via Koopman VAEs.\" ICLR 2024.\n\nOther related work:\n\n[5] Coletta, Andrea, et al. \"On the constrained time-series generation problem.\" NeurIPS 2023.\n\n[6] Jeon, Jinsung, et al. \"GT-GAN: General purpose time series synthesis with generative adversarial networks.\" NeurIPS 2022", "was_revised": true, "final_justification": "I appreciate the authors' response that address my primary concerns. In light of their thorough rebuttal, I have raised my score to 4.\n\nWhile the reason for omitting the ImagenTime results on the Cont.-FID metric remains unclear, I have two recommendations for the final version. The authors should incorporate the work on  SDformer[1], ImagenTime [2], and KoVAE [3] into the related work section and discuss key difference with Koopman approach in KoVAE. Furthermore, the new results, which are currently only in the rebuttal, should be moved into the main paper. Their inclusion is important as they effectively position the authors' method as being on par with or competitive against other state-of-the-art techniques.\n\nReferences\n\n[1] \"SDformer: Similarity-driven discrete transformer for time series generation.\" NeurIPS 2024\n\n[2] \"Utilizing image transforms and diffusion models for generative modeling of short and long time series.\" NeurIPS 2024\n\n[3] \"Generative modeling of regular and irregular time series data via Koopman VAEs.\" ICLR 2024.", "rebuttal": "We thank the reviewer for the constructive feedback. We address the main concerns and questions below.\n\n(note: we revised our architecture after submission, which led to improvements on high-dimensional Energy and fMRI in Table 1. Please see **Table A** for reviewer k6rf. All new experiments are done with this version. We are happy to elaborate on the details during the discussion period.)\n\n> W1. “All empirical studies are restricted to fixed 24 timesteps. The authors do not investigate scalability to longer or variable-length sequences.” “Important recent baselines are omitted, making it hard to judge the method’s relative performance.”\n\nTo address the concerns jointly, we ran additional experiments on more challenging setups: (1) modeling long time series and (2) learning on irregularly sampled time series and then generating full time series, comparing with recent SoTA ImagenTime [1] and Koopman VAE (KoVAE) [2] suggested by the reviewer.\n\n(note: while the reviewer also suggested SDFormer [3] and FIDE [4], we could not find code for [3], and official code of [4] did not contain sufficient information for reproducibility, preventing us from experimenting. Still, we will cite and discuss [1-4] in our next revision and also [5, 6] as suggested by the reviewer.)\n\nFor long time series, we used FRED-MD and NN5 Daily with 728- and 791-length data [7]. Their modest sizes allowed testing within the rebuttal period. Our models are designed similarly to **Table A**, equipped with time-delay observables used in operator methods [9]. We use marginal, classification, and prediction metrics [1].\n\n**Table B.** Long time series. Baselines are from [1], except for ImagenTime which was reproduced using official code.\n| | FRED-MD | | | NN5 Daily | | |\n|--|--|--|--|--|--|--|\n| | marg. ↓ | class. ↑ | pred. ↓ | marg. ↓ | class. ↑ | pred. ↓ |\n| Ours | **0.0145** | **1.547** | 0.0261 | **0.00453** | **1.07** | 0.427 |\n| ImagenTime | 0.0292 | 1.154 | 0.0279 | 0.00745 | 0.935 | 0.496 |\n| LS4 | 0.022 | 0.544 | 0.037 | 0.007 | 0.636 | **0.241** |\n| SaShiMi | 0.048 | 0.001 | 0.232 | 0.020 \t| 0.045 | 0.849 |\n| SDEGAN | 0.084 | 0.501 | 0.677 | 0.085 | 0.085 | 1.01 |\n| TimeGAN | 0.081 | 0.029 | 0.058 | 0.040 | 0.001 | 1.34 |\n| Latent ODE | 0.042 | 0.327 | **0.013** | 0.107 | 0.000 | 1.04 |\n| ODE2VAE | 0.122 | 0.028 | 0.567 \t| 0.211 | 0.001 | 1.19 |\n| GP-VAE | 0.152 | 0.016 \t| 2.05 | 0.117 | 0.002 | 1.17 |\n| RNN-VAE | 0.132 | 0.036 | 1.47 | 0.137 | 0.000 | 0.967 |\n\nThe results are in **Table B**. Our method shows strong results, achieving the best metric in 4 out of 6 cases and outperforming the previous best (ImagenTime, LS4) in 5 out of 6 cases.\n\nFor irregular time series, we used 3 datasets derived from Google stock data by omitting 30%, 50%, and 70% of observations [2]. The task is generating full signal. Our models are designed similarly to **Table A**. To handle missing values at training, a simple idea of interpolating observed values worked well. Based on the code of [2], we use discriminative metric.\n\n**Table C.** Irregular time series. 30%, 50%, and 70% dropped. Baselines are from [2] except for KoVAE which was reproduced using official code.\n| | Stocks, no drop | 30% | 50% | 70% | | |\n|--|--|--|--|--|--|--|\n| disc. ↓ | | | | |\n| Ours | **0.0150** | **0.1109**| **0.0130**| **0.0628**|\n| KoVAE | 0.0203 | 0.1336 | 0.0899| 0.0740|\n| GT-GAN | 0.077| 0.251 | 0.265 | 0.230 |\n| TimeGAN | 0.102 | 0.411 | 0.477 | 0.485 |\n| RCGAN | 0.196 | 0.436 \t| 0.478 | 0.381 |\n| C-RNN-GAN | 0.399 | 0.500 \t| 0.500 | 0.500 |\n| T-Forcing | 0.226 | 0.305 | 0.308 | 0.365 |\n| P-Forcing | 0.257 | 0.341 | 0.388 | 0.317 |\n| GP-VAE | 0.152 | 0.016 \t| 2.05 | 0.117 |\n| RNN-VAE | 0.132 | 0.036 | 1.47 | 0.137 |\n\nThe results are in **Table C**. Our method achieves the best metric in all 4 cases.\n\n> W2. “The paper does not convincingly articulate the practical advantages over standard flow-matching or diffusion. Because the operator is not modeled explicitly, potential benefits such as interpretability through mode analysis remain out of reach, and the added complexity may not translate into clear gains.”\n\nWhile interpretation of eigenfunctions is not direct due to our end-to-end design, we believe that a closely related benefit of operator methods is retained to some degree: incorporation of physics-based prior knowledge. Thanks to linearity of tensor network, this can be done via explicit spectral regularization. We test this using a challenging physics-constrained problem from [2].\n\nWe consider a nonlinear pendulum, governed by an ODE of angular displacement $\\theta$ from an equilibrium [2]: $\\frac{d^2}{dt^2}\\theta + 9.8\\sin\\theta = 0, \\dot{\\theta}(0)=0$.\nAs the pendulum is stable and conservative, the physical knowledge is that eigenvalues of underlying operator have non-positive real part and some have zero real.\n\nSince our model is governed by a linear product of matrix hidden states $Ow(x^n)$ (Eq. 3.12), we can implement this knowledge in an end-to-end manner.\nIn KoVAE [2], this stability constraint is built-in using a loss on matrix $A$ for hidden states $z_{t+1} = Az_t$, specifically $|s_1 - 1|^2 + |s_2 - 1|^2$ where $(s_1, s_2)$ are the two largest eigenvalues of $A$. We compare our model against [2] using the same spectral loss for each matrix hidden state’s two largest eigenvalues.\nWe use correlational score for eval. [9]\n\n**Table D.** Physics-informed nonlinear pendulum modeling.\n|  | w/o stability loss | w/ stability loss |\n|--|--|--|\n| corr. ↓ | | |\n| Ours | 2.298e-3 | 2.320e-4 |\n| KoVAE | 2.164e-3 | 6.954e-4 |\n\nThe results are in **Table D**. Without physics-based loss, ours performs on par with KoVAE. When the physics-based loss is incorporated, we observe a clear improvement, outperforming KoVAE which uses the same stability constraint. This supports incorporating physical knowledge into our model is possible, resulting in a better match with the true data-generating process.\n\n> W3. “A detailed complexity analysis is missing.”\n\nWe conducted detailed measurements of wall-clock training and inference times, as well as real and complex parameter counts, for the models in Table 1. The time measurements are done on a single NVIDIA RTX A6000 GPU with 48GB and Intel Xeon Gold 6330 CPU @ 2.00GHz. The results in **Tables E, F, G** show that our models for Table 1 has reasonable computational complexity, and their generation is actually faster than the previous SoTA diffusion model (Diffusion-TS).\n\n**Table E.** Wall-clock training times.\n| | Sines | Stocks | ETTh | MuJoCo | Energy | fMRI |\n|--|--|--|--|--|--|--|\n| Ours | 18min | 15min | 27min | 60min | 44min | 35min |\n| Diffusion-TS | 16min | 15min | 30min | 25min | 60min | 48min |\n\n**Table F.** Wall-clock sampling times per 1,000 samples.\n| | Sines | Stocks | ETTh | MuJoCo | Energy | fMRI |\n|--|--|--|--|--|--|--|\n| Ours | 30sec | 30sec | 29sec | 109sec | 95sec | 109sec |\n| Diffusion-TS | 66sec | 66sec | 67sec | 133sec | 202sec | 262sec |\n\n**Table G.** Numbers of parameters. The parameter counts of baselines are from Table 7 of [10].\n| | Sines | Stocks | Energy |\n|--|--|--|--|\n| Ours | 315,136 | 315,200 | 1,088,160 |\n| Diffusion-TS | 232,177 | 291,318 | 1,135,144 |\n| Diffwave | 533,592 | 599,448 | 1,337,752 |\n| Cot-GAN | 40,133 \t| 52,675 | 601,539 |\n| TimeGAN | 34,026 | 48,775 | 1,043,179 |\n| TimeVAE | 97,525 \t| 104,412 | 677,418 |\n\n---\n\n[1] Utilizing image transforms and diffusion models for generative modeling of short and long time series, NeurIPS 2024\n\n[2] Generative modeling of regular and irregular time series data via Koopman VAEs, ICLR 2024\n\n[3] SDformer: Similarity-driven discrete transformer for time series generation, NeurIPS 2024\n\n[4] FIDE: Frequency-inflated conditional diffusion model for extreme-aware time series generation, NeurIPS 2024\n\n[5] On the constrained time-series generation problem, NeurIPS 2023\n\n[6] GT-GAN: General purpose time series synthesis with generative adversarial networks, NeurIPS 2022\n\n[7] Monash time series forecasting archive, NeurIPS 2021\n\n[8] Time-delay observables for Koopman: Theory and applications, SIAM 2020\n\n[9] Conditional Sig-Wasserstein GANs for time series generation, arXiv 2020\n\n[10] Diffusion-TS: Interpretable diffusion for general time series generation, ICLR 2024"}]

Author rebuttal:

Meta-review: The paper proposes a sequence modeling framework Spectral Mean Flow (SMF), which combines an operator‑theoretic view of Hidden Markov Models, a neural network implementation of spectral decomposition of the underlying operator, and a time‑dependent MMD gradient flow that is trained with the Flow‑Matching objective. Empirically, the method achieves good performance on several time sequence modeling tasks with reasonable time efficiency. There have been quite some discussions, additional related works, and new experiments added during the discussion period and the authors shall incorporate them into the final version.
Decision: Accept (poster)