Data-driven spectral decomposition and forecasting of ergodic dynamical systems

1.1. Background and motivation

• Ergodicity enables nonparametric dimension reduction and forecasting from long time series without prior knowledge of the governing equations, including when systems are unknown, infeasible to simulate, or complex.Long-time averages equal expectations under the invariant measure, allowing time-series data to sample the full phase space.Source: §1.1 · p. 1
• Koopman and Perron-Frobenius operators provide complementary linear representations of nonlinear dynamics for spectral analysis, coherent-set identification, mode decomposition, dimension reduction, and nonparametric prediction.The paper uses diffusion maps to learn smooth orthonormal bases for Galerkin and spectral methods involving these operators.Source: §1.1 · p. 2
• Continuous spectra, dense eigenfrequencies, high-dimensional ambient spaces, and non-intrinsic standard Koopman modes create challenges requiring regularization, geometry-adapted bases, and modality-consistent mode representations.For weak-mixing systems, nonconstant Koopman eigenfunctions may be absent; for irrational torus flows, dense frequencies can produce poorly conditioned schemes, while ambient-space bases may scale exponentially with dimension.Source: §1.1 · pp. 3–4

1.2. Contributions of this work

• The paper constructs diffusion-maps orthonormal bases for Koopman spectral and Galerkin methods, vector-field reconstruction, and forecasting from single time series, including partial observations satisfying a Takens delay-embedding condition.The basis is supported on the data manifold, converges under mild assumptions to a complete L2 basis, scales linearly with ambient dimension for computation, and supports Koopman–Perron-Frobenius conversion by matrix transpose.Source: §1.2 · p. 4
• For pure point spectra, the method uses H1 Galerkin regularization and diffusion geometry to select smooth Koopman eigenfunctions, generate commuting vector-field components, and obtain modality-natural tensorial modes.Delay-coordinate diffusion maps connect Laplace-Beltrami and Koopman operators, while factor maps produce projectible simple harmonic oscillations and noise-robust Galerkin schemes.Source: §1.2 · p. 5
• Koopman eigenfunctions forecast arbitrary L2 observables and probability measures as uncoupled harmonic oscillators, while time change extends decomposition and forecasting to complex spectra with coupled, time-dependent dynamics.The time-change transformation preserves orbits, uses the ambient vector-field norm without governing equations, and can recover pure point spectra in special cases.Source: §1.2 · pp. 5–6

• The setting is a compact, connected, orientable smooth manifold with a smooth invertible flow preserving a probability measure, observed at fixed intervals through a smooth embedding whose induced metric is extrinsic.If the observation is not an embedding, delay-coordinate maps can provide one with high probability.Source: §2 · p. 6
• The Koopman group is unitary and strongly continuous on L2, with a skew-adjoint generator extending differentiation along the flow; measure preservation gives the vector field the Liouville property.The generator acts as a derivation on a dense invariant subalgebra, while Koopman operators distribute over products.Source: §2 · p. 7
• Koopman eigenvalues are purely imaginary, distinct eigenfunctions are orthogonal, and each eigenfunction evolves along trajectories as a simple harmonic oscillator with intrinsic frequency.Ergodicity also makes the zero eigenvalue simple and normalizes eigenfunctions to the unit circle, yielding circle-rotation factor maps.Source: §2 · pp. 7–8
• Ergodicity makes temporal averages converge to spatial averages almost everywhere, allowing L2 inner products and related Sobolev quantities to be estimated by Monte Carlo sums from sampled time series.For nonergodic systems, the same sums converge on the ergodic component containing the initial point.Source: §2 · p. 8
• For pure-point systems, eigenfunctions form an orthonormal basis with additive eigenvalues and multiplicative eigenfunctions, enabling a Koopman Fourier representation in which observables evolve as uncoupled oscillators and observations decompose into Koopman modes.With l rationally independent basic frequencies, all spectral frequencies are integer combinations of those generators, although for l ≥ 2 they are dense and non-isolated.Source: §2 · pp. 8–9
• Mixing systems have a zero eigenvalue plus continuous spectrum, so the paper restricts initial analysis to pure-point systems and later uses time-change regularization to seek generators more amenable to numerical eigendecomposition.The main continuous-time formulation is retained because it provides skew-adjointness, the Leibniz rule, data-space reconstruction, and coercivity of a regularized generator.Source: §2 · p. 9

