On the difficulty of learning chaotic dynamics with RNNs
Authors: Jonas Magdy Mikhaeil
Zahra Monfared
Daniel Durstewitz
Venue: NeurIPS 2022 Accept (2022)
Keywords: Recurrent neural networks
Dynamical systems
Attractors
Time series analysis
Chaos
Exploding and vanishing gradient problem
Teacher forcing

Abstract: Recurrent neural networks (RNNs) are wide-spread machine learning tools for modeling sequential and time series data. They are notoriously hard to train because their loss gradients backpropagated in time tend to saturate or diverge during training. This is known as the exploding and vanishing gradient problem. Previous solutions to this issue either built on rather complicated, purpose-engineered architectures with gated memory buffers, or - more recently - imposed constraints that ensure convergence to a fixed point or restrict (the eigenspectrum of) the recurrence matrix. Such constraints, however, convey severe limitations on the expressivity of the RNN. Essential intrinsic dynamics such as multistability or chaos are disabled. This is inherently at disaccord with the chaotic nature of many, if not most, time series encountered in nature and society. It is particularly problematic in scientific applications where one aims to reconstruct the underlying dynamical system. 
Here we offer a comprehensive theoretical treatment of this problem by relating the loss gradients during RNN training to the Lyapunov spectrum of RNN-generated orbits. We mathematically prove that RNNs producing stable equilibrium or cyclic behavior have bounded gradients, whereas the gradients of RNNs with chaotic dynamics always diverge. 
Based on these analyses and insights we suggest ways of how to optimize the training process on chaotic data according to the system's Lyapunov spectrum, regardless of the employed RNN architecture. 

On the difficulty of learning chaotic dynamics with RNNs
Jonas M. Mikhaeil Zahra Monfared 1.4" and Daniel Durstewitz!25
j. mi mikhaiil Sataaeel@@co www.oolllaa....ommmmmoom edu, {zahra. ( monfared, daniel.. www.duuitttttteeeo 0zi mrizii-mnnheimi de
Department Theoretical Neuroscience, Central Institute Mental Health, Medical Faculty Mannheim.
Heidelberg University, Mannheim, Germany
Faculty of Physics and Astronomy, Heidelberg University, Heidelberg, Germany
IInterdisciplinaryy Center for Scientific Computing, Heidelberg University
'Department Mathematics & Informatics and Cluster of Excellence STRUCTURES, Heidelberg University,
Hcidelberg, Germany
"These authors contributed equally
Abstract
Recurrent neural networks (RNNs) are wide-spread machine learning tools for modeling sequential and time series data. They notoriously hard to train be- cause their loss gradients backpropagated in time tend to saturate or diverge during training. This is known the exploding and vanishing gradient problem. Previous solutions to this issue either built on rather complicated, purpose-engineere archi- tectures with gated memory buffers, or more recently imposed constraints that ensure convergence to fixed point or restrict (the eigenspectrum of) the recurrence matrix. Such constraints, however, convey severe limitations on the expressivity of the RNN. Essential intrinsic dynamics such multistability or chaos are disabled. This inherently at disaccord with the chaotic nature many, not most, time series encountered in nature and society. It is particularly problematic in scientific applications where one aims to reconstruct the underlying dynamical system. Here offer comprehensive theoretical treatment of this problem by relating the loss gradients during RNN training to the Lyapunov spectrum of RNN- -NN-greerated orbits. We mathematically prove that RNNs producing stable equilibrium cyclic behavior have bounded gradients, whereas the gradients of RNNs with chaotic dynamics always diverge. Based these analyses and insights suggest ways how to optimize the training process on chaotic data according to the system's 's Lyapunov spectrum, regardless of the employed RNN architecture.
Introduction
Recurrent neural networks (RNNs) are widely used across various fields in engineering and science for learning sequential tasks or modeling and predicting time series [55]. Yet, they struggle when long-term temporal dependencies, very slow, or hugely varying time scales are involved 6 34 53]
73
80
Time series sequential data with such properties are, however, very common in fields
like climate physics [84], neuroscience [22] 75), ecology [87], or language processing [T]. Training
RNNs on such data is hard because the loss gradients backpropagated in time easily saturate or diverge in this process. This commonly referred to as the exploding and vanishing gradient problem
(EVGP)
34 66].
One solution the EVGP based on specifically designed RNN architectures with gating mecha- nisms, such as long short-term memory (LSTM) 35 or gated recurrent units (GRU) [10]. These architectures allow states at earlier time steps to more easily influence activity much later through
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

