Learning interpretable control inputs and dynamics underlying animal locomotion
Authors: Thomas Soares Mullen
Marine Schimel
Guillaume Hennequin
Christian K. Machens
Michael Orger
Adrien Jouary
Venue: ICLR 2024 poster (2024)
Keywords: computational neuroscience
interpretable dynamics
motor control
animal behavior
dynamical systems
system identification
unsupervised learning
zebrafish
Primary area: applications to neuroscience & cognitive science
TL;DR: We proposed a novel approach to modeling time series of behavior observations by combining two existing methods in order to learn a reduced and interpretable model of behavioral dynamics

Abstract: A central objective in neuroscience is to understand how the brain orchestrates movement. Recent advances in automated tracking technologies have made it possible to document behavior with unprecedented temporal resolution and scale, generating rich datasets which can be exploited to gain insights into the neural control of movement. One common approach is to identify stereotypical motor primitives using cluster analysis. However, this categorical description can limit our ability to model the effect of more continuous control schemes. Here we take a control theoretic approach to behavioral modeling and argue that movements can be understood as the output of a controlled dynamical system. Previously, models of movement dynamics, trained solely on behavioral data, have been effective in reproducing observed features of neural activity. These models addressed specific scenarios where animals were trained to execute particular movements upon receiving a prompt. In this study, we extend this approach to analyze the full natural locomotor repertoire of an animal: the zebrafish larva. Our findings demonstrate that this repertoire can be effectively generated through a sparse control signal driving a latent Recurrent Neural Network (RNN). Our model's learned latent space preserves key kinematic features and disentangles different categories of movements. To further interpret the latent dynamics, we used balanced model reduction to yield a simplified model. Collectively, our methods serve as a case study for interpretable system identification, and offer a novel framework for understanding neural activity in relation to movement.

Published as a conference paper at ICLR 2024
LEARNING INTERPRETABLE CONTROL INPUTS AND DYNAMICS UNDERLYING ANIMAL LOCOMOTION
Thomas Soares Mullen', Marine Schimel Guillaume
Hennequin', Suullauee Christian K. Machens
Michael B. Orger & Jouary
Department of Neuroscience, Champalimaud Foundation
2 Department of Engineering, University of Cambridge
E fadr j mi chael (rgger)rresearch fchampal imaud. org
ABSTRACT
A central objective in neuroscience is understand how the brain orchestrates movement. Recent advances in automated tracking technologies have made it possible document behavior with unprecedented temporal resolution and scale, generating rich datasets that can be exploited to gain insights into the neural con- trol of movement. One common approach is to identify stereotypical motor prim- itives using cluster analysis. However, this categorical description can limit our ability to model the effect of more continuous control schemes. Here, we take a control-theoretie approach to behavioral modeling and argue that movements be understood as the output of a controlled dynamical system. Previously, models of movement dynamics, trained solely behavioral data, have been effective in reproducing observed features of neural activity. These models addressed specific scenarios where animals were trained execute particular movements upon re- ceiving a prompt. In this study, we extend this approach to analyze the full natural locomotor repertoire of an animal: the zebrafish larva. Our findings demonstrate that this repertoire can be effectively generated through a sparse control signal driving a latent Recurrent Neural Network (RNN). Our model's learned latent space preserves key kinematic features and disentangles different categories of movements. To further interpret the latent dynamics, we used balanced model reduction yield simplified model. Lastly, we demonstrate the flexibility of our model by successfully applying it to the study of continuous locomotion another organism, C. elegans. Collectively, interpretable system identification, and offer a novel framework for understanding neural activity in relation movement.
methods serve case study for
INTRODUCTION
Over the past decade, several techniques have emerged to facilitate the high-throughput, automated tracking of animal behavior (Pereira et al.. 2020). This capability is crucial for various applications, from understanding genetics and neural activity (Pereira et al., 2020; Krakauer et al., 2017) to testing the effects of new pharmaceutical compounds (Wiltschko et al., 2020).
A commonly used approach to quantifying animal behavior has been to segment postural time series into a series of discrete actions. This approach has been used effectively represent the behavioral repertoires of species such as Drosophila (Berman et al., 2014), zebrafish (Marques et al., 2018) 8), and mice (Wiltschko et al., 2015). However, this discretization can obscure continuous variations in motor patterns and does not provide insights directly into the mechanisms that generate these actions. Instead of compartmentalizing behaviors into distinct categories, we suggest a novel approach: treat these observable movements as the output of underlying latent dynamical system related to motor function. This perspective reframes the study of motor actions as a system identification problem, aiming to uncover the latent control signals and dynamics that govern animal's full range of movements.
"These authors contributed equally this work.

---

Published as a conference paper at ICLR 2024
Prior work (Sussillo et al., 2015; Saxena et al., 2022) has shown that RNNs trained to replicate mus- cle activation in simple tasks, such as reaching towards a target. closely mimic features of concurrent brain activity. We aim to extend this line of work to natural behaviors, where the animal's actions are not constrained to a specific task, and the cues for movement may be unobserved or internally generated. To achieve this, we propose fit latent dynamical model that captures the movement repertoire of an animal, and, at the same time, to infer the unknown control signals that instigate or modify actions.
Here, our goal is thus to identify both the latent control signals and the underlying dynamics that make up the complete locomotor repertoire of an animal. We use iLQR ILQR-VAE (Schimel et al., 2022) to Iearn the latent RNN and infer the control signals from time series of behavioral observations. While we find that more expressive and overparametrized RNNs are easier to learn and provide better fits of the data, even when they are considerably higher-dimensions (Sedler et al., 2022), we additionally propose scheme to reduce our models to a minimal form, that is easier to dissect and relate to behavior. More specifically, we turn our focus to linear models of the data, and use balanced model reduction (Moore, 1981) to recover the minimal dynamical system that can still capture the system's input-output response.
than the ground truth
In the following, we first apply this framework to a system where rich behavioral data sets are avail- able: the zebrafish. This species is well suited to test our approach because its movement repertoire has already been systematically characterized (Marques et al., 2018 8; Mearns et al., 2020; Johnson al., 2020), offering a useful benchmark. We train family of RNNS and show that sparse trol signals activated at the onset of actions are sufficient for reproducing the observed postural sequences. These sparse control signals form a compressed representation that encompasses both discrete and continuous descriptions of movement. In the latent space, the trajectories associated with different movements are efficiently "untangled". Notably, this representation also yields spatial navigation information that was not explicitly provided in the training data. Finally, we show that the learned high-dimensional dynamics can be systematically reduced to a low-order, interpretable linear model that reveals the principal modes underlying behavior. Next, we apply the method to modeling continuous locomotion in C. elegans. There, sparse inferred control signals successfully match external stimuli applied to the animal, and model reduction reveals different patterns of crawl- ing and coiling.
We believe our investigation provides novel framework for studying movements, that will facilitate the understanding of neural activity that generates them, as well as aiding in the investigation of perturbations induced by pharmaceutical optogenetic treatments. Furthermore, our work serves as a valuable case study for interpretable system identification.
2 RELATED WORK
Below, we review existing methods that have been proposed to uncover the dynamics underlying behavior, before detailing how our approach differs from previous work.
Switching Linear Dynamical Systems In the context of behavioral data modeling, switching Linear Dynamical Systems (LDS) models have been particularly popular, and are at the core of some of the most widely used algorithms for behavioral segmentation (Datta, 2019). These models assume the existence of a hidden, discrete state that parameterizes the data-generating process, with different dynamics being used depending on the value of the hidden state. MoSeq (Datta, 2019) employs an autoregressive generative process coupled with sticky Hidden Markov Model. An alternative strategy employs adaptive segmentation algorithms that rely on statistical model testing instead of explicit state transitions (Costa al., 2019). Such models have been used successfully segment behavioral time series in various contexts (Johnson et al., 2016; Batty et al., 2019), and relating switches in behavior to concurrent neural recordings has yielded exciting insights into the neural control of behavior (Markowitz et al., 2023). However, potential limitation of switching LDS models is their reliance on multiple underlying dynamical systems, which complicates the interpretation of the learned dynamics.
Input-driven latent dynamics Another way of modeling behavioral dynamics while accounting for potential discontinuities due e.g. to descending commands from the brain, is to use input- latent dynamical systems models. However, sinultaneously learning system's dynamics and any
2

---