3.1. Intrinsic dimension reduction coordinates

• Koopman eigenfunctions define projection coordinates whose composite map is a diffeomorphism onto an m-torus, with each coordinate evolving as a simple harmonic oscillator and remaining consistent across observation modalities.The individual projections map onto the unit circle, while rationally independent frequencies yield an ergodic torus rotation; compatible generating eigenfunctions can parameterize data from different sensors.Source: §3.1 · pp. 10–11

3.2. Vector field decomposition

• For smooth generating eigenfunctions, the dynamical vector field decomposes into mutually commuting, measure-preserving, non-ergodic components whose composition recovers the full evolution.The components are observation-map-independent and their pushforwards reconstruct tangent vector fields in data space that transform naturally under changes of observation modality.Source: §3.2 · pp. 11–13

3.3. Nonparametric forecasting

• Koopman-basis expansions provide nonparametric forecasts of evolving probability densities, observable expectations, and variances through diagonal harmonic evolution of the spectral coefficients.Forecasts can also start from a single ambient-space observation using out-of-sample eigenfunction extension, including observations off an imperfect model manifold.Source: §3.3 · pp. 13–14

3.4. Irrational flow on the 2-torus

• In the irrational 2-torus flow, low-roughness generating eigenfunctions produce unit-circle coordinates with frequencies 1 and α, and the projected dynamics are simple harmonic oscillators.The reconstructed commuting components are explicitly depicted for the standard R3 torus embedding in Fig. 1.Source: §3.4 · p. 14, Figure 1 · p. 15
• The torus example shows that eigenfunctions can have arbitrarily small frequencies but arbitrarily large Dirichlet energy, motivating diffusion-based suppression of highly oscillatory modes.This spectral pathology can impair numerical conditioning and means that identifying slow observables also requires controlling smoothness.Source: §3.4 · p. 14
• The forecast observable is the first component F^1 of the observation map, evaluated under a circular-Gaussian density with location (0, 0) and concentration 30.The time-dependent expectation and standard deviation are computed relative to Haar measure.Source: §3.4 · p. 15, Figure 2 · p. 16
• The forecast moments can be evaluated analytically using modified Bessel functions and the expansion of F^1 in the Koopman eigenfunction basis.The relevant Bessel function is I_n, the modified Bessel function of order n.Source: §3.4 · p. 15
• The statistical forecast produces separate mean and standard-deviation forecasts for x1, shown in Fig. 2.The initial density uses frequencies (1, 301/2) and the stated circular-Gaussian parameters.Source: §3.4 · p. 15, Figure 2 · p. 16

4.1. Choice of Galerkin approximation space

• The conformal metric h = gσ2/m is observation-map invariant within conformal classes, has uniform invariant-measure volume, and yields low-roughness basis functions with enhanced resolution where sampling density is high.Its Laplace–Beltrami eigenfunctions minimize expected roughness, while rescaled eigenfunctions provide a well-conditioned H1 basis with equal nonconstant Dirichlet energies.Source: §4.1 · pp. 16–18
• The h-based basis is more robust to sampling-density fluctuations than the ambient-metric alternative, although analogous error estimates are not established for correlated ergodic samples.The same Galerkin scheme can also use the weighted Laplacian basis, but the h-based construction is reported to behave more stably when sampling density varies strongly.Source: §4.1 · pp. 17–18

4.2. Data-driven orthonormal basis

• Diffusion maps construct a compact ergodic Markov operator whose eigenvalues and eigenfunctions converge to those of the heat operator and Laplace–Beltrami operator for h, yielding the data-driven basis.For suitable bandwidth scaling, the discrete eigenvectors converge to basis functions and −log ˆκi/ϵ converges to ηi; the weighted discrete inner product also approaches the L2(M, µ) inner product.Source: §4.2 · pp. 19–20

4.3. Spectral Galerkin method

• The regularized generator Lε = v + ε∆h is solved in weak form on H1(M, µ), where diffusion suppresses highly oscillatory eigenfunctions and skew-symmetry of v ensures well-posedness.The finite-dimensional approximation uses the data-driven basis, discrete Dirichlet and mass matrices, and second-order time differences to estimate vector-field action from trajectory data.Source: §4.3 · pp. 20–21
• The Galerkin discretization solves a generalized eigenvalue problem in a basis scaled to the weak problem’s H1 regularity, keeping the highest-order form equal to the identity.This scaling avoids the unbounded condition-number growth expected in the unscaled eigenfunction basis.Source: §4.3 · pp. 21–22
• Approximation errors vanish under sampling, diffusion-map, finite-difference, and Galerkin limits, but diffusion regularization remains difficult for non-H1 Koopman eigenfunctions or continuous spectra.Non-orthogonal regularized eigenfunctions require inverse transforms using a generally nonidentity Gram matrix.Source: §4.3 · pp. 22–23

4.4. Spectral properties of the regularized generator

• When the dynamical and Laplace-Beltrami operators commute, regularized eigenfunctions are orthogonal and exact for the dynamics; otherwise diffusion suppresses oscillatory modes and can produce problematic spectra.For weak-mixing systems, diffusion may yield point spectrum whose eigenfunctions have questionable physical meaning, motivating time change rather than diffusion alone,.Source: §4.4 · pp. 23–24

4.5. Applications to variable-speed flows on the 2-torus

• Using only time-ordered data, the method detected the variable-speed torus’s integrable structure, with generating eigenfunctions showing simple harmonic oscillations and variation purely along angles θ1 and θ2.The system mixes two timescales and exhibits a density contrast of approximately 34 relative to the Haar measure.Source: §4.5 · p. 25
• The reconstructed commuting vector fields achieved 0.2% root-mean-square error using 3721 basis functions, while 121 functions yielded 5.9% error.The components separately describe flows along the θ1 and θ2 directions, and their expansion coefficients are strongly sparse.Source: §4.5 · pp. 25–26
• The nonparametric density forecast’s means closely matched a 10,000-sample perfect-model ensemble, while standard deviations were generally accurate but overestimated during low-uncertainty periods.Finite-l truncation may limit accuracy when the evolving density develops small-scale structure and broad expansion bandwidth.Source: §4.5 · p. 27, Figure 6 · p. 27

5.1. Relationship between Koopman and Laplace-Beltrami operators in delay-coordinate space

• Delay-coordinate diffusion maps produce a Koopman-commuting operator, making the Galerkin approximation space efficient for computing Koopman eigenfunctions in pure-point systems.The limiting diffusion eigenfunctions converge to Koopman eigenfunctions and provide maximal timescale separation because each evolves at a single frequency.Source: §5 · p. 27, §5.1 · pp. 29, 31
• For systems with nonzero Lyapunov exponents, increasing delays can produce a nonsmooth limiting metric and bias diffusion eigenfunctions toward stable directions, potentially hindering Koopman approximation.This limitation motivates schemes that approximate Koopman eigenfunctions without delay-coordinate maps.Source: §5.1 · pp. 28, 29, 31
• In pure-point systems with smooth Koopman eigenfunctions, the infinite-delay metric is smooth, flow-invariant, flat, and has a volume form with uniform density relative to the invariant measure.These properties make the associated vector fields Killing fields and ensure commutation with the Laplace-Beltrami operator.Source: §5.1 · p. 29
• 800 delays make the Koopman eigenfunction spectra significantly sparser, with ζ1 projecting almost entirely onto the leading two-dimensional diffusion eigenspace.The delay-coordinate dataset has ambient dimension 2400 but intrinsic dimension 2, reflecting the torus structure.Source: §5.1 · p. 30, Figure 8 · p. 30

5.2. Koopman eigenfunctions from noisy data

• Delay-coordinate maps can denoise i.i.d. observational noise because the noise-induced distance bias cancels in diffusion-map normalization as the number of delays grows, enabling consistent Laplace-Beltrami approximation.Theoretically, arbitrarily large noise variance can be removed with sufficiently many delays, while finite samples constrain feasible noise levels through s ≪ N.Source: §5.2 · pp. 31–32
• The modified noisy-data algorithms recover generating Koopman eigenfunctions under moderate and strong Gaussian noise, matching noise-free frequencies within 3 significant figures and outperforming EDMD in the strong-noise case.The approach and EDMD perform comparably under moderate noise, whereas the proposed method yields higher-quality numerical eigenfunctions under strong noise.Source: §5.2 · p. 32

6.1. Time change in dynamical systems

• Time change multiplies a vector field by a positive smooth function, preserving its orbits while changing the time parametrization and invariant measure.The transformed flow remains ergodic, whereas suitable smooth time changes can convert an ergodic flow into a mixing one.Source: §6.1 · pp. 32, 34

6.2. Dimension reduction and nonparametric forecasting with time change

• The time-change approach regularizes mixing by replacing simple harmonic oscillators with variable-frequency oscillators, enabling Koopman-based dimension reduction, vector-field decomposition, and forecasting.Time change preserves ergodicity but can alter mixing; the resulting vector-field components generally do not commute.Source: §6.2 · pp. 34–35, §6.1 · p. 34
• For systems related to pure-point-spectrum dynamics by time change, Koopman eigenfunctions retain constant-amplitude projections with nonconstant frequency, while the combined eigenfunction map remains an embedding with projectible dynamics.Forecasting uses finite-dimensional truncations because the non-multiplicative operator prevents a generally closed-form solution; raw forecasts perform well at short to medium lead times but develop longer-lead biases.Source: §6.2 · pp. 35–36

6.3. Data-driven implementation

• The data-driven method uses the observable vector-field norm ξ, estimated by finite differences, as an accessible time-change function and constructs eigenfunctions through a ξ-adjusted Laplace–Beltrami basis and advection–diffusion operator.It can recover the true time change in special cases and applies modified vector-field decomposition and forecasting algorithms accounting for the nonunitary Fourier transform.Source: §6.3 · pp. 36–37

6.4.1. Mixing flow on the 3-torus

• In the 3-torus experiment, the speed variations exceeded a factor of 10, producing phase-modulated wavetrains and a highly nonuniform invariant density; the reconstructed vector-field decomposition had 1% error in R6.Without time change, regularized eigenfunctions showed near-periodic or strongly amplitude-modulated patterns rather than coherent oscillator behavior.Source: §6.4.1 · p. 38, Figure 14 · p. 40
• In the mixing 3-torus experiment, time change recovered generating eigenfunctions with frequencies 0.5587, 0.2498, and 0.7902, agreeing within three significant figures with the irrational-flow frequencies.Their real eigenvalue parts were approximately 5 × 10^-4, indicating minimal diffusion effects after time change.Source: §6.4.1 · p. 38, Figure 13 · p. 39
• In the mixing 3-torus example, the recovered vector fields V1, V2, and V3 generate flows along the θ2, θ1, and θ3 directions, respectively.Source: §6.4.1 · p. 40
• For the mixing 3-torus, the nonparametric forecast tracks ensemble means and standard deviations before the equilibrium relaxation time but later develops spurious oscillations.Finite deterministic oscillators cannot reproduce the arbitrarily fine density scales required for correct long-time relaxation without diffusion; increasing bandwidth may suppress fixed-time aliasing but cannot remove eventual failure.Source: §6.4.1 · p. 40, Figure 15 · p. 41

6.5. Ergodic flow on the 2-torus with a fixed point

• Time change markedly improves diffusion-map eigenfunctions for the fixed-point torus, where unmodified methods fail because trajectory sampling becomes singular near the fixed point.The induced degenerate metric strengthens heat-kernel smoothing near the vanishing-speed fixed point, producing a well-behaved basis.Source: §6.5 · pp. 41–42
• With ε = 0.02, the fixed-point system yields phase-modulated generating eigenfunctions with frequencies {0.735, 0.165}, energies {1.54, 2.42}, and time-dependent scale separation.The faster mode is transverse to the flow, while the slower mode has level sets aligned with dynamical orbits; lagged-state EDMD produces noisier or missing analogues.Source: §6.5 · p. 42
• The fixed-point torus vector-field decomposition reconstructs the full field with 7% time-averaged error at bandwidth l = 30, separating distinct directional flow components.V1 is primarily directed along θ2 except near θ1 = θ2, while V2 has substantial magnitude in that band and opposing inflow/outflow regions.Source: §6.5 · pp. 42–43
• Over [0, 10], time-change forecasting of the fixed-point torus closely matches ensemble means and reasonably estimates standard deviations, though uncertainty errors are larger.The approach also improves the diffusion-map basis and regularizes the dynamics through an orbit-equivalent system with improved spectral properties.Source: §6.5 · p. 43

Core framework

• The framework learns a smooth complete L2 basis from time-ordered data, enabling Koopman-based dimension reduction and mode decomposition alongside Perron-Frobenius forecasting of densities and observables.The basis connects the Koopman and Perron-Frobenius representations through standard linear algebra operations and supports equation-free forecasting.Source: §7 · p. 45

Pure-point systems

• For pure-point-spectrum systems, Koopman eigenfunctions yield smooth projectible coordinates, commuting vector-field decompositions, improved delay-coordinate Galerkin methods, and uncoupled harmonic-oscillator dynamics.These methods were demonstrated under multiple timescales, strongly nonuniform sampling density, and substantial i.i.d. observational noise.Source: §7 · p. 45

Delay coordinates

• The paper establishes that diffusion maps on delay-coordinate data approximate Laplace-Beltrami operators linked to Koopman operators, explaining timescale separation in diffusion coordinates for pure-point systems.The connection also improves the efficiency and noise robustness of Galerkin eigenvalue computations.Source: §7 · p. 45

Spectral regularization

• Adding small basis-tailored diffusion regularizes Koopman eigenvalues and eigenfunctions by removing highly oscillatory, high-Dirichlet-energy modes, including when irrational flows lack isolated generator eigenvalues.This regularization targets smoothness for the observation modality and is important even in simple pure-point systems.Source: §7 · p. 45

Mixing systems

• For mixing systems, time-change regularization rescales the generator by the ambient-data vector-field norm, preserves orbits while altering time parametrization, and supports vector-field decomposition and nonparametric forecasting.The approach produces non-commuting vector fields and coupled oscillators, and was demonstrated on mixing flows on the 3-torus and a fixed-point flow on the 2-torus where diffusion maps otherwise became ill-conditioned.Source: §7 · p. 45

Future research

• Future work should extend the methods beyond toroidal systems to mixed spectra and nonzero-Lyapunov systems, broaden time-change transformations, and replace inefficient global L2 bases with multiscale bases for localized densities.These directions address both the paper’s dynamical-system scope and the limitations of global spectral expansions.Source: §7 · p. 46

Algorithm 1: Data-driven orthonormal basis

• Algorithm 1 constructs a data-driven orthonormal L2(M, µ) basis using diffusion maps, kernel density estimation, bandwidth tuning, and a symmetric eigenproblem.For large datasets, nearest-neighbor truncation and the similarity-transformed matrix S improve computational efficiency and stability, with small accuracy loss from truncation.Source: Appendix A · pp. 46–47

Implementation considerations

