Source-linked AI summary
Variational approach for learning Markov processes from time series data
Hao Wu, Frank Noé
TL;DR
Extracting dynamical models and their main characteristics from time-series data is a recurring problem in science and engineering. This paper introduces VAMP, a variational approach for learning from time series, and reports that its framework leads to an optimal choice of hyperparameters.
Problem
Extracting dynamical models and their main characteristics from time-series data is a recurring problem in science and engineering.
Source: §1 · p. 2Method
VAMP uses Koopman-operator singular components to define data-calculable scores for optimizing feature mappings and Markovian models.
Source: §1 · p. 5, §2.1 · p. 7Results
The framework leads to an optimal choice of hyperparameters.
Source: §1 · p. 6Takeaways & Limitations
VAMP-r and VAMP-E scores support model optimization and cross-validation for hyperparameter tuning and model selection.
Source: §1 · p. 6, §6 · p. 34The authors identify a major challenge in applying VAMP to real-world applications but do not specify it in the supplied passage.
Source: §6 · p. 34
Abstract
from arXiv · showhide
Inference, prediction and control of complex dynamical systems from time series is important in many areas, including financial markets, power grid management, climate and weather modeling, or molecular dynamics. The analysis of such highly nonlinear dynamical systems is facilitated by the fact that we can often find a (generally nonlinear) transformation of the system coordinates to features in which the dynamics can be excellently approximated by a linear Markovian model. Moreover, the large number of system variables often change collectively on large time- and length-scales, facilitating a low-dimensional analysis in feature space. In this paper, we introduce a variational approach for Markov processes (VAMP) that allows us to find optimal feature mappings and optimal Markovian models of the dynamics from given time series data. The key insight is that the best linear model can be obtained from the top singular components of the Koopman operator. This leads to the definition of a family of score functions called VAMP-r which can be calculated from data, and can be employed to optimize a Markovian model. In addition, based on the relationship between the variational scores and approximation errors of Koopman operators, we propose a new VAMP-E score, which can be applied to cross-validation for hyper-parameter optimization and model selection in VAMP. VAMP is valid for both reversible and nonreversible processes and for stationary and non-stationary processes or realizations.
Mathematics Subject Classification (2010) 37M10 · 37L65 · 47N30 · · 1 Introduction
The section presents VAMP as a framework for choosing features, Markov models, and hyperparameters from time-series data, including for irreversible and non-stationary processes. It derives scores from Koopman-operator singular components to optimize models and support cross-validation.
1 Introduction
- Regression error cannot select the form or number of basis functions, whereas variational scores can; optimal choices for f and g remain unclear at fixed dimension or fixed data amount.The text identifies f(x)=g(x)=1 and K=1 as an uninformative model in this context.Source: §1 · pp. 4–5
- VAMP optimizes models for arbitrary Markov processes, including irreversible and non-stationary dynamics, using the Koopman operator’s singular-value decomposition.This addresses the limited usefulness of eigenvalue decompositions for time-irreversible and non-stationary processes.Source: §1 · p. 5
- The Koopman-operator approximation error is minimized when the model’s feature maps are its top left and right singular functions.The introduction uses this result to motivate variational scores for selecting features and model parameters.Source: §1 · p. 5
- VAMP-r scores measure similarity to the true singular functions, and maximizing them yields optimal model parameters, algorithmically equivalent to CCA on featurized time-lagged data.The scores are defined for r = 1, 2, … .Source: §1 · p. 6
- The CCA-based approach can learn feature transformations with nonlinear function approximators, including deep neural networks.The introduction also connects the VAMP-2 score to approximation error relative to the true Koopman operator.Source: §1 · p. 6
- The negative Koopman-model approximation error, computable up to a constant, defines the VAMP-E score for cross-validation.The introduction proposes optimizing VAMP-E in cross-validation.Source: §1 · p. 6
- The framework leads to an optimal choice of hyperparameters.Source: §1 · p. 6
2 Theory
This section develops an SVD-based variational theory for optimal finite-rank Koopman models, including a data-based score for estimating leading components and relationships between model error and transition-density error. It specifies the assumptions and reversible-process connections that delimit the results, with optimality applying at a fixed lag time.
2.1 Koopman analysis of dynamical systems and its singular value decomposition
- With dim(f), dim(g) ≤ k, the linear model (1) attaining the smallest modeling error in Hilbert-Schmidt norm is given by f=(ψ1,…,ψk)ᵀ, g=(φ1,…,φk)ᵀ, and K=diag(σ1,…,σk).The first singular component is (σ1, φ1, ψ1) = (1, 1, 1).Source: §2.1 · pp. 7–8
- The theory requires a Hilbert-Schmidt Koopman operator for its SVD, and the paper says its conclusions and methods do not apply to deterministic systems in usual cases.The theorem is otherwise described as universal for Markov processes.Source: §2.1 · p. 8
- The Hilbert-Schmidt error between Koopman operators equals a weighted L2 error of transition densities, giving the operator approximation error a transition-density interpretation.The illustrative system’s first four singular components yield a relative Koopman approximation error of 6.6%.Source: §2.1 · pp. 8–9
- In the illustrative system, the second left and right singular functions indicate metastable states, while the third and fourth provide finer dynamical detail.A rank-4 approximation also predicts the state-distribution evolution with small error.Source: §2.1 · p. 9
2.2 Variational principle for Markov processes
- The VAMP variational principle makes the k dominant Koopman singular components the maximizers of VAMP-r, providing a data-driven route to estimating them.The score can be evaluated on time-series transition pairs; for r = 1, its objective uses cross-time feature products.Source: §2.2 · pp. 11–12
2.3 Comparison with related analysis approaches
- For a time-reversible stationary Markov process with ρ0 = ρ1, Koopman SVD is equivalent to eigenvalue decomposition, generalizing reversible conformation-dynamics variational analysis; VAMP-1 maximizes Rayleigh trace, while VAMP-2 maximizes kinetic variance.Source: §2.3 · p. 12
- Singular functions can embed kinetic distances in low-dimensional coordinates, and dominant singular functions can identify coherent sets in nonstationary Markov processes [25].The passages describe coherent sets as generalizations of metastable states.Source: §2.3 · p. 12
- SVD-based analysis overcomes the variational and projection-error limitations for irreversible processes and yields optimal finite-rank approximate models.The cited passages contrast it with the lack of a variational principle for approximate eigenfunctions and the difficulty evaluating projection errors onto invariant subspaces.Source: §2.3 · p. 13
- Because singular functions depend on lag time, model (4) is optimal only for a fixed τ; multi-lag optimization and error analysis are left for future work.The passages state that the lag-time dependence is a major limitation of the SVD-based approach.Source: §2.3 · p. 13
3 Estimation algorithms
The section develops estimation algorithms for optimal dynamical models from time-series data using basis-function representations. It presents feature TCCA, which finds the best linear model in a given feature space.
3 Estimation algorithms
- The feature functions f and g are linear combinations of basis functions, with U and V combining basis sets to approximate k singular components.The formulation allows different basis sets for f and g.Source: §3 · p. 13
- Despite the linear Ansatz, the resulting feature functions can be strongly nonlinear in the system state x.The authors state that this does not restrict the generality of representable functions.Source: §3 · p. 14
- The estimation section addresses optimizing U and V, optimizing the feature functions, and assessing the resulting dynamical model.Source: §3 · p. 14
- Covariance and time-lagged covariance matrices of basis functions are estimated from trajectory data, with multiple trajectories handled by averaging.Source: §3 · p. 14
3.1 Feature TCCA: finding the best linear model in a given feature space
- Feature TCCA finds the best linear model in a given feature space by applying linear CCA to time-lagged basis data.The method solves for the optimal parameter matrices U and V when the basis functions are known.Source: §3.1 · pp. 14–15
2. Perform the truncated SVD
The section describes truncated-SVD feature TCCA, nonlinear optimization of its basis functions, and an error-based VAMP-E score. It also reports how these methods estimate Koopman singular components and relate to EDMD.
2. Perform the truncated SVD
- The truncated SVD selects k corresponding left and right singular vectors and outputs estimates of singular values and singular functions for the Koopman operator.The resulting model uses Kii, fi = u_i^T χ0, and gi = v_i^T χ1.Source: §2 · pp. 15–16
- With the numerical modifications, feature TCCA singular values are bounded by 1, and the first estimated singular component is exactly (K11, f1, g1) = (1, 1, 1), even with statistical noise.Source: §2 · p. 16
- Feature TCCA is a generalization of EDMD: with equal basis sets and full-rank covariance matrices, its output is equivalent to a least-squares regression solution.When χ0 = χ1, that regression is identical to the linear model of EDMD.Source: §2 · pp. 16–17
3.2 Nonlinear TCCA: optimizing the basis functions
- Nonlinear TCCA optimizes basis parameters w by maximizing fixed-w VAMP-r, then estimates components with feature TCCA; candidate bases include neural networks and decision trees.The optimization can avoid explicitly computing U and V; the estimates generally depend on r, and the authors suggest r = 2.Source: §3.2 · pp. 17–18
- Nonlinear TCCA gives more accurate singular-function and singular-value estimates than feature TCCA, although both underestimate dominant singular values as the variational principle predicts.The estimates are compared in Fig. 3.Source: §3.2 · p. 19
- Unlike EDMD-DL, nonlinear TCCA excludes uninformative zero-regression-error models without extra feature constraints and targets Koopman approximation error; EDMD-DL does not directly handle stochastic Koopman operators.The paper describes nonlinear TCCA as similar to EDMD-DL, which optimizes feature transformations by minimizing regression error.Source: §3.2 · p. 20
3.3 Error analysis
- The approximation error separates into an unknown constant Hilbert-Schmidt-norm term and a model-dependent term RE estimable from data; maximizing RE is equivalent to maximizing R2 in both TCCA methods.RE is named the VAMP-E score and can serve as an alternative to VAMP-r scores; the scores behave differently for hyper-parameter optimization.Source: §3.3 · pp. 20–21
4 Model validation
Section 4 addresses VAMP model validation and selection, including VAMP-E-based cross-validation for hyperparameters and rank and Chapman–Kolmogorov testing for lag times. It also states that an unspecified set of factors dominates the Markov-process dynamics.
4.1 Cross-validation for hyper-parameter optimization
- Cross-validation partitions trajectories into J folds, fits each candidate hyperparameter set on training data, scores estimated singular components on held-out data, and selects the set with maximal MCV(θ).Hyperparameters can include basis-function number and form or neural-network architecture and connectivity; models of different types may be compared.Source: §4.1 · pp. 22–23
- A test-set VAMP-r subspace score cannot directly use the training singular functions because they are generally not orthonormal on test data, and it fixes k and may be numerically unstable.Source: §4.1 · p. 23
- VAMP-E validation scores Koopman-operator approximation error on test data, permits selecting k as a hyperparameter, and avoids matrix inverses for stable computation.Feature TCCA and VAMP-2-based nonlinear TCCA maximize the training-set VAMP-E score.Source: §4.1 · pp. 23–24
- A smaller basis set (m = 13) poorly approximates singular functions, while m = 250 overfits and produces an estimated singular value larger than the true value.The paper attributes the latter to statistical noise.Source: §4.1 · pp. 24–25
- Cross-validation and exact VAMP-E scores both reach their maximum at m = 33, identifying the preferred basis-function count in Example 3.The training-set average score continues to rise with m; Fig. 4a shows the cross-validated scores and exact score.Source: §4.1 · p. 24, Figure 4 · p. 25
4.2 Chapman-Kolmogorov test for choice of lag times
- VAMP scores cannot compare models at different lag times, so the Chapman–Kolmogorov test assesses whether multi-lag observable predictions agree and favors τ below the timescale of interest.The paper applies the test to estimated leading singular functions and requires approximate agreement across multiple observables and lag times.Source: §4.2 · pp. 25–27
- An unspecified set of factors is said to dominate the Markov-process dynamics.The passage does not identify the factors.Source: §4.2 · p. 27
5 Numerical examples
The numerical examples show VAMP selecting a basis size that accurately captures double-gyre dynamics, while stochastic Lorenz singular functions reveal coarse states and support transformation invariance. The estimated models also predict long-time evolution and pass lag-time validation in the double-gyre example.
5.1 Double-gyre system
- Diffusion leaves the gyre’s two subdomains almost invariant, and singular components reveal these sets alongside rotational kinetics within them.Source: §5.1 · p. 27
- Five-fold VAMP-E cross-validation selected m=37 basis functions, producing accurate leading singular-component estimates; smaller m underfit and larger m showed statistical-noise effects.The figure caption reports that both the largest test-set MCV and the exact VAMP-E score were achieved at m=37.Source: §5.1 · p. 28, Figure 5 · p. 29
- The estimated Koopman operator predicts long-time state distributions, and Chapman–Kolmogorov tests support τ=2 as a suitable lag time.Source: §5.1 · p. 28, Figure 6 · p. 30, Figure 7 · p. 30
5.2 Stochastic Lorenz system
- Singular-function patterns coarse-grain stochastic Lorenz dynamics into four macrostates—the inner and outer basins of both attractor lobes—and ψ1’s sign boundary matches the almost-invariant-set boundary.Source: §5.2 · p. 31
- Projected leading singular components from a higher-dimensional observable are almost the same as those computed directly from three-dimensional data, illustrating VAMP’s transformation invariance.Source: §5.2 · p. 31, Figure 9 · p. 33
6 Conclusion … A.1 Definition of empirical distributions
The conclusion presents VAMP as a framework for evaluating and optimizing linear Markov models, while the appendix defines notation and states consistency results for empirical covariance estimates. It also identifies high-dimensional optimization and probability preservation as open challenges.
6 Conclusion
- VAMP unifies analysis of linearized coarse-grained Markov models, with VAMP-r and VAMP-E scores evaluating accuracy and data-driven methods optimizing finite-data, finite-dimension models.The methods include feature TCCA, nonlinear TCCA, and VAMP-E-based cross-validation.Source: §6 · p. 34
- For high-dimensional systems, the paper identifies the curse of dimensionality and effective variational optimization as challenges, and points to VAMPnet and tensor-decomposition approximations.VAMPnet approximates singular components with deep neural networks.Source: §6 · p. 34
- Operator-error optimization may yield models without nonnegative transition densities, leaving probability-preserving modeling open; applications to metastable states, coherent sets [18], and dominant cycles are future work.Source: §6 · p. 34
Appendix
- The appendix defines transition-density and scalar-product notation and uses N(·|c, σ²) for the normal density with mean c and variance σ².Source: Appendix · p. 35
A.1 Definition of empirical distributions
- For S independent trajectories, the estimates of C00, C11, and C01 in (16–18) are unbiased and consistent as S → ∞.The passage frames this result using the Markov propagator Pt defined in (63).Source: §A.1 · p. 35
- Similar conclusions are stated for trajectories with different initial conditions and different lengths.Source: §A.1 · p. 35
A.2 Proof of Theorem 1 … B Markov propagators
These appendices establish the best-rank-k Koopman approximation and describe how its errors and limitations relate to transition densities. They also state sufficient conditions for Theorem 1, a noncompactness result for deterministic systems, and a Markov-propagator formulation.
A.2 Proof of Theorem 1
- The model-derived approximate Koopman operator is the best rank-k approximation to Kτ in Hilbert-Schmidt norm, by the generalized Eckart-Young theorem.The proof also concludes that (σ1, φ1, ψ1) = (1, 1, 1) [14]; the approximation theorem is attributed to.Source: §A.2 · p. 36
A.3 Transition densities deduced from Koopman operators
- The operator error between the approximate and true Koopman operators can be represented by the error between their transition densities.Source: §A.3 · p. 37
- The approximate transition density is normalized but may be negative, so the resulting approximate operators and densities are not guaranteed to define valid probabilistic models.They can still be used for quantitative analysis of Markov processes.Source: §A.3 · p. 37
A.4 Sufficient conditions for Theorem 1
- The stated sufficient conditions for Theorem 1 include a finite state space or a state space Rd with ρ0(x) and ρ1(y) positive for all states.Source: §A.4 · p. 38
A.5 Koopman operators of deterministic systems
- For deterministic systems, the Koopman operator is not compact on L2 and therefore is not Hilbert-Schmidt.Source: §A.5 · p. 39
B Markov propagators
- The Markov-propagator formulation uses p_t(x)=P(x_t=x) as the density of x_t and applies the Koopman-operator SVD with stated normalizations.Source: B Markov propagators · p. 39
C Proof of the variational principle … E.2 Feature TCCA of projected Koopman operators
The appendices justify the variational principle and show that feature TCCA recovers optimal singular components of a projected Koopman operator, while reversible processes admit an eigenfunction variational principle.
C Proof of the variational principle
- The proof rewrites the optimization problem in singular coordinates and applies the Cauchy-Schwarz inequality to establish its variational bound.The supplied proof passages give the rewritten problem and inequality-based argument, but omit the intervening formulas.Source: C Proof of the variational principle · p. 40
- The bound is attained when the first k rows of D0 and D1 are identity matrices, proving the variational principle.The proof uses the constraints D0ᵀD0 = I and D1ᵀD1 = I.Source: C Proof of the variational principle · p. 41
D Variational principle of reversible Markov processes
- Under time reversibility and stationarity, the variational principle for r ≥ 1 is maximized by eigenfunctions corresponding to the largest eigenvalues.The reversible-process result follows because the Koopman operator’s eigendecomposition is equivalent to its SVD.Source: D Variational principle of reversible… · p. 41
E.1 Correctness of feature TCCA
- Feature TCCA solves optimization problem (19) under its orthonormality constraints.The appendix identifies this result as the correctness of the feature TCCA algorithm.Source: §E.1 · pp. 41–42
- The solution uses the first k left and right singular vectors of C under the stated orthonormality constraints.The passages describe the singular-vector choice as the equality case in the optimization argument.Source: §E.1 · p. 42
E.2 Feature TCCA of projected Koopman operators
- Ignoring statistical noise, feature TCCA returns the singular components of the projected Koopman operator Kproj.The algorithm is interpreted variationally for the Koopman operator projected onto the feature subspaces.Source: §E.2 · pp. 42–43
- The optimality of the estimation result is invariant for any choice of r ≥ 1.This invariance is stated for the feature TCCA estimation result in the projected-operator setting.Source: §E.2 · p. 43
E.3 An example of nonlinear TCCA … H.1 Proof of (28)
The appendices illustrate nonlinear TCCA score optimization, describe basis de-correlation and parameter optimization, and outline the stated VAMP–EDMD and VAMP-E proof setups.
E.3 An example of nonlinear TCCA
- Across 9,901 equally spaced w values in [0.01, 1], maximal R1 and R2 scores occur at w=0.3157 and w=0.7069, respectively.The score was evaluated for r = 1 and 2.Source: §E.3 · p. 44
F.1 De-correlation of basis functions
- PCA de-correlates basis functions by centering them, forming covariance matrices, and retaining positive eigenvalues above ϵ0 and absolute values of negative eigenvalues.The applications set ϵ0 = 10^-10.Source: §F.1 · p. 44
- The estimated leading singular value satisfies σ1 ≤ 1, and including constant basis functions makes 1 the largest singular value.Source: §F.1 · p. 46
F.2 Parameter optimization in nonlinear TCCA
- For high-dimensional w, nonlinear TCCA uses gradient descent; for large datasets it can use random-subset stochastic gradient descent.The examples also solve the optimization by direct search.Source: §F.2 · p. 46
- When C00 or C11 is rank deficient, nonlinear TCCA can be numerically singular; gradient methods regularize Rr(w) with ε selected by cross-validation.For direct search, de-correlation of basis functions addresses the issue.Source: §F.2 · pp. 46–47
G Relationship between VAMP and EDMD
- The VAMP–EDMD relationship is examined through the eigenvalue problem for K̂τ, assuming χ0 = χ1 = χ and invertible C00, consistent with EDMD spectral approximation theory.Source: G Relationship between VAMP and EDMD · p. 47
H.1 Proof of (28)
- The proof of (28) invokes the Hilbert-Schmidt inner product and the definition of the Hilbert-Schmidt norm.Source: §H.1 · p. 49
H.2 Relationship between VAMP-2 and VAMP-E … K.2 Double-gyre system
The appendices relate VAMP-E optimization to TCCA, extend the variational principle to subspaces, and provide computation and numerical-example details.
H.2 Relationship between VAMP-2 and VAMP-E
- Feature TCCA maximizes the first term in (116), and nonlinear TCCA also maximizes VAMP-E because its first term is zero while its second is maximized over w.The passages state that maximizing the first term yields VAMP-E for feature TCCA; for nonlinear TCCA, the first term vanishes at the optimal model.Source: §H.2 · pp. 49–50
I Subspace variational principle
- The subspace score Rspace_r relaxes VAMP-r by measuring consistency between feature subspaces and dominant singular spaces, with equality when the subspaces match.The equality condition specifies matching spans for the feature functions and the corresponding singular functions.Source: I Subspace variational principle · p. 50
- Rspace_r is invariant to invertible linear transformations of f and g, but at maximal k it is independent of estimated singular components, leaving k undetermined by cross-validation.The passage also notes that V_k may be singular, making the validation score unreliable.Source: I Subspace variational principle · pp. 50–51
J Computation of ˆKn
- Substituting f = U⊤χ0 and g = V⊤χ1 into (123) yields (37), as stated in the cited passage.Source: J Computation of ˆKn · p. 51
K.1 One-dimensional system
- The one-dimensional system is discretized into 2,000 uniform bins over [−20, 20], and the appendix describes computing its stationary distribution and singular components.For Examples 2 and 3, nonlinear TCCA’s smoothing parameter w is optimized with a golden-section search.Source: §K.1 · pp. 52–53
K.2 Double-gyre system
- The double-gyre example uses Euler–Maruyama with step size Δ = 0.02 and a uniform partition into 1,250 bins across [0, 2] × [0, 1].The appendix says Fig. 6’s transition density at lag time nτ is computed from estimated singular components.Source: §K.2 · p. 54