Published as a conference paper at ICLR 2024
unobserved input drive is very challenging, and often degenerate problem. Two methods, LFADS VAE (Schimel et al., 2022), have recently been proposed and shown to give promising results various synthetic and real datasets. These two methods for- mulate the learning problem as a variational autoencoder (VAE; Kingma & Welling, 2013) and approximate the posterior distribution over unobserved inputs with Gaussian distribution. They differ critically in the way the posterior over inputs is computed: in LFADS, the posterior mean is amortized using bidirectional recurrent neural network (bi-RNN), while iLQR-VAE uses iLQR (Li & Todorov, 2004), nonlinear control algorithm, to obtain the inputs most likely to have gen- erated the observations given the model parameters. Additionally, LFADS has only been used with choice of autoregressive prior over the inputs, while iLQQ- has been shown to work with spurse prior. While there have been very few examples of work modeling behavioral observations using such input-driven, latent dynamical models, we note that Wimalasena et al. (2022) used Au- toLFADS (Keshtkaran et al., 2022) to model musele activity in mice and monkeys. Importantly however, they made the assumption that the muscles were driven by continuous, autoregressive in- puts. The approach we present below differs from this by enforeing a sparse prior over the inputs, which encourages autonomy in the dynamics, and leads to a more interpretable control signal.
(Pandarinath et al., 2018) and -LQR-VAE
3 METHODS
3.1 SYSTEM IDENTIFICATION USING ILQR- -R--VAE
In this work, we modeled time series of behavioral observations o as the output of a latent dynamical system whose state # is being driven by set of control inputs u, i.c. as
(A) 2++1
fz
(R) U
(1)
(R) N(Cz, (R) E).
(2)
In Equation (1). o E F. is the temporal evolution of the behavior in the k-th set of observations, (k)
(R) Zt E I is latent representation of that behavior, and u, E Ih' are the control inputs driving
the latent dynamies. Note that we used inputs to set the initial state of the dynamics in the latent space, that is, = f0 uo. 0). Therefore, given set of parameters 0, the input sequence fully determines the latent trajectory z(u) {21, 22, .2. ......,,,,, and thus the corresponding output. this model, inputs therefore represent the true latent variables of interest. As behavioral datasets often consist of real-valued observations, we assumed Gaussian observation noise throughout our work, although other noise models could be used in other settings.
Learning the dynamics of a system driven by unobserved inputs is a challenging system identifica- tion problem. We tackled it using the recently proposed iLQR- VAE method (Schimel et al., 2022, see Appendix A) which optimizes a lower bound on the log marginal likelihood w.r.t. the model parameters 0 (which include any parameters in fo(.). C, Σ, and the prior parameters; see below). To facilitate identifiability, here we assumed that the dynamics were driven by sparse set of inputs. Specifically, we followed Schimel et al. (2022) and used a Student-t prior distribution (see details in Appendix 1. Because of its heavy-tailedness, this choice of prior is consistent with the model using small inputs most of the time, whilst still accommodating the occasional large input. It is thus well suited to inferring large sporadic inputs that may trigger e.g. switches in behavior, while encouraging the model to capture most behavioral segments via strong, near- aua-mouuomonous dynamics.
In iLQR- VAE, the posteriors over inputs are obtained using iLQR for each set of contiguous ob- servations; inference can be performed on any new time series not seen during training, even its length is different from the training data. After training, the model described above thus pro- vides us with a) a set of parameters 0 that can be dissected further, and b) a posterior distribution over the inputs, for every set of observations. Assuming the behavior can be generated with sparse, low-dimensional inputs, the posterior over inputs po (u (o) yields a compressed representation of the high-dimensional time series, which is typically more interpretable and easier visualize than the original data.
3.2 MINIMAL INTERPRETABLE DYNAMICS USING BALANCED MODEL REDUCTION
Alongside nonlinear RNNs, we also fitted the data using linear dynamics, i. modeled fo Equa-
tion (1) as
fo(2)
(K) uu (b)))
Az (k) + Bug (k)
(3)
3

---

Published as a conference paper at ICLR 2024
The linear model is characterized by its parameters (A E R= B E RTX", C ") However, the smallest dimension n required to accurately model the dynamics is unknown. One way of determining the appropriate size of the dynamical system involves fitting models of various sizes and comparing their generalization performance on test data. Here, we instead propose to fit a high-dimensional, potentially overparametrized linear model, and use balanced truncation (Moore, 1981; Antoulas, 2005) post-hoc to obtain reduced-order model that still accurately captures the transfer function of the original model.
Balanced truncation is performed in two stages; in a first step, the original dynamics (A, B, C) are transformed into an equivalent balanced system (A, B, C), i.c. dynamical system with the same transfer function as the original model, but in which each mode is controllable as it is observable. In a second stage, this balanced realization of the original system is truncated to a desired order r < n. Balancing the system in step relies on computing its observability and controllability Gramians. Those are positive-definie definite matrices defined respectively as
Wo-/ C'
(4)
We -
e Al BB" A (5)
They capture the information about the state space directions in which the system is most sensitive to inputs, and those in which it most likely to elicit outputs. More precisely, the top eigenmodes of W. correspond to the most observable directions of the system, i.e. directions which will lead to most variance in the output when the network is initialized along them. The top eigenmodes of W. are the most controllable directions of the linear system, i.e correspond to directions along which the dynamics can be steered at a minimal input energy cost (Kao & Hennequin, 2019). Importantly, given (A, B, C), W. and We can be computed in closed form as the solutions to two dual Lyapunov equations. Finding a balanced realization corresponds to finding a similarity transformation T such that, in the transformed system characterized by parameters {Á= = P-IAT, B = T-1B, = CT}, the resulting Gramians W. = T W.T and W. T W 17" become equal to each other (balanced) and diagonal. In other words,
W. W 2 = diag(oi, 00.0000.............._.... .2....,.
(6)
In Equation (6), o,, = 1, 1,.... n are the ordered, positive Hankel singular values (HSVs) of the system. These play the same role for the dynamical system as singular values do for constant, finite- size matrices. Thus, in the same way as the decay of the singular values can be used to choose the dimensionality cutoff' in principal component analysis (PCA), evaluating the HSV spectrum provides a principled way of deciding the truncation dimension r. Once r has been set, truncation of the system is achieved by selecting the first r X r block of Á, A the first r rows of B, B and the first r columns of C, C(r), The dynamics of the reduced dynamical system can then be summarized as follows:
c+1 fo(E(( u (1 = A
(k)
B (8) u
The removal of dimensions that are uncontrollable and unobservable provides simple, yet effective, way reduce the latent dynamics of the system to minimal form that can be interpreted and analyzed with more ease.
As a proof of concept, we validated the balanced model truncation's ability to recover the true dynamics of a system using synthetic example. In Appendix C, we generated data from set of lightly damped pendulums coupled by spring and driven by sparse forcing inputs. We modeled the pendulum dynamics as linear, and fitted 20-dimensional latent linear dynamical system to the data. The spectrum of the Hankel singular values exhibited a drop after modes corresponding to the true dimensionality of the data. Finally, truncating the dynamical system to 4 dimensions successfully revealed the modes corresponding to the underlying synthetic dataset (see Figure S1). This validation emphasizes how combining the iLQR- LLQR-VA with balanced model reduction can yield a minimal latent dynamical system that accurately describes the data.
RESULTS
4.1 MODELING ZEBRAFISH BIIIHHOOOOR LOCOMOTION
Next, we turned to modeling a dataset of behavioral observations in larval zebrafish. This case study highlights how modeling approach enables to (i) gain insights into the latent space underlying
4

---

Published as a conference paper at ICLR 2024
this complex behavior, and (ii) extract minimal model of the underlying dynamies, which is more amenable to further analyses.
4.1.1
ZEBRAFISH BEHAVIORAL DATASET
A.
B.
C.
-M
-m
Tail angle [mod
20 Latent dynamic
M
Ann An
Ann
Mr
0.0
segments
Mit Motur command inpurs
Thil kinematic
50 ms
iLQR Inference
Time [116]
Figure 1: Overview of the model setup for zebrafish Locomotion (A) Tail angle tracked along 7 tail segments (top). Example swim bout captured by displacement of tail angles along tail segments through time displays wave propagating along the tail. (B) Illustration of the model setup. We learn a generative model from inputs to observations, and use iLQR- VAE to perform inference in this model. (C) Illustration of the behavioral repertoire, one movement is shown for each category.
The larval zebrafish is a model organism studied in thousands of biological research labs, owing to its optical transparency, rapid developmental cycle, and the extensive genetic library available (Mrinalini et al., 2022). Zebrafish larvac swim using discrete episodes of propulsion, known as swim bouts, which typically last around 200ms. The body shape of the zebrafish can be characterized by measuring the bending angle as a function of distance along the length of its body (Figure IA). We compiled a dataset of 30800 tail angle times series classified according the known repertoire of zebrafish larva (Marques et al. (2018), details in Appendix G).
4.1.2 LOW-DIMENSIONAL, SPARSE CONTROL OF LARVAL ZEBRAFISUI LOCOMOTION
To characterize zebrafish locomotion, we sought to learn dynamical systems that were able to gen- erate its locomotor repertoire (see Figure 1B).
Model selection To find the models that best described the data, we varied the dimensions of the control signal and of the latent size, as well as the RNN architecture. We considered both MGU RNNs (Zhou et al., 2016; details in Appendix B.2), and Linear Dynamical Systems (LDS). We evaluated models based on both (i) their ability to reconstruct the data (fraction R² variance explained by the mean of the posterior predictive distribution over observations), and (ii) the sparsity of the inferred control signals. We quantified temporal sparsity by the ratio of the L L'-ro to the norm. Therefore, a sparsity score of corresponds to an input signal that non-zero time point, resembling an impulse function. We found that all trained models could generate accurate reconstructions of our dataset (R2 > 0.94; Table S1). The output signal predicted by the dynamical system's 's flow was smoother than the behavioral observation, i.c. our method acted as a denoising autoencoder. The high reconstruction fidelity can be attributed to the inherently low-dimensional nature of behavioral recordings. Yet, critical concern was ensuring the control signal captured the genuine dynamies rather than simply mirroring zebrafish posture. Thus, we also evaluated models on the sparsity of control. We found that matching the observations using very sparse control was challenging, and possible only for the largest and most expressive models. Both an LDS and an MGU network of latent size 120 driven by = 10 inputs could accurately reconstruct the zebrafish behavior while using sparse inputs (Figure 2A). As we detail below, however, the MGU network was able to compress the behavioral information more effectively. For dissection of the latent controls driving the system, we thus turned our focus to this network.
single
Sparse control signal and latent state representation Our assumption regarding the sparsity of the control signal was grounded in biological observations. Indeed, sparse electrical activation of brainstem neurons projecting to the spinal cord has been shown to be sufficient to elicit forward locomotion escape (Severi et al., 2014; Xu et al., 2021). In addition, Markov et al. (2021) estab- lished that zebrafish bouts are modified by unexpected visual feedback only after 220ms, with the intial portion of the bout proceeding in ballistic fashion after the onset.
5

---

Published as a conference paper at ICLR 2024
A.
Tall ingle
MGU
[mdd]
Tail kng
100
200
LDS
Time [ms]]
100
200
0
100
200
100
200
Nomsalzze contol
Time [mss]
Time [ms)]
Time [ms)]
0.5
D.
energ
Posture PCA
MGU latent PCA
Time
100
Classificaiion
securncy [%]
50
50
100
PC3
PC
PC
Time (mss)
possure time LFADS
posture LUS
shuffled MGU
Figure 2: Sparse inferred inputs generate range of swim bouts. (A) Example bout reconstruc- tion (R² 0.97 between the observation O and the reconstruction Ô) with a sparse input driving the latent trajectories for MGU and LDS (B) Mcan + standard deviation across dataset of the normalized Lo-morm of aligned on the first peak. MGU RNN (black) learns a more transient input to generate movements across the dataset compared LDS (gray). The alignment the first peak the same is B., C. and D. (C) Performance of linear classifier to predict the movement category using latent posture data (see Appendix D.3) (D) PCA projections of the posture time series (left) and latent state (right) color-coded according to bout category.
Here, the model's multivariate control impulse arose right before the onset of movement (which started at different time points in different trials), thus setting the initial conditions for the latent dynamical system (Figure 2B). After this point, the bout was then generated by quasi-autonomous dynamics, before decaying back a fixed point. As shown in Figure 2A, the MGU could reconstruct bouts using a sparser control compared to the LDS.
Most of the information generate movement was contained in this diw-dimnnsionll impulse. Indeed, restricting the input to its initial peak was still sufficient to reconstruct bouts with R² = 0.82 across the test dataset. Although a single time-step impulse lasting for single time bin (1.42ms) might appear unrealistically fast from biological standpoint, our model was robust to changes in the duration the control signal up to 10ms (Figure S2a). During this impulse signal, we found that, rather than using one-hot encoding of the movement category, the different control units were continuously modulated different ratios to generate specific bout types (see Figure S3). The state-space position at the onset of movement acts compressed representation of the forthcom- ing actions. Here, five principal components explained 95% of the variance at movement onset. This representation not only segregates known categories of movements but also retains continuous variations in kinematic features like frequency and amplitude (Figure S5). We found that the latent representation was robust across random initialization and splits of the training data, even those that excluded bout categories (see Figure S2B). Following this initial state, the movement unfolds by following the learned flow field. In the absence of an additional control signal, the trajectory state space should therefore be untangled, with similar positions in state space leading to similar patterns in the near future (Russo et al., 2018). We tested this in the models, by measuring how well we could decode the category of movement from a single snapshot of the latent state (Figure 2C). In contrast to postural trajectories, which were highly tangled, the low-tangling of the latent state space allowed for high classification accuracy from that state. Surprisingly, the state space of the MGU 40ms after the control peak provided a higher classification accuracy even when compared with a linear classifier trained on the full time series of tail movements Figure 2C. This separation between

---

Published as a conference paper at ICLR 2024
A.
B.
18
Tail angle [rad]]
300
0.2
Frequency (Hb)
200
Reduged Fall
Inpor
0.0
H1V
Toll mil
100
[
-0.2
-200
100
Time [m0e]
Time
Time [me]]
100
D.
Index
Sum of2 modes
t 16.5ts
f- DEL
53 Stns
Forward
Tall Angle
0.5
0.5
Tall Angle
40
W
-0.3
-0.5
Input [sal
Froguancy (Hd)
100
t
Tall Ande
Thans swriim
Tal Angle
20
Topus Boaj
Time
Time [me]]
300
Time (1mte)
100
300 Time (mm]
-80
60
[ms)
Figure 3: Linear model reduction used to obtain compositional structures of behavioral reper- toire. (A) Hankel singular values of the system partitioned at r = 18 (red vertical line). (B) Input signal projected to full scale network A (top) and reduced network Á (r) (bottom) able to recon- struct swim R² = 0.818. (C) Eigen spectrum of balanced truncated dynamics (r)()) compared to the full state dynamics X(A). (D) Top two eigenmodes that capture the most variance are con- volved through the input signal construct forward swim (top) and turn swim with R 0.76 and R² = 0.74, respectively. The oscillatory and decay mode are highlighted by a square and triangle in the truncated dynamics spectrum (right).
categories of movement cannot be attributed a difference in the dimensionality of the latent space (see Appendix D.3). It demonstrates the quality of the representation of the movement within the latent space. Spreading of target variables in the neural state space has been shown be associated with robustness to noise (Russo et al., 2018; Susman et al., 2021).
To benchmark our method, we compared with the non-autonomous LFADS (Pandarinath et al., 2018) with inferred inputs and an autoregressive prior AR(1) for the control signal (see Appendix E). The method successfully reconstructed the postural observation (R² = 0.94). However, in this predominantly input-driven. Indeed, we could regress
regime, the latent dynamical system projection of the control signal to display a correlation of r= 0.73 with the postural time series (see Figure S6). This leakage from the input data resulted in less informative dynamical system, with a worse classification accuracy of the movement category from the latent space Figure 2C. This result suggests that for low-dimensional behavioral observation, sparse input priors better suited at learning a meaningful dynamical system.
Control of spatial navigation The control signal accounts for the influence of sensory and decision-making areas. As such, this signal is expected to encode information relative to naviga- tional landmarks such as the relative position of a prey that the fish want to capture. To achieve its goal, the fish has to solve inverse problem: determining what control signal will propel its body in given position in space.
We propose that straightforward mapping exists between the control impulse and the ensuing fish trajectory. This proposal is non-trivial as our training dataset contained only posture information. Nonetheless, Figure S4 demonstrates that spatial displacement be linearly predicted from the 120-dimensional initial latent state. The control signal and dynamic inferred by our method therefore allow for simple sensorimotor coupling.
7

---

Published as a conference paper at ICLR 2024
4.1.3 DISSECTING THE DYNAMICS
In the nonlinear dynamical system described above, input sparsity aids in dissecting the control mechanisms. Here, we transition to a linear dynamical system (LDS) to gain insights into the dy- namics. Although L.DS models did not match the performance of their nonlinear counterparts, their dynamics could be more easily interpreted. This illustrates a common tradcoff between the use of more expressive models which can lead to more compressed representations, versus simpler mod- els which can be more easily interpreted. We note that both types of models are accommodated by our framework, and can be switched between depending on the focus of the application.
Using -LQR-VAE, we fitted an LDS model using 120 latent dimensions and 10 control inputs. We truncated the system to the minimal dimension after which the decay in Hankel Singular Values started plateauing (see Figure 3). Here, this corresponded to 18 eigenvalues, and 10 distinct eigen- modes, corresponding to temporal decays of up to 100 ms and frequencies up to 45 Hz ( Figure 3c). This reduced model could accurately reproduce the output across range of frequencies, given the sparse forcing inputs given to the original model (Figure 3B). Note that, as we show in Table SI, was not possible to directly train models of such small size that were both sparse and yielded good reconstruction accuracy. In particular, when directly training linear models of smaller sizes, we were unable to find solutions that relied on sparse control inputs. Thus, model reduction was necessary step to reach the smaller size models we needed to further dissect the dynamics. As we can see in Figure 3C, reducing the model corresponded to removing many eigenvalues corresponding high frequencies and short time scales while adding cigenvalues with longer time scales and lower fre- quencies. The reduced system modes are therefore closer to the temporal and spectral components critical to generate the movement.
With our reduced model in hand, we were then able decompose postural time series into the contri- butions of the individual modes. We focused on forward movements and tumns, shown in Figure 3D (respectively top and bottom). Interestingly, turn movements were generated by the superposition of oscillatory (corresponding to complex conjugate eigenvalues) and non-oscillatory modes (with real- valued cigenvalue), while forward movements primarily relied on oscillatory modes (see Figure S7). This aligns with the findings from Huang et al. (2013), who identified populations of reticulospinal neurons essential for converting forward swims into turns. These neurons could be implicated in triggering the non-oscillatory modes observe. This highlights how our model may facilitate the dissection of neural circuits responsible for generating movements, and help investigate how neural ablations or perturbations may relate specific dynamic modes.
4.2 MODELING LOCOMOTION OF C. ELEGANS
To map the relevance of inferred control signals driving behavior requires controlled activation of ground truth external stimulus. Here, we switched animal models, using a dataset from Broekmans et al. (2016) to study the input-driven dynamics of the C. elegans response to 100ms aversive heat shock to the head. In response, C. elegans exhibited a robust escape sequence of behavioral maneu- vers: reverse crawl, high-angle Q-turn to change heading direction, the forward crawl (Figure 4B). Unlike the larval zebrafish, the continuous locomotion of the C. elegans is captured as mixture of principal modes, termed eigen worms (EW), to capture its 2D body posture ( FFigure 4A).
We trained a LDS of network size = 60 driven by m 5 inputs which could accurately recover the continuous dynamics (Figure 4C). Across all 90 trials, a salient input peak was aligned on the heat shock activating the reverse crawl (Figure 4D). Moreover, there is secondary peak driving the Q-turn changing the heading direction away from the stimuli. While driving inputs were less sparse due to the continuous locomotive nature of the C. elegans locomotion, large inputs were recruited in the transition between behavioral maneuvers.
The choice of LDS model allowed for balanced model truncation, (Figure S8), where we interpreted how the task-related inputs influenced the dynamics of behavioral switches in the escape sequence. We projected the three slowest decaying eigenmodes of the reduced model into the EW space. Each mode contribution reflected specific behavioral maneuvers, shown in Figure S8C. The frequency of modes that were recruited in forward and reversal crawls aligned with results found in Costa et al. (2023). The reduced LDS provides an alternative approach to segment behavioral states in continuous locomotion and infer switches between behavioral modes.
8

---

Published as a conference paper at ICLR 2024
C.
-j
ag
Forward
Reverse
2 fum
Forward
EW 2
EW
EW
x
a4
B.
Bestabock
Revers craid
Qumm
D.
Normalized inpart (a:0
0.5
Forwaard crawe
10
Time 20
25
Figure 4: iLQR- captures C. elegans continuous locomotion and infers heat shock-induced behavioral transitions. (A) C. elegans 2D body posture was expressed as a combination of eigen-worm (EW) postures (95% explained variance). (B) Time-lapse of C. elegans body posture following an aversive heat shock that triggers behavioral sequence of reverse crawling, high angled Q-tum, and forward crawling away from the stimuli. (C) (From top to bottom) example trial of the postures around the heat shock (red marker). The inferred driving inputs align with the heat shock 10 s, with secondary burst of inputs transitioning the reversal mode reconstruction of the posture from the latent dynamics Ô matches the observation O. (D) Mean = standard deviation across all escape trials of the normalized L2-norm of the control is aligned to the external heat shock.
U-turn. The
5 DISCUSSION
In this work, we introduced a novel framework to learn the dynamics underlying the generation of movement from unstructured postural datasets. Our approach simultaneously identifies the control signals that drive actions and the dynamical systems that generate them. Our empirical findings revealed that temporally sparse control signals are sufficient to replicate the complete locomotor repertoire of larval zebrafish. In dataset of C. elegans continuous crawling recordings, could infer control signals that were aligned with external perturbations. While high-dimensional expressive dynamics were necessary to model the data accurately, proposed to balanced truncation to extract a minimal, interpretable linear dynamical system model of the data. In both modelled species, our approach yielded insights consistent with prior findings.
model
While balanced truncation is a well-established method in the control literature, it had not yet been combined with data-driven approaches reduction, one must have access to all the parameters of a controlled dynamical system. If the inputs the system are unobserved (as is often the case behavioral datasets), thus necessary to take the approach we did, and propose constraints on the inputs that allow to infer the most likely inputs have generated the data.
learning behavioral dynamics. Indeed, to apply balanced
A limitation of our work is that we relied on two separate dynamical systems to be able to extract interpretable inputs and interpretable dynamics from the data. This could be handled in a more unified manner in the future, for instance by obtaining the LDS by linearization (Sussillo & Barak, 2013), or by co-training the linear and MGU models, as in Smith et (2021).
Overall, separating out the contributions of intrinsic dynamics and control signals that shape an- imal movements can be helpful for interpreting neural recordings, and provide insight into how brains control movement. Our model could be used to estimate how either control signals, or out- put dynamics, change in different behavioral contexts, or pathological states, and make mechanistic predictions testable with brain imaging studies. In subsequent research, adapted to other animal behavior, such as limb-based locomotion in flies or mice. Alternative con- alternative to temporal sparsity.
framework could be
trol priors could also be investigated

---

Published as a conference paper at ICLR 2024
ACKNOWLEDGMENTS
We would like thank II Memming Park for insightful discussions the manuscript. T.S.M re- ceived funding from the European Union's 's Horizon 2020 research and innovation program under the Marie Skłodowska Currrr grant agreement 813457. C.K.M. acknowledges support from the Simons Collaboration on the Global Brain (543009); C.K.M and A.J. received support from Fundação para CTTTT//CBBIL-OTT/337777777---A-
a Ciência Tecnologia
M.B.O received support from the European Research Council (ERC NEUROFISH 773012) and Volkswagen Stiftung "Life?" initiative.
was supported by a Fellowship from the La Caixa Foun- dation (LCF/BQ/PR20/1 17700007) M.S. was funded by EPSRC DTP studentship. This work was performed using resources provided by the Cambridge Service for Data Driven Discovery (CSD3) operated by the University of Cambridge Research Computing Service (www.csd3.cg csd3.cam.ac.uk) pro- vided by Dell EMC and Intel using Tier-2 funding from the Engineering and Physical Sciences Research Council (capital grant EPPT022159/1/ and DiRAC funding from the Science and Tech- nology Facilities Council www.dirac.ac.oo Laborde for their help collecting the zebrafish dataset. Fish husbandry was supported by the Champalimaud Foundation Fish Facility, which received support from the research infrastructure CONGENTO, co-financed by the Lisboa Regional Operational Programme (Lisboa2020), under the PORTUGAL 2020 Partnership Agreement, through the European Regional Development Fund (ERDF), and the Fundação para a Ciência Tecnologia (Portugal) under the project LISBOA-01- 0145-FEDER-022170.
We thank Pedro Tomás Martins Silva and Alexandre
10

---

Published as a conference paper at ICLR 2024
REFERENCES
Athanasios C Antoulas. Approximation of large-scale dynamical systems. SIAM, 2005.
Eleanor Batty, Matthew Whiteway, Shreya Saxena, Dan Biderman, Taiga Abe, Simon Musall, Winthrop Gillis, Jeffrey Markowitz, Anne Churchland, John P Cunningham, et al. BehaveNet: nonlinear embedding and Bayesian neural decoding of behavioral videos. Advances in Neural Information Processing Systems, 32, 2019.
Gordon Berman, Daniel M Choi, William Bialek, and Joshua W Shaevitz. Mapping the stereotyped behaviour of freely moving fruit flies. Journal of The Royal Society Interface, 11(99):20 1066722 2014.
Onno D. Broekmans, Jarlath B. Rodgers, William S. Ryu, and Greg J. Stephens. Data from: Re- solving coiled shapes reveals new reorientation behaviors in c. elegans, 2016. URL https: / //datadryyd.cc wttpoo//ttt//aaataeaeeeeeeeeeeee
1066:111111111111 dryad. 0m6p.
Antonio C Costa, Tosif Ahamed, and Greg Stephens. Adaptive, locally linear models of complex dynamics. Proceedings of the National Academy of Sciences, 116(5):101-11500, 2019.
Antonio C. Costa, Tosif Ahamed, David Jordan, and Greg Stephens. Maximally predictive states: From partial observations to long timescales. Chaos: An Interdisciplinary Journal of Nonlinear Science, 33(2):02313,, February 2023. doi: 10.1063/5.0129398. URL https: /doi .rr/ 10. 0010663//66 0129398.
Sandeep Robert Datta. Q & A: Understanding the composition of behavior. BMC Biology, 17 (1), May 2019, doi: 10.1 (0.1188611911-1990000 s1 s12915-00 01 9- 9-0663-3.
URL https:/ http:://dd....gg////
Joel Heck and Fathi M. Salem. Simplified minimal gated unit variations for recurrent neural networks. 20/7 IEEE 60th International Midwest Symposium on Circuits and Systems (MWS- CAS). pp. 1593-1596, 2017. URL https: api. swmanticscholar........
control turning behaviors in zebrafish. Current Biology, http:p/://.gggeennnncccooooo.... org/ Corg/oorpsss:
14334496.
Kuo-Hua Huang, Misha B Ahrens, Timothy W Dunn, and Florian Engert. Spinal projection neurons 2013.
Matthew Johnson, David K Duvenaud, Alex Wiltschko, Ryan Adams, and Sandeep Datta. Composing graphical models with neural networks for structured representations and fast infer- ence. Advances in neural information processing systems, 29, 2016
Robert Evan Johnson, Scott Linderman, Thomas Panier, Caroline Lei Wee, Erin Song, Kris- tian Joseph Herrera, Andrew Miller, and Florian Engert. Probabilistic models larval zebrafish behavior reveal structure on many scales. Current Biology, 30(1):70-82.e4 January 2020. doi: 10.1016/j.cub.2019.11 11.02.. URL https:tt//// //doi./ooiio// org/ 10. 1016/ / j. cub. 2019 1.11... 11. 0266.
Ta-Chu Kao and Guillaume Hennequin. Neuroscience out of control: control-theoretic perspec- tives on neural circuit dynamics. Current Opinion in Neurobiology, 58:122-129, October 2019. doi: 10. 00.1016j.comb.2119.000 09. 001.
URL https: /tts:://ooiiie ortp:0..0001011111611 conb. 2019.
Mohammad Reza Keshtkaran, Andrew R. Sedler, Raeed H. Chowdhury, Raghav Tandon, Diya Bas- rai, Sarah L. Nguyen, Hansem Sohn, Mehrdad Jazayeri, Lee E. Miller, and Chethan Pandar- inath. A large-scale neural network training framework for generalized estimation of single- trial population dynamics. bioRxiv, 2022. doi: 10. 10.1100/2001.01.13342000 / //www. biorxiv. wwwpo//////o/////o/oooooooeooe/e hatpm//eaaaarr//2eme/ ///1///22222000///222122 11.13.1...... 426570.
URL https:
Diederik Kingma and Max Welling. Auto-encoding variational Bayes.
arXiv preprint
arXiv: 14,661114000 2013.
Diederik P. Kingma, Tim Salimans, and Max Welling. Variational dropout and the local reparame- terization trick, 2015.
11

---

Published as a conference paper at ICLR 2024
John W Krakauer, Asif A Ghazanfar, Alex Gomez- Malcolm A MacIver, and David Poeppel. Neuroscience needs behavior: correcting a reductionist bias. Neuron, 93(3):480-499) 2017.
A. Laub, M. Heath, C. Paige, and R. Ward. Computation of system balancing transformations and other applications of simultaneous diagonalization algorithms. IEEE Transactions Automatic Control, ( February 1987. doi: 10.1 1099/a// 1909///...98..111111 104549. URL https: /doi. org/10, 1110/ttaeee 1987. 1104549.
Weiwei Li and Emanuel Todorov. Iterative linear quadratic regulator design for nonlinear biological movement systems. In Helder Araújo, Alves Vieira, José Braz, Bruno Encarnação, and Marina Carvalho (eds.), ICINCO (I), pp. 222-229. INSTICC Press, 2004. ISBN 972-8865-12444 URL http:///bl/ uni memi-trier... de/db////// conf www.ccccciininnn0/..... htm #LiT04.
Daniil A Markov, Luigi Petrucco, Andreas M Kist, and Ruben Portugues. A cerebellar internal model calibrates a feedback controller involved in sensorimotor control. Nature Communications, 12(1):6694, 2021.
Jeffrey E Markowitz, Winthrop F Gillis, Maya Jay, Jeffrey Wood, Ryley W Harris, Robert Cieszkowski, Rebecca Scott, David Brann, Dorothy Koveal, Tomasz Kula, et al. Spontaneous behaviour structured by reinforcement without explicit reward. Nature, 61477946):00--117, 2023.
João C. Marques, Simone Lackner, Rita Félix, and Michael B. Orger. Structure of the zebrafish locomotor repertoire revealed with unsupervised behavioral clustering. Current Biology, 28(2): 181-195.e5, January 2018. doi: 10.1016j..ub.2011.2..000 URL https://doi. wwg/oo.gg/o.ooo /016/j. cub. 2017 12.12.... 002.
Dunean Mearns, Joseph C. Donovan, António Fernandes, Julia Semmelhack, and Herwig Baier. Deconstructing hunting behavior reveals tightly coupled stimulus-respone loop. Current Biology, 30(1):54-69.e9 January 2020. doi: 10.10066j.cub.2199..... URL httpt:/// org/ htg/10.. 10. 101.10166/// (/.0066jj-. cub. 2019. 11 .022.
B. Moore. Principal component analysis in linear systems: Controllability, observability, and model reduction. IEEE Transactions on Automatic Control, 26(1 1673), February 1981, doi: 10.1
1ac.1191.1 100566.. URL https: http//o.oroopo 10. 1109 /tac.taccco. 1981. 1102568.
Mrinalini. T Tamilanban, V Naveen Kumar, and K Manasa. Zebrafish-the neurobehavioural model in trend. Neuroscience, 2022.
Chethan Pandarinath, Daniel J. O'Shea, Jasmine Collins, Rafal Jozefowicz, Sergey D. Stavisky, Jonathan C. Kao, Eric M. Trautmann, Matthew T. Kaufman, Stephen I. Ryu, Leigh Hochberg, Jaimie M. Henderson, Krishna V. Shenoy, I.. F. Abbott, and David Sussillo. Inferring single-trial neural population dynamies using sequential auto-encoders. Nature Methods, 15(10)88058155 September 2018. doi: 10. 10.10338/441522-18-00099 s41592- 018 1008800100 0109-9. 9.
URL http:::/tc...ogg////po
Talmo D Pereira. Joshua Shaevitz, and Mala Murthy. Quantifying behavior to understand the brain. Nature neuroscience, 23(12):1537-1549 2020.
Maithra Raghu, Justin Gilmer, Jason Yosinski, and Jascha Sohl-Dickstein. Sveca: Singular vector canonical correlation analysis for deep learning dynamics and interpretability. Advances in neural information processing systems, 30, 2017.
Abigail Russo, Sean R. Bittner, Sean M. Perkins, Jeffrey S. Seely, Brian M. London, Antonio H. Lara, Andrew Miri, Najja J. Marshall, Adam Kohn, Thomas M. Jessell, Laurence F. Abbott, John Cunningham, and Mark M. Churchland. Motor cortex embeds muscle-like commands in an untangled population response. Newron, 97(4):953396666 February 2018. doi: 10. 1016// 10t.on.211..0..000 URL https: /ttp:///oiiieee.t http//ii.ogg/.//////// wtt.://.././..///..... neuron. 2018 01.004.
Shreya Saxena, Abigail A Russo, John Cunningham, and Mark M Churchland. Motor cortex activity across movement speeds is predicted by network-level strategies for generating muscle activity. eLife, 11:e67620, may 2022. ISSN 2050-084X. doi: 10.7554/eLife URL https:/ doi doi.org/10. 7554 4/ / fe. 67620.
12

---

Published as a conference paper at ICLR 2024
Marine Schimel, Ta-Chu Kao, Kristopher T Jensen, and Guillaume Hennequin. iLQR VAE: control- based learning of input-driven dynamics with applications to neural data. In International Confer- ence Learning Representations, 2022. URL https: http:///ppeeeeee.eeeeeeeeeeeeeee WRODLDHaAi..
Andrew R Sedler, Christopher Versteeg, and Chethan Pandarinath. Expressive architectures enhance interpretability of dynamics-based neural population models. arXiv preprint arXiv:2212.03771, 2022.
Kristen E. Severi, Ruben Portugues, João C. Marques, Donald M. O'Malley, Michael B. Orger, and Florian Engert. Neural control and modulation of swimming speed the larval zebrafish. Neuron, 83(3):692-707, August 2014. doi: 10. http://.eeroo.211...00. 10.1016/ wtw.:o//...uuonn..... www.on...n..... 2014 16.003....
URL https: /oi.. org/
Jimmy Smith, Scott Linderman, and David Sussillo. Reverse engineering recurrent neural networks with jacobian switching linear dynamical systems. Advances in Neural Information Processing Systems, 11116600066713 2021.
Lee Susman, Francesca Mastrogiuseppe, Naama Brenner, and Omri Barak. Quality of internal rep- resentation shapes learning performance feedback neural networks. Physical Review Research, 3(1):013176, 2021.
David Sussillo and Omri Barak. Opening the black box: Low-dimensional dynamics in high- dimensional recurrent neural networks. Neural Computation, 25(3):626-649) 649, March 2013. doi: 10.1 2,1116///cooaa.04400 URL httph:://doooo www.orgggoo.ooo....oo 100 11 62 httr//ncoo___0000
David Sussillo, Mark M Churchland, Matthew T Kaufman, and Krishna Shenoy. neural network that finds naturalistic solution for the production of musele activity. Nature Neuroscience, 18 (7):1025-100 1033, June 2015. doi: 10. 10.1008/mn..00 URL https: doi. wrw/oorgg......o...ooooooo 1038/ 4042.
Alexander B Wiltschko, Matthew Johnson, Giuliano Iurilli, Ralph E Peterson, Jesse M Katon, Stan Pashkovski, Victoria E Abraira, Ryan P Adams, and Sandeep Robert Datta. Mapping sub-second structure in mouse behavior. Neuron,
2015.
Alexander B Wiltschko, Tatsuya Tsukahara, Ayman Zeine, Rockwell Anyoha, Winthrop Gillis, Jeffrey E Markowitz, Ralph E Peterson, Jesse Katon, Matthew J Johnson, and Sandeep Robert Datta. Revealing the structure of pharmacobehavioral space through motion sequencing. Nature neuroscience, 23(1 (13(1)):43331111 1443, 2020.
Lahiru N Wimalasena, Jonas F Braun, Mohammad Reza Keshtkaran, David Hofmann, Juan Álvaro Gallego, Cristiano Alessandro, Matthew Tresch, Lee Miller, and Chethan Pandarinath. Esti- mating muscle activation from emg using deep learning-based dynamical systems models. Jour- nal of neural engineering, 19(3):036013, 2022.
Lulu Xu, Na N. Guan, Chun-Xiao Huang, Yunfeng Hua, and Jianren Song. A neuronal circuit that generates the temporal motor sequence for the defensive response in zebrafish larvae. Current Biology, 31(15):3343-3357.c
August 2021. doi: 10. 10.1016/j.uub.0021.06 202..00000.22.00.54..... URL https:
//doi. www.oor./0.............1 1016/
j. cub. 2021 054.
Guo-Bing Zhou, Jianxin Wu, Chen-Lin Zhang, and Zhi-Hua Zhou. Minimal gated unit for recurrent neural networks. International Journal of Automation and Computing, 13(3):226-25) 2344, 2016.
13

---

Published as a conference paper at ICLR 2024
Appendix
A SUMMARY OF THE ILQR-VA METHOD
iLQR-VAE is recently proposed method to learn input-driven dynamical systems from data. In this section, we provide a brief description of the various components of the methods. More details can be found in Schimel al. (2022). LLQR-VAE models time series observations o as arising from an input-driven nonlinear latent dynamical system. The generative models considering in iL.QR VAE have the following form
latent state observations 0 (o)))(z)).
21+1 = fe( Ue,
(S1)
(S2)
where z denote latent trajectories, u the inputs driving those trajectories, and are model parame- ters. Inputs in iLQR- VAE are used capture process noise, ics, and to model control inputs driving the dynamics. The parameters are learned by optimizing the evidence lower-bound (EL BBO) (Kingma Welling, 2013), a lower bound on the log- of the observations given the model parameters, (o). Computing the ELBO requires introducing an approximate posterior distribution q(ulo). The ELBO objective then takes the following form
the initial condition of the dynam-
((0,0,0)
Σ
E go (a)
Σ log Po (0 (z))
(z) + logpo(u;) log (
(S3)
1o OO)
(S4)
with respect to both 0 and d. In iLQR-VAE, the approximate posterior over inputs for a set of observations (*)) defined as Gaussian distribution, centered at the maximum of the true posterior over inputs,
" (x) argmax log (x)
(S5)
This set of inputs can be obtained the solution of an optimization problem, which solved here using the iLQR algorithm (Li & Todorov, 2004). iLQR-VA PAE then uses a shared set parameters to model the posterior covariance. Those consist of of spatial covariance parameters E,, and set of temporal covariance parameters Σ.). The full covariance is then defined as the Kronecker product of those terms.
Given this parameterization of the approximate posterior, iL-QR-VAE then optimizes the ELBO by drawing samples from (()))(k)) and using the reparameterization trick (Kingma et al., 2015) to obtain gradients.
B DETAILS OF THE GENERATIVE MODEL
B. PRIOR DISTRIBUTION
We used a Student-t distribution as our prior over inputs. This took the following form
[[ + m)/2 T [2/2] //
(v4m)
pa (u+)
+ - - u
(S6)
where diag(ss...... 5.). The Student-t prior Equation (S6) arises from the following con- struction; every point in time, the scale of the input variance d,, is sampled from x distribution with y degrees of freedom, i.e Gg ~ X (v) where larger values of y shift the distribution towards a larger mean. In second step, inputs scaled by ous ie u AO(0, S). The advantage of this construction is that promotes the use of small inputs most of the time, whilst still accommodating the occasional large input.
sampled from Gaussian distribution whose variance S
14

---

Published as a conference paper at ICLR 2024
B.2 DYNAMICS
To fit the larval zebrafish data, we selected Minimal Gated Unit (MGU) RNN presented by Heck & Salem (2017). The dynamics are defined as,
ft - a (
(S7)
â. (UU) 1 + Wu:
(S8)
Ze = (1 ft) 21-1 + ft &
(S9)
where the non-linearity functions 0 and & were chosen to be tanh functions.
C DETAILS OF THE PENDULUM EXAMPLE
We generated observations from synthetic system consisting of two driven pendulums. The equa- tions describing the autonomous dynamics of the system are as follows:
mři -g/ g/ sin 01 + k02 fāz
(S10)
2 -9/ L0 + k02 f0}
(S11)
m8² -9/ / sin 02 k0ı f6z
(S12)
g/ L02 + k0 fÚ2,
(S13)
where 81 and 82 are the angles of the pendulums relative to their equilibrium position. In Equa- tion (S10), m 1kg is the pendulum mass, g 9.8ms 2 is the gravitational acceleration term, L = Im the length of the pendulum, the pendulums, and f the damping term. As the angles of the pendulums remained small, ap- proximated the dynamics as linear.
160Nm the strength of the coupling spring in-between
We modeled the effect of inputs forcing switches in the pendulum motion by integrating the system's position from random initial conditions, and stitching together several trajectories. In both our training and test sets, the dynamies were driven by either I or 2 inputs. We sampled observations at time intervals of dt = 0.06s. Note that even though we only had two pendulums, the differential equation underlying the dynamics of order two, such that modeling the dynamics of this system requires 4-dimensional linear dynamical system.
To recover the coupling constant from the frequencies learned by the model, we used the fact that the two oscillatory modes have frequencies of wp and W1 v w + 2k: where wo = g/ L. As such, (w) w22
C.I RECOVERING THE DYNAMICS OF DRIVEN COUPLED PENDULUMS
To illustrate the use of our approach, we first validated its ability to recover the true dynamics of a system synthetic example. We generated data from toy dynamical system consisting of a set of two lightly damped pendulums, coupled by a spring, and driven by sparse forcing inputs (Figure $1). We modeled the pendulums linear systems, as they remained in the small angle regime throughout the experiment. The state trajectories of these coupled pendulums (which were 4-dimensional) were then projected onto lower-dimensional (2D) space and corrupted with Gaussian noise give rise to a set of observations. The equations governing the system's 's motion are described in Equation (S10).
We fitted the observations using a 20-dimensional latent linear dynamical system. Computing the balanced realization of this model revealed the existence of 4 dominant modes in the dynamics (see Figure S1B). Reducing the dynamics eigenvalues, corresponding to two oscillatory modes (Figure S1C). From these eigenvalues, we were able to recover the coupling constant of the underlying synthetic system highly accurately (see details in Appendix C), with the true coupling constant & = 160 Nm falling within the error of our prediction, kpred = (160 + # 1) across random initializations and data splits. This highlights how combining iLQR- VAE with balanced model reduction allows to recover the minimal linear dynamical system that can describe the data, and to extract valuable information about the underlying dynamical system.
dimensions thus allowed us to recover dominant
15

---

Published as a conference paper at ICLR 2024
A.
B.
Hankel Singular Values
Fit
Octpar
Time
10
15
20
Dimensions
R(AA)
Figure S1: System Identification (A) Latent coupled oscillator driven by unknown sparse forcing inputs. The control signal u. is fed coupled pendulums parameterized by 2. iL.QR- RARE learns from observations O and can reconstruct Ô after training. Here, blue and orange are used to denote the two input channels of the model, and the two observed outputs of the model. (B) Spectrum of the Hankel singular values. (C) Ground truth eigenvalues of the dynamics (blue), and eigenvalues of the truncated model (orange).
D MODEL COMPARISON
D.1 MODEL SELECTION
The linear dynamical system and MGU RNN were benchmarked on the zebrafish dataset. We eval- uated each model's performance reconstructing the tail dynamics and determining the sparsity of the control signal required to generate the dynamics, see Table S1. Across all combinations of network types and input dimensions, all models achieved high R² scores in the validation dataset. However, the sparsity of the control signal varied, depending on the scale of the network and input dimensions throughout the parameters space.
Table S1: Benchmark model parameters (mean + SEM)
m
MGU
Linear
R²
|
/
60
5
0.974
21.7
60
10
0.975
13.9
90
5
0.976
7.88
90
10
0.976
9.43
90
15
0.976
4.59
120
5
0.976
8.26
120
10
0.983 + 0.000
2.29 0.05
120
15
0.986 + 0.000
4.71 + 0.03
120
20
0.974
3.78
120
10
0.971 + 0.000
2.75 + 0.03
80
10
0.979 + 0.00
3.06 + 0.08
60
10
0.985 0.00
3,97 0.10
40
10
0,986 + 0.00
5.76 + 0.20
20
10
0.987
15.84
18
10
0.985 + 0.002
23.5 + 6,2
15
10
0.985 + 0.000
28.7 5.2
16

---

Published as a conference paper at ICLR 2024
D.2 MODEL ROBUSTNESS
We evaluated the robustness of the model fit across random initializations. To do so, we apply Singular Vector Canonical Correlation Analysis (SVCCA) the Principal Component Analysis (PCA) projection of the initial condition of the RNNs, at movement onset Raghu et al. (2017). SVCCA allows to compare representations by measuring their similarity a way that is independent to rotations of the data. SVCCA is based on Canonical Correlation Analysis (CCA), but includes an additional step during which the dimensionality of the whole data is reduced, greatly speeding up the analysis. In Figure S2, we applied SVCCA to dataset of size K X P where / was the number of points in the dataset, and K the dimension of the PCA projection of the initial condition.
A
0:
Input
u.
0.5
0.5
0.0
0.0
0..0
0.0
108
200
100
200
1001
2000
0
100
2001
1.00
200
2
1.5
Tail angle
[rad.]
0.5
0.0
100
200
100
200
100
200
0
100
200
100
200
B
Time [ms]
D
1.00
4.75
0.986
Simalarity
0.95
4.50
0.985
0.8
us stulfed
0.90
4.25
CCC valuee
0.6
Models
0.85
4.00
0.983
0.80
Sparzds
3.75
0.4
0.75
3.50
0.70
0.981
0.2
3.25
0.980
2
0.60
3,00
Models
0
15
3.75
CCA coof
R.
Figure S2: Robustness of the learned model. Smoothing the input signal to generate movements using a Gaussian kernel a {1.00, 1.68, 2.83, 4.76, 8.00) had small effect on reconstructing the original movement. (b) Using SVCCA to compare the initial latent state representation of two independently initialized models (blue) and control with shuffled labels for the initial state of the same model (orange). (c) The similarity matrix of the top 3 principal components of the initial state across 5 independently initialized models.
(d) A comparison of sparsity and reconstruction for these models when trained the full bout dataset (F) and the dataset excluding specific swim type during training (E).
D.3 CLASSIFICATION MOVEMENT CATEGORY FROM THE LATENT SPACE
Here we detail the method to assess the ability of the latent representation to untangle movement from different categories. At every time step we trained and tested a linear classifier using a 80%/20% split of the latent state trajectory at that time. The of the MGU and LDS were projected into their respective 20 first principal components as reg- ularization. We only considered movement presenting positive first half beat so that rectification would not be required to identify the movement category. With categories and a balanced dataset the chance level was at 9%.
120 dimensional latent state
17

---

Published as a conference paper at ICLR 2024
A.
positive first half-beat
negaláve first ball-bed
B.
Categary
Mean control activation
PC2
control unit
-0.
-10
Figure S3: MGU control peak for different category of movement. (A) Average activation during the first peak for each category of tail movement (B) 2D projection of the 10-dimensional peak in control signal. The control continuously parameterize the movement category
D.4 LINEAR PREDICTION OF TRAJECTORY FOLLOWING INITIAL CONTROL IMPULSE
We measured if the initial control input used to drive the network was predictive of navigational landmarks. By aligning input signals by their initial peak and defining the initial state the network, fit a linear regression of the initial 120-dimensional state of the network to each navigation feature, Figure S4. This was further constrained by T.1NNorm regularization.
r
150
10r
100
100
Prodiced distange [mon]
Proflcted vireatiinn 8
50
Predlcttd pow
-50
-50
-100
-100
Tate direction. (1))
-50
50
100
-150
Thte yow tum, O 1"]
0
2
True distano, [mm]
Figure S4: Linear decoding of navigation features following the first control peak of the net- work. Schematies demonstrating important movement features to navigate in the environment are swim direction 0, swim yaw o, and swim distance r. Linear regression with LI constraint of the initial network state (state reached after the peak control signal) to predict swim trajectory direction (R2 = 0.90), turn yaw (R² = 0.93), and swim distance (R² 0.72).
D.5 LINEAR PREDICTION OF TAIL MOVEMENT FEATURES
We extracted the tail frequency and the maximum tail amplitude from dataset. Using a similar method to Appendix D.4, we could predict the tail beat frequency (R² 0.85) and maximum tail bend amplitude (R² 0.96).
18

---

Published as a conference paper at ICLR 2024
Tail beat frequency
Extremum tail angle
70
Predicted frequency [[1x]
Predicted angle [rud.]
40
30
40
50
True frequency [Hz]
tail angle [rad.]
Figure S5: Linear decoding of tail oscillations features following the first control peak of the network. Linear regression with L1-Norm constraint of the initial network state (state reached after the peak control signal) to predict the tail beat frequency (R² = 0.85) and maximum tail beat amplitude (R² 0.96).
E IMPLEMENTATION OF LFADS
The non-autonomous LFADS model with inferred inputs is composed of RNNs in the encoder, controller and generator. The encoder is constituted of bidirectional RNNs. The controller RNN generates the control signal from the encoder and past factors. To facilitate comparison with iLQR- VAE, we used similar architecture where the generator RNN hidden state directly read into the output. As in Pandarinath et al. (2018), we used an autoregressive prior AR(1) for the control signal. During training, we implemented linear KI. warmup from 0 during the first 100 epochs. See Table S2 for important hyperparameter settings.
Table S2: LFADS Hyperparameters
RNN type
GRU
Encoder Dimension
64
Controller Dimension
32
Control Dimension
10
Generator Dimension
120
AR(I) tau
10
AR(1) noise variance
0.1
Dropout rate
0.3
Batch size
75
Learning rate
3e-4
IMPLEMENTATION OF BALANCED MODEL REDUCTION
F.1 STATE CONTROLLABILITY
OBSERVABILITY
The controllability of the LDS is determined by, both, the coupling of the state dynamics A, and, the number control inputs in the LDS m. The controllability matrix C measures the behavior in which the input matrix propagates through the state dynamics A, defined as,
C(AB) B
AB
A
(S14)
where C (mmmn)) If the control reach all directions in the A subspace, the LDS is com- pletely controllable, thus the will have column rank n. In such case, any state in the network from initial state eventually be reached by arbitrary input signal {u}.
In contrast, the observability matrix assesses whether states within the LDS can realized from its measured observable O measures how much information be retrieved from the
19

---

Published as a conference paper at ICLR 2024
2.5
O
0.0
-2.5
2.5
Wu(c) u(O
0.0
-2.5
2.5
0.0 AM
2.5
Olt)
0.0
-2.5
50
100
50
100
50
100
50
100
Time [ms)]
Time [me]]
Time [ms]]
Time [mss]
Figure S6: Example of LFADS AR(1) Each columns corresponds to an example movements. By regressing the 10-dimensional control signal to the last tail segments across the dataset found projection of the input Wt(t) presenting a correlation of r=0.73 with the last tail segments of the input movement Ott).
readout matrix C' about the state dynamics A, defined by,
C CA
O(C:A) =
(S15)
CC
where O E ()) The LDS is observable providing O has full rank otherwise regions in the state space cannot be reached via observations o. Observability important in order to control the LDS.
F.2 IDENTIFICATION OF BALANCED MODEL TRANSFORMATION
In order to find balanced basis transform, T for the LDS, we employed the framework from Laub al. (1987). Due to the LDS state dynamics being constrained to be asympotically stable, both, W. and W. are positive definite. Therefore, we can compute the Cholesky factors of the Gramians,
W. L
(S16)
W. = L. L.D
T
(S17)
where L. and Lo are the lower triangles of the respective factorization. Using the product of the two Cholesky factors, L L.. we computed the singular value composition,
L. UAV UAV
(S18)
where U E R"xn and V E R"*a are right and left basis, respectively, and A E Rax* is a diagonal matrix where A " corresponding the Hankel singular values. From this decomposition, the balancing transformation can be obtained, yielding,
T = VV
(S19)
T A L
(S20)
This basis is then projected on the LDS, yielding a system to that Equation (7).
F.3 COMPARISON OF REDUCED MODELS
We varied the order of truncation, r, and compared the balanced truncated model performance in reconstructing the data. We found that the LDS truncated to r 18 could reconstruct most move- ments, with the exception of escape swims.
20

---

Published as a conference paper at ICLR 2024
Table S3: Reduced models reconstruction
R
120
0.970
60
0,960
40
0.904
20
0,770
10
0.119
F.4 CONTRIBUTION OF EIGEN MODES BALANCED TRUNCATED TURNS AND FORWARD MOVEMENTS
Using our LDS truncated to 1818, we can decompose the readout of the dynamical system as the sum of the contribution of each mode. When quantify the contribution of each mode to the final tail movement we found that overall routine turns are most differentiated from slow2 forward movement in their use of the mode corresponding a real eigenvalue with 7=-17ms similar the examples in Figure 3D.
Mode comtibution [ax)]
10.0
Slow 2
7.5
Routine Turn
5.0
2.5
0.0
T 92 ms
r= .75 ms
I= ms
I- ms
t ms
20 ms
ms Te
t 17 ms f=
ms
I el ms
f= Mh
f= li
2=2
f= 38 ll
== S He
Hi
f= 29
f=0
Thil aingl
[cad]
100
100
100
00
100)
100
100
100 Time [ss)
100
100 Time [es)
Time [ms)
Time [ms]
Time (ms)
Time [ms]
Time [us)]
Time (ns
Timee [ms]]
Time [ms]
Figure S7: (Top) Distribution of the - of each mode across routine turn and slow2 move- ments. (Bottom) Projection of impulse response into the readout for each of the 10 modes (8 com- plex conjugates eigenvalues and real eigenvalues).
G ZEBRAFISH DATASET
Fish Care Adult fish were maintained at 25 C 14:10 hour light cycle following standard methods. Embryos were collected and larvae were raised at 28 in E3 embryo medium (5 mM NaCI, 0.17 mM KCI, 0.33 mM CaC12 and 0.33 mM MgSO4) in groups of 30. The wild-type line Tübingen (Tu) was used. Sexual differentiation occurs at later larval stage, therefore the sex of the animals cannot be reported. All experimental procedures were approved by the Champalimaud Foundation Ethics Committee and the Portuguese Direcção Geral Veterinária, and were performed according to the European Directive 2010/63/EU.
Behavioral experiments The dataset was compiled from observations from over 100 freely- swimming larvae at 6 or days post fertilization. The larvae were exposed to various ethologically relevant visual and auditory stimuli, including moving gratings of different speeds and orientations, and transitions between light and dark environments. We monitored the larvae's movements by tracking seven landmarks along their tails The oscillations of the tail, driven by activity in the spinal cord, were characterized by measuring the angle of each segment relative to the fish trunk. Using methods from Marques et al. (2018), we segmented the data into swim bouts which we classified into 11 distinct kinematic categories: approach swim, slow 1, slow 2, burst swim, turn, high angle turn, routine turn, shadow avoidance, bend, long latency C start, short latency C start (see Fig- IC and Supplementary Video 1). We excluded capture swim categories since the dataset did not
21

---

Published as a conference paper at ICLR 2024
contain prey capture episodes. For each category, the dataset contained both movements presenting a positive first half beat as in Figure IC or a negative first half beat. The final dataset was balanced and contained 30800 tail angle time series lasting 230ms.
H DETAILS OF C. ELEGANS
H.1 C. ELEGANS DATASET
Behavioral experiments We used the dataset from Broekmans et al. (2016), composed of two experiments. First, foraging dataset where 12 worms were recorded at 16 Hz for 35 min. The other set consisted of 98 worms recorded at 20 Hz for 30 s. Each worm was subjected to a 100 ms shock on the head with an infrared laser seconds after acquisition and tracked for an additional 20 s post-zap. The frame rate was down sampled to 16Hz, and recordings were segmented into 12.5 s trials. To train the LDS model we used these 12.5 trials. To validate reconstruction the dynamics and inferred inputs, we use the entire 30 trials for each worm in the escape dataset.
Behavioral tracking C. elegans behavior was acquired by tracking 100 segments along the body. At each segment, the tangential angles along the centerline calculated. Four principal compo- nents explain over 95% of the variation, and each mode was referred to as an eigen worm (EW), shown in Figure EWI and EW2 map the sine-cosine pair of the body posture, different com- binations reflect wave along the body at different phases, important for forward and reverse crawl maneuvers. EW3 represents the cear-consttnt curvature along the central half of the body line, this mode is mainly recruited for turns. Finally, EW4 denotes the bend in the head and tail of the posture.
H.2 REDUCED MODEL EIGEN MODES CAPTURE BEHAVIORAL STATES IN ESCAPE SEQUENCES
We applied balanced model reduction to enhance our understanding of the LDS model's s dynamics. In Figure S8A, compressed the LDS model, size 60, to 20 with balanced model truncation. The reduced model was able to reconstruct the EW dynamics with R² 0.702. The model with 20 states covered wide spectrum of frequencies and timescales comparable those in the full LDS model, as shown in the eigen spectrum Figure S8B. Interestingly, the three slowest decaying modes mapped the reverse, Q-turn, and forward crawl, see Figure S8C.
22

---

Published as a conference paper at ICLR 2024
Heat shock
1=-10s f=066 Hz
D.
2.0
0 a. a
Reverse
11SV [A.U]
1.0
== s0.5
0.5
Qawm
f= 0.1 Hz
0.0
20
40
60
Index
0.5
Forward
=m 0.6 Hz
B.
101
Full Reduced
Time
3
Forward
Reverse Q turn
Forward
Frequency [12
0.
-5
-10
-1.0
-0.5
0.0
10
15 Time [s]
20
25
[s]
Figure S8: Eigen-modes from balanced model reduction reflect behavioural states in escape sequence. (A) Hankel singular values of the LDS truncated at r 20 (red vertical line). (B) Eigen spectrum of the reduced model (open marker) compared the full LDS (black marker). (C) From top to bottom, are the projections of an impulse response for the three slowest decaying eigen modes in the reduced model. (D) The individual contribution the three slowest eigen modes capture the highest variance across at different behavioral states within the escape sequence. Moving from top to bottom, each modes corresponds to the reverse crawl, escape Q-turn, and forward crawl, respectively. Their contribution is reflected in the in the full reconstruction (middle-bottom pancl). The driving-input to generate mode contribution (bottom).
23


Reviews:
[{"review_id": "R93l60qQLv", "reviewer": "Reviewer_AsPq", "summary": "This paper uses an unsupervised method of discovering latent control signals from behavioural sequences, to discover behavioural motifs. The main methodological contribution is using a model reduction technique to gain interpretability insights into the inferred models. Overall I think this is a good paper.", "questions": "See weaknesses", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "This paper analyses various behavioural sequences and infers intuitive behavioural motifs.\n\nWhile the core method is not new, the interpretability aspect is, and the insights are interesting.\n\nThe paper is well written, and clearly presented.", "weaknesses": "The paper could do with a proper model comparison, i.e. vs moseq.\n\nThe Balanced truncation bit needs to be better explained with some intuition. \n\nEqn 5 should have $ \\tilde{W}_o = \\tilde{W}_c = \\dots$\n\nFig 1 A: what does blue and orange correspond to? Different sequences?\n\nSome of the Figures seem to be mis-referenced, e.g., “As shown in Figure 3B, the MGU could reconstruct bouts using a sparser control compared to the LDS. \"\n\nAll the figures could do with more explanations in their captions.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "", "extra_scores": {}}, {"review_id": "Yyt9I9gCUL", "reviewer": "Reviewer_9AdG", "summary": "The paper proposes a framework to identify the latent control signals and the underlying dynamics that make up the natural behavior of zebrafish. It utilizes iLQR-VAE for learning the parameters and inferring the control signals. The paper trained recurrent neural network (RNN) models and linear dynamical systems (LDS) to reproduce the postural sequence with sparse control signals. The underlying dynamics of those models are explored and how they relate to the behaviors is studied. The paper further demonstrates the model order reduction on the LDS model and relates the reduced mode to the observed behaviors.", "questions": "Please see above weaknesses.", "limitations": "", "rating": 5, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "The paper looks at a new approach for studying naturalistic behaviors, and defines the problem in a clear way.\n\nThe model reduction on the LDS model gives new insights into modeling behavior.", "weaknesses": "The paper utilizes iLQR-VAE for learning the parameters and inferring the control signals. It would be helpful if the authors stated clearly in the text what is the difference between the proposed model and the previous model. It seems like the differences are minimal if any, in which case, this is a good application paper, however, it doesn't introduce many novel elements from a modeling perspective. While the model reduction technique for linear models is not usually applied to this field of behavioral modeling, it is a very well-known concept and does not offer novelty from a methodological perspective. There are interesting takeaways from a neuroscience perspective, but these would have to be validated more thoroughly and may be more suitable for a different venue.\n\nThe kinematics of the zebrafish are somewhat simple; however, the naturalistic behaviors of other animals (mice and monkeys, for example) are much more complex. Will the model still be able to capture these dynamics? One real-world dataset may not be able to provide enough insight about the generalizability of the model for this question.\n\nThere is no comparison with other models presented.\n\n\nMinor:\n1. Section 4.1, mismatched figure labels. Figure 1C->Figure 1B; Figure 1D -> Figure 1C\n2. Figure 3 Caption: missing space: ‘udriving’ -> ‘u driving’\n3. More explanations on Figure 3D in the main text.\n4. Page 15: ‘all states withing the LDS’ -> ‘all states within the LDS’\n5. Page 15: ‘cannot be reach via observations’ -> ‘cannot be reached via observations’\n6. Page 16: ‘Therefore, we can the compute the’ -> ‘Therefore, we can compute the’", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "", "extra_scores": {}}, {"review_id": "fZ03jhdEdY", "reviewer": "Reviewer_325F", "summary": "In this paper the authors propose a control framework for estimating movement dynamics from behavioural/postural observations. Movements are considered as generated by an unknown dynamical system controlled through sparse input signals.\nThey employ two alternative methods to fit the latent dynamics: \n- One that uses the iLQR-VAE from Schimel et al., 2022 that approximates the observed dynamical system through an RNN.\n- And a linear model of “large dimensionality” that is consequently reduced to a lower dimensional model using balanced truncation.\nThey demonstrate their method on a toy example of a linearly approximated system comprising two pendulums, and on a behavioural dataset of zerbafish movements. \nThey are able to further provide insights into the identified dynamics by dissecting the spectrum of the reduced linear model.", "questions": "- In the toy example with the pendulum, as I understand, the authors created observations by simulating the linear approximation of the system (assuming small angles). I wonder how  the proposed framework would perform if the authors created the observations for the same system parameters with the nonlinear dynamics, and for parameter sets that result in increasingly larger angles. I think in this toy example it is interesting to demonstrate the robustness of the linear-reduced order framework.\n\n- In Figure 3 why does the classification accuracy decrease with time (observation time?)?\n\n- In Figure 4C and associated main text, the authors mention that the reduced model added eigenmodes with larger timescales, while it neglected the small timescale modes of the the non-reduced system. While the second part of the previous sentence is expected, I am not sure how should I understand the first one. Can you provide some intuition. Also related to this, when I first looked at Figure 4C I thought that the large the large time scale values of the reduced and non-reduced model overlap, therefore the light grey ones are not visible. Can you probably make the circles of the reduced system non-opaque or non-filled to make the plot clearer?\n\n- I think it would be helpful if the authors mention in the supplement how they fit the linear model of their framework, and provide a brief description of the iLQR-VAE framework that is a crucial component of the proposed approach.\n\n- Model Selection C.1. section is missing from the supplement.\n\n- There is a typo in the subscripts in Eq. 5.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "- Applying balanced model reduction together with dynamical and control inference is an interesting contribution to the existing system identification literature.\n- The linear approximation through the model reduction allows for interesting insights into the approximated dynamics (like in Figure 4 for analysing the eigenmodes active during different movements), that would be otherwise more cumbersome to perform for an RNN.", "weaknesses": "- The dissection of the identified dynamics through the reduced linear model requires the non-reduced linear system to already approximate accurately the observed system. I wonder whether the authors could provide a systematic analysis on the robustness of their framework when fit on the nonlinear toy model.\n- As the authors already mention in their discussion, the framework they propose requires the fitting of two dynamical systems: an RNN that captures accurately the observed behavioural trajectories, and an linear system that provides interpretability.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "", "extra_scores": {}}, {"review_id": "1iN8K6HZyB", "reviewer": "Reviewer_j5DC", "summary": "This manuscript provides a framework for modeling the behavior of the larval zebrafish using a control theoretic and dynamical systems framework. They model the behavior of a system as arising from a dynamical system driven by sparse inputs with a students-t distributed prior. They assume behavior is generated from the dynamical system by either a linear or RNN model, and use iQLR-VAE to infer model parameters and the sparse inputs. In each case they use a model truncation approach to reduce the number of dimensions considered. They begin by validating the approach on a double pendulum oscillator system, showing they could reproduce the number and values of eigenvalues. Then they analyze a dataset of larval zebrafish behavior, which is organized into discrete bouts and thus a natural test case for this approach. They find they need high capacity models (120 latent dimensions, 10 control inputs) to fit the behavior. They show that the resulting representation is useful, in that the dynamical representations provides a more disentangled input that can be better decoded from the top 5 PCs compared to the postural representation or even the full time series in some cases. By using model reduction, they show the dynamics underlying behavior can be decomposed  into separate modes.", "questions": "* Does the number of control inputs match the number of bout types people have reported in the literature? \n\n* There is no ground truth in behavior for the driving inputs and thus there is no way to validate that they are correct, however with an appropriate experiment (e.g. applying a perturbation or stimulus at a specific time) perhaps they could be useful. \n\n* Generally, more concretely linking any of the observations to real neural or behavioral variables would improve the manuscript's impact.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "* The formulation is well described and I believe novel, and an addition to other methods fitting dynamical systems to neural and muscular systems, but here in the case of sparsely driven behavior. \n\n* The result about better disentangling of behavior is interesting, and the spatial decoding result as well. \n\n* The model reduction approach to achieve small models is interesting.", "weaknesses": "* There is a lack of benchmarks comparisons to other techniques in the literature such as LFADs. Even though some suggest continuous inputs for their control scheme, could they be used to model the behavioral dynamics in Figure 3C for instance. \n\n* I found the applications a bit limited. It wasn’t obvious from this analysis that there were fundamentally new types of experiments that were enabled by this approach. \n\n* Overall, I think the lack of benchmarks combined with the lack of novel applications was a major sticking point for me and I view this manuscript as very borderline. I would be happy to reconsider the manuscript (which was very well written) if these were added or addressed. I see this as just 1 compelling figure short of a reasonable submission.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (poster)