• The algorithms are conceptual rather than fully optimized implementations, and the appendix recommends fixed choices knn = 8 and ϵ_l = 2^l because the procedure is not highly sensitive to these parameters.The stated bandwidth grid uses l from −30 through 10 in increments of 0.1.Source: Appendix A · pp. 46–47

Algorithm 2: Generating frequencies and eigenfunctions

• Algorithm 2 computes Koopman generating frequencies and eigenfunctions by solving a regularized generalized eigenproblem, normalizing eigenfunctions, and selecting the lowest-energy rationally independent frequencies.The number of basic frequencies may be set to the dimension estimated by Algorithm 1.Source: Appendix A · pp. 47–48

Algorithm 3: Vector field decomposition

• Algorithm 3 reconstructs vector fields by expanding products of Koopman eigenfunctions, forming their weighted Gram matrix, and applying frequency-indexed diagonal operators to the observation map.The procedure uses multi-index products with bounded spectral order and reconstructs each vector field as a data-space matrix.Source: Appendix A · p. 48

Algorithm 4: Nonparametric forecast

• Algorithm 4 forecasts probability densities and vector-valued observables by expanding them in Koopman products and advancing each density coefficient as e^(-iω_jt).The forecast reconstructs densities from the time-evolved coefficients after projecting the initial density with the inverse Gram matrix.Source: Appendix A · pp. 48–49

Metric identity

• The constructed metric satisfies ¯g(v_i, v_j) = B_ij because the vector fields are nowhere-vanishing and linearly independent.The calculation uses the defining relations and Theorem 5; Fourier-coefficient decay justifies the required limit and summation interchanges.Source: Appendix B · pp. 49–50

Isometry

• The dynamical map Φ_n,t* is an isometry of ¯g, as verified by checking preservation of inner products on the vector-field basis.The verification applies Theorem 5(iv) to ¯g(v_i|_a, v_j|_a) and its pushforward under Φ_n,t*.Source: Appendix B · p. 50

Flatness

• The metric ¯g is globally flat because the v_i form globally defined, linearly independent, mutually commuting Killing vector fields.Each v_i generates a one-parameter group of isometries, and the standard local flatness result extends globally.Source: Appendix B · p. 50

Volume form

• The volume ratio dvol¯g/dµ is globally constant because ¯g has constant matrix B in vector-field coordinates and the v_i have zero µ-divergence.Locally, the ratio is √det B/γ(θ), while vanishing divergence makes γ constant in every coordinate direction.Source: Appendix B · pp. 50–51

Asymptotic setting

• For pure-point ergodic systems with smooth Koopman eigenfunctions, the noisy diffusion-map procedure is analyzed as ε→0, s→∞, and N(s)≫s.The conclusions hold for almost every training starting state, almost surely in the noise, and uniformly over states.Source: Appendix C · p. 51

Noise in delay space

• With non-overlapping delay sequences, noisy squared distances have expectation ∥Xi−Xj∥2+2R2, while their noise contribution converges almost surely to 2R2.The fraction of overlapping delay pairs vanishes as s/N(s)→0.Source: Appendix C · pp. 51–52

Kernel normalization

• Fixed-bandwidth kernels acquire a constant multiplicative bias e−2R2/ϵ that cancels between numerator and denominator, whereas variable-bandwidth normalization would retain an index-dependent bias.This cancellation supports the noisy-data convergence of the normalized diffusion-map quantities.Source: Appendix C · p. 52

Sampling-density estimation

• In the large-data, infinite-delay limit, the density estimate τϵ is constant up to O(ϵ2) and estimates the true sampling density up to proportionality.The limiting delay-space metric has uniform volume form relative to the invariant measure, producing the constant-density result.Source: Appendix C · p. 52

Averaging-operator approximation

• With sufficiently many delays, i.i.d. noise induces an O(ϵ2) bias in approximating Pϵf, matching the intrinsic diffusion-map discretization error.Thus, (I−P̃f_i)/ϵ approximates Δḡf(ai) with the same accuracy as the noiseless operator approximation.Source: Appendix C · p. 53