{'fraction_non_alphanumeric': 0.040191047202886035, 'fraction_numerical': 0.015141507037244043, 'mean_word_length': 4.537985368598762, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Neural Dynamic Mode Decomposition for End-to-End Modeling of Nonlinear Dynamics

Koopman spectral analysis has attracted attention for understanding nonlinear dynamical systems by which we can analyze nonlinear dynamics with a linear regime by lifting observations using a nonlinear function. For analysis, we need to find an appropriate lift function. Although several methods have been proposed for estimating a lift function based on neural networks, the existing methods train neural networks without spectral analysis. In this paper, we propose neural dynamic mode decomposition, in which neural networks are trained such that the forecast error is minimized when the dynamics is modeled based on spectral decomposition in the lifted space. With our proposed method, the forecast error is backpropagated through the neural networks and the spectral decomposition, enabling end-to-end learning of Koopman spectral analysis. When information is available on the frequencies or the growth rates of the dynamics, the proposed method can exploit it as regularizers for training. We also propose an extension of our approach when observations are influenced by exogenous control time-series. Our experiments demonstrate the effectiveness of our proposed method in terms of eigenvalue estimation and forecast performance.


Introduction
The data-driven analysis of nonlinear dynamics is important in a wide variety of fields, such as sociology [10], epidemiology [24], neuroscience [4], physics [3], finance [26], biology [6], and engineering [22]. Recently, an approach based on the Koopman operator [18,27] has garnered attention. With this approach, a nonlinear dynamical system is lifted to the corresponding linear one in a possibly infinite-dimensional space by a nonlinear function, by which various methods that were developed for analyzing and controlling linear dynamical systems can be used for nonlinear ones. For example, such properties of nonlinear dynamical systems as frequencies and growth rates can be identified with the eigenvalues of a Koopman operator in the lifted space.
One representative data-driven method based on the Koopman operator is dynamic mode decomposition (DMD) [33,32,20], which analyzes the properties of nonlinear dynamics based on the eigenvalues and eigenvectors of a low-rank approximation of a state transition matrix. However, since DMD assumes that observations are obtained in a lifted space, an appropriate lifted space must be defined. To automatically find a lifted space, several methods have been proposed that use neural networks for modeling lift functions [35,25,42,21,29]. They train neural network-based autoencoders and linear dynamics in the lifted space, where neither a low-rank approximation nor eigen decomposition is used. After the autoencoders are trained, the dynamics is analyzed by the eigen decomposition of the estimated state transition matrix as post processing. Therefore, eigenvalues and eigenvectors, which contain important information for characterizing the dynamics, are not directly used for training the autoencoders.
In this paper, we propose neural dynamic mode decomposition (NDMD), which trains autoencoders such that the forecast error is minimized when the dynamics is analyzed by DMD in lifted space. With NDMD, a forecast is made based on the following procedure. First, observations are encoded by a neural network-based encoder in a lifted space, Section 4.2, we present NDMD's training procedure. In Section 4.3, we propose an extension of NDMD when the auxiliary information on the eigenvalues of dynamics is available. In Section 4.4, we propose another extension of NDMD for learning a lift function when observations are influenced by exogenous control time-series.

Neural dynamic mode decomposition
Let x t ∈ R M be an observation vector at timestep t, which is obtained by a set of observablesh := [h 1 , . . . , h M ] as x t =h(w t ). Since the space spanned byh is not necessarily large as mentioned in the previous section, we consider a K-dimensional lifted space, where observation vectors are mapped onto the lifted space by a nonlinear function, where ψ t ∈ R K is the encoded vector at timestep t, and f is modeled by a neural network called an encoder. That is, we expect that the space spanned by composition f •h is sufficiently large (in the sense of [36]) for the underlying dynamics behind the data. Thus, in the lifted space, encoded vectors are assumed to follow a discrete linear dynamics: where A ∈ R K×K is a linear transition matrix. With the proposed method, infinite-dimensional Koopman operator K is approximated by low-dimensional linear transition matrix A.
Our model obtains the linear dynamics by applying DMD on the encoded vectors. We define two matrices, Ψ 1 and Ψ 2 , which are constructed from the encoded vectors: where Ψ 1 consists of the encoded vectors at timesteps τ = {τ 1 , · · · , τ S }, and Ψ 2 consists of the encoded vectors at their next timesteps. Timesteps τ can be non-consecutive, and unordered.
The transition matrix can be estimated by the least squares method: where † denotes the pseudo-inverse. Our model considers the low-rank approximation of the transition matrix, A = U AU, as it did with DMD. Here, U ∈ R N ×R is obtained by the singular value decomposition (SVD) of Ψ 1 : where denotes the transpose, Σ ∈ R R×R , V ∈ R S×R , and R is the reduced rank. Then an estimate of the low-ranked transition matrix is given:Ã Using the low-rank approximation, we can reduce the noise in the encoded vectors for modeling the dynamics.
Let the columns of Y ∈ C R×R be the eigenvectors ofÃ, where Λ ∈ C R×R is a diagonal matrix containing corresponding eigenvalues λ k :Ã The encoded vector at timestep t is forecasted using the eigenvalues and the eigenvectors: where is the DMD modes, and α = Φ † ψ τ1 ∈ C R . Using eigen decomposition, we can efficiently obtain the estimation at any timestep since the power of diagonal matrix Λ t−1 is calculated by the power of its diagonal elements. In addition, when auxiliary information on the eigenvalues and/or eigenvectors is available, we can use it as a regularizer for training our model as described in Section 4.3. In Eq. (10), the encoded vectors are decomposed into R different dynamics, where the dynamics of the rth component is defined by eigenvalue λ r . Randomly sample S timesteps τ = {τ 1 , · · · , τ S } from {1, · · · , T − 1}.

11:
Calculate loss L = 1 |t| t∈t x t −x t 2 , and its gradient. 12: Update parameters Θ using the loss and the gradient. 13: end while The observation vector at timestep t is forecasted by mapping forecasted encoded vectorψ t into the observation space: where g is the decoder modeled by a neural network.
The number of eigenvalues estimated by our model is at most K, which is the dimensionality of the lifted space. We can easily limit the number of eigenvalues for obtaining sparse but essential dynamics by tuning the number of output units of the encoder, as well as the number of input units of the decoder, even if the dimensionality of the observation space is high.

Training
Suppose that we are given time-series data of length T , [x 1 , · · · , x T ]. The parameters in the encoder and decoder neural networks, Θ, are estimated by minimizing the expected forecast error of our model in the observation space: where E denotes the expectation. In particular, we minimize the objective function using a stochastic gradient descent based method with batch size S, as shown in Algorithm 1.

Neural DMD with auxiliary information
The real part of the eigenvalues represents the growth and decay of the dynamics, which diverges when real(λ k ) > 0 and converges when real(λ k ) < 0. The imaginary part of the eigenvalues represents the frequency of the dynamics, . When auxiliary information on the dynamics is available, we can use it for training our model by adding regularizer η in the objective function in Eq. (13): where β > 0 is a hyperparameter, η is a regularization function, and E is auxiliary information. For example, suppose we are given true eigenvalues E = λ * = {λ * r } R r=1 as auxiliary information. The regularizer can be modeled such that estimated eigenvalues Λ in Eq. (9) in our model become similar to the true ones: where the first term is the absolute distance between true eigenvalue λ * r and estimation λ r , which is closest to the true one, and the second term is the absolute distance between estimated eigenvalue λ r and true λ r , which is closest to the estimated one. We can incorporate other types of auxiliary information as regularizers. For example, when only frequency information is available, a regularizer on the imaginary parts of the eigenvalues can be used. When there is a limit cycle dynamics, we can use a regularizer that turns the real part of an eigenvalue to zero.

Neural DMD with control
We propose an extension of neural DMD for analyzing dynamics using time-series data with control, where we observe the exogenous variables that affect the dynamics.
Let z t ∈ R D be the observed exogenous control vector at timestep t. We encode z t by nonlinear function c: where ξ t ∈ R N is the encoded exogenous vector at timestep t, and c is modeled by a neural network. Our model assumes the following linear dynamics in the lifted space: where the encoded observation vector at the next timestep is determined by the encoded observation vector at the previous timestep and the encoded exogenous vector.
Our model obtains linear dynamics in the lifted space based on DMD with control [30]. We define the following three matrices, Ψ 1 , Ψ 2 , and Ξ 1 , from the encoded observation and exogenous vectors: Note that Ψ 1 , Ψ 2 , and Ξ 1 should be consecutive, unlike Eq. (5). Let Ω be a concatenation of Ψ 1 and Ξ 1 : We take the SVD of Ω: where U ∈ R (K+N )×P ,Σ ∈ R P ×P ,V ∈ R P ×(T −1) , and P is the reduced rank. We obtain estimates of A and B by least squares using the SVD result:Â where LetÛ ∈ R K×R ,Σ ∈ R R×R , andV ∈ R (T −1)×R be the SVD of Ψ 2 with reduced rank R: Then the low-rank approximation of A is given: The dynamics of Eq. (17) is characterized by the eigenvalues and eigenvectors ofÃ: Algorithm 2 Training procedures of the proposed method with control. Input: Observation time-series data x 1 , · · · , x T , exogenous time-series data z 1 , · · · , z T , batchsize S Output: Estimated neural network parametersΘ 1: Initialize parameters Θ. 2: while End condition is satisfied do 3: Randomly sample timestep τ ∈ {1, · · · , T − S}.

9:
Forecast encoded observation vectorsψ t by Eq. (26) for t ∈ t. x s −x s 2 , and its gradient. 12: Update parameters Θ using the loss and the gradient. 13: end while With Eq. (17), the encoded observation vector at timestep t is given: Since A is approximated by A ≈ ΦΛΦ † using the eigen decomposition result, we can forecast the encoded observation vector at timestep t:ψ where and α = Φ † ψ 1 . The observation vector at timestep t is forecasted byx t = g(ψ t ) by decoder g in the same way as with Eq. (12). We estimate the parameters of the encoder and decoder neural networks by minimizing the expected error in Eq. (13). Algorithm 2 shows the estimation procedure using a stochastic gradient descent-based method.

Eigenvalue estimation
We evaluated the proposed NDMD described in Section 4.1 using synthetic data, which were generated by the following procedure. First, the dynamics in a two-dimensional lifted space was obtained with linear transition matrix, A = [0.9, −0.5; 0.4, 0.9]. Then the dynamics in a ten-dimensional observation space was calculated by a nonlinear transformation using Gaussian processes (GPs) [31] with the RBF kernel. The time-series in the observation space are shown in Figure 2(a, top). The first 70 timesteps of the data were used for training, and the following 10 were used for validation.
We compared the eigenvalue estimation performance of the proposed method (NDMD) with DMD [36], a DMD with rank-two approximation (DMD2), and a sparse DMD [14] (SDMD). With the proposed method, we used four-layered feed-forward neural networks with 256 hidden units for the encoder and decoder. The dimensionality of the lifted space was two. We optimized the parameters using Adam [16] with a learning rate of 10 −3 , a dropout rate of 0.1 [34], and a batch size of 128. We implemented the proposed method with Tensorflow [1]. The validation data were used for early stopping, and the maximum number of epochs was 1,000. With DMD2, we selected two eigenvalues based on the singular values. With SDMD, we tuned the best hyperparameter from {10 −3 , 10 −2 , 10 −1 , 1, 10 1 , 10 2 , 10 3 }.

Eigenvalue estimation with auxiliary information
We evaluated the proposed NDMD with auxiliary information described in Section 4.3 using synthetic data. The dynamics in a four-dimensional lifted space was obtained with linear transition matrix A = [0.9, −0.5, 0, 0; 0.4, 0.9, 0, 0; 0, 0, 0.8, −0.5; 0, 0, 0.6, 0.8], and the dynamics in a ten-dimensional observation space was calculated using GPs with an RBF kernel. The observation time-series data are shown in Figure 2(b, top). The first 70 timesteps of the data were used for training, the next 10 were used for validation, and the remaining 20 were used for testing.
We calculated the four true eigenvalues from A, and used them for auxiliary information. Figure 2(b, bottom) shows the result. The NDMD without auxiliary information (NDMD) failed to estimate two eigenvalues in this case. On the other hand, the NDMD with auxiliary information (NDMDr) appropriately estimated all four eigenvalues. The mean squared error for forecasting the test data was 0.146 with NDMD and 0.074 with NDMDr. This result indicates that the proposed method can improve the forecasting performance using auxiliary information.

Eigenvalue estimation with control
We evaluated the proposed NDMD with control described in Section 4.4 using synthetic data. The data were generated by the following procedure. First, a one-dimensional exogenous time-series was obtained by a standard Gaussian distribution. Second, the dynamics in a two-dimensional lifted space was obtained with linear transition matrices A = [0.9, −0.5; 0.4, 0.9], and B = [1; 0]. Finally, we calculated the dynamics in a ten-dimensional observation space by nonlinear transformation using GPs with the RBF kernel. The observation time-series data are shown in Figure 2(c, top). The first 140 timesteps of the data were used for training, and the next 20 were used for validation.   Figure 2(c,bottom) shows the estimated eigenvalues. The neural DMD with control (NDMDc) successfully estimated the true eigenvalues. The NDMD without control (NDMD) failed to estimate them since it ignores the influence of the exogenous time-series data. The DMD with control (DMDc) estimated ten eigenvalues, only two of thewhich were close to the true ones. The estimated eigenvalues by the DMD with control with rank-two approximation (DMDc2) were different from the true ones.

Fluid flow forecast
We evaluated our proposed NDMD with fluid flow forecasts. We used the cylinder wake data shown in Figure 3. The original data consisted of the vorticity in 199×499 field for 151 timesteps. We randomly subsampled {0.01%, 0.03%, 0.1%} of the field for the observation time-series data, and forecasted them. We used 70% of the timesteps for training, 10% for validation, and the remainder for testing.
We compared NDMD with the following nine methods: DMD, EDMD, KDMD(RBF), KDMD(Poly), AEAR, LKIS, AR, NN, and LSTM. With NDMD, we used four-layered neural networks with 256 hidden units and the 256 dimensional lifted space, for the encoder and decoder. The rank of the linear transition in the lifted space was determined by removing the elements that had relative singular values lower than 0.1%. DMD is the standard DMD without neural networks. EDMD is the extended DMD [39], where observation vectors were augmented by the product of two elements in the original observation vector. KDMD is the kernel DMD [40], where we used the RBF kernel and the polynomial kernel with a degree of three. AEAR is a combination of neural network-based autoencoders and linear autoregressive models [25]. It forecasted the value at the next timestep using an encoder, a linear transition matrix, and a decoder. The observation vectors were first encoded by the encoder. Next, the encoded vectors were transformed to those at the next timestep by the linear transition matrix. Then, the transformed encoded vectors were decoded in the observation space. The neural networks were trained by minimizing the mean squared error for the next timestep forecast. LKIS is learning Koopman invariant subspaces [35], which also consisted of an encoder, a linear transition matrix, and a decoder. LKIS was trained by minimizing the sum of the reconstruction error in the observation space, and the next timestep forecasting error in the lifted space. AR is the linear autoregressive model. NN is a neural network that forecasts the next timestep observation vectors. We used four-layered neural networks with 256 hidden units. LSTM is a long short-term memory recurrent neural network [12], where 256 hidden units were used. AR, NN, and LSTM were trained by minimizing the mean squared error for the next timestep forecast.
The test mean squared error averaged over ten experiments with different random subsamplings is shown in Table 1.
The proposed NDMD achieved the lowest error. As the number of observation fields increased, the error decreased. The DMD failed to produce a forecast since it could not extract the nonlinear dynamics. Although EDMD and KDMD were outperformed DMD by incorporating higher order information, they were worse than NDMD since the predefined nonlinear functions were not appropriate for the dynamics. On the other hand, the proposed method learned a proper lifted space from the data by neural networks. The errors by existing neural network-based methods, i.e., AEAR, LKIS, NN, and LSTM, were higher than NDMD. The reason is that NDMD was trained based on long-term forecasting performance using the eigenvalues, and the noise was removed by SVD. The forecasting values by AR diverged, and the error was very high, especially when there were few observations.

Conclusion
We proposed a neural network-based method for analyzing nonlinear dynamics by the Koopman operator. With the proposed method, we estimated a lifted space with a linear dynamics such that the forecast error is minimized when the dynamics is modeled based on spectral decomposition in the lifted space. We demonstrated the effectiveness of the proposed method for estimating the eigenvalues of the nonlinear dynamics and forecasting future values compared with existing DMD-based and neural network-based methods.
Although our results are encouraging, we must extend our approach in several directions, First, we will apply it for controlling the dynamics. Since the proposed method can estimate lifted linear dynamics, we can control the nonlinear dynamics using techniques developed in optimum control theory for linear dynamics [17,11,28,23], such as linear quadratic control [13] and pole placement [5]. Second, we want to use other types of neural networks, which include convolutional neural networks [19] and transformers [38], for extracting rich information from data. Third, we will use our proposed end-to-end learning approach for supervised learning tasks, such as time-series regression and classification [8]. Finally, we plan to evaluate our method with various types of auxiliary information, e.g., a regularizer on dynamic modes.