kind of protected memory buffer, thus alleviating the EVGP by structural design. In practice, such models need to be backed up by further techniques like gradient clipping to keep the gradients in check [67]. 'he relatively complex architectural design of these networks impedes their mathe-
matical analysis and requires reverse engineering after training [57 63
80]. Partly to forego
these complications, a variety of other solutions has been proposed recently, imposing restrictions on the recurrence matrix to bound the gradients [4] 9). or enforcing global stability by design or regularization
19 49]. Often these procedures dramatically curtail the expressivity of the RNN 45 66 [8]; in particular, they rule out chaotic dynamics (for reasons discussed further below).
This is at odds with the plethora of chaotic phenomena in nature, engineering, and society. Chaotic dynamics are commonplace, almost default any complex physical or biological system. This
includes scientific areas as diverse as neuroscience [17 93 physiology [46], geophysics (8T), climate
systems [88), astrophysics (52), ecology [14). chemical reactions (2). cell (65] or population (59)
biology. Chaotic phenomena are also crucial for the understanding of societal and epidemiological
processes, such as the spread of diseases [58] or in economics [20]. They are further relevant in
purely technical contexts such as electrical engineering [1 [3 or laser optics
+ They have even been suggested to play an up to now largely neglected, but potentially very significant role in speech recognition
77 and natural language processing [66.. Hence, in almost any practical setting, chaotic phenomena abound. They cannot, in general, be ignored when devising RNN training algorithms.
Here we offer comprehensive theoretical treatment of the relation between RNN dynamics and the behavior of the loss gradients during training. We find a close connection between an RNN's loss gradients and the largest Lyapunov exponent of its freely generated orbits. We mathematically prove that RNNs producing stable fixed point or cyclic behavior have bounded gradients. Crucially, however, the loss gradients of RNNs producing chaotic dynamics always diverge. Hence, the chaotic nature of many time series data induces principle problem, and, despite significant efforts in the past to solve the EVGP, training RNNs on such data remains an open issue. We illustrate the implications of our theory for RNN training on several simulated and empirical chaotic time series, and adapt the idea of sparsely forced Back- Pack-Propagation Through Time (BPTT) as simple yet effective remedy that enables to learn the underlying dynamics despite exploding gradients.
2 Related works
Exploding and vanishing gradients. While 'classical' remedies of the EVGP 6 purpose-tailored architectures with gating mechanisms, which safeguard information flow across longer temporal distances (TO 35), the focus has recently shifted to simpler RNNs that address the EVGP by restricting the recurrence matrix to be orthogonal 31
67 rest on
32 37 unitary 44, antisymmetric 9.. or by ensuring globally stable fixed point solutions [8]
for example through co-trained Lyapunov functions [49]. However, all these approaches impose strong limitations on the dynamical repertoire of the RNN, enforcing global convergence to fixed points or simple cycles [ In doing so, they drastically reduce the expressiveness of these models 45
66). To address this problem, Erichson et al. [19] somewhat relaxed the constraints on the recurrence matrix by introducing a skew- symmetric decomposition combined with Lipschitz condition on the activation function. Another recent approach discretizes oscillator ordinary differential equations (ODEs) to arrive at a stable system of coupled [73] or independent (74 oscillators which increase the RNN's expressiveness while bounding its gradients. By design (and as acknowledged by the authors), neither of these architectures is capable of producing chaotic dynamics, however, as the underlying ODEs do not allow for exponential divergence of close-by trajectories (a prerequisite for chaos). Given these often principle limitations of parametrically or dynamically strongly constrained models, a fruitful direction may be to modify the training process itself, e.g. through modified or auxiliary loss functions (80 85). or special procedures for parameter updating (39 or loss truncation [6]] (95.. Our empirical evaluation will follow up on such ideas, but also highlight that simple loss truncation, windowing, or architectural solutions like LSTMs are not sufficient.
Learning dynamical systems. Surprisingly disconnected from the work on the EVGP and learning long-term dependencies, huge and long-standing literature deals with training RNNs on nonlinear dynamical systems (DS) 770 86 89]. including chaotic systems like the famous Lorenz equations (56) or chaotic turbulence in fluid dynamics
Of those, methods based on reservoir computing 69 are special in that they start with a large complex dynamics-enabling repertoire to begin with
We make this point more formal in Appx. A.1.6
2

---

for which a linear mapping onto the observations in feedback loop with the reservoir is learned (see [8] for issues associated with this strategy in the context of DS reconstruction). Teacher forcing see also (8) is one of the earliest techniques introduced to keep RNN trajectories on track while training. The idea behind TF is to simply replace RNN states by observations when
(TF: 13 70 95,,
available, thereby also effectively cutting off the gradients. TFF essentially derives from ideas in dynamical control theory, and adaptive schemes that increasingly hand over control to the RNN throughout training have been devised ( 2 S A related technique with applications in control theory "multiple shooting" [92] Here the whole observed time series is chopped into chunks, and for each chunk of trajectory new initial condition is estimated. Explicit constraints ensure continuity between the separate trajectory bits during optimization. State space models and the Expectation- Maximization algorithm became popular particularly in the 90es for uncovering the latent dynamics underlying a set of time series observations (23). and remain an important tool until today [16 [5]]
Most recently, approaches based variational inference and the reparameterization trick [48] like
sequential variational autoencoders (SVAE), gained in popularity for DS approximation
33 51
"Deterministic" RNNs (i.e., with latent states not treated random variables), like conventional LSTMs [90]. remain top choices for DS reconstruction, however.
Although connections between DS ideas and loss gradients have been drawn early on [6), so far only particular scenarios (like fixed point attractors) have been considered. Closest to our work is recent work by Schmidt et al. [80), where non-divergence of loss gradients is established when RNNs converge to fixed points or cycles. However, this was done only for the particular class of piecewise-linear RNNs (PLRNNs), more restrictive conditions for cycles were imposed than assumed here, and most importantly the chaotic case on which we focus here was not considered. Recent studies [18 9T] also point out the general connections between Lyapunov exponents and loss gradients that we develop in sect. 5.T. but do not provide any in-depth theoretical treatment, proofs, or empirical evaluation of methods alleviate exploding gradients in training, as we do here. Thus, systematic theoretical framework that relates RNN dynamics more generally, and across a range of different RNN architectures, to the behavior of its training gradients, still lacking so far.
3 Theoretical analysis: Relation between RNN dynamics and loss gradients
In our analysis, we will cover all major types of asymptotic dynamics (fixed points, cycles, chaos, and qurii--priodiccity). and mathematically investigate their implications for the loss gradients. We will do this for all major classes of RNNs, including standard RNNs with largely arbitrary activation function, LSTMs, GRUs, and PLRNNs. The next section will first develop and illustrate the basic intuition behind the relations between RNN dynamics and loss gradients.
3.1 Preliminaries: RNN dynamics and loss gradients
Formally, all popular RNN architectures, including LSTMs, GRUs, or PLRNNs, are discrete time DS, defined by Markkvvvaa recursive prescription for the temporal evolution of the latent states 2, E M of the general form
(1) where E R is the input at time t and 0 are RNN parameters. Map Fş may be instantiated by any the common RNN architectures: For instance, for standard RNNs we have F0 s:) fW 2;-1 Bs, + h), where is an element- wise activation function like tanh or a rectified linear unit (ReLU), W a connection matrix, matrix B weighs the inputs, and h is the usual bias term
2, = Fo(z1-1, 8t),
(see sect. A.1.3 A. for the definition of other RNN models explored here).
Assuming we start at some initial value 2 E MM and given a sequence of external inputs = (8t). can recursively rewrite eq.
ZT = (F), =:
(s.).
In DS theory, we characterize the long-term behavior of such sequences by its spectrum of Lyapunov exponents. The Lyapunov exponents estimate the exponential growth rates in different local directions of the system's 's state space, and the largest Lyapunov exponent gives the dominant exponential behavior. Let us denote the system's ' Jacobian at time by
(Z
Oz
Jr θz.-1
(3)
dz,-1
3

---

Then, the maximum Lyapunov exponent along an RNN trajectory {21, 22, ,2,,"" "r) is defined as
T-2 log
Xmaz = lim T-000
1 J
(4)
-=0
where || || denotes the spectral norm (or any subordinate norm) of a matrix. Amar < 0 nearby trajectories will ultimately converge to a fixed point or cycle, while for Xmaz > 0 (a necessary condition for chaos) initially nearby trajectories will exponentially separate, i.e. we will have divergence along one (or more) directions in state space. This accounts for the sensitive dependence on initial conditions in chaotic systems.
Now let ((0) be some loss function employed for RNN training that decomposes in time as L = E Lt- Suppose fancy BPTT our training algorithm (similar derivations could be performed for Real Time Recurrent Learning [RTRL]) we recursively develop the loss gradients w.r.t some RNN parameter in time (i.e., across layers of the RNN unrolled in time)
OC
T
OL.
Oh, Oz. 0 2r
with
(5)
00
00
00
Oz. dz, 00
t=1
and
azi
dz
0zi-1
0z,+1
Oz.
0z,-1
0z,-2
Ozr
t-リー1
Ōz1-k
t-p-1
II
J
(6)
Ōz,
k=0
k=0
where a+ denotes the immediate derivative. Now observe that the behavior of the loss gradients
crucially depends on the product series of Jacobians in eqn. 6 If the maximum absolute eigenvalues
1/T
of the Jacobians J, will, in the geometric mean, be larger than 1 (i.e.,
> 1).
1/T
gradients will explode as T + 0o, while they will saturate if II Jr
< 1. Thus, the key
point to note is that the same terms that occur in the definition of the Lyapunov spectrum, eqn. + resurface in the loss gradients, eqn. 55 & (6. This accounts for the tight links between system dynamics and gradients.
3.2 Fixed points and eyclic dynamics
Let us start by considering the simplest types of limit dynamics that can occur in RNNs (or any discrete-time DS): fixed points and cycles. In fact, by far most of the literature on global stability in RNNs and on loss gradients focused on just fixed points 9 19 49]. with only few authors who recently started to also connect eyclic behavior to loss gradients 73 80]. Recall that fixed point of
a recursive map 2e = F(2(11) is defined as a point z" for which we have x" = F(x*)) ) Likewise. a k-cycle (k > 1) is a set of temporally consecutive periodic points Pe := {zr,, Ztg: ................................ Zt,
(Ze,) ((z.,)..... 1 (2,,), }} that we obtain from recursive application of the map such that each of the cyclic points z,, P& is fixed point of the k times iterated map F* (with i being the smallest positive integer for which this holds). To simplify the subsequent treatment, we will collectively refer to fixed points and cycles as k-cycles (k: 2 1). Further recall that a fixed point or k-cycle is called stable if the maximum absolute eigenvalue of the Jacobian evaluated at that point is smaller than 1, neutrally stable if exactly 1, and unstable otherwise. Although the results we develop in this and the following sections will hold more widely, we will restrict our attention to recursive maps Fg from the class of RNNs R- (standardRNN, LSTM, GRU, PLRNN} ) (see Appx. A.1 for details).
Based on the observations made in the previous sections we can state the following theorem that links RNN dynamics and loss gradients:
Theorem 1. Consider an RNN Fe E R parameterized by 0, and assume that it converges to a stable fixed point or k-cycle I'k (k: 2 1) with Br, as its basin of attraction. Then for every 21 Br, (i)
the Jacobian Øzy exponentially vanishes as T oc; for I' the tangent vectors Gzz. and thus
dzi
From here on will suppress the explicit dependence external inputs s; notation-wis see Remark 2

---

the gradient of the loss function, and (ii) for the PLRNN 27 both T + OC.
0Lr, will be bounded from above, i.e. will not diverge for T + oo;
86
dzr and
@C will remain bounded for every 21 Br, as
Proof. (i) Assume that Fa is stable k-cycle (k: 2 1) denoted by ... = {
Ts {z1, #2, Z1 - (1- zt- -1- Zish - (k- 1)--1
(7)
Then, the largest Lyapunov exponent of &2 is given by
Ar lim In J J ... J22 6- t
h-1
lim
In
(8)
jk
By assumption of stability of T& we have Ar, radius), which implies
<
and also p( [Il=0 Je-n_s) < (the spectral
Jemnss ) O.
lim J7 J2
lim 5-000
(9)
Now suppose that O2, is an orbit of the map eqn. and I's have the same largest Lyapunov exponent, we have
①I converging to l'ks i.c. 21 E Br. Since Oz.
1 T T
Ao = lim
In JT Jr-1 ... J2ll = Ar < 0,
(10)
and hence for z1 E Br
0zy
lim T+000
= lim T
J Jr-1 ... J2ll = O.
(11)
0zi
(fi) & (iii) See Appx. A.2.1
Remark 1. The result of Theorem part (i) will be generally true for any furst-order-Markovinn
recursive map 4 but the conclusions in part (ti) may hinge on its specific definition.
Remark None of the results above and throughou sect. 3 require the dynamics to autonomous, the theory applies whether there is external input not. in fact, mathematically, non-autonomons (externally forced) systems can always be rewritten as autonomous dynamical systems
3 71
see Appx. A. for details.
The results above ensure that loss gradients will not diverge (explode) as T oc in RNNs that are "well-hehaved" in the sense that they converge to fixed point or cycle (see Fig. [T). This is a generalization of the results given Theorem
Schmidt et al. (80), where this was shown only a)
for the specific class of PLRNNs and b) for specific constraints imposed on the eigenvalue spectrum of the RNN's Jacobians which were relaxed in our theorem above.
While our treatment above is centered on the "exploding-gradienns" case, various architectural modifications or regularization techniques ensure that gradients do not vanish either, i.e. remain bounded from below as well. This was established, for instance, in Schmidt et al. (80) for PLRNNs using manifold attractor regularization'. In Appx. A.2. we show that the results from Theorem from Schmidt et al. [80] on doubly bounded (from below and above) loss gradients can indeed be extended to the more general covered by Theorem above.
3.3 Chaotic dynamics
We will now consider the all-important chaotic case. Let F be a recursive map and Oz, {21, Z2, 23. ... } be an orbit of F. The orbit is chaotic if (i) it is not asymptotically periodic and (ii) has at least one positive Lyapunov exponent (26 60). If the system's 's invariant set is bounded, condition (ii) is considered a standard signature of chaos, as in this case two nearby orbits separate exponentially fast, but at the same time their mutual separation cannot go to infinity so that there are also folds. The following theorem states the sufficient condition for exploding gradients:
5

---

(a)
(b)
[| |0ar (0sr//0ill
(0a Wall
10-7
10a8
10"
10
||
10-
10
10F
10
10
10
10
state separation T
slate seperation T
Figure 1: Illustration the exploding gradient problem when training RNNs on dynamical systems. Jacobians
(a) decay away across time separation T for RNN trained simple cycle (ef. Thm. []. but (b) quickly shoot through the roof when training performed chaotic time series (Lorenz system; cf. Thm. 2 Note the doubly-logarithnn scale these graphs.
Theorem Suppose that an RNN Fe R. (parameterized by 0) has chaotic attractor I with Br- as its basin of attraction. Then, for almost every orbit with z1 E Br., (i) the Jacobians connecting temporally distal states zr and 2₁ (T > 1), B2r will exponentially explode for T - xc, and (ii) the tangent vector
Us.'
Szr and so the gradients of the loss function,
8Cr will diverge as T + co.
09
Proof. Let the RNN Fe E R have chaotic orbit denoted by I" { =i. Z7. Then, denoting by J7. the Jacobian at 2 E I", the largest Lyapunov exponent I is given by
lim T
J7 J
(12)
T
Since I* is chaotic, so À > 0. Hence, from
12 it is concluded that
lim 2-000
Jill
lim +00
00, TSL.
(13)
Now, according to Oseledec's multiplicative ergodic Theorem, almost all the points in the basin of attraction of 1* have the same largest Lyapunov exponent A. Thus, (3) holds for almost every z1 E Br..
(it) See Appx. A.2.2
Remark The first part of Theorem holds for all Sesst-oder-Maarkon recursive maps 4 Note özT.
that for LSTMs,
0z, (z == (h, )) denotes the full Jacobian of both hidden and cell states.
We collect some further mathematical results and remarks related Theorem in Appx.
Hence, the essential result is that for all popular RNNs R. and activation functions, loss gradients will inevitably diverge if the RNN latent states converge to chaotic attractor (as illustrated in Fig. Ip).
3.4 Quasi-periodicity
Quasi-periodicity is long-term behavior which occurs on a torus and, superficially, bears some similarity to chaos in the sense that, strictly speaking, orbits also aperiodic. That is, + oc, trajectories will never close up with themselves. Moreover, every trajectory becomes arbitrarily close to any point on the torus, that is, it is dense. One important difference between quasi-pert quasi-periodic and chaotic systems is, however, that in quasi-periodic system, as time passes, two close initial conditions are linearly diverging, while in a chaotic system the divergence is exponential.
Theorem 3. Assume that an RNN Fo E R (parameterized by 0) has quasi-periodic attractor with Bp as its basin of attraction. Then, for every E B
VO<E<1 3T. > 1 s.t.
VT 2 To =
(1 (1- <
(1+6)
(14)
0z
Proof. See Appx.
A.2.3

---

According to Theorem for every orbit converging to a quasi-periodie attractor, the Jacobians Ozr 8z may diverge or vanish as 7 + oo, but this will not occur exponentially fast as T" - oc. Thus, even for bounded non-chaotic RNNs we may sometimes stumble into the problem of diverging gradients. Although this may be a less common scenario, we point out it may occur if we train RNNs on real data from oscillatory systems with incommensurate frequencies, as for instance encountered in electronic engineering,
In Appx. A.3.2 we have collected further mathematical results on the connection between RNN dynamics and loss gradients that hold regardless of the RNN's limiting behavior.
4 Empirical evaluation
Our theoretical results imply that chaotic time series pose a principle challenge for RNN training that cannot easily be circumvented through specifically designed architectures, constraints, or regulariza- tion criteria. If the underlying DS we aim to capture is chaotic, loss gradients propagated back in time will inevitably explode. Hence we need to curtail gradients in an ideal way. The issue arises especially in scientific ML where time series from chaotic systems are ubiquitous and the aim is to reconstruct the generating DS with its limiting behavior. Thus, our exposition will focus on this area.
4.1 Training on systems with exploding gradients by sparse teacher forcing
To illustrate the connections between theory and RNN training, we revive the old idea of TF 995 as a mechanism for truncating error gradients and keeping model-generated trajectories on track while training. However, we would like do this such that important information about the system dynamics does not get lost, for which Lyapunov theory offers some guidance. Specifically, we should not force the system back onto the true trajectory all most of the time (as in "classical TF"), but should effectively "re-calibrate" it only at certain time points chosen wisely according to the system's 's local divergence rates. This procedure will be referred to as sparsely forced BPTT in the following. Assume we want to train an RNN with hidden states z, E IR and linear (or affine) output layer on a time-series 1 222." Dr) generated by chaotic system Fl The linear output layer át Bzg B E R. M maps the RNN hidden states into the observation space. This allows us to modify the original TF procedure by constructing a control series {~1, $2, ... ,rr from the observations by "inverting" the linear output mapping
z BB=BBB) -0)-B''.
B xt.
(15)
The idea is to supply this control signal only sparsely, separated by the learning interval r between consecutive forcings. Hence, defining T {nr + nen as the set of all time points at which we force the RNN onto the 'true' values, the RNN updates can be written as
(2,) ift
21+1
(16)
(NN(2))
else
This forcing is applied after calculation of the loss, such that C, ||3, Bzel2 irrespective of whether t is in T or not (and of course it is applied only during training, not at test time!). Replacing hidden states 2, with their teacher-foreed signals ε, simply breaks divergence between true and predicted trajectories at time points E TT and also cuts off the Jacobians by breaking the temporal contingency (for details see Appx. A.7) The learning interval 7 hence controls how many time steps are included in the gradient calculation and has to be chosen with care such as to balance the effects of exploding gradients vs. those of losing relevant time scales and long-tern term dependencies. While it is general wisdom that an optimal batch size will facilitate training the point here is thus much
'Note that in DS reconstruction one usually considers the data observations (unsupervised problem). "To ensure invertibility, one could add a regularizer AI to B' B in eqn. 1.5 as in ridge regression, but we did not find this necessary in any of our examples.
'It also reminiscent of truncated BPTT. but with the all-important differences that we suggest 1) a theoretically informed choice of the optimal 'truncation length' (forcing interval) and 2) specific procedure for replacing current latent states by control values. As shown in sect. 4.2 & 4.3 both these aspects are indeed crucial to avoid diverging gradients and trajectories whilst not loosing relevant longer time scales.
7

---

(a)
10-
stap
10
Did
PLRNN LSTM
0.0
0
Fjea
100
200
Opeol
100
200
learning interval
learning interval
(b)
1.0-
15
st.pp
Frr
10
0.5
RNN PLRNN
LSTM
0.0
Fyowt
?gend
learning interval
learning interval 5
Figure 2: Overlap in attractor geometry (Dstap. lower = better) and dimension wise comparison of power- spectra (Dy. lower better) against learning interval 7 for (a) the Lorenz and (b) the chaotically forced Dufting oscillator. Continuous lines sparsely forced BPTT. Dashed lines classical BPTT with gradient clipping. Prediction time indicated vertically in black.
more specific: Ideally should be chosen in accordance with the system's Lyapunov spectrum, for instance based on the predictability time
In 2
Tpred
(17)
Amax
Various open-source packages exist for calculating the maximal Lyapunov exponent Amax from empirical time series data (c.g., Julia: DynamicalSystemss..] [12]. C++: TISEAN (30), potentially after delay-embe the data (see Appx. A.4 more details). Note that this needs to be done the empirical data and only once before RNN training, as Amax an invariant characteristic of the empirical system that
aiming reconstruct by RNN (i.e., after successful training the RNN should have the same invariant properties as the underlying DS). We also emphasize that such simple recipe for addressing the exploding gradient problem is based modifying the training routine, and is thus in principle applicable to any model architecture.
4.2 Example 1: Lorenz system and externally forced Duffing oscillator in chaotic regime
Let us illustrate these ideas on two classical textbook examples of chaotic DS, the chaotic Lorenz attractor as an autonomous system, and the chaotically forced Duffing oscillator as an example with explicit external input (see Appx. A.4 for details). Trajectories were repeatedly drawn from these systems, on which we trained a PLRNN, vanilla RNN with tanh activation function, and LSTM by stochastic gradient descent (SGD) to minimize the MSE loss between predicted and actual observations. As optimizer we used Adam [47 from PyTorch (68] with a learning rate of 0.001. For all models, training proceeded solely by sparsely forced BPTT and did not employ gradient clipping or any other technique that may interfere with optimal loss truncation.
forced
In nonlinear DS reconstruction, we are mainly interested in reproducing invariant properties of the underlying system such as the attractor geometry (or topology [78 (2) or the frequency composition (i.e., time-indee ieme-independent properties), while measures like ahead prediction errors are less meaningful especially on chaotic time series [50|96]. Thus, in evaluating training performance, here we follow Koppe et al. [50] in using Kullback-Leible divergence Detsp to quantify the agreement between observed and generated probability distributions across state-space to asses the overlap attractor geometry (Appx. A.5). Moreover, we calculate the dimension-wis Hellinger distance Dnr between power spectra quantify the temporal agreement of the observed and generated time-series (Appx.
A.5).

---

Fig. 2 shows the dependence of the reconstruction quality on the learning interval for all RNN architectures on (a) the Lorenz and (b) the externlly forced Duffing system. Fig. X provides particular examples of reconstructions for chosen too small, too large, or about right.
(a)
(b)
2
0
0
0.0
2.5
-2
Figure 3: Lorenz attractor (blue) and example reconstructions by an LSTM (orange) trained with learning
interval (a) chosen too small (T = 5), (b) chosen optimally (T = 30), and (c) chosen too large (r = 200). Sec
Fig. L
for vanilla RNN example.
For all models we find system-dependent range for the optimal learning interval that agrees well with the predictability time defined in eqn. [7]. where estimates for the maximal Lyapunov exponent were taken from the literature [24] 72). As reference, dashed lines represent the reconstruction performance for all architectures when trained with classical BPTT and gradient clipping. The training procedure control-signal, gradients were normalized to prior to each parameter update (see Fig. T3 for a more systematic evaluation of different clipping procedures and thresholds). As evidenced by the much worse performance, gradient clipping does not effectively address the EVGP, even for LSTMs. As further shown in Fig. 10 in Appx. A.6.4 using the optimal resetting the initial condition for each chunk to either zero or the last forward-iterated state Fe (z) (instead of its control value 2,) equally destroys performance. This suggests that neither mere gradient normalization nor simple windowing are sufficient, but will wipe out essential information about the dynamics.
the same for spansely forced BPTT, except that instead of supplying
any other) window length but
In Appx. A.6 we collect further results on the chaotic Rössler attractor (Fig. & high-dimensional Mackey-Gila Glass equations (Fig. R. and the Lorenz attractor with partial observations (Fig. D
4.3 Example 2: Chaotic weather data
())
(b)
-0.3
16
(log d(e))
Dapp
Du
-1.5
10-4
PERSS USEM
n.-
time-steps
Tyee
learning interval
200
Teed
learning interval
1001
200
Figure 4: (a) maximal Lyapunov exponent was determined as slope of average log-divergence of nearest neighbors in embedding space (n = embedding dimension). (b) Reconstruction quality assessed by attractor overlap (lower better) and dimension-wise comparison of power-spectra (Dz, lower = better). Black vertical lines Tpeed-
As for an empirical example, we trained all RNNs (vanilla RNN, PLRNN, LSTM) on temperature time series recorded at the Weather Station the Max Planck Institute for Biogeochemistry in Jenal Germany To expose the chaotic behavior and obtain a robust estimate of the maximal Lyapunov exponent, trends and yearly cycles were removed, and nonlinear noise-reduction
performed [3]
Appx. A.4p. The maximal Lyapunov exponent was determined with the TISEAN package [30 as
shown Figure (a). The value obtained is in close agreement with the literature (Amaz R 0.017
62
Figure shows that also for these empirical data the optimal training interval agrees well with the predictability time, eqn. (17,, for all trained RNNs. Furthermore, as was the case for the DS
benchmarks, gradient clipping not able to satisfactorily tackle the EVGP, even when paired with

---

architectures like LSTMs explicitly designed for alleviating this problem. Similar results are reported for another real-world dataset, electroencephalogram (EEG) recordings, in Appx. A.6.55
Discussion and conclusions
In this paper proved that RNN dynamics and loss gradients are intimately related for all major types of RNNs and activation functions. If the RNN is "well behaved" in the sense that its dynamics converges to fixed point or cycle, loss gradients will remain bounded, and established remedies (35] 80] can be used to refrain them from vanishing. However, if the dynamics are chaotic, gradients will always explode. This constitutes principle problem in RNN training that cannot easily be mastered through architectural design or gradient clipping. This is because to avoid exploding gradients while training on time series from chaotic systems, one either needs to constrain the RNN so much that chaotic behavior is completely disabled to begin with (i.e., ultimately by forcing all Lyapunov exponents be smaller equal zero), implying very poor fit to such data. Or one needs to be a bit more lenient and thereby allow for the possibility of exploding loss gradients (as LSTMs or PLRNNs in fact do). This problem is furthermore practically highly relevant, as most time series we encounter in nature, and many from man--made systems as well, are inherently chaotic.
While we do not offer a full solution to this problem here, suggest it might be tackled in training by taking a system's 's local divergence rates as measured through the Lyapunov spectrum into account. Hence, rather than conquering the EVGP by structural design specific constraints regularization terms, we recommend to put the focus more on the training process itself. We illustrated this point empirically using sparsely forced BPTT, a training technique that pulls trajectories back on track at times determined by the maximal Lyapunov exponent. Doing so leads to optimal reconstruction results for a variety of simulated and real-world benchmarks, regardless of the specific RNN architecture employed in training.
As noted in sect. exponents from data. Some background knowledge, as provided in classical textbooks (c.g. Ch. 5 in (42). ), may be required for properly reading the output from these packages: Essentially, one would be looking be for a linear scaling region as in Figs. 44 & [2h, ignoring both the initial noise transient well as the plateau caused by reaching the full attractor extent. If unsure about the exact value, a moderate amount of jittering around the estimated mean value may help (see Appx. Fig. 15 A further interesting direction for improvement might be to regulate the forcing interval through an annealing procedure 2). for instance starting at 7 = and ramping up to 7 Tpred throughout training, similar as in adaptive schemes short-term behavior right, and then challenge the system more and more for longer time spans until the predictability time is reached.
fairly standard packages are available for computing maximal Lyapunov
5 The idea here would be to first get the
We stress that our goal here above all was to provide a mathematically grounded perspective on the problem, with the empirical section focused on elucidating the practical implications of the theoretical results. We believe that more thorough theoretical understanding important and needed for guiding future research into more powerful training procedures that avoid exploding gradients without compromising expressiveness. In our application examples, we developed the case from the perspective of scientific machine learning, which by now is a broad area in its own right with huge societal relevance (e.g., climate or epidemiological time series), and where the reconstruction of geometrical or topological (invariant) properties is important, beyond mere prediction. Nevertheless, we believe that our theoretical results will also have implications for other domains, like NLP (36). While scientific time series problems traditionally have been extensively considered from a DS perspective (e.g.. (42) ) much more groundwork is needed, however, in areas like NLP, where, for instance, it may not even be immediately clear how to best define a Lyapunov spectrum.
All code from this paper is available at http://s/g.tttlooocco/////tttata
Acknowledgements
This work was funded by the German Research Foundation (DFG) under Germany's Excellence Strategy EXC-2181 390900948 (STRUCTURES), and through grant Du 354/10-1 to DD.
We thank of the referees for pointing this out.
10

---

References
[1] Abarbanel. Predicting the Fuuure: Completing Models of Observed Complex Systems. Understanding Complex Systems. Springer- 2013.
[2] H. D. I. Abarbanel, P. J. Rozdeba, and S. Shirman. "Machine Learning: Deepest Learning as Statistical Data Assimilation Problems". In: Neural Computation 30.8 (2018), pp. 2025-2055.
[3] K. T. Alligood, T. D. Sauer, and J. A. Yorke. Chaos: An Introduction to Dynamical Systems. Springer, New York, NY, 1996.
[4] M. Arjovsky. A. Shah, and Y. Bengio. "Unitary Evolution Recurrent Neural Networks". In: Proceedings of the 33rd International Conference on International Conference on Machine Learning. Vol. 48. ICML' 16. New York, NY, USA: JMLR, 2016, pp. 1120-1128.
[5] S. Bengio, O. Vinyals. N. Jaitly, and N. Shazeer. "Scheduled sampling for sequence prediction with recurrent Neural networks". In: Proceedings of the 28th International Conference on Neural Information Processing Systems Volume 1. NIPS' 15. Cambridge, MA, USA: MIT Press, 2015. pp. 1171-1179.
[6] Y. Bengio, P. Simard, and P. Frasconi. "Learning long-term dependencies with gradient descent is difficult". In: IEEE Transactions on Neural Networks 5.2 (1994), pp. 157-166.
[7] B. P. Bezruchko and D. A. Smirnov. Extracting Knowledge From Time Series: An Introduction to Nonlinear Empirical Modeling. Springer Series in Synergetics. Berlin Heidelberg: Springer-
Verlag, 2010. DO1: 10. 1007/0009991 0977/998330 642 12601-7
[8] M. Brenner, F. Hess, J. M. Mikhaeil, L. F. Bereska, Z. Monfared, P.-C. Kuo, and D. Durstewitz. "Tractable Dendritic NNNS for Reconstructing Nonlinear Dynamical Systems". In: Proceedings of the 39th International Conference Machine Learning. Vol. 162. Proceedings of Machine Learning Research. PMLR, 2022, pp. 2292-2320.
[9] B. Chang, M. Chen, E. Haber, and F. H. Chi. "AntisymmetricRNN: A Dynamical System View on Recurrent Neural Networks". In: Proceedings of the International Conference on Learning Representations (ICLR). 2019.
[10] K. Cho, B. van Merriënboer, D. Bahdanau, and Y. Bengio. "On the Properties of Neural Machine Translation: Encoder- Approaches". In: Proceedings of SSST-8, Eighth Workshop Syntax, Semantics and Structure in Statistical Translation. Association for Computational Linguistics, 2014, pp. 103-111.
[11] K. Cho, B. Merriënboer, C. Gulcehre, F. Bougares, H. Schwenk, and Y. Bengio. "Learning phrase representations using RNN encoder-decoder for statistical machine translation". In: Proceedings of the 2014 Conference on Empirical Methods in Natural Language Processing (EMNLP). 2014.
[12] Datseris. "DynamicalSSstems...i
Julia software library for chaos and nonlinear dynamics". In: Journal of Open Source Software 3.23 (2018), p. 598.
[13] K. Doya. "Bifurcations in the learning of recurrent neural networks". In: Proceedings of the IEEE International Symposium Circuits and Systems. Vol. 6. 1992, pp. 2777-2780.
[14] J. Duarte, C. Januário, N. Martins, and J. Sardanyés. "Quantifying chaos for ecological stoichiometry". In: Chaos: An Interdisciplinary Journal of Nonlinear Science 20.3 (Sept. 2010), p. 033105.
[15] G. Duffing. Forced vibrations with variable natural frequencies and their technical significance. F. Vieweg und Sohn, Braunschweig, 1918.
[16] D. Durstewitz. Advanced Data Analysis in Neuroscience: Integrating Statistical and Com- putational Models. Bernstein Series in Computational Neuroscience. Springer International Publishing, 2017.
[17] D. Durstewitz and T. Gabriel. "Dynamical Basis of Irregular Spiking in NMDA-Driven Preffronta Cortex Neurons". In: Cerebral Cortex 17.4 (2007), pp. 894-908.
[18] R. Engelken, F. Wolf, and I.. F. Abbott. "Lyapunov spectra of chaotic recurrent neural net- works" In: arXiv:2006. 02427 (2020).
[19] N. B. Erichson, O. Azencot, A. Queiruga, L. Hodgkinson, and M. W. Mahoney. "Lipschitz Recurrent Neural Networks". In: Proceedings of the International Conference on Learning Representations (ICLR). 2021.
11

---

[20] M. Faggini. "Chaotic time series analysis in economics: Balance and perspectives". In: Chaos: An Interdisciplinary Journal of Nonlinear Science 24.4 (2014). Publisher: American Institute of Physics. p- 042101.
[21] R. J. Field and L. Györgyi. Chaos in Chemistry and Biochemistry. World Scientific, 1993.
[22] S. Fusi, W. F. Asaad, E. K. Miller, and X.-J Wang. "A neural circuit model of flexible sensorimotor mapping: learning and forgetting multiple timescales". In: Neuron 54.2 (2007), pp. 319-333.
[23] Z. Ghahramani and S. Roweis. "Learning Nonlinear Dynamical Systems Using an EM Algo- rithm".. In: Advances in Neural Information Processing Systems. Vol. 11. MIT Press, 1999.
[24] W. Gilpin. "Chaos as an interpretable benchmark for forecasting and data-driven modelling". In: Proceedings of the Thirty-fitth Conference on Neural Information Processing Systems Datasets and Benchmarks Track. 2021.
[25] L. Glass and M. C. Mackey, "Pathological conditions resulting from instabilities in phys- iological control systems". In: Annals of the New York Academy of Sciences 316 (1979), pp. 214-235.
[26] P. A. Glendinning and D. J. W. Simpson. "A constructive approach to robust chaos using invariant manifolds and expanding cones". In: Discrete & Continuous Dynamical Systems 41.7
(2021), pp. 3367-3387.
[27] A. L. Goldberger, L. A. N. Amaral, Glass, J. M. Hausdorff, P. C. Ivanov, R. G. Mark, J. E. Mietus, G. B. Moody, C.-K. Peng, and H. E. Stanley. "PhysioBank, PhysioToolkit, and PhysioNet: Components of New Research Resource for Complex Physiologic Signals". In: Circulation 101. .33 (2000), e215-e220
[28] I. Goodfellow, Y. Bengio, and A. Courville. Deep Learning. MIT Press, 2016.
[29] A. Graves, G. Wayne, M. Reynolds, and et al. "Hybrid computing using neural network with dynamic external memory". In: Nature 538 (2016), pp. 471 -476.
[30] R. Hegger, H. Kantz, and T. Schreiber. "Practical implementation of nonlinear time series methods: The TISEAN package".. In: Chaos: Interdisciplinary Journal of Nonlinear Science 9.2 (1999), pp. 413-435.
[31] K. Helfrich, D. Willmott, and Q. Ye. "Orthogonal Recurrent Neural Networks with Scaled Cayley Transform". In: Proceedings of the International Conference on Machine Learning. PMLR, 2018, pp. 1969-1978.
[32] M. Henaff, A. Szlam, and Y. LeCun. "Recurrent Orthogonal Networks and Tasks". In: Proceedings of the International Conference Machine Learning. PMLR, 2016, pp. 2034-2042.
[33] Hernandez, A. K. Moretti, Z. Wei, S. Saxena, Cunningham, and Paninski. "Nonlinear Evolution via Spatially-Dependdee Linear Dynamics for Electrophysiology and Calcium Data". In: arXiv preprini 18Xiv:1811.00 1.02459 (2020).
[34] S. Hochreiter. Untersuchungen zu dynamischen neuronalen Netzen. Diploma thesis. Institut für Informatik, Technische Universität München. 1991.
[35] S. Hochreiter and J. Schmidhuber. "Long Short-Term Memory". In: Neural Computation 9.8
(1997), pp. 1785-1780.
[36] K. Inoue, S. Ohara, Y. Kuniyoshi, and K. Nakajima. "Transient Chaos in BERT". In: arXiv:2106. 101.28006 (2021).
[37] L. Jing, C. Gulcehre, J. Peurifoy, Y. Shen, M. Tegmark, M. Soljacic, and Y. Bengio. "Gated Orthogonal Recurrent Units: On Learning to Forget". In: Neural Computation 31.4 (2019),
pp. 765-783.
[38] A. Kag and V. Saligrama. "Time Adaptive Recurrent Neural Network". In: Proceedings of the IEEE/CVE Conference on Computer Vision and Pattern Recognition (CVPR). IEEE, 2021, pp. 15144-15153.
[39] A. Kag and V. Saligrama. "Training Recurrent Neural Networks via Forward Propagation Through Time". In: Proceedings of the 38th International Conference on Machine Learning. Vol. 139. Proceedings of Machine Learning Research. PMLR, 2021, pp. 5189-5200.
12

---

[40] A. Kag, Z. Zhang, and V. Saligrama. "RNNs Incrementally Evolving on an Equilibrium Mani- fold: A Panacea for Vanishing and Exploding Gradients?" In: Proceedings of the International Conference on Learning Representations. 2020.
[41] L. Kamdjeu Kengne, J. R. Mboupda Pone, and H. B. Fotsin. "On the dynamics of chaotic circuits based on memristive diode-bridge with variable symmetry: A case study". In: Chaos, Solitons & Fractals 145 (2021), p- 110795.
[42] H. Kantz and T. Schreiber. Nonlinear Time Series Analysis. 2nd ed. Cambridge University Press, 2003.
[43] H. Kantz, T. Schreiber, I. Hoffmann, T. Buzug, G. Pfister, L. G. Flepp. J. Simonet, R. Badii, and E. Brun. "Nonlinear noise reduction: A case study on experimental data". In: Physical Review E 48 (1993), pp. 1529-1538.
[44] M. B. Kennel, R. Brown, and H. D. I. Abarbanel. "Determining embedding dimension for phase-space reconstruction using a geometrical construction". In: Phys. Rev. A 45 (1992).
[45] G. Kerg, K. Goyette, M. Puelma Touzel, G. Gidel, E. Vorontsov, Y. Bengio, and G. Lajoie. "Non-normal Recurrent Neural Network (nnRNN): learning long time dependencies while im- proving expressivity with transient dynamics". In: Advances in Neural Information Processing Systems. Vol. 32. Curran Associates, Inc., 2019.
[46] M. Kesmia, S. Boughaba, and S. Jacquir. "Control of continuous dynamical systems modeling physiological states", In: Chaos, Solitons & Fractals 136 (2020), p. 109805.
[47] D. P. Kingma and J. Ba. "Adam: A Method for Stochastic Optimization". In: 3rd International Conference on Learning Representations. Ed. by Y. Bengio and Y. LeCun. ICLR, 2015.
[48] D. P. Kingma and M. Welling. "Auto-Encoding Variational Bayes". In: Proceedings of the 2nd
International Conference on Learning Representations. 2014.
[49] J. Z. Kolter and G. Manek. "Learning Stable Deep Dynamics Models", In: Advances in Neural Information Processing Systems. Vol. 32. Curran Associates, Inc., 2019.
[50] G. Koppe, H. Toutounji, P. Kirsch, S. Lis, and D. Durstewitz. "Identifying nonlinear dynamical systems via generative recurrent neural networks with applications to fMRI". In: PLOS Computational Biology 15.8 (2019), pp. 1-35.
[51] D. Kramer, P. L. Bommer, C. Tombolini, G. Koppe, and D. Durstewitz. "Reconstructing Nonlinear Dynamical Systems from Multi-Modal Time Series". In: Proceedings of the 39th International Conference on Machine Learning. Vol. 162. PMLR, 2022.
[52] J. Laskar and P. Robutel. "The chaotic obliquity of the planets". In: Nature 361.6413 (1993), pp. 608 6122. (Visited on 09/10/2021).
[53] S. Li. W. Li. C. Cook, C. Zhu, and Y. Gao. "Independently Recurrent Neural Network (In- dRNN); Building A Longer and Deeper RNN", In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition. 2018, pp. 5457-5466.
[54] L. B. Kovachki, K. Azizzadeneshelli B. liu, K. Bhattacharya, A. Stuart, and A, Anandku- mar. "Fourier Neural Operator for Parametric Partial Differential Equations". In: Proceedings of the International Conference on Learning Representations. 2021.
[55] Z. C. Lipton, J. Berkowitz, and C. Elkan. "A Critical Review of Recurrent Neural Networks for Sequence Learning". In: arXiv: 00019 [cs] (2015).
[56] E. N. Lorenz. "Deterministic Nonperiodic Flow". In: Journal of the Atmospheric Sciences 20.2 (1963), pp. 130-141.
[57] N. Maheswaranathan, A. Williams, M. Golub, S. Ganguli, and D. Sussillo. "Reverse engi- neering recurrent networks for sentiment classification reveals line attractor dynamies". In: Advances in Neural Information Processing Systems. Vol. 32. Curran Associates, Inc., 2019.
[58] S. Mangiarotti, M. Peyre, Y. Zhang, M. Huc, F. Roger, and Y. Kerr. "Chaos theory applied to the outbreak of COVID-19: ancillary approach to decision making in pandemic context". In: Epidemiology and Infection 148 (2020), e95.
[59] R. M. May. "Chaos and the Dynamics of Biological Populations". In: Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences 413.1844 (1987), pp. 27-44.
[60] J. Meiss. Differential Dynamical Systems. Society for Industrial and Applied Mathematics,
2007.
13

---

[61] J. Menick, E. Elsen. U. Evci, S. Osindero, K. Simonyan, and A. Graves. "Practical Real Time Recurrent Learning with Sparse Approximation", In: Proceedings of the International Conference on Learning Representations. 2021.
[62] H. Millán, B. Ghanbarian- -laviihh,,, and Garcia-Fornaris "Nonlinear dynamics of mean daily temperature and dewpoint time series Babolsar, Iran, 1961-2005". In: Almospheric Research. Clouds, Aerosols and Radiation 98.1 (2010), pp. 89-101. DOT: 10.1016/] atw.sres..................o 2010
06.001 001001
[63] 7. Monfared and D. Durstewitz. "Existence of n-cycles and border-collision bifurcations in piecewise-linear continuous maps with applications to recurrent neural networks". In: Nonlinear Dynamics 101.2 (2020), pp. 1037-1052.
[64] 7. Monfared and D. Durstewitz. "Transformation of ReILU-based recurrent neural networks from discrete-time to continuous-time" In: Proceedings of the 37th International Conference on Machine Learning. PMLR, 2020, pp. 6999-7009.
[65]
I.. Olsen and Degn. "Chaos in an enzyme reaction". In: Nature 267.5607 (1977), pp. 177-
178.
[66] E. Orhan and X. Pitkow. "Improved memory in recurrent neural networks with sequential non-normal dynamies". In: Proceedings of the International Conference on Learning Repre- sentations. ICLR, 2020.
[67] R. Pascanu, T. Mikolov, and Y. Bengio. "On the Difficulty of Training Recurrent Neural Networks". In: Proceedings of the 30th International Conference on International Conference on Machine Learning. Vol. 28, ICML' 13. JMLR, 2013, 111-111--11--3188
[68] Paszke, S. Gross, S. Chintala, G. Chanan, E. Yang, Z. Vito, Z. Lin, A. Desmaison, Antiga, and A. Lerer. "Automatic differentiation in PyTorch". In: Proceedings of the 31st Conference on Neural Information Processing Systems. NIPS, 2017.
[69] J. Pathak, B. Hunt, M. Girvan, Z. Lu, and E. Ott. "Model-Free Prediction of Large Spatiotem- porally Chaotic Systems from Data: A Reservoir Computing Approach". In: Phys. Rev. Lett. 120 (2018), p. 024102. DO1: 10. 1103 Phwysheyeseettttt.meeeee 120. 024102
[70]
B. Pearlmutter. Dynamic recurrent neural networks. Carnegie Mellon University. 1990.
[71] L. Perko. Differential Equations and Dynamical Systems. Vol. 7. Springer, New York, NY, 2001.
[72] M. T. Rosenstein, J. J. Collins, and C. J. De Luca. "A practical method for calculating largest Lyapunov exponents from small data sets". In: Physica D: Nonlinear Phenomena 65.1 (1993), pp. 117-134. DOT: 10. 1016/0167771111 2789 (93) 90009-P
[73] T. K. Rusch and S. Mishra. "Coupled Oscillatory Recurrent Neural Network (coRNN): An ac- curate and (gradient) stable architecture for learning long time dependencies". In: Proceedings of the International Conference on Learning Representations. ICLR, 2021.
[74] T. K. Rusch and S. Mishra. "UnICORNN: A recurrent model for learning very long time dependencies". In: Proceedings of the 38th International Conference on Machine Learning. Vol. 139. PMLR, 2021, pp. 9168-911911 9178.
[75] E. Russo and D. Durstewitz. "Cell assemblies at multiple time scales with arbitrary lag constellations". In: elife 6 (2017). Ed. by M. Howard, e19428.
[76] E. Rössler. "An equation for continuous chaos". In: Physics Letters 57.5 (1976), pp. 397- 398.
[77] S. Sabanal and M. Nakagawa. "The fractal properties of vocal sounds and their application in the speech recognition model". " In: Chaos, Solitons Fractals' 7.11 (1996), pp. 1825-1843.
[78] T. Sauer, J. A. Yorke, and M. Casdagli. "Embedoogy.. In: Journal of Statistical Physics 65.3 (1991), pp. 579-616.
[79] G. Schalk, D. McFarland, T. Hinterberger, N. Birbaumer, and J. R. Wolpaw. "BCI2000: general-purpose brain-computer interface (BCI) system". In: IEEE transactions on bio-medical e engineering 51.6 (2004), pp. 1034-1043.
[80] D. Schmidt, G. Koppe, Z. Monfared, M. Beutelspacher, and D. Durstewitz. "Identifying nonlinear dynamical systems with multiple time scales and long-range dependencies". In: Proceedings of the International Conference Learning Representations. ICLR, 2021.
14

---

[81] B. Sivakumar. "Chaos theory in geophysics: past, present and future".. In: Chaos, Solitons & Fractals. Fractals in Geophysies 19.2 (2004), pp. 441-462.
[82] F. Takens. "Detecting strange attractors in turbulence". In: Dynamical Systems and Turbulence, Warwick 1980. Springer Berlin Heidelberg, 1981, pp. 366-381.
[83] R. Tchitnga, B. A. Mezatio, T. Fonzin Fozin, R. Kengne, P. H. Louodop Fotso, and A. Fomethe.
"A novel hyperchaotic three-component oscillator operating at high frequency", In: Chaos, Solitons Fractals 118 (2019), pp. 166-180.
[84] D. J. Thomson. "Time series analysis of Holocene climate data". Philosophical Transactions of the Royal Society of London. Series A, Mathematical and Physical Sciences 330.1615 (1990), pp. 601-616.
[85] T. H. Trinh, A. M. Dai, M.-T. Luong, and Q. V. Le. "Learning Longer-term Dependencies in RNNs with Auxiliary Losses". In: arXiv:1803.001444 (2018).
[86] A. P. Trischler and G. M. D'Eleuterio. "Synthesis of recurrent neural networks for dynamical system simulation". In: Neural Networks 80 (2016), pp. 67-78.
[87] P. Turchin and A. D. Taylor. "Complex Dynamics in ecological Time Series". In: Ecology 73.1 (1992), pp. 289-305.
[88] E. Tziperman, H. Scher, S. Zebiak, and M. Cane. "Controlling Spatiotemporal Chaos in a Realistic El Niño Prediction Model". In: Physical Review Letters 79 (1997), pp. 1034-1037.
[89] P. R. Vlachas, G. Arampatzis, C. Uhler, and P. Koumoutsakos. "Learning the Effective Dynam- ics of Complex Multiscale Systems". In: arXiv:2006. 13431 [nlin, physics:physics] (2020).
[90] P. R. Vlachas, W. Byeon, Z. Y. Wan, T. P. Sapsis, and P. Koumoutsakos. "Data-driven fore- casting of high-dimensional chaotic systems with long short-term memory networks". In: Pro- ceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 474.2213 (2018), p. 20170844.
[91] R. Vogt. M. Puelma Touzel, E. Shlizerman, and Lajoie. "On Lyapunov Exponents for RNNs: Understanding Information Propagation Using Dynamical Systems Tools". In: Frontiers in Applied Mathematics and Statistics 8 (2022).
[92] H. U. Voss, J. Timmer, and J. Kurths. "Nonlinear dynamical system identification from uncer- tain and indirect measurements". In: Proceedings of the International Journal of Bifurcation and Chaos 14.06 (2004), pp. 1905-1933.
[93] C. van Vreeswijk and H. Sompolinsky. "Chaos in Neuronal Networks with Balanced Excitatory and Inhibitory Activity". " In: Science 274.5293 (1996), pp. 1724-1726.
[94] J. H. M. Wedderburn. Lectures on Matrices. New York: American mathematical society, New York: Dover Publications, 1964.
[95] R. J. Williams and D. Zipser. "A Learning Algorithm for Continually Running Fully Recurrent
Neural Networks". In: Neural Computation 1.2 (1989), pp. 270-280.
[96] S. N. Wood. "Statistical inference for noisy nonlinear ecological dynamic systems". In: Nature 466.7310 (2010).
[97] Zhang, D. Liu, and Wang. Controlling Chaos. Springer, London, 2009.
15

---

1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's
contributions and scope? [Yes]
(b) Did you describe the limitations of your work? [Yes] See Sec.
S
(c) Did you discuss any potential negative societal impacts of your work? [N/A]
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [Yes] We provide comprehensive proofs of all presented theorems in the Appendix.
(b) Did you include complete proofs of all theoretical results? [Yes]
3. If you ran experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main experi- mental results (either in the supplemental material or as URL)? [Yes] The code will be made public after publication.
(b) Did you specify all the training details (e.g., data splits, byperparameters, how they were chosen)? [Yes]
(c) Did you report error bars (e.g., with respect the random seed after running experi- ments multiple times)? [Yes]
(d) Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, cloud provider)? [N/A]
4. If you using existing assets (e.g., code, data, models) curating/releasing assets...
(a) If your work uses existing assets, did you cite the creators? [Yes]
(b) Did you mention the license of the assets? [Yes]
(c) Did you include any new assets either in the supplemental material or as URL? [No]
(d) Did you discuss whether and how consent was obtained from people whose data you're using/curating? [N/A]
(e) Did you discuss whether the data you using/curating contains personally identifiable information or offensive content? [N/A]
5. If you used crowdsoureing or conducted research with human subjects...
(a) Did you include the full text of instructions given participants and screenshots, applicable? [N/A]
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent participant compensation? [N/]]
16


Reviews:
[{"review_id": "qkrpWC52I_v", "reviewer": "Reviewer_Dedw", "summary": "This paper concerns the issue of using RNNs to learn chaotic dynamics, a phenomenon exhibited in many important applications. To address the exploding/vanishing gradient problem (EVGP) when training RNNs, other works have typically proposed specific RNN architectures, however, these architectural choices rule out chaotic dynamics. This work performs theoretical analysis to formalize these ideas by relating the loss gradients of RNNs to the RNN dynamics using Lyapunov exponents. A result is presented for RNNs that converge to a stable fixed point or cycle showing that the loss gradients of these models will not diverge. A result for RNNs converging to chaotic dynamics shows that these models' loss gradients will inevitably diverge. A result for quasi-periodicity shows that there exist cases for bounded non-chaotic RNNs whose loss gradients can diverge.\n\nGiven these theoretical results, the work then performs an empirical evaluation to show the practical implications of the theoretical results. The main idea is that since we cannot both control EVGP with architectural choices and learn chaotic dynamics, we need to use a different approach to control EVGP. This paper suggests the use of a form of teacher forcing where at sparse time points the \"true\" state is computed given the observation and weights at that timestep. This has the effect of truncating the loss gradients.  The particular sparse time points are determined by the system's Lyapunov spectrum using the predictability time. Examples are provided using a Lorenz system, Duffing oscillator and weather data, showing that this approach performs better than backpropagation through time BPTT with gradient clipping as well as truncated BPTT where the initial states are set to zero and the truncated timepoints are chosen using the predictability time. ", "questions": "1. If I am not mistaken, it seems to me that \"sparsely forced BPTT\" is basically a form of truncated BPTT where the truncation points are determined by the Lyapunov exponents, and the initial state that kicks off the truncated segment is the computed \"true\" state (using the TF approach) from the last step of the previous segment. It seems it would be helpful to the reader to point this different viewpoint out and perhaps highlight any practical differences. Section A.6.4 starts to make this connection but it seems a more thorough discussion in the main text would be helpful to the reader. Please let me know if I am misunderstanding something though.\n\n2. The baselines used for the empirical results are: BPTT with gradient clipping and truncated BPTT where the initial state of each truncated segment is set to zero and the truncation points are determined by the Lyapunov exponents. This latter baseline shows that just truncating the segments by  $\\tau_{pred}$ is not sufficient and using the teacher forcing approach is important. But this raises the question of the necessity of using $\\tau_{pred}$ for the windowing? I.e., how would this approach do using evenly (or randomly) spaced truncation points but using the TF method to initialize the initial state? It is currently difficult for the reader to judge if using $\\tau_{pred}$ is actually important. Another baseline, (though I would rank this as less important than the previously mentioned one), might be to learn the initial state as opposed to setting it to zero, similar to the \"multiple shooting\" ideas mentioned in the related work. I think these additional baselines, particularly the first one, would help strengthen the claims related to the proposed training method.\n\n3. How practical is computing the predictability time, $\\tau_{pred}$, for general systems? It is noted that the maximal Lyapunov exponents for the Lorenz and Duffing systems were taken from the literature and this exponent was computed using the TISEAN package for the weather data. It seems a discussion of computing this exponent is warranted to give the reader a better sense of what is required to get this value. The Vogt et al. 2022 paper that is cited seems to propose a method to estimate the Lyapunov spectrum. Perhaps the usefulness or limitations of that method could be discussed?\n\n", "limitations": "The authors acknowledge they have not proposed a full solution to the main problem, but rather have helped to formalize the problem and offer potential directions. ", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Originality: While other works (cited in the related work section) have explored the connections between Lyapunov exponents and loss gradients before, this appears to be the first work to present formal theory for this connection in general RNN architectures, particularly for the chaotic case.  While the \"sparsely forced BPTT\" appears to be a specific variant of truncated BPTT (though see Question #1 related to this below), the teacher forcing approach seems useful and the idea to base the truncation points on the Lyapunov spectrum appears to be novel. \n\nQuality:  The theoretical results appear to be sound (though I did not go as closely through every line of the additional proofs in the appendix). While the authors acknowledge they have not fully solved the difficult problem of controlling EVGP for problems with chaotic dynamics, the theoretical results point to attempting to address this with training methods as opposed to architectural constraints and the empirical results help to illustrate the proposed training method.  However, it seems some lingering questions related to the effectiveness of using the predictability time to choose the truncation points could be resolved with a few additional baselines (see Question #2 below). It is also currently unclear how practical using this predictability time is for more general problems (See question #3 below)?\n\nClarity: The paper is well-written and clear.\n\nSignificance: I think this paper is useful in formalizing the problem of attempting to control the EVGP while also attempting to learn chaotic dynamics with RNNs. Given the prevalence of chaotic systems, this is an important problem. The point that this problem cannot be addressed through architectural constraints is important and can hopefully lead to future training methods that can successfully address this.\n\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "nSHfZ3X4A6i", "reviewer": "Reviewer_94xf", "summary": "This paper studies the problem of training RNNs to learn chaotic dynamics. First, the paper draws analytic connections between exploding and vanishing gradients (due to the repeated product of Jacobian terms when computing gradients), and dynamical systems characterized by Lyapunov exponents. The paper proves that in order for the gradients and loss of an RNN to be well behaved, the RNN must converge to a limit cycle or fixed point. The consequence of this is that it will be difficult to use gradient-based optimization to train an RNN to learn chaotic dynamics.\n\nThe paper proposes a new solution to address this problem: the idea of resetting network dynamics periodically to a \"best-fit\" state given the readout weights. They call this method \"sparsely forced\" BPTT, and compare it to other methods of stabilizing RNN training (gradient clipping, and windowing the time series into distinct sections). The paper also provides a heuristic for selecting how often to force/reset the state, based on the maximum Lyapunov exponent of the underlying DS.", "questions": "Rather than use a fixed learning interval (tau), would there be any benefit to annealing tau over time? A natural choice in my mind would be to start small and slowly grow tau, and use early stopping with a validation set to determine when the network has started to diverge? Curious if the authors think this is reasonable.", "limitations": "yes", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths\n1. The paper is well written, clearly motivated, and has an appropriate amount of context (with additional detail present in the supplement). I found it a pleasure to read.\n2. The problem the paper addresses is significant and, to my knowledge, their proposed solution is novel.\n3. I really like the simplicity of sparsely forced BPTT.\n\nWeaknesses\n1. For the empirical evaluation, the authors mention they compare against a baseline of BPTT with gradient clipping as well as with a simple windowing scheme (resetting the state to zero). I think these comparisons should be more fleshed out in the main text. For gradient clipping, what values of gradient clipping were tried? Was it clipping by value or shrinking the norm? For windowing, what window lengths were tried? Does the sparsely forced BPTT outperform all possible choices for those hyperparameters?\n\nMinor suggestions\n- I don't have a concrete suggestion here, but I feel like Figures 1 and 2 emphasize the wrong information. As drawn, they emphasize the difference between the three architectures (that's what the eye gets drawn to first), but the important comparison is between the solid lines (the proposed method) and the dashed lines (the baseline method). I wonder if there is someway of regenerating those figures that puts the emphasis on the comparison between ST-BPTT and gradient clipping.\n- Line 310: Instead of saying that the value of the exponent is in close agreement with previous work, why not just say what the value is from that previous work?\n- I wonder if there is a simple figure or toy example that can demonstrate Theorem 1 and 2. For example, if one were to train a toy RNN on a chaotic system, can you show how the gradient terms diverge?", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "jFrW-KDbCuO", "reviewer": "Reviewer_3EGL", "summary": "The hidden states of recurrent neural networks are analyzed for both chaotic and non-chaotic time series and proved to suffer from exploding gradients in case of approximating the former type of dynamics. Furthermore, an efficient method to still train RNNs on chaotic time series is proposed and exhaustively validated on synthetic and real-world data.", "questions": "- In line 121 an explanation of $h$ could not be found. How does the author's formulation relate to the conventional $\\mathbf{z_t} = f((\\mathbf{W_s}\\mathbf{s_t} + \\mathbf{B_s}\\mathbf{s_t}) + (\\mathbf{W_z}\\mathbf{z_{t-1}} + \\mathbf{B_z}\\mathbf{z_{t-1}}))$ notation, where two different weight matrices $\\mathbf{W_s}$ and $\\mathbf{W_z}$ are used for the input and past hidden state mapping, respectively?\n- In equation (17) the authors suggest to choose the teacher forcing frequency of the hidden state as $\\tau_{\\text{pred}} = \\frac{\\ln2}{\\lambda_{\\text{max}}}$. Hence, it depends on the largest Lyapunov exponent $\\lambda_{\\text{max}}$ (which likely changes during training). Does this mean, $\\tau_{\\text{pred}}$ should be assessed and updated regularly during training, based on a repeatedly computed $\\lambda_{\\text{max}}$? I might be wrong with my assumption of a changing $\\lambda_{\\text{max}}$ during training, when it is calculated on the training data and not on $\\mathbf{z_t}$ (I understood it to be the latter, though). A brief explanation of how $\\lambda_{\\text{max}}$ was computed (using the training data or the RNN's hidden state) would be highly appreciated.\n- Figure 1 indicates that LSTMs behave more consistent with the proposed theory than RNNs. Did the RNNs have more difficulties in learning the processes in general, which shows up in the $D_{stsp}$ and $D_{H}$ analysis, or how do the authors explain this deviation between RNNs and LSTMs? A plot like Figure 2 for RNNs would be appealing and informative.", "limitations": "Limitations are addressed by the authors.", "rating": 9, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "### Originality\nStrengths:\n- Computing the largest Lyapunov exponent $\\lambda_{\\text{max}}$ of the RNN orbit appears to be an appealing and highly informative measure about the RNN's nature.\n\n### Quality:\nStrengths:\n- The problem of modeling chaotic time series and the proposed analysis are of high value for the ML community.\n- Theoretical proofs are excellently verified empirically on a broad variety of synthetic and real-world processes with varying difficulty, suggesting a fairly simple but efficient method of how to train RNNs on chaotic time series.\n\nWeaknesses:\n- The authors might include [this work](https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.120.024102), where an ESN has been shown to successfully predict chaos. \n\n### Clarity and significance:\nStrengths:\n- According to my assessment, the analysis of exploding gradients when trianing an RNN on chaotic time series is of very high relevance for the ML community and thus constitutes a high significance of this work.\n- The relevance and significance of this work is further amplified by the proposed effective method to train RNNs on chaotic time series.\n- Despite the complex nature of the discussed problem, the paper is written and structured particularly clear and applicable.\n\nWeaknesses:\n- In line 172, does $\\rho$ relate to the spectral radius? Might be helpful for the readers to specify.\n- In line 173, please add \"equation\" before (1) to a facilitate application.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "PMPYyDrLJbR", "reviewer": "Reviewer_okuc", "summary": "This paper demonstrates that the loss gradients are inherently linked to RNN dynamics. When RNN dynamics converge to a fixed point or cycle, the paper shows that loss gradients will remain bounded. When RNN dynamics are chaotic, gradients will always explode and cannot be addressed by architectural designs or gradient clipping. The paper proposes one way of addressing the problem on chaotic data, by tracking the Lyapunov spectrum during training and uses teacher forcing when local divergence rates becomes too large. ", "questions": "- It is unclear to me why other forms of constraints and regularization cannot be used to alleviate the exploding gradient problem, which the sparsely forced BPTT can. I'd appreciate a more intuitive explanation from the authors.", "limitations": "Yes", "rating": 7, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "### Strengths:\n- This paper is decently well written and has a good overview of related works.\n- The paper shows that the exploding and vanishing gradient problem in RNNs can be easily addressed when the RNN is well-behaved, but will necessarily be problematic when modeling chaotic dynamics with RNNs. This principled view into the problem and the conditions on which they occur appear to be significant. \n- The paper proposes sparsely forced BPTT as a way to enable the learning of chaotic dynamics.\n\n### Weaknesses\n- It's not clear how applicable the sparsely forced BPTT method would be in real-world applications, as it would be expensive to calculate the Lyapunov spectrum throughout training.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept