The manifold structure of limb coordination in walking Drosophila
Abstract
Terrestrial locomotion requires animals to coordinate their limb movements to efficiently traverse their environment. While previous studies in hexapods have reported that limb coordination patterns can vary substantially, the structure of this variability is not yet well understood. Here, we characterized the symmetric and asymmetric components of variation in walking kinematics in the genetic model organism Drosophila. We found that Drosophila use a single continuum of coordination patterns without evidence for preferred configurations. Spontaneous symmetric variability was associated with modulation of a single control parameter—stance duration—while asymmetric variability consisted of small, limbspecific modulations along multiple dimensions of the underlying symmetric pattern. Commands that modulated walking speed, originating from artificial neural activation or from the visual system, evoked modulations consistent with spontaneous behavior. Our findings suggest that Drosophila employ a lowdimensional control architecture, which provides a framework for understanding the neural circuits that regulate hexapod legged locomotion.
Introduction
Legged locomotion requires flexible coordination of limbs in varied, changing environments. This coordination has been studied in a variety of model organisms, from camels and rhinoceros to ferrets, rats, beetles, and ants (Alexander and Jayes, 1983; Büschges et al., 2008). In hexapods, many studies have noted considerable variability in limb movements (Ayali et al., 2015; Mendes et al., 2013; Pereira et al., 2019; Strauß and Heisenberg, 1990; Wosnitza et al., 2013; Zollikofer, 1994), but the structure and origin of this variability remain unclear. A comprehensive characterization of this variability is required to understand locomotor control in insects (Krakauer et al., 2017). To understand how insects coordinate their limbs during terrestrial locomotion, we investigated the variability in walking patterns in the fruit fly Drosophila, where we could record behavior in many flies and manipulate neural signals with genetic tools.
In many animals, variability in walking patterns takes the form of distinct locomotor gaits (Alexander, 1989; Alexander and Jayes, 1980; Alexander and Jayes, 1983; Srinivasan and Ruina, 2006; Thorstensson and Roberthson, 1987). In the study of hexapod locomotion, it is unclear how the observed variability in walking coordination relates to potential preferred gaits or to forms of continuous variability. Although some studies report distinct gaits in insects (Ayali et al., 2015; Bender et al., 2011; Burrows, 1996; Graham, 1972; Mendes et al., 2013; Pereira et al., 2019), many note a high degree of variability in limb coordination, with intermediate patterns intermixed with movements that resemble canonical gaits (Ayali et al., 2015; Mendes et al., 2013; Pereira et al., 2019; Strauß and Heisenberg, 1990; Szczecinski et al., 2018; Wosnitza et al., 2013; Zollikofer, 1994). Drosophila have been reported to use different coordination patterns at low and high walking speeds, with significant variability around these canonical patterns (Mendes et al., 2013; Pereira et al., 2019; Strauß and Heisenberg, 1990; Szczecinski et al., 2018; Wosnitza et al., 2013). The nature of this variability has not been fully investigated, as these studies focused primarily on small numbers of individual bouts of walking.
When walking speed is varied, it represents a form of symmetric variability in walking patterns. These walking patterns are themselves often antisymmetric, with legs moving out of phase, but we refer to modulations of these patterns as symmetric if they are equal on the two sides of the body. This symmetry maintains heading. Locomotion also contains asymmetric variability, in which modulations are unequal on the two sides of the body. Such asymmetry is required to change heading. In hexapods, measurements of limb kinematics (Cruse et al., 2009; Domenici et al., 1998; Dürr and Ebeling, 2005; Franklin et al., 1981; Frantsevich and Mokrushov, 1980; Graham, 1972; Gruhn et al., 2009; Mu and Ritzmann, 2005; Strauß and Heisenberg, 1990; Zollikofer, 1994; Zolotov et al., 1975) and dynamics (Jindrich and Full, 1999) have described small, limbspecific modulations of movement associated with turning. As in the study of symmetric variability, these studies focused primarily on individual bouts of locomotion, and most employed tethered preparations. To understand how limb movements are modulated during turns, and how symmetric and asymmetric forms of variability are related to one another, it is necessary to observe full distributions of behavior in intact and unrestrained animals.
In this study, we generated large datasets of locomotor behavior using a simple method for robustly tracking the body and limb movements of freelywalking Drosophila. By analyzing hundreds of thousands of steps, we found that the locomotor variability in flies exists on a single continuous manifold, without evidence for discrete preferred patterns. Symmetric variability could be described by changes to a single global parameter governing walking speed. In contrast, turning involves small, asymmetric, and limbspecific modulations that are precisely timed relative to the underlying symmetric pattern. When slowing was evoked by optogenetic activation of command neurons or by visual stimuli, the resulting symmetric modulations of limb movement were consistent with the manifold present in spontaneous walking. The structure of the variability in spontaneous and evoked coordination patterns suggests a simple model for limb coupling during locomotion.
Results
A simple automated method robustly identifies limb positions
We constructed a planar arena in which freelywalking flies were illuminated from above and could be filmed from below in silhouette by a highspeed camera at 150 frames per second (Figure 1A). Within each frame, flies were first located and then rotated into a common orientation (see Materials and methods). From the centroid positions and heading angles, we computed the fly’s forward velocity (parallel to its heading), its lateral velocity (perpendicular to its heading), and its yaw velocity, which we denote as ${v}_{\parallel}$, ${v}_{\perp}$, and ${v}_{r}$, respectively (Figure 1B). Exploiting the stereotyped anatomy of flies, we then used a simple linear model to identify the position of each limb in each frame (see Materials and methods, Figure 1—figure supplement 1, Table 1, and Videos 1–3). To ensure that correlations between limb poses in our training set did not systematically bias predicted limb positions, we required the model to predict limb positions from image regions limited to each limb’s approximate range of motion (Figure 1—figure supplement 1). We trained our model on a set of 5000 handannotated images, and used 10fold crossvalidation to test for overfitting. This crossvalidation showed that the mean error between predicted and handannotated limb positions was 0.15 mm (3.5 pixels) (Figure 1—figure supplement 1). Importantly, the accuracy of this simple linear model in our experimental context was comparable to that achieved by published deep neural network methods (Mathis et al., 2018; Pereira et al., 2019) (Figure 1—figure supplement 1).
Measurements of spontaneous wildtype walking are consistent with previous studies
Many previous studies have investigated the statistics of walking in Drosophila (Berman et al., 2014; Bidaye et al., 2014; Branson et al., 2009; Geurten et al., 2014; Kain et al., 2013; Katsov and Clandinin, 2008; Katsov et al., 2017; Martin, 2004; Mendes et al., 2013; Strauß and Heisenberg, 1990; Wosnitza et al., 2013). As the flies walked in the arena, their heading, lateral, and forward velocities oscillated continuously due to the periodicity of limb forces (Figure 1—figure supplement 2). Since we were interested in longertimescale modulations of the path, we smoothed the centroid kinematics to eliminate this oscillation (Figure 1—figure supplement 2, see Materials and methods) (Katsov and Clandinin, 2008; Katsov et al., 2017). Consistent with prior studies, we observe forward velocities between −1.3 and 30.4 mm/s (2.5th and 97.5th percentiles), with peaks in the distribution of forward walking speeds at 0 mm/s and ~17.5 mm/s (Figure 1Ci). Variation across studies in the measured distributions of walking speeds may be due to variation in experimental conditions that affect the motivational state of the fly, such as temperature and hunger (Chadha and Cook, 2014; Crill et al., 1996; SotoPadilla et al., 2018).
Our large dataset enabled us to characterize the full distributions of centroid kinematics, rather than measuring only a small number of individual examples. The distributions of lateral and angular velocities were roughly symmetric, illustrating that we do not observe populationlevel handedness in Drosophila walking (Figure 1C). The joint distribution of the forward velocity ${v}_{\parallel}$ and the yaw velocity ${v}_{r}$ showed that flies turn at a broad range of yaw rates across many forward speeds (Figure 1Di–ii). Crabwalking, in which the fly’s heading vector is no longer tangent to its path, predominately occurs during highyawrate turns at slow to moderate forward velocities (Figure 1Di–iii). These characteristics of aggregate kinematic behavior are consistent with those reported in previous studies (Strauß and Heisenberg, 1990; Geurten et al., 2014; Katsov and Clandinin, 2008; Katsov et al., 2017).
Our measurements of limb kinematic parameters were also broadly consistent with previous studies. Importantly, we restricted this and subsequent analyses of limb coordination to locomotion by excluding flies moving below 0.5 mm/s. This criterion excluded flies that were stopped or grooming, behaviors that involve interesting patterns of limb coordination (Seeds et al., 2014), but lie outside of the focus of this study. In our data, step length in the stationary frame of the camera increased roughly linearly with forward walking speed (Mendes et al., 2013) (Figure 1Ei). Across walking speeds, swing duration increased, but this modulation is small in comparison to the change in stance duration (Figure 1Eii–iii). Stance duration was approximately inversely proportional to forward walking speed (${\tau}_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}~{{v}_{\left\right}}^{1.025},{R}^{2}=0.59$) (see Materials and methods) (Figure 1Eiii). Therefore, our data were broadly consistent with previous observations that walking speed differences are dominated by changes in stance duration, while swing duration remains relatively constant (Wilson, 1966; Wosnitza et al., 2013; Dürr et al., 2018; Mendes et al., 2013).
Swingstance patterns change in a speeddependent manner
In hexapod locomotion, studies have described three canonical gaits (Figure 2A–B) (Bender et al., 2011; Burrows, 1996; Collins and Stewart, 1993; Graham, 1972; Spirito and Mushrush, 1979; Wilson, 1966). In tripod gait, two groups of three limbs swing simultaneously, creating two alternating configurations of swing and stance. In tetrapod gait, three groups of two limbs swing in sequence, while in wave gait, each limb swings individually. In general, distinct gaits can be distinguished by discontinuous changes in one or more parameters of the coordination pattern (Alexander, 1989; Srinivasan and Ruina, 2006). However, defining gaits in terms of discontinuous transitions can be problematic in the presence of measurement noise or behavioral variability, which can blur sharp transitions. Instead, one may characterize the distribution of walking coordination patterns, where distinct gaits correspond to multiple modes in the distribution (Hoyt and Taylor, 1981).
To examine the structure of limb coordination in flies, we first investigated the structure of the swingstance patterns used at different forward walking speeds (Figure 2C). We determined whether a limb was in swing or stance phase by measuring its frametoframe displacement in the camera frame. If a limb position moved more than roughly our error in positional measurement (~0.13 mm in the camera frame) then we considered the limb to be in swing phase (see Materials and methods). Movements below this threshold were classified as stance phase. Thus, for each flyframe, we generated a sixdigit binary vector describing the swingstance configuration of the fly’s six limbs.
These individual configurations can be categorized by how many limbs are in stance phase. Previous studies have reported an enrichment of the threefoot stance category at faster forward walking speeds and a corresponding enrichment of the fourfoot stance category at slower forward walking speeds (Mendes et al., 2013; Pereira et al., 2019; Strauß and Heisenberg, 1990; Wosnitza et al., 2013). These studies suggested that the proportion of feet in stance at different speeds were consistent with flies using tetrapod gaits at lower speeds and a tripod gait at higher speeds. Our data also showed speeddependent enrichments of five and fourfoot stance categories at slow speeds, and enrichments of the threefoot stance category for flies walking at faster speeds (Figure 2C).
To understand the structure of these different stance categories, we investigated the specific configurations of limbs in each of the three, four, and fivefoot stance categories. For each of the canonical gaits, one would expect the constituent discrete limb configurations to occur in roughly equal proportions (Figure 2B). In the case of the threefoot stance category, the two configurations that comprise the canonical tripod gait occur in roughly equal proportions at all forward walking speeds (Figure 2D). In contrast, different configurations of the fivefoot and fourfoot stance categories do not occur with equal abundance across walking speeds (Figure 2D). Interestingly, the fourfoot stance configurations that occur with greatest frequency are those that are a single swingstance transition away from a canonical tripod gait.
A large fraction of four and fivefoot stance configurations are transient
To investigate the dynamics of stance configurations in our data, we examined the dwell times in all individual stance configurations within each canonical gait (Figure 2E). We defined the dwell time as the number of contiguous frames in the configuration. If each canonical gait existed in the data, one would expect that the dwell time in each stance configuration associated with the gait would approximately match the average swing duration of the limbs. Because previous studies reported that Drosophila preferentially use tetrapod and tripod gaits at slow and fast speeds respectively, we partitioned the data into thirds based on forward walking speed (Mendes et al., 2013; Pereira et al., 2019; Wosnitza et al., 2013). At slow speeds (v_{ }≤ 10.2 mm/s), while many instances of canonical configurations were transient (one frame or ~6.7 ms), there existed longer durations for all tripod, tetrapod, and wave gait configurations (Figure 2E). At higher speeds, virtually all tetrapod and wave gait configurations were transient, lasting only one frame.
Next, we found the most frequent forward velocities within each stance category, which were 7 mm/s, 13 mm/s, and 24 mm/s for 5foot, 4foot, and 3foot stance categories, respectively (Figure 2D). If each canonical configuration represented a distinct, preferred gait, then we might expect these walking speeds to occur with greater frequency than nearby walking speeds, as is observed in horses (Hoyt and Taylor, 1981). Yet, no peaks in the distribution of forward velocities exist in these locations, suggesting that this phenomenon is not present in flies (Figure 1Ci).
The abundance of each stance category peaks twice per singlelimb stride
Each canonical gait contains a different number of discrete limb configurations and thus a different number of transitions between these configurations (Figure 2B). In canonical gaits, swinging limbs move simultaneously, but deviations in step timing can transiently change the number of feet down. For tripod gait, changes in the number of feet down would occur twice per cycle, corresponding to the two transitions between canonical stance configurations. For a tetrapod gait, we would expect a threecycle associated with the three transitions between the three tetrapod stance configurations. Correspondingly, a wave gait would produce a sixcycle corresponding to the transitions between each of the independently swinging limbs. Previous studies have noted that deviations in step timing can produce intermediate noncanonical stance configurations (Mendes et al., 2013). Yet, even if these intermediate configurations are present, the periodicity of the relative abundance of three, four and fivefoot stance categories as a function of limb phase still provides information about potential underlying gait patterns.
To calculate the phase of each limb at each point in time, we decomposed the time series of the limb position in the direction parallel to the fly’s body axis into amplitude and phase components using the discretetime analytic signal method (Figure 2—figure supplement 1) (see Materials and methods) (Boashash and Reilly, 1992; Gabor, 1946; Lamb and Stöckl, 2014; Marple, 1999; Vakman and Vaĭnshteĭn, 1977). This transformation gives an estimate of the phase of each limb at each timepoint and has previously been applied to analyze cockroach limb oscillations (CouzinFuchs et al., 2015). We then visualized the number of feet down as a function of each limb’s phase to count the number of transitions per step cycle (Figure 2F, Figure 2—figure supplement 2). If canonical wave and tetrapod gaits occurred with significant frequency in our data, we would expect to see walking speeds where the number of feet down oscillates as a 6 or 3cycle. Yet, the frequency of all limb categories peaked twice per cycle at all forward walking speeds. This twocycle persists even when the threshold that distinguishes between swing and stance is halved, a manipulation that accentuates any effects of measurement noise. The persistent twocycle in the data is inconsistent with the existence of tetrapod or wave gaits in a large fraction of measured trajectories. However, if tetrapod or wave gaits accounted for only a small fraction of the walking data, their characteristic periodicity might not be visible in this analysis.
Contralateral limbs are antiphase across all walking speeds
To investigate variability in the coordination among limbs, we examined the full distributions of pairwise relative limb phases (Figure 3A–B). Within the set of canonical gaits, each has a characteristic set of phase offsets among the six limbs (Collins and Stewart, 1993). The most probable relative phase configuration in our data was not consistent with any single canonically defined gait (Figure 3—figure supplement 1). Like the canonical tripod gait, all contralateral pairwise relative phases maintain a mean relative phase of a halfcycle ($\u27e8\mathrm{\Delta}\varphi \u27e9\approx 0.5$) (Wosnitza et al., 2013). However, unlike a perfect tripod gait, ipsilateral midfore and hindmid pairings are not on average in antiphase ($\u27e8\mathrm{\Delta}\varphi \u27e9\approx 0.4$), while hindfore pairings are also not in phase ($\u27e8\mathrm{\Delta}\varphi \u27e9\approx 0.85$), consistent with previous observations (Wosnitza et al., 2013). The observed pairwise relative phases are also inconsistent with expected phasing for canonical tetrapod and wave gaits (Figure 3—figure supplement 1).
To understand how phase relationships depend on forward walking speeds, we estimated the distributions of all pairwise limb relative phases conditioned on forward walking speed (Figure 3—figure supplement 2). Across all forward velocities, these distributions were unimodal, allowing us to characterize them by their first and second circular moments. Strikingly, at all forward walking speeds, all three contralateral limb pairings maintain a constant mean phase offset of a halfcycle (Figure 3C). Mean ipsilateral limb relative phases, in contrast, are not constant across walking speeds, and approach a tripodlike coupling at the fastest forward walking speeds (Figure 3C).
Notably, all pairwise relative phases exhibit a monotonic decrease in angular deviation as forward walking speed increases (Figure 3D). One interpretation for this decrease in variance is that interlimb coupling is weaker at low speeds than at fast speeds (Aminzare and Holmes, 2018; Aminzare et al., 2018). This possibility is consistent with the observation that lowspeed walking is more sensitive to sensory information and correspondingly more variable (Berendes et al., 2016; Isakov et al., 2016). A second interpretation is that some of the variability results from the characteristics of phases for waveforms with duty cycles that deviate from onehalf (CouzinFuchs et al., 2015; Lamb and Stöckl, 2014; Vakman and Vaĭnshteĭn, 1977) (Figure 3—figure supplement 3). When measuring the relative phases between limbs, each limb sweeps out one halfcycle while in stance and another halfcycle while in swing. When the limbs have asymmetric duty cycles, they sweep out different fractions of a cycle per unit time when in swing versus stance. Thus, slow walking, which is empirically characterized by longer stances and approximately the same swing durations as fast walking, will necessarily have greater variance in relative phasing. Both explanations likely contribute to the observation of greater variance in limb relative phases at low speed.
The joint distribution of limb relative phases is unimodal
To determine the preferred phase offsets used by the control circuits that govern walking, we examined the joint distributions of interlimb relative phases (Figure 3E). If Drosophila used multiple distinct gaits, we would expect to see multiple peaks in this probability distribution (Figure 3—figure supplement 1). For instance, in this view, tripod gaits exist in one region, tetrapod gaits exist in two, and wave gaits exist in a fourth. Yet, the observed distribution has only a single peak, suggesting that the spontaneously walking fruit fly uses a single continuum of limb relative phases (Figure 3E). This observation holds across all forward walking speeds; in particular, there was only a single mode even at the slowest walking speeds (Figure 3—figure supplements 2 and 3). The most abundant relative phasing is distinct from the peaks predicted for any canonical gait; however, it is closest to that of a canonical tripod gait. There were no peaks in the distribution at locations corresponding to canonical tetrapod or wave gait configurations. However, the tails of this distribution extend to the location of both canonical tetrapod peaks, which allowed us to identify individual trajectories with tetrapodlike characteristics.
Consistent with previous studies, we could extract individual trajectories that are similar to canonical tripod and tetrapod gaits by selecting trajectories with the expected duty cycles (Figure 3F–G) (Mendes et al., 2013; Pereira et al., 2019). Interestingly, the tetrapodlike examples appear to have a characteristic crossbody offset in step timing, so that the limbs that should swing together in the canonical tetrapod gait are actually slightly offset in time. This offset maintains approximate antiphase between contralateral pairs of limbs, rather than the expected phase offset of the canonical tetrapod gait (Figure 3F). In addition, we observed a previously unidentified noncanonical coordination pattern in which fore and hindlimbs swing together while each midlimb swings alone (Figure 3F). This pattern is distinct from all previously described canonical gaits, but might be expected when considering a single continuum of coordination patterns. In particular, this noncanonical configuration was predicted from an early continuum view of limb coordination (Wendler, 1964; Wilson, 1966).
Templatematching suggests little evidence for preferred canonical tetrapod and wave gaits
We next sought to systematically investigate how closely our data matched the configurations of relative limb phasing specified by each of the canonical gaits. To do so, we defined a templatematching coherence metric based on a simple model for networks of coupled oscillators (see Materials and methods). The coherence falls between zero and one, with unity corresponding to an exact match to the template phase configuration. The tripod coherence distribution is sharply peaked around ~0.9, indicating a close template match, while the distributions for the other canonical gaits were broader and peaked at lower coherences (Figure 3H). If the fly preferentially used different gaits at different forward walking speeds, then we would expect to see a peak in the conditional coherence distribution for a given canonical gait at some characteristic forward velocity. Yet, there were not prominent peaks in our data, and all coherences yielded low average values at slow walking speeds (Figure 3I). At those slower speeds, where mean coherences were low, the tetrapod coherences became more likely to be the highest coherence (Figure 3J). This is consistent with the observation of tetrapodlike coordination patterns at low speeds (Figure 3F–G), and with the increased spread of the distribution of relative phases at low speeds (Figure 3D,E, Figure 3—figure supplements 2 and 3; Wosnitza et al., 2013). However, even when the largest coherence was not tripod, these nontripod coherences were lower than when the tripod coherence was largest (Figure 3—figure supplement 3). Since the tetrapod coherence of a perfect tripod gait is 3^{–1/2} ~ 0.6, this coherence metric may not show strong contrasts under some deformations in coordination. Therefore, under these analyses, we do not observe strong evidence for the preferred use of canonical tetrapod and wave gaits.
A single lowdimensional manifold describes limb coordination during walking
To directly visualize the structure of the continuum of coordination patterns indicated by our phase analysis, we first applied principal component analysis to our limb kinematic data (Figure 4—figure supplement 1). The principal components of these data occur in degenerate, phaseshifted pairs. The projection into the first two principal components provides some information about the phasic structure of our data. However, variations in walking frequency are not wellrepresented by this projection. In particular, the first two principal components emphasize a single frequency, so that frequencies that are both higher and lower are mapped to the same radius in the plane. Furthermore, the projection into the first three principal components does not provide greater insight. Therefore, principal component analysis provides some information about the structure of the data, but does not permit easy visualization of its full structure in two or three dimensions.
We therefore applied Uniform Manifold Approximation and Projection (UMAP) to generate a lowdimensional embedding of our data (McInnes et al., 2018). UMAP minimizes the crossentropy between topological representations of the local distance metric on some highdimensional data manifold and on a lowdimensional embedding manifold. Among nonlinear dimensionality reduction techniques, UMAP is ideal for this application because it attempts to preserve the local and global topology of the highdimensional data, in contrast to algorithms like tSNE, which is designed to emphasize local structure at the expense of global structure (Becht et al., 2019; McInnes et al., 2018).
Using UMAP, we generated lowdimensional representations of 10^{5} randomlysampled segments of limb positional data, each with a halfwindow length of 100 ms (15 frames). This sampling thus had a dimensionality of 31 frames by 6 limbs by two spatial coordinates, or 372 total dimensions. Our data was represented by a bellshaped manifold embedded in threedimensional space (Figure 4). The linear dimension along the axis of the bell is parameterized by the mean stepping frequency (Figure 4A, Figure 4—figure supplement 2), while the cyclic coordinate is defined by a single global phase (Figure 4B, Figure 4—figure supplement 2). Antiphase between contralateral limbs is preserved at all phases of the global oscillator (Figure 4B, Figure 4—figure supplement 3). Embedding individual trajectories illustrates a clockwise rotation during forward walking when viewed from the mouth of the bell (Figure 4C). The number of feet in stance at the central timepoint of the sampled segments oscillates as a twocycle as a function of the cyclic coordinate (Figure 4D). Consistent with the absence of multiple preferred configuration patterns, the density of points within this manifold along its axial dimension is unimodal (Figure 4—figure supplement 2). If one halves or doubles the segment length, the axial extent of the manifold is compressed or extended, intuitively corresponding to one’s ability to extract frequency information from the time series, following the frequencytime uncertainty principle (Figure 4—figure supplement 4). Thus, this manifold learning analysis of our data reveals a structure that exactly corresponds to the intuitive representation of a system of coupled oscillators with a shared frequency that varies between segments.
To test how this UMAP visualization would act on a dataset known to contain distinct gaits, we generated a synthetic dataset in which all three canonical gaits were present, each at a continuum of stepping frequencies (see Materials and methods). To generate an embedding comparison that provided evidence for the absence of distinct gaits, one would need to have a model for the structure of gait transitions. However, there is not a clear null model for such transitions. In the absence of a null model, this analysis considers only pure canonical gaits. Even when canonical tetrapod and wave gaits each accounted for only 1% of the segments, this embedding clearly showed distinct manifolds for the canonical gaits (Figure 4—figure supplement 5). Therefore, although this simulation does not provide direct evidence for the absence of multiple distinct gaits, it provides additional evidence that the manifold structure of our data is consistent with a single continuum, since data containing discrete gaits could have yielded distinct manifolds for each.
Varying one parameter in a simple model generates key attributes of observed walking patterns
Our experimental data suggest that walking speed modulation in Drosophila may be consistent with a simple control circuit with speedindependent coupling. We tested this possibility by constructing a minimal, rulebased generative model for phase dynamics (Figure 5A–B) (see Materials and methods). Conceptually, this model is related to rulebased control algorithms for walking robots, but it differs from most such models in that it generates phase dynamics rather than target limb placement positions (Kwak and McGhee, 1989; Song and Waldron, 1987; Wettergreen and Thorpe, 1992). In this model, we varied the stance duration while keeping the interlimb coupling and swing duration fixed, thus generating posteriortoanterior metachronal waves along each side of the body (Figure 5A). The interval between metachronal waves was governed by the stance duration, which is the single modulated parameter in the model. The two sides of the body were coupled such that they were in antiphase, independent of stance duration.
This simple construction recapitulated several key observations in our data. First, it can generate canonical tripod and wave gaits, as well as a continuum of intermediate coordination patterns (Figure 5B). Importantly, because the phase difference across the body is always one halfcycle, this model cannot produce the canonical tetrapod gait, but it can produce the tetrapodlike pattern with a characteristic offset in swings observed in our experimental data. With an appropriate choice of stance duration, it also generated the observed noncanonical limb pattern in which fore and hindlimbs swing simultaneously while midlimbs swing individually. With this model, the number of feet in stance phase varied as a function of the forward walking speed, as observed in experimental data (Figure 5C). The model also reproduced the measured twocycle in the number of feet down as a function of limb phase (Figure 5D, Figure 5—figure supplement 1). This twocycle occurs in the model because each ipsilateral posteriortoanterior wave is effectively a single event, and they are always exactly out of phase. The mean contralateral relative phases are constant across forward walking speeds, while the mean ipsilateral relative phases vary in a speeddependent manner, consistent with our experiments (Figure 5E). Finally, this model recapitulates the decrease in variance in the relative phases with increasing forward speed observed in our data (Figure 5F). Consistent with experimental observations, embedding the modelgenerated data using UMAP produces a single manifold (Figure 5G). Overall, this simple model with speedindependent limb coupling and a single control parameter qualitatively captures major characteristics of the symmetric variability in Drosophila spontaneous walking coordination.
Turns are achieved through asymmetric, segmentspecific modulations of limb movements
Our simple model shows that symmetric variability in limb movements can be described using a single parameter to control speed. Turning is inherently a form of asymmetric variation, as it requires that the outside legs of the fly traverse a greater distance than the inside legs. This path length differential could be achieved by modulating any combination of three parameters of limb kinematics: the frequency of stepping, the length of steps, and the direction of steps. We therefore sought to determine which of these parameters are modulated during turning. To exclude the effect of forward speed modulation on these parameters, we restricted our analysis to forward speeds between 15 and 20 mm/s, the range with the greatest quantity of turning data.
On average, the swing durations of the inside mid and hindlimbs decrease as a function of turning rate, while their stance durations increase (Figure 6A–B). The swing and stance durations of the remaining limbs are oppositely modulated (Figure 6A–B). At the highest yaw rates, the net modulation of stepping frequency is about 25%. The step length of the inside forelimb is not significantly modulated during turning (Figure 6C, Figure 6—figure supplement 1). The step lengths of the outside limbs increase with yaw rate while step lengths of the inside midlimb and hindlimb decrease with increasing turning rate. Overall, the observed limbspecific modulations of step frequency and step length are consistent with reports in other hexapods (Franklin et al., 1981; Frantsevich and Mokrushov, 1980; Graham, 1972; Gruhn et al., 2016; Jindrich and Full, 1999; Strauß and Heisenberg, 1990; Zollikofer, 1994; Zolotov et al., 1975).
The remaining kinematic parameter is the direction of limb movement during stance, which can also change the total path traversed by each limb (Figure 6D, Figure 6—figure supplement 1).
Our measurements of the two forelimbs show nearly identical modulations of stepping frequency and modest differences in step length modulations. The remaining path length differential between forelimbs is achieved through modulation of the stance direction of the inside forelimb, which is modulated more dramatically than that of other limbs, by up to about 45°. This suggests that the fly achieves the required difference in forelimb path length by directing forelimb motions away from, rather than along, its path. This is consistent with measurements in cockroaches (Jindrich and Full, 1999; Mu and Ritzmann, 2005) and in tethered stick insects (Gruhn et al., 2009). Thus, unlike forward walking speed changes, which are dominated by global modulations of a single parameter, turns are achieved through asymmetric and limbspecific modulations of multiple parameters of limb movement.
Turns are aligned to preferred phases of the limb oscillator
Next, we sought to understand whether symmetric and asymmetric patterns of limb movement are coupled. In particular, we examined whether flies execute turns at all phases of their walking pattern or if turns occur at particular limb configurations. To visualize the gross placement of limbs during turning, we estimated the spatial distributions of limb positions over all walking behaviors (Figure 7A) and at yaw rate extrema (Figure 7B). These distributions demonstrate that turns occur at a preferred spatial configuration of the walking oscillator. To understand how turns are timed relative to the underlying straightwalking coordination pattern, we examined the distributions of the instantaneous phases of each of the limbs at yaw rate extrema (see Materials and methods). Turns were preferentially aligned to a particular phase for each leg, which corresponds to a preferred phase of the walking oscillator (Figure 7C–D). In particular, there exists a statistically significant difference between the distributions of phases for all timepoints and for yaw extrema (p<10^{−5} for all limbs, see Materials and methods, Table 4). Thus, turns consist of preciselytimed modulations of the underlying walking coordination pattern, which couple the asymmetric variability to the symmetric coordination patterns.
Perturbations of forward walking speed modulate stance duration
In the analyses presented so far, we have characterized the structure of variability in spontaneous walking behavior. We next sought to test whether the structure of evoked behavior was similar to spontaneous behavior. Our experiments and model both suggest that stance duration controls forward walking speed during spontaneous behavior in Drosophila. Since flies also transiently modulate their walking speed in response to stimuli, we wondered whether this transient modulation was along the same dimension as the spontaneous behavior. For instance, while spontaneous walking speed is regulated primarily by stance duration, transient modulations of walking speed could instead be modulated by changing the step length, swing duration, or both. To investigate these possibilities, we generated transient slowing in walking flies both by directly activating a set of command neurons and by presenting visual stimuli that induced slowing.
Previous work identified a set of command neurons in a pathway named moonwalker, which regulates the speed and direction of walking (Bidaye et al., 2014). When these neurons are tonically activated, the fly walks backwards. We expressed Chrimson in these neurons and then stimulated freelywalking flies with 8 ms pulses of red light (see Materials and methods, Table 1) (Klapoetke et al., 2014). Shorttimescale activation of moonwalker neurons caused a transient reduction in forward walking speed that lasted ~100 ms (Figure 8A). This slowing is consistent with the suppression of forward walking elicited by activating a subset of the moonwalker neurons (Bidaye et al., 2014). This reduction in forward walking speed was accompanied by a transient increase in stance duration, the control parameter for speed regulation in our coordination model (Figure 8B). In individual trials, all limbs increased their stance duration over hundreds of milliseconds in response to the brief neural activation (Figure 8C). When visualized using UMAP, moonwalker activation corresponded to a deflection along the stepping frequency dimension of the manifold, consistent with a modulation of stance duration (Figure 8D). Thus, symmetrically activating a small subset of command neurons shifts the fly’s limb coordination along the same manifold as spontaneous symmetric variation.
Last, we sought to test whether more naturalistic slowing commands also modulate stance duration. Previous studies in tethered and freelywalking flies have shown that a variety of visual stimuli evoke slowing (Creamer et al., 2018; Katsov and Clandinin, 2008; Silies et al., 2013). We presented rigidly translating random dot patterns above freelywalking flies. In response to this visual stimulus, flies slowed transiently, consistent with previous studies (Figure 8E). As with the moonwalker manipulation, visual stimulation produced a transient increase in stance durations (Figure 8F–G). The visual stimulusbased perturbations also induced shifts along the frequency dimension of the UMAP manifold, corresponding to the change in stance duration (Figure 8H). Therefore, slowing evoked by visual stimuli generated transient changes along the same manifold as the spontaneous symmetric modulations that control walking speed.
Discussion
In this study, we comprehensively characterized the variability in Drosophila limb coordination during spontaneous walking behaviors. To fully sample the space of limb coordination patterns, we developed an automated method for video annotation (Figure 1). We investigated the distribution of limb coordination patterns across forward walking speeds in freelymoving flies (Figures 2 and 3). Our experimental data lie on a single manifold, consistent with a continuum of coordination patterns (Figure 4). A simple model with fixed coupling could reproduce the observed symmetric patterns of walking at a wide range of speeds by varying only stance duration (Figure 5). In contrast, asymmetric variability associated with turning is characterized by small, preciselytimed, and limbspecific modulations of this underlying straightwalking coordination pattern (Figures 6 and 7). Evoked changes in walking speed generated changes in limb coordination similar to the structure of spontaneous modulations in walking speed (Figure 8).
A single continuum of leg coordination patterns suggests a simple neural control circuit
A variety of computational models have been developed to describe patterns of hexapod locomotion (Aminzare and Holmes, 2018; Aminzare et al., 2018; Collins and Stewart, 1993; Collins and Stewart, 1994; Cruse, 1979; Cruse, 1980; Cruse, 1990; Cruse et al., 1995; Delcomyn, 1999; Graham, 1977; Ijspeert, 2008). A central component of these models is the nature of interlimb coupling, which in turn determines the predicted limb coordination patterns. Our simple model suggests that one might expect to find mutual inhibitory coupling between contralateral neuropil of the ventral nerve cord (VNC) and posterior to anterior inhibitory coupling between ipsilateral neuropil. Importantly, neither contralateral nor ipsilateral couplings need vary with walking speed. If a single continuum can describe fly walking, then models of the circuit controlling walking may be made substantially simpler, since they need not account for multiple distinct coordination patterns. Instead, a simpler control circuit need only be able to vary stance duration while maintaining ipsilateral and contralateral coupling. This suggests the intriguing possibility that walking Drosophila may use a qualitatively different locomotor control circuit than previously studied animals in which distinct gaits have been observed (Alexander, 1989; Alexander and Jayes, 1980; Alexander and Jayes, 1983; Brett and Sutherland, 1965; Rayner et al., 1986; Srinivasan and Ruina, 2006).
With a simple, oneparameter model, we recapitulated several major elements of the symmetric variability in spontaneous walking (Figure 5). However, there were aspects of our experimental data that were not captured by our simple model. Though the ipsilateral relative phases exhibited qualitatively similar speeddependent shifts, the magnitude of these shifts is larger in our model than in our data, with greater deviations at slow forward walking speeds. There are several potential explanations for this larger range. First, our experimental data included turning, which our model does not consider. Second, the relative timing of ipsilateral swings was less variable in our model than in our experimental data. Third, each swing event in the model was of constant duration while our experimentallymeasured swings were more variable in duration overall and shorter on average at slow forward walking speeds. Taken together, these differences may account for the deviations in ipsilateral phasing. These differences do not diminish the utility of this simple model, since its purpose is to illustrate the sufficiency of a single set of couplings to describe the symmetric variability in spontaneous limb coordination.
Optimality and expected gait transitions
Prior studies have suggested that the interlimb coordination patterns used by insects reflect an optimization against physical constraints (Nishii, 2000; Szczecinski et al., 2018). Theoretical work has suggested that continua of coordination patterns, like those identified here in Drosophila, minimize the energetic cost of transport when ground reaction forces are balanced across limbs (Nishii, 2000). In contrast, distinct gaits are optimal when nonuniform ground reaction forces are necessary to balance the body (Nishii, 2000). The magnitude and robustness of mechanical stability during walking may also influence the presence or absence of distinct gaits (Szczecinski et al., 2018). In addition to these mechanical considerations, we suggest that the simple control of the observed continuum of coordination patterns might make it preferable in animals with small circuits for controlling limbs.
Turning limb kinematics represent small modifications of straight walking
The asymmetric variability in limb kinematics associated with turning represents a minor modification of the coordination patterns used during forward walking, with ~25% modulations of limb movement parameters at the highest yaw rates. Detailed neuromechanical and phasereduced models in cockroaches (Kukillaya and Holmes, 2009; Kukillaya et al., 2009; Proctor and Holmes, 2018; Proctor et al., 2010) and simplified bipedal models (Proctor and Holmes, 2008) have shown that ~20% modulations in forelimb placement parameters can generate 90° degree turns within three strides. These theoretical findings are consistent with experimental results in cockroaches, which show that limb kinematics and dynamics during turning represent a minor modification of the symmetric coordination pattern (Jindrich and Full, 1999). The small magnitude of these modulations suggests a conceptual parallel between turning during walking and during flight. In flying Drosophila, detailed kinematic and dynamical measurements and modeling have shown that yaw is regulated through small modulations of wing rotation angle while gross wingbeat movements remain mostly unaffected (Dickinson and Muijres, 2016; Fry et al., 2003; Lindsay et al., 2017; Muijres et al., 2015).
Distinct control mechanisms imply distinct circuit bases for yaw and speed regulation
Mechanically, turning could be achieved through asymmetric modulations of limb movement that do not affect the coordination between ipsilateral limbs (Proctor and Holmes, 2008). In this case, the command signals regulating turning could be unilateral but not segmentspecific. However, swing duration, stance duration, step length, and step direction are all modulated on an individuallimb basis during turning. This suggests that the motor output of individual limbs in freelywalking flies is regulated during turning at the hemisegmental level by control signals, a conclusion which is consistent with electrophysiological measurements and pharmacological manipulations in dissected stick insects (Gruhn et al., 2016). In contrast to the limbspecific modulation of many kinematic degrees of freedom during turning, forward speed changes are associated primarily with modulations in stance duration that are comparable across all limbs (Figure 6). The fact that limb movements are modulated along different dimensions for yaw and speed control suggests differentiation between the circuit bases of symmetric and asymmetric variability beyond the simple requirement of asymmetric limb movement in turning. However, as turns are precisely timed relative to the walking coordination pattern (Figure 7), the control of these symmetric and asymmetric modulations must be neurally or mechanically coupled.
Differentiated roles for the posterior and anterior neuropils of the VNC
Our detailed characterization of walking suggests that, surprisingly, the observed complexity in symmetric interlimb coordination can be reproduced with fixed, unidirectional posteriortoanterior coupling between limbs. This model implies that the most posterior of the limb neuropils could have a differentiated role in walking speed regulation, as swing events originate in these posterior segments. Anatomical evidence supports this differentiated role. In particular, the moonwalker ascending neurons receive dendritic input almost exclusively from the metathoracic neuromere, with axons projecting to more anterior neuropils in the thorax and to the subesophageal ganglion (Bidaye et al., 2014). Their activation causes the fly to slow its walking (Bidaye et al., 2014), and they appear most active when all six limbs are in stance phase (Chen et al., 2018). In contrast, our kinematic measurements suggest that the forelimbs, and thus the anterior neuropil, play a distinct role in spontaneous turning, a result consistent with dynamical measurements in cockroaches (Jindrich and Full, 1999). Symmetry arguments suggest that walking speed and turning must be regulated through different coordinated patterns of neural activity. The difference in the roles of anterior and posterior neuropils in the two locomotor behaviors supports a hypothesis that these behaviors employ distinct circuitry.
Separation of action selection and motor implementation in insects
Despite the phasic relationships between limbs during locomotion and grooming, a variety of studies have demonstrated that tonic activation of descending neurons (DNs) suffices to generate complex, structured motor patterns in the fruit fly (Bidaye et al., 2014; Cande et al., 2018; Seeds et al., 2014). Interestingly, previous electrophysiological recordings in orthopteran insects showed that neurons projecting from the insect brain to the thoracic ganglia are tonically active in all regions except for those emanating from the subesophageal ganglion (Heinrich, 2002). This organization suggests that the control of motor behaviors in insects physically segregates action selection and the lowlevel implementation of motor behaviors into the brain and VNC respectively. The manifold structure of limb coordination is compatible with tonic descending inputs that modulate limb movements, while local circuits in the VNC coordinate limbs.
These experiments have shown that variability in walking coordination patterns exists on a continuous manifold. In Drosophila, explorations of VNC circuits are increasingly made possible by automated behavioral tracking (Branson et al., 2009; Mathis et al., 2018; Pereira et al., 2019; Williamson et al., 2018), by preparations to visualize neural activity in the VNC (Chen et al., 2018; Tuthill and Wilson, 2016), and by genetic tools to target neurons in the VNC (Cande et al., 2018; Mamiya et al., 2018; Namiki et al., 2018) and to manipulate neural activity (Luo et al., 2018). The manifold structure of limb coordination in Drosophila provides a framework for dissecting the circuits that regulate walking.
Materials and methods
Fly strains and husbandry
Request a detailed protocolFemale Drosophila melanogaster used in experiments were grown at 25°C on a 12 hr/12 hr light/dark cycle. Wildtype flies were staged on CO_{2 }12–24 hr after eclosion, and run after 24 hr of staging. Flies used in optogenetic experiments were grown at 25°C and staged 12–24 hr after eclosion on CO_{2}. When staged, flies used in optogenetic experiments were transferred to food supplemented with alltransretinal (ATR) following previous protocols (de Vries and Clandinin, 2013). Flies remained on ATRsupplemented food for 4 days prior to behavioral experiments. In all cases, flies were grown at near 50% relative humidity, and experiments were performed at 50% relative humidity.
Experimental setup
Request a detailed protocolBehavioral experiments were performed in an illuminated 5 cm diameter circular planar arena consisting of two plates of glass separated by 2.5 mm. The top glass plate was coated with RainX (wax) to prevent flies from walking on this surface during experiments. Above the arena, we mounted a diffusing screen and a 530 nm green Luxeon SP01G4 LED (Quadica Developments Inc, Lethbridge, AB, Canada) that provided background illumination for the arena (~1 µW/mm^{2}). Flies were recorded from below at 150 fps using a Point Grey Flea3 FL3U313Y3MC camera (FLIR Systems, Wilsonville, OR). The camera was positioned such that its field of view (2.2 × 2.75 cm) was approximately centered within the fly arena. All experiments were performed at 34°C to increase the frequency of fly walking (SotoPadilla et al., 2018).
During walking experiments, groups of 12–15 female flies were loaded into the arena. We allowed the flies to acclimate to the arena for 20 min prior to image acquisition. We then recorded the spontaneous activity of flies in the arena for 1.1 hr. Data from eight separate recordings, with a total of 114 flies, were merged to generate the aggregate wild type, freewalking dataset analyzed in this manuscript.
Feature extraction
Request a detailed protocolTo track freelywalking flies in an open arena, we decomposed the tracking process into distinct steps of centroid and footfall tracking. Our body tracking method relies on the roughly stereotyped size of fruit flies. To track the body, we binarized camera frames to separate out the darker fly pixels from the lighter background pixels. In this binarized camera frame, we identified body ellipsoids restricting our analysis to regions with roughly the correct area and eccentricity (Bradski, 2000). After aligning these ellipsoids, we then extracted stereotyped images of individual flies oriented such that the fly was aligned across all frames. We smoothed estimates of velocities by zerophase filtering with a Gaussian kernel with a standard deviation of 20 ms (three frames), which eliminated this oscillation (Figure 1—figure supplement 2) (Katsov and Clandinin, 2008; Katsov et al., 2017). Due to the height of the arena and the temporal resolution of our camera, we were, in most circumstances, able to distinguish individual flies, even when they were in close proximity. To end trajectories when fly bodies touched and became difficult to distinguish with our algorithm, we included a threshold on the maximum ellipsoid size.
Footfall tracking was performed using linear regression. To train our footfall prediction model, we manually labeled 5,000 randomly sampled oriented fly frames output by our body tracking algorithm. Exploiting the bilateral symmetry of the fly, we augmented our training data by appending the mirror symmetric images and corresponding manual labels. We used principal component analysis to reduce the dimensionality of the image space from the full number of pixels (19,881) to 1000 (Turk and Pentland, 1991). This reduction in the number of model parameters was important because we wished to train our model with a limited set of handannotated frames. We then trained 12 separate linear regression models to independently predict each footfall location variable (Pedregosa et al., 2011). To ensure that correlations between the positions of limbs did not influence the predictions made for an individual limb, we masked input variables to exclude portions of the image outside of the target limb’s range of motion prior to training (Figure 1—figure supplement 1). We fit our model with 5000 frames of annotated data, and used 10fold crossvalidation to test for overfitting. This showed that the mean Euclidean distance between predicted and hand annotated limb positions was 0.15 mm (3.5 pixels), or about 10% of stride length (Figure 1—figure supplement 1). This error distance is comparable to that reported for distal limb positions by published deep neural network methods for pose estimation (Mathis et al., 2018; Pereira et al., 2019).
Fit of functional relationship between forward velocity and stance duration (Figure 1)
Request a detailed protocolTo compare our measurements of the relationship between forward velocity and stance duration with previous studies (Mendes et al., 2013; Szczecinski et al., 2018; Wosnitza et al., 2013), we fit a power law of the form:
where ${v}_{0}=1$mm/s is a fixed scaling factor used to nondimensionalize the velocity (hence $a$ has units of milliseconds). We fit this power law using the nonlinear leastsquares solver built into MATLAB in the fit function, giving $a=932.8$ ms and $b$= –1.025 with a coefficient of determination of ${R}^{2}=0.59$.
Swing and stance determination (Figures 1–4, 6 and 8)
Request a detailed protocolMovement in the camera frame was used to determine swing and stance in each of the six limbs. Frames were categorized as stance if the smoothed instantaneous velocity of the limb was less than 20 mm/s (3.1 pixels per frame), approximately matching our error in positional measurement. Frames above this threshold were categorized as swing. Smoothing of limb velocities was performed using a fiveframe moving average filter. Both thresholds were chosen manually in order to give a reasonable representation in the presence of limb position noise.
Estimation of limb oscillator phase (Figures 2, 3, 5 and 7)
Request a detailed protocolTo estimate the phases of the limbs, we used the discretetime analytic signal method. Briefly, this method uses the Hilbert transform to construct the harmonic conjugate of a realvalued signal, and then estimates the phase as the argument of the resulting complexanalytic function. In continuous time, the Hilbert transform of a real signal $s\left(t\right)$ is related to the Fourier transform as
which allows it to be defined in discretetime by replacing the continuoustime Fourier transform by its discretetime analog. Using the Hilbert transform, one may define the analytic signal as:
From the analytic signal, the instantaneous phase is defined as
where the argument function is defined such that the instantaneous phase is a continuous function of time (Vakman and Vaĭnshteĭn, 1977; Boashash and Reilly, 1992; Gabor, 1946). We computed the analytic representations of the positions of each limb in the direction parallel to the fly’s body axis using the hilbert function in Matlab (Mathworks, Natick, MA, USA). This function implements the fastFouriertransformbased algorithm introduced in Marple (1999). The zeropoint of the resulting phase estimate corresponds to the maximally posterior point of the limb cycle. The discretetime analytic signal method was previously used to estimate the phases of cockroach limb oscillations (CouzinFuchs et al., 2015), where it was shown that this method produced comparable results to the phase estimation algorithm introduced in Revzen and Guckenheimer (2008). We smoothed the resulting phase estimates using a thirdorder SavitzkyGolay smoothing filter with a window length of 15 frames, and estimated instantaneous frequencies using a corresponding differentiating filter (Savitzky and Golay, 1964). To calculate the relative phase differences between limbs, instantaneous phase measurements were subtracted from each other at each frame and expressed in fractions of a limb cycle.
Coherence analyses (Figure 3)
Request a detailed protocolTo measure how well our data matched the configurations of relative limb phasing specified by each of the canonical gaits, we defined a metric of coherence based upon the Kuramoto model for networks of coupled oscillators (Acebrón et al., 2005). For each canonical gait, we defined a template of relative phases ${\left\{{\psi}_{k}\right\}}_{k=1}^{6}$. The template phases for each canonical gait are listed in Table 2. Then, given the set of instantaneous phases of the six limbs ${\left\{{\varphi}_{k}\left(t\right)\right\}}_{k=1}^{6}$ we defined for each template pattern the global phase $\mathrm{\Phi}\left(t\right)$ and coherence $r\left(t\right)$ by the relation
which corresponds exactly to the definition of the order parameters of the Kuramoto model with the addition of the phase template. The resulting coherence $r\left(t\right)$ ranges from a value of zero for asynchrony to a value of one for perfect alignment with the phase template.
Dimensionality reduction with UMAP (Figures 4, 5 and 8)
Request a detailed protocolTo generate lowdimensional representations of our data, we applied the nonlinear dimensionality reduction algorithm UMAP (McInnes et al., 2018). Using UMAP, we generated lowdimensional representations of 10^{5} randomly sampled segments of limb positional data, each with a halfwindow length of 100 ms. Given our sampling rate of 150 frames per second, each segment was thus an element of a 372dimensional vector space. Before embedding, each time series was meansubtracted on a persegment basis, and each timepoint of each variable was standardized across trajectories to have zero mean and unit variance. This standardization does not substantially affect the resulting embedding, since all flies are roughly the same size, with the same mean limb positions. We note also that comparable results may be obtained with segments of length 100 ms or 400 ms rather than 200 ms (Figure 4—figure supplement 4). In Figure 4B, we colored the UMAP embedding by meansubtracted limb position, a raw form of the input data.
Synthetic canonical gait data (Figure 3—figure supplement 1; Figure 4—figure supplement 5)
Request a detailed protocolTo generate synthetic data matching canonical gait patterns, we adapted a model originally designed to generate quadruped gaits (Righetti and Ijspeert, 2008). Briefly, this model is composed of six coupled nonlinear oscillators, each with two degrees of freedom. We denote the state of the i^{th} oscillator as $\left\{{x}_{i},{y}_{i}\right\}$, where ${x}_{i}$ is the position of the limb in the direction parallel to the body axis, and ${y}_{i}$ is a variable which allows for the incorporation of feedback and defines whether the limb is in swing (${y}_{i}<0$) or in stance phase (${y}_{i}>0$). In this model, uneven duty ratios are implemented by introducing the statedependent frequency
where the stance frequency is given in terms of the swing frequency ${\omega}_{\mathrm{s}\mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}$ and the duty factor $\beta $ as
Then, defining ${r}_{i}^{2}={x}_{i}^{2}+{y}_{i}^{2}$, the dynamical equations of the model are given as
where $\alpha $ is a positive parameter governing the strength of selfinteractions, $\mu $ is a positive parameter governing the radius of the limit cycle, and ${k}_{ij}$ is a coupling matrix. This system of nonlinear ordinary differential equations admits stable limit cycles for a broad range of coupling matrices (Righetti and Ijspeert, 2008). By adjusting the coupling matrix and the duty cycle appropriately, this model can produce synthetic limb traces with the duty cycle and phasing of the four canonical gaits. We integrated these equations using Matlab’s ode45 mediumorder adaptive timestep RungeKutta integrator. For each canonical gait we generated trajectories with swing durations iid from a uniform distribution corresponding to stepping frequencies between 5 and 12.5 Hz for a tripod gait. To correct for the fact that this model (with the coupling matrices chosen) is timereversal symmetric, we adjust limb ordering posthoc to enforce posteriortoanterior propagation of metachronal waves. Finally, we add white Gaussian noise with a signaltonoise ratio of 12.5 to the synthetic limb data to model the noise in our experimental data.
Numerical modeling (Figure 5)
Request a detailed protocolTo build a generative model for phase dynamics during straightwalking, we make use of the fact that the crossbody phasing is constant, and first construct a simple rulebased discretetime algorithm to generate metachronal waves along one side of the body. Following the empirical observation that the swing duration is approximately constant across all forward walking speeds, we fix the swing duration ${\tau}_{\mathrm{s}\mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}$ and have as a free parameter the stance duration ${\tau}_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}$. As the swing and stance durations are dimensionful quantities, we define a timestep $\mathrm{\Delta}t$ such that the ratios $\mathrm{\Delta}t{\tau}_{\mathrm{s}\mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}^{1}$ and $\mathrm{\Delta}t{\tau}_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}^{1}$ are dimensionless. We denote the phases of the fore, mid, and hindlimbs at timestep $t$ as ${\theta}_{1}\left(t\right)$, ${\theta}_{2}\left(t\right)$, and ${\theta}_{3}\left(t\right)$, respectively. Defining phases modulo $2\pi $, we represent swings by values less than $\pi $, and stances by values greater than $\pi $. We fix a desired total number of timesteps $N$, and an initial condition ${\left\{{\theta}_{i}\left(0\right)\right\}}_{i=1}^{3}$.
Since the metachronal waves propagate from the posterior limb to the anterior limb, we model the phase dynamics of the posterior limb as being independent from those of other limbs; at each timestep its phase is simply incremented by $\mathrm{\Delta}t\pi {\tau}_{\mathrm{s}\mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}^{1}$ during swing and by $\mathrm{\Delta}t\pi {\tau}_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}^{1}$ during stance. We now must choose a mechanism to ensure that the forelimb and hindlimb swing after their respective posterior neighbors. A simple means to do this is to reduce the rate of a given limb’s phase advance during swing phase if its anterior neighbor is also in swing phase. We make the arbitrary but simple choice that the rate of phase advance is halved. These choices lead to the feedforward rulebased algorithm to generate the phase dynamics described in pseudocode in Algorithm 1. Conceptually, this is similar to classical rulebased algorithms for generating robot limb placement targets (Wettergreen and Thorpe, 1992; Kwak and McGhee, 1989; Song and Waldron, 1987).
Algorithm 1 Rulebased generation of metachronal waves. 
Parameters $\mathrm{\Delta}t$, ${\tau}_{\mathrm{s}\mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}$, ${\tau}_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}$, $N$, ${\left\{{\theta}_{i}\left(0\right)\right\}}_{i=1}^{3}$ 
for $t=1,\text{}2,\text{}3,\text{}\dots ,\text{}N$do 
if ${\theta}_{3}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ then 
$\begin{array}{c}{\theta}_{3}\leftarrow {\theta}_{3}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{swing}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
else 
$\begin{array}{c}{\theta}_{3}\leftarrow {\theta}_{3}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{stance}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
end if 
if ${\theta}_{2}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ and ${\theta}_{3}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ then 
$\begin{array}{c}{\theta}_{2}\leftarrow {\theta}_{2}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}{2}^{1}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{swing}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
else if ${\theta}_{2}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ then 
$\begin{array}{c}{\theta}_{2}\leftarrow {\theta}_{2}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{swing}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
else 
$\begin{array}{c}{\theta}_{2}\leftarrow {\theta}_{2}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{stance}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
end if 
if ${\theta}_{1}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ and ${\theta}_{2}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ then 
$\begin{array}{c}{\theta}_{1}\leftarrow {\theta}_{1}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}{2}^{1}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{swing}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
else if ${\theta}_{2}\left(t1\right)<\pi \text{}\text{mod}\text{}2\pi$ then 
$\begin{array}{c}{\theta}_{1}\leftarrow {\theta}_{1}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{swing}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
else 
$\begin{array}{c}{\theta}_{1}\leftarrow {\theta}_{1}(t1)\phantom{\rule{thinmathspace}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\mathrm{\Delta}t\phantom{\rule{thinmathspace}{0ex}}\pi \phantom{\rule{thinmathspace}{0ex}}{\tau}_{stance}^{1}\phantom{\rule{thinmathspace}{0ex}}\text{mod}\phantom{\rule{thinmathspace}{0ex}}2\pi \end{array}$ 
end if 
end for 
We now note that Algorithm 1 simply implements a forward Euler method for integrating the following system of ordinary differential equations:
where we have implemented the discretetime rules using indicator functions for swing and stance phase, defined as
and
respectively. As shown in Figure 5, this simple model suffices to generate metachronal waves that appear qualitatively similar to those observed empirically. We note that this system of piecewiselinear ODEs could be replaced with a system that uses continuous, differentiable approximations for the indicator functions. However, this would introduce additional parameters governing the sharpness of the approximation. As the piecewisecontinuous system is wellbehaved but for the discontinuities, it may be numerically integrated with satisfactory accuracy using a fixedstep loworder RungeKutta method. Therefore, we chose to integrate the discontinuous rulebased model.
To generate sixlegged phase dynamics, we coupled together two such generative models for metachronal waves. Since crossbody antiphase is conserved across forward walking speeds, our model is chiefly concerned with the generation of metachronal waves. The nature of the crossbody coupling is not central to our model; rather, it is a convenience used to generate phase dynamics for all the limbs in a feedforward manner. With this in mind, we chose to enforce crossbody antiphase by adding additional divisive terms depending on the sine of the contralateral phase differences. In doing so, we introduce an additional parameter $0<\alpha <1$ governing the strength of the crossbody coupling. In Figure 5, a value of $\alpha =\frac{1}{8}$ was used, though we find empirically that varying the value of $\alpha $ between $\frac{1}{8}$ and $\frac{1}{2}$ does not produce qualitative changes in the behavior of the model. Then, denoting phases of the left forelimb, left midlimb, left hindlimb, right forelimb, right midlimb, and right hindlimb as ${\theta}_{1}$ through ${\theta}_{6}$, respectively, the dynamical equations of the sixlimb model are
We integrated this system of ODEs with the parameters given in Table 3 using Heun’s method, a secondorder explicit RungeKutta integrator (Ascher and Petzold, 1998) with a timestep of 0.025 ms, as it provides a reasonable balance of computational efficiency and numerical stability.
Optogenetic manipulations (Figure 8)
Request a detailed protocolFor optogenetic experiments, we augmented the arena used in freewalking experiments by mounting a DLP Lightcrafter 4500 projector (Texas Instruments, Dallas, TX) below the walking surface. During experiments, flies were loaded into the arena and allowed to acclimate for 20 min prior to recording. For each experiment, 12–15 females were allowed to walk freely within the arena and were recorded from below for 1.1 hr. Throughout the recording, at 4 s intervals, 8 ms fullfield flashes of ~ 0.05 mW/mm^{2} red light (peak 624 nm) were projected into the arena in order to optogenetically activate moonwalker neurons of the flies. In both experimental and random trigger control conditions, we selected trajectories where the fly was walking > 5 mm/s on average in the 100 ms prior to stimulation. Additionally, in order to exclude stopping flies, we removed trajectories in which the flies’ forward velocity was less than 0.1 mm/s at any point post stimulus. Data from seven separate recordings, with a total of 96 flies, were merged to generate the aggregate dataset analyzed in this manuscript.
Visual stimuli (Figure 8)
Request a detailed protocolVisual stimuli were designed using Matlab and its associated Psychophysics Toolbox (Brainard, 1997; Kleiner et al., 2007; Pelli, 1997). They were projected onto a virtual flat torus masked with a disk placed on a flat screen 29.5 mm above the freewalking arena using a DLP Lightcrafter projector (Texas Instruments, Dallas, TX). The spatial resolution of the screens was ~ 0.3° and the screen update rate was 180 Hz. All visual stimuli used green light (peak 520 nm) with mean luminance of 100 cd/m^{2}. The video data were aligned to the stimulus by presenting periodic flashes measured with a photodiode. All stimuli were presented in 50 ms bouts separated by 950 ms interleaves. The stimulus presented consisted of an array of 15°diameter white and black circular dots, placed on a meangrey background at locations iid from a uniform distribution with a density of 10%. Overlap between dots was handled such that the contrast of the stimulus did not increase (two white dots add to white, two black dots sum to black, and a white and black dot sum to gray). The resulting image was translated rigidly in a randomly selected cardinal direction during stimulation, and was static during the interleave. Three velocities of translation, 720, 960, and 1440 ^{o}/s, were grouped and used as the visual stimulus condition in Figure 8. Trajectories were selected in a manner identical to the optogenetic experiments described above. For these experiments, the flies were backlit with infrared light (peak 850 nm) at an intensity of ~ 2 µW/mm^{2}. Data from three separate recordings, with a total of 44 flies, were merged to generate the aggregate dataset analyzed in this manuscript.
Statistics
All statistical analysis was conducted using Matlab 9.4 (The MathWorks, Natick, MA). Throughout this manuscript, presented error bars are bootstrapped 95% confidence intervals obtained using the biascorrected and accelerated percentile method (Efron, 1987) unless otherwise indicated. In Figure 1, error bars for both body and limb variables were generated from bootstraps over videos (eight videos total from eight experiments, each with approximately 1/8 of the total data). In Figures 3, 6 and 7, error bars for analyses of relative limb phases and coherences were also generated by bootstrapping over videos. In Figure 8, unlike other analyses, the confidence intervals were generated by bootstrapping over individual selected fly walking trajectories. In Figure 8, distributions of stance durations were compared using a twosample KolmogorovSmirnov test. Throughout this manuscript, circular probability density functions of limb relative phases were estimated using kernel methods (Batschelet, 1981; Fisher, 1989). In Figure 7, to test against the null hypothesis that the distribution of phases at yaw extrema is indistinguishable from that over all instants, we performed twosample Monte Carlo resampling tests (Dwass, 1957; Gandy, 2009) using the Kuiper Vstatistic (Batschelet, 1981; Kuiper, 1960). We used 10^{5} permutations, and computed ClopperPearson confidence intervals for the resulting pvalues (Clopper and Pearson, 1934; Dwass, 1957; Gandy, 2009) (Table 4). As the value of the test statistic for the null distribution did not exceed the observed test statistic in 10^{5} permutations, this test only allows us to bound the true pvalue from above.
Data and code availability
Request a detailed protocolThe datasets analyzed in this work are available from the Dryad Digital Repository: https://doi.org/10.5061/dryad.3p9h20r. Matlab code to integrate all models, reproduce all statistical analyses, and generate all figure panels is available at https://github.com/ClarkLabCode/GaitPaperCode (DeAngelis et al., 2019; copy archived at https://github.com/elifesciencespublications/GaitPaperCode).
Data availability
The datasets analyzed in this work are available from the Dryad Digital Repository: https://doi.org/10.5061/dryad.3p9h20r. Matlab code to reproduce all statistical analyses and figure panels is available at https://github.com/ClarkLabCode/GaitPaperCode (copy archived at https://github.com/elifesciencespublications/GaitPaperCode).

Dryad Digital RepositoryData from: The manifold structure of limb coordination in walking Drosophila.https://doi.org/10.5061/dryad.3p9h20r
References

The Kuramoto model: a simple paradigm for synchronization phenomenaReviews of Modern Physics 77:137–185.https://doi.org/10.1103/RevModPhys.77.137

Optimization and gaits in the locomotion of vertebratesPhysiological Reviews 69:1199–1227.https://doi.org/10.1152/physrev.1989.69.4.1199

Fourier analysis of forces exerted in walking and runningJournal of Biomechanics 13:383–390.https://doi.org/10.1016/00219290(80)900196

A dynamic similarity hypothesis for the gaits of quadrupedal mammalsJournal of Zoology 201:135–152.https://doi.org/10.1111/j.14697998.1983.tb04266.x

Gait transitions in a phase oscillator model of an insect central pattern generatorSIAM Journal on Applied Dynamical Systems 17:626–671.https://doi.org/10.1137/17M1125571

BookComputer Methods for Ordinary Differential Equations and DifferentialAlgebraic EquationsSIAM: Society for Industrial and Applied Mathematics.

The comparative investigation of the stick insect and cockroach models in the study of insect locomotionCurrent Opinion in Insect Science 12:1–10.https://doi.org/10.1016/j.cois.2015.07.004

Dimensionality reduction for visualizing singlecell data using UMAPNature Biotechnology 37:38–44.https://doi.org/10.1038/nbt.4314

Kinematic and behavioral evidence for a distinction between trotting and ambling gaits in the cockroach blaberus discoidalisJournal of Experimental Biology 214:2057–2064.https://doi.org/10.1242/jeb.056481

Speeddependent interplay between local patterngenerating activity and sensory signals during walking in DrosophilaThe Journal of Experimental Biology 219:3781–3793.https://doi.org/10.1242/jeb.146720

Mapping the stereotyped behaviour of freely moving fruit fliesJournal of the Royal Society Interface 11:20140672.https://doi.org/10.1098/rsif.2014.0672

Time Frequency Signal Analysis: Methods and Applications163–168, Algorithms for TimeFrequency Signal Analysis, Time Frequency Signal Analysis: Methods and Applications, Longman Cheshire, Melbourne.

Highthroughput ethomics in large groups of DrosophilaNature Methods 6:451–457.https://doi.org/10.1038/nmeth.1328

Respiratory Metabolism of Pumpkinseed ( Lepomis gibbosus ) in Relation to Swimming SpeedJournal of the Fisheries Research Board of Canada 22:405–409.https://doi.org/10.1139/f65039

Organizing network action for locomotion: insights from studying insect walkingBrain Research Reviews 57:162–171.https://doi.org/10.1016/j.brainresrev.2007.06.028

Hexapodal gaits and coupled nonlinear oscillator modelsBiological Cybernetics 68:287–298.https://doi.org/10.1007/BF00201854

A grouptheoretic approach to rings of coupled biological oscillatorsBiological Cybernetics 71:95–103.https://doi.org/10.1007/BF00197312

Intersegmental coupling and recovery from perturbations in freely running cockroachesJournal of Experimental Biology 218:285–297.https://doi.org/10.1242/jeb.112805

A new model describing the coordination pattern of the legs of a walking stick insectBiological Cybernetics 32:107–113.https://doi.org/10.1007/BF00337442

A quantitative model of walking incorporating central and peripheral influences II. the connections between the different legsBiological Cybernetics 37:137–144.

What mechanisms coordinate leg movement in walking arthropods?Trends in Neurosciences 13:15–21.https://doi.org/10.1016/01662236(90)90057H

A modular artificial neural net for controlling a sixlegged walking systemBiological Cybernetics 72:421–430.https://doi.org/10.1007/BF00201417

Tight turns in stick insectsJournal of Comparative Physiology A 195:299–309.https://doi.org/10.1007/s0035900804063

Optogenetic stimulation of escape behavior in Drosophila melanogasterJournal of Visualized Experiments 50192.https://doi.org/10.3791/50192

Walking robots and the central peripheral of locomotion in insectsAutonomous Robots 7:259–270.

The aerodynamics and control of free flight manoeuvres in DrosophilaPhilosophical Transactions of the Royal Society B: Biological Sciences 371:20150388.https://doi.org/10.1098/rstb.2015.0388

Curve walking in freely moving crayfish (Procambarus clarkii)The Journal of Experimental Biology 201:1315–1329.

Motor flexibility in insects: adaptive coordination of limbs in locomotion and nearrange explorationBehavioral Ecology and Sociobiology 72:15.https://doi.org/10.1007/s0026501724123

The behavioural transition from straight to curve walking: kinetics of leg movement parameters and the initiation of turningJournal of Experimental Biology 208:2237–2252.https://doi.org/10.1242/jeb.01637

Modified randomization tests for nonparametric hypothesesThe Annals of Mathematical Statistics 28:181–187.https://doi.org/10.1214/aoms/1177707045

Better bootstrap confidence intervalsJournal of the American Statistical Association 82:171–185.https://doi.org/10.1080/01621459.1987.10478410

Smoothing a sample of circular dataJournal of Structural Geology 11:775–778.https://doi.org/10.1016/01918141(89)900126

Rotational locomotion by the cockroach blattella germanicaJournal of Insect Physiology 27:249–255.https://doi.org/10.1016/00221910(81)900585

Turning and righting inGeotrupes (Coleoptera, scarabaeidae)Journal of Comparative Physiology ? A 136:279–289.https://doi.org/10.1007/BF00657348

Sequential implementation of monte carlo tests with uniformly bounded resampling riskJournal of the American Statistical Association 104:1504–1511.https://doi.org/10.1198/jasa.2009.tm08368

Saccadic body turns in walking DrosophilaFrontiers in Behavioral Neuroscience 8:1–10.https://doi.org/10.3389/fnbeh.2014.00365

A behavioural analysis of the temporal organisation of walking movements in the 1st instar and adult stick insect (Carausius morosus)Journal of Comparative Physiology 81:23–52.https://doi.org/10.1007/BF00693548

The effect of amputation and leg restraint on the free walking coordination of the stick insectCarausius morosusJournal of Comparative Physiology ? A 116:91–116.https://doi.org/10.1007/BF00605519

Straight walking and turning on a slippery surfaceJournal of Experimental Biology 212:194–209.https://doi.org/10.1242/jeb.018317

Impact of descending brain neurons on the control of stridulation, walking, and flight in orthopteraMicroscopy Research and Technique 56:292–301.https://doi.org/10.1002/jemt.10033

Recovery of locomotion after injury in Drosophila depends on proprioception and de bivort, BThe Journal of Experimental Biology 219:1760–1771.https://doi.org/10.1242/jeb.133652

Manylegged maneuverability: dynamics of turning in hexapodsThe Journal of Experimental Biology 202:1603–1623.

Legtracking and automated behavioural classification in DrosophilaNature Communications 4:1910–1918.https://doi.org/10.1038/ncomms2908

Independent optical excitation of distinct neural populationsNature Methods 11:338–346.https://doi.org/10.1038/nmeth.2836

Tests concerning random points on a circleIndagationes Mathematicae 63:38–47.https://doi.org/10.1016/S13857258(60)500060

Neuromechanical models for insect locomotion: stability, maneuverability, and proprioceptive feedbackChaos: An Interdisciplinary Journal of Nonlinear Science 19:026107.https://doi.org/10.1063/1.3141306

A model for insect locomotion in the horizontal plane: feedforward activation of fast muscles, stability, and robustnessJournal of Theoretical Biology 261:210–226.https://doi.org/10.1016/j.jtbi.2009.07.036

Rulebased motion coordination for a hexapod walking machineAdvanced Robotics 4:263–282.https://doi.org/10.1163/156855390X00297

Computing the discretetime "analytic" signal via FFTIEEE Transactions on Signal Processing 47:2600–2603.https://doi.org/10.1109/78.782222

DeepLabCut: markerless pose estimation of userdefined body parts with deep learningNature Neuroscience 21:1281–1289.https://doi.org/10.1038/s415930180209y

Kinematics and motor activity during tethered walking and turning in the cockroach, blaberus discoidalisJournal of Comparative Physiology A 191:1037–1054.https://doi.org/10.1007/s003590050029x

Body saccades of Drosophila consist of stereotyped banked turnsJournal of Experimental Biology 218:864–875.https://doi.org/10.1242/jeb.114280

Legged insects select the optimal locomotor pattern based on the energetic costBiological Cybernetics 83:435–442.https://doi.org/10.1007/s004220000175

SoftwareScikitlearnMachine Learning in Python.

Fast animal pose estimation using deep neural networksNature Methods 16:117–125.https://doi.org/10.1038/s4159201802345

A phasereduced neuromechanical model for insect locomotion: feedforward stability and proprioceptive feedbackPhilosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 368:5087–5104.https://doi.org/10.1098/rsta.2010.0134

Steering by transient destabilization in piecewiseholonomic models of legged locomotionRegular and Chaotic Dynamics 13:267–282.https://doi.org/10.1134/S1560354708040047

Estimating the phase of synchronized oscillatorsPhysical Review E 78:51907.https://doi.org/10.1103/PhysRevE.78.051907

ConferencePattern generators with sensory feedback for the control of quadruped locomotionIEEE International Conference on Robotics and Automation.

Smoothing and differentiation of data by simplified least squares proceduresAnalytical Chemistry 36:1627–1639.https://doi.org/10.1021/ac60214a047

An analytical approach for gait study and its applications on wave gaitsThe International Journal of Robotics Research 6:60–71.https://doi.org/10.1177/027836498700600205

Thermosensory perception regulates speed of movement in response to temperature changes in Drosophila melanogasterThe Journal of Experimental Biology 221:jeb174151.https://doi.org/10.1242/jeb.174151

Interlimb coordination during slow walking in the cockroach: I. effects of substrate alterationsThe Journal of Experimental Biology 78:233–243.

Coordination of legs during straight walking and turning in Drosophila melanogasterJournal of Comparative Physiology A 167:403–412.https://doi.org/10.1007/BF00192575

Static stability predicts the continuum of interleg coordination patterns in DrosophilaThe Journal of Experimental Biology 221:jeb189142.https://doi.org/10.1242/jeb.189142

Adaptations to changing speed in human locomotion: speed of transition between walking and runningActa Physiologica Scandinavica 131:211–214.https://doi.org/10.1111/j.17481716.1987.tb08228.x

Eigenfaces for recognitionJournal of Cognitive Neuroscience 3:71–86.https://doi.org/10.1162/jocn.1991.3.1.71

Amplitude, phase, frequency—fundamental concepts of oscillation theorySoviet Physics Uspekhi 20:1002–1016.https://doi.org/10.1070/PU1977v020n12ABEH005479

Laufen und stehen der stabheuschrecke carausius morosus: sinnesborstenfelder in den beingelenken als glieder von regelkreisenZeitschrift Für Vergleichende Physiologie 48:198–250.https://doi.org/10.1007/BF00297860

ConferenceGait generation for legged robotsIn IEEE International Conference on Intelligent Robots and Systems.https://doi.org/10.1109/IROS.1992.594568

Insect walkingAnnual Review of Entomology 11:103–122.https://doi.org/10.1146/annurev.en.11.010166.000535

Interleg coordination in the control of walking speed in DrosophilaJournal of Experimental Biology 216:480–491.https://doi.org/10.1242/jeb.078139

Stepping patterns in ants  INFLUENCE of speed and curvatureThe Journal of Experimental Biology 192:95–106.

Kinematik Der Phototaktischen Drehung bei der HonigbieneApis mellifera LJournal of Comparative Physiology 97:339–353.https://doi.org/10.1007/BF00631970
Decision letter

Ronald L CalabreseSenior and Reviewing Editor; Emory University, United States
In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.
Thank you for submitting your article "The manifold structure of limb coordination in walking Drosophila" for consideration by eLife. Your article has been reviewed by 3 peer reviewers, and the evaluation has been overseen by Ronald Calabrese as the Senior and Reviewing Editor. The reviewers have opted to remain anonymous.
The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.
Summary:
This work investigates the repertoire of forward locomotion dynamics in the fruit fly (D. melanogaster) using a combination of limb tip and centroid tracking, dimensionality reduction and mathematical modeling, and neuronal perturbations. The key claims are that the speed of forward locomotion are controlled using a single control parameter, and that the traditional division of walking patterns into alternating tripod, tetrapod, and wave gaits is not apparent in spontaneous walking and is, accordingly, better understood through the manifold framework proposed here. Lastly, the authors perform both optogenetic and sensory perturbations during spontaneous locomotion, finding that the flies altered their movements along the trajectory of the manifold described.
Essential revisions:
There is much to interest the fly community and the greater neuroscience community in the thorough analysis presented and the interesting manifold structure described.
Several concerns detailed in the expert reviews should be addressed.
1) The UMAP embedding (dimensionality reduction) is wellsuited for preserving topological features of a data set requiring warping of global metrics in order to preserve local topological structure. PCA may do a better job of describing the full dynamic range of the variables than the UMAP does. The authors should directly compare PCA to the UMAP manifold and should be very clear about what features of the data UMAP is preserving.
2) The effect of error on gait definitions (e.g. those in Figure 3) is not explored/discussed. The authors make stance/swing definitions based on the smoothed instantaneous velocity of the limb (with the tracking error threshold chosen as the separation threshold). While this is a sensible approach given the type data obtained here there is no analysis or discussion about how error in this value might affect the results found.
3) Previous studies have shown that limb coordination patterns in Drosophila are highly variable and that while there are changes with speed that the overall shifts are more continuous than observed in many other animals. The authors would be better served by not emphasizing how this work shows that flies do not exhibit discrete gaits and emphasizing more the structure of the manifold discovered and how it is relevant to locomotion.
4) The paper would be substantially improved if the coupled phase oscillator model were motivated and presented more clearly.
5) A common feature in the gait transition literature is hysteresis – the transition occurs at a different speed going from slowto fast than fasttoslow. Either adapting the analysis to include this potential effect (or explaining why it is not relevant here) is a necessary addition to the manuscript.
6) The spontaneous turning data from your other eLife submission and its relationship to the manifold described would enhance this paper and further take the emphasis off the current discrete gaits controversy. This is not a requirement but a suggestion.
Reviewer #1:
Note: I am also a reviewer of your other eLife submission, so several of my comments will mirror those for that paper, due to the similarity of the methods used.
The manuscript by DeAngelis and colleagues investigates the repertoire of forward locomotion dynamics in the fruit fly (D. melanogaster) using a combination of limb tip and centroid tracking, dimensionality reduction, mathematical modeling, and neuronal perturbations. The key claims are that the speed of forward locomotion are controlled using a single control parameter, and that the traditional division of walking patterns into alternating tripod, tetrapod, and wave gaits is not apparent in spontaneous walking (and is, accordingly, better understood through the manifold framework proposed here). Lastly, the authors perform both optogenetic and sensory perturbations during spontaneous locomotion, finding that the flies altered their movements along the trajectory of the found manifold.
Overall, I thought that this manuscript was creatively performed and would be of general interest to the community. I do have some serious reservations about the findings, however, and would like to see some further analysis of the data before recommending publication. Specifically, I have questions about the authors' choice of embedding to find the manifold, the effects of footfall error on their results and interpretation of their results in light of gait hysteresis. If these concerns are resolved satisfactorily, however, I envision that this paper will be wellsuited to publication in eLife.
1) After reading and digesting the paper (as well as the cosubmitted paper), I remain a bit perplexed why the authors decided to use UMAP for the dimensionality reduction technique upon which many of their results lie. UMAP is wellsuited for preserving topological features of a data set requiring warping of global metrics in order to preserve local topological structure. Despite this, the authors only use UMAP as a global coordinate, treating it like one often observed PCA (which does the opposite, preserving global metrics while warping local ones), pointing at directions where yaw or lateral or forward velocities alter. In fact, the results from Figure 2—figure supplement 1 in ZavatoneVeth et al. show that PCA actually does a better job of describing the full dynamic range of the variables than the UMAP does. As a result, it also would be interesting to see what PCA provides for the kinematic variables.
I don't make these points as merely a technicality, however. As many of the results in the paper come from this embedding, it is important to actually note what the structure of this manifold is. Moreover, the density within the manifold is not shown in any of the plots as well. Are there islands of high density with valleys in between? Or is It smoothly distributed? I think that more care is needed to more fullycharacterize the manifold that is the one of the key mathematical descriptions in this manuscript (and is in fact in the title!).
2) Perhaps my largest concern with the paper, however, is the effect of error on gait definitions (e.g. those in Figure 3). The authors make stance/swing definitions based on the smoothed instantaneous velocity of the limb (with the tracking error threshold chosen as the separation threshold). While this is a sensible approach given the type data that they obtain here, as far as I saw in this manuscript at any rate, I did not see any discussion about how error in this value might affect the found results. This is rather important, as many of the disagreements with the previous literature about the discreteness of individual gaits rely on these measurements. I think that it would be important to see how the results change if the threshold for swing/stance was changed within some reasonable limits, or if swing/stance was assigned using a confidence (using, say, a logistic function of the limb speed). Specifically, I would be curious to see if the 2wave vs. 3wave and 6wave results could be altered by changing this threshold.
That being said, the values for seem that far off. The results from Mendes et al., 2013, for example, are built upon a much more accurate measurement of swing/stance determination, since the footfalls are directly measured using frustrated total internal reflection (although the authors should please correct me if they feel this assessment is incorrect). The results in Figure 1 here are similar to those found in Figure 2 of that paper, so no major tweaking is likely needed, but some assessment of the effect of potential error is important. Relatedly, it might be useful to see what happens to the "tripod index" as defined in the Medes paper with the data collected here. Are the results replicated?
3) Related to both the previous comments, I think that in some senses, this paper is offering a potentiallyfalse dichotomy between discrete and continuous gaits. For example, although the fly may have some preferred forward velocities or gaits that are peaks in the distribution, there may be nonzero valleys between them that, although rare, are important parts of the overall manifold that describes locomotion. For example, in Figure 4—figure supplement 1, the separate regions are potentially the result of the fact that there is no intermediate data. Because UMAP explicitly exaggerates topological differences, my guess is that the inclusion of a relatively small number of intermediate gaits will link the discrete rings together into something that looks more like the real data. I also worry that the phase error potentiallyintroduced from the uncertainty of the stance/swing transition time might affect the phase difference peaks in Figure 3 and Figure 3—figure supplements 1 and 2.
4) A common feature in the gait transition literature is hysteresis – the transition occurs at a different speed going from slowto fast than fasttoslow. All of the analyses here fundamentally assume that gait is a singlevalued function of speed, however, which may cause for some misinterpretations of the data, since all spread of the distribution would just be assigned to a single mean value. Either adapting the analysis to include this potential effect (or explaining why it is not relevant here) is a necessary addition to the manuscript.
Reviewer #2:
In this paper, DeAngelis and colleagues investigate the temporal organization of the fly walking gait. They use a simple tracking algorithm to measure the body orientation and distal leg position of flies walking in a flat, featureless environment. They confirm previous observations that stance duration changes as a function of walking speed. They then extracted swing and stance phases from each leg and used these patterns to classify different walking gaits. They find that flies typically have three feet in stance (tripod gait), but they observe that flies occasionally (and transiently) have 4 or 5 feet in stance, particularly at low walking speeds. Contralateral limbs reliably swing in antiphase, while ipsilateral limb coordination vary a bit more. A lowdimensional embedding visualization of limb coordinates resulted in cyclic (global oscillator phase) and one linear coordinate (stepping frequency). Overall, this suite of analyses suggest that there is no clear differentiation of walking behavior into distinct gaits (e.g., tripod vs. tetrapod). The authors then construct a model that can predict the continuum of walking behavior by tuning of a single parameter. Finally, they find that perturbing walking behavior through activation of the Moonwalker Descending Neuron or visual stimulation alters stance duration, and thus induces a shift along the length of the walking manifold.
Overall, the paper is clearly organized and wellwritten. Although the subject of the fly walking gait is already welltrodden, this is certainly the most exhaustive and definitive treatment to date. I think it is a very strong paper and have only a handful of suggestions to improve it prior to publication.
1) My biggest point of confusion was the coupled phase oscillator model. The paper would be improved if this model were motivated and presented more clearly. The equations are not intuitive at first glance (e.g., the 1/sin terms, the discontinuous indicator functions); is there an easy interpretation of them? And does the specific form of the model matter? How is it related to the other coupled oscillator models in the literature? I would have liked to see a sensitivity analysis that provides some intuition for how the model works and which parameters matter. For example, τ_{stance} is the relevant parameter that the authors tune to get different interlimb coordination patterns, but this is not clearly explained (especially in the main text). Finally, how does this model contribute to our understanding or intuition of the control strategy used by the fly to coordinate walking? The authors should include a deeper discussion of the motivation for and interpretation of this model.
Reviewer #3:
This paper is notable for collecting a very comprehensive dataset of kinematic gait patterns in Drosophila over a wide range of spontaneous and evoked speed changes. While the optogenetic and visually induced changes are nice as supporting evidence, the core of the paper comes from a detailed and thorough gait analysis. As such the paper has some very good insights when connected to modern dimensionality reduction and some simple models of oscillations, but needs to distinguish itself carefully from the long history of papers that also do gait analysis especially in Drosophila.
1) The structure of the argument is that animals in general are thought to have distinct gaits. The authors cite a number of studies on larger legged animals for this, but in those cases walking is an out of phase oscillation of gravitational potential energy and kinetic energy (the inverted pendulum model) and running is inphase oscillations (the springloaded inverted pendulum). The transition in gaits in many animals from insects to horses corresponds to the transition between these underlying dynamics If this transition is unlikely (although possible?) for Drosophila given their size, it is not clear why we would expect distinct gaits even if larger animals have them.
More critical to the conclusions drawn about Drosophila, previous studies (Wosnitza et al., 2012, and Mendes et al., 2013) My reading of these previous studies is that they show that limb coordination patterns in Drosophila are highly variable and that while there were changes with speed that the overall shifts are more continuous than observed in many other animals. These are cited in the paper but the relationship of their conclusions to this paper are only discussed superficially. For example Wosnitza et al., 2012 conclude that their findings "imply that Drosophila's walking behavior is more flexible than previously thought (Strauss and Heisenberg, 1990): there are no clearly separable gaits and, more specifically, the neural controller producing interleg coordination is not restricted to a fixed tripod pattern." This is very close to the initial conclusions of this paper.
So while some data from those studies do fit tetrapod and tripod gaits (as do some in this study), the conclusion that Drosophila use a variety kinematic patterns that encompass the canonical gaits is already somewhat established. Therefore the hypothesis of distinct gaits in this paper seems like a bit of a strawman. The authors do make some nice analyses of these data and show the point perhaps more convincingly than in previous work. However, prior to the manifold section, the initial presentation of the data did not leave me thinking that I had learned something new. I recognize that the authors may not agree with interpretations in previous works, but the onus is on the authors to clearly articulate what cannot be concluded from the previous work that their work in turn justifies.
More importantly this whole argument of distinct gaits takes away from what for me was the real impact of the paper. The central advance of the paper is the finding that kinematics patterns across a large range of speeds can collapse down onto a manifold with a single parameter governing speed variation. The authors are still making an advance here, but the motivation of the introduction and the conclusion of the first sentence of the Discussion highlight the continuum vs. discrete nature of Drosophila gaits. The emphasis should be about the manifold.
2) My concern about the manifold discussion is how this approach is different from the many other dimensionality reductions that have been done on gait. What really sets it apart from the tSNE analyses of behavior that show all walking lying in a low dimensional cluster (Berman et al., 2014) and the reduced phase and oscillator models of stick insects and cockroaches (e.g. CouzinFuchs et al., 2015). I think there is novelty here but it must be articulated more in context.
https://doi.org/10.7554/eLife.46409.sa1Author response
Essential revisions:
There is much to interest the fly community and the greater neuroscience community in the thorough analysis presented and the interesting manifold structure described.
Several concerns detailed in the expert reviews should be addressed.
1) The UMAP embedding (dimensionality reduction) is wellsuited for preserving topological features of a data set requiring warping of global metrics in order to preserve local topological structure. PCA may do a better job of describing the full dynamic range of the variables than the UMAP does. The authors should directly compare PCA to the UMAP manifold and should be very clear about what features of the data UMAP is preserving.
Thank you for this suggestion. We initially considered using PCA to visualize our locomotion dataset. Unfortunately, PCA does not provide a particularly clarifying description of the overall structure of variability in limb movement patterns. We have added a supplementary figure (Figure 4—figure supplement 1) to illustrate this point. The principal components of our limb kinematic data occur in degenerate, approximately phaseshifted pairs. Therefore, though the projection of these data into the first two principal components provides some information about their phasic structure, PCA does not allow for easy visualization of the full structure of variability in two or three dimensions. Specifically, the relationship between the projection and the frequency of stepping is multivalued.
We would like to highlight that the principal component analysis presented in the companion manuscript was performed on the body kinematic variables. Due to longtimescale correlations among those coordinates, the body kinematic data are approximately threedimensional, even when segments of trajectories are embedded. Therefore, projection into two dimensions using PCA yields more informative results than the application of PCA to the higherdimensional limb data. We have adjusted the prose to address the difficulties encountered with PCA and tSNE, another commonly used dimensionalityreduction technique.
We acknowledge that UMAP is a relatively new dimensionality reduction technique. To further clarify the features of the data that are preserved by UMAP, we have added a supplementary figure (Figure 4—figure supplement 2) including the additional visualizations requested by Reviewer #1. We have also expanded the prose to more clearly describe the optimization performed by the UMAP algorithm. We have added citations to several other publications, notably (Becht et al., 2019), in which it was applied to reveal the global topology of highdimensional gene expression data, consistent with its usage in this manuscript.
2) The effect of error on gait definitions (e.g. those in Figure 3) is not explored/discussed. The authors make stance/swing definitions based on the smoothed instantaneous velocity of the limb (with the tracking error threshold chosen as the separation threshold). While this is a sensible approach given the type data obtained here there is no analysis or discussion about how error in this value might affect the results found.
Thank you for this suggestion. We have taken two steps to address the sensitivity of the swing/stance analyses. First, the analyses presented in Figure 3 are independent of our determination of swing and stance, as the instantaneous phase measurements used throughout the figure are estimated directly from raw limb position time series. We have modified the prose to state this more clearly. Second, in Figure 2, the analyses presented do depend on swing/stance determinations. If one believed that measurement noise strongly affected these analyses, one would expect that lowering the separation threshold would accentuate the effect. To address this concern, we replicated Figure 2 using a separation threshold of 10 pixels per frame rather than 20 pixels per frame as in the main figure. As shown in Author response image 1, the twocycle in the instantaneous number of feet in stance is preserved even under halving of the separation threshold. Thus, we believe this analysis is relatively insensitive to measurement noise in our swingstance definition. We have updated the prose in the manuscript to describe this sensitivity analysis.
Following this comment, we also wanted to test whether our phasebased analyses were sensitive to noise in the limb measurements. To test the noise sensitivity of the discretetime analytic signal method we used for phase estimation, we simulated a set of trajectories using our canonical gait model. As this phase estimation method is known to perform best with waveforms with duty cycle equal to 1/2 (Boashash and Reilly, 1992; Lamb and Stöckl, 2014), we chose to simulate canonical right tetrapod gaits, which have a duty cycle equal to 2/3, so that we tested the algorithm where it might not work as well. We then corrupted the resulting simulated limb position data with additive white Gaussian noise to yield signaltonoise ratios (std of signal/std of noise) of 12.5 (matching Figure 3—figure supplement 1 and Figure 4—figure supplement 4; also roughly matching our estimate of the experimental signal to noise) and 1.25. Following the protocol applied to our experimental data, we then estimated the phase from these noisy traces and smoothed the resulting estimate using a thirdorder SavitzkyGolay filter. As shown in Author response image 2, the spread of the relative phase distribution increased with decreasing SNR, but the location and number of modes were not qualitatively altered. Therefore, our phase estimation method appears to be relatively robust to noise.
3) Previous studies have shown that limb coordination patterns in Drosophila are highly variable and that while there are changes with speed that the overall shifts are more continuous than observed in many other animals. The authors would be better served by not emphasizing how this work shows that flies do not exhibit discrete gaits and emphasizing more the structure of the manifold discovered and how it is relevant to locomotion.
We agree that previous studies have presented evidence that limb coordination patterns in Drosophila are highly variable, and believe that the primary contribution of our study is to characterize the structure of that variability. We have adjusted the prose throughout the manuscript to better reflect this understanding. In the literature, the structure of variability in walking patterns in many animals has most often been studied using the framework of distinct gaits, so we believe it is important to address the question of gaits. Indeed, this framework has persisted in the descriptions (Aminzare and Holmes, 2018; Aminzare et al., 2018; Mendes et al., 2013; Pereira et al., 2018; Strauss and Heisenberg, 1990) and analyses (Berman et al., 2014; Mendes et al., 2013; Pereira et al., 2018) used in the Drosophila locomotion literature, despite the observations of variable coordination patterns. We have therefore retained some analysis of and discussion of gaits, while simultaneously shifting the primary focus of the manuscript to describing the variability in coordination patterns. Adding results about spontaneous turning (Essential Revision 6) has also decreased the focus on gaits.
4) The paper would be substantially improved if the coupled phase oscillator model were motivated and presented more clearly.
Thank you for this suggestion. We have updated the main text and Materials and methods to clarify and motivate the description of the coupled phase oscillator model.
We have updated the prose to clarify the rulebased nature of the model, and its relation to conceptually similar rulebased algorithms for robot limb placement, such as (Kwak and McGhee, 1989; Song and Waldron, 1987; Wettergreen and Thorpe, 1992). To that end, we now present an algorithmic description in pseudocode. Additionally, we regret the lack of clarity surrounding the various model parameters. Parameters other than the stance duration are not critical to the model; we have updated the prose to clarify this. We therefore did not include analyses of the sensitivity to parameters other than the stance duration.
We have also updated the prose to better highlight the intuition generated by the model regarding the fly’s walking control strategy. Briefly, in this framework a single descending command conveying the desired stance duration would suffice to control the forward speed of the fly. This single control parameter permits a simpler control architecture than models in which interlimb coupling is actively varied with walking speed. On each side of the fly, a posteriortoanterior set of limb swings can be thought of as a single event. The two sides of the body are then coupled to each other in a manner that maintains a halfcycle difference in phase between contralateral limbs.
5) A common feature in the gait transition literature is hysteresis – the transition occurs at a different speed going from slowto fast than fasttoslow. Either adapting the analysis to include this potential effect (or explaining why it is not relevant here) is a necessary addition to the manuscript.
We agree that, based upon the literature, it is important to address the possibility of hysteresis in our experimental data. As we do not observe any evidence of multiple fixed points in our limb coordination data, it is difficult to directly address the question of hysteresis. To address the possibility that multiple preferred limb coordination patterns exist at certain speeds, we have added a supplementary figure (Figure 3—figure supplement 2) showing the distributions of pairwise relative phases conditioned with high resolution on forward velocity. At all forward walking speeds, these distributions are unimodal, which leads us to believe that our dataset does not include substantial hysteretic effects.
6) The spontaneous turning data from your other eLife submission and its relationship to the manifold described would enhance this paper and further take the emphasis off the current discrete gaits controversy. This is not a requirement but a suggestion.
Thank you for this suggestion. We agree that incorporating the spontaneous turning data originally presented in our rejected submission would enhance this paper, and have done so.
Reviewer #1:
Note: I am also a reviewer of your other eLife submission, so several of my comments will mirror those for that paper, due to the similarity of the methods used.
The manuscript by DeAngelis and colleagues investigates the repertoire of forward locomotion dynamics in the fruit fly (D. melanogaster) using a combination of limb tip and centroid tracking, dimensionality reduction, mathematical modeling, and neuronal perturbations. The key claims are that the speed of forward locomotion are controlled using a single control parameter, and that the traditional division of walking patterns into alternating tripod, tetrapod, and wave gaits is not apparent in spontaneous walking (and is, accordingly, better understood through the manifold framework proposed here). Lastly, the authors perform both optogenetic and sensory perturbations during spontaneous locomotion, finding that the flies altered their movements along the trajectory of the found manifold.
Overall, I thought that this manuscript was creatively performed and would be of general interest to the community. I do have some serious reservations about the findings, however, and would like to see some further analysis of the data before recommending publication. Specifically, I have questions about the authors' choice of embedding to find the manifold, the effects of footfall error on their results and interpretation of their results in light of gait hysteresis. If these concerns are resolved satisfactorily, however, I envision that this paper will be wellsuited to publication in eLife.
1) After reading and digesting the paper (as well as the cosubmitted paper), I remain a bit perplexed why the authors decided to use UMAP for the dimensionality reduction technique upon which many of their results lie. UMAP is wellsuited for preserving topological features of a data set requiring warping of global metrics in order to preserve local topological structure. Despite this, the authors only use UMAP as a global coordinate, treating it like one often observed PCA (which does the opposite, preserving global metrics while warping local ones), pointing at directions where yaw or lateral or forward velocities alter. In fact, the results from Figure 2—figure supplement 1 in ZavatoneVeth et al. show that PCA actually does a better job of describing the full dynamic range of the variables than the UMAP does. As a result, it also would be interesting to see what PCA provides for the kinematic variables.
I don't make these points as merely a technicality, however. As many of the results in the paper come from this embedding, it is important to actually note what the structure of this manifold is. Moreover, the density within the manifold is not shown in any of the plots as well. Are there islands of high density with valleys in between? Or is It smoothly distributed? I think that more care is needed to more fullycharacterize the manifold that is the one of the key mathematical descriptions in this manuscript (and is in fact in the title!).
Our primary response to this comment is under Essential Revision 1. To address the concerns regarding the specific structure of the manifold, we have added a supplementary figure (Figure 4—figure supplement 2) showing the parameterization of the manifold with cylindrical coordinates. Using this projection, we visualized the density of points within the manifold. The distribution along the phase dimension is roughly uniform. As one would expect given the nonuniform distribution of forward speeds, the density along the frequency dimension of the manifold is not uniform, though it changes smoothly. The projection into cylindrical coordinates also allows us to generate more easily interpretable visualizations relating the UMAP dimensions to limb phases and frequencies, as well as a clear visualization of how the radius of the manifold is related to the axial dimension.
2) Perhaps my largest concern with the paper, however, is the effect of error on gait definitions (e.g. those in Figure 3). The authors make stance/swing definitions based on the smoothed instantaneous velocity of the limb (with the tracking error threshold chosen as the separation threshold). While this is a sensible approach given the type data that they obtain here, as far as I saw in this manuscript at any rate, I did not see any discussion about how error in this value might affect the found results. This is rather important, as many of the disagreements with the previous literature about the discreteness of individual gaits rely on these measurements. I think that it would be important to see how the results change if the threshold for swing/stance was changed within some reasonable limits, or if swing/stance was assigned using a confidence (using, say, a logistic function of the limb speed). Specifically, I would be curious to see if the 2wave vs. 3wave and 6wave results could be altered by changing this threshold.
That being said, the values for seem that far off. The results from Mendes et al., 2013, for example, are built upon a much more accurate measurement of swing/stance determination, since the footfalls are directly measured using frustrated total internal reflection (although the authors should please correct me if they feel this assessment is incorrect). The results in Figure 1 here are similar to those found in Figure 2 of that paper, so no major tweaking is likely needed, but some assessment of the effect of potential error is important. Relatedly, it might be useful to see what happens to the "tripod index" as defined in the Medes paper with the data collected here. Are the results replicated?
Our primary response is in the comments to Essential Revision 2. Specifically, the results regarding the 2, 3, and 6 wave results were preserved under halving of the separation threshold, suggesting that choice of threshold does not strongly influence this finding. Additionally, we have reproduced in Author response image 3 the gait index analysis from (Mendes et al., 2013) on the individual example trajectories presented in Figure 3FG; that analysis is shown in Author response image 3. Consistent with Mendes et al., 2013, we are capable of extracting individual trajectories with tripodlike and tetrapodlike characteristics under this measure. This, in combination with other presented analyses (Figures 12, specifically Figure 2CD), suggests that our experimental data is largely consistent with data collected in previous studies of Drosophila walking. Importantly, the tetrapodlike patterns identified by this index are consistent with our smoothly varying generative model.
3) Related to both the previous comments, I think that in some senses, this paper is offering a potentiallyfalse dichotomy between discrete and continuous gaits. For example, although the fly may have some preferred forward velocities or gaits that are peaks in the distribution, there may be nonzero valleys between them that, although rare, are important parts of the overall manifold that describes locomotion. For example, in Figure 4—figure supplement 1, the separate regions are potentially the result of the fact that there is no intermediate data. Because UMAP explicitly exaggerates topological differences, my guess is that the inclusion of a relatively small number of intermediate gaits will link the discrete rings together into something that looks more like the real data. I also worry that the phase error potentiallyintroduced from the uncertainty of the stance/swing transition time might affect the phase difference peaks in Figure 3 and Figure 3—figure supplements 1 and 2.
Our primary response is in the comments to Essential Revision 3. Regarding the specific analysis presented in Figure 4—figure supplement 1, we would like to note that we do not view the UMAP embedding as a direct line of evidence for the absence of distinct gaits. Rather, it is an ancillary line of evidence that is consistent with the existence of the continuum suggested by the analyses presented in Figures 2 and 3. To generate an embedding comparison that provides direct evidence for the absence of distinct gaits, it would be necessary to have a model for the structure of gait transitions. However, there does not exist a priori a clear null model for intermediate coordination patterns, nor does our experimental data suggest such a model. In the absence of such a null model, Figure 4—figure supplement 1 presents the case in which transitions do not occur within the embedded segments of coordination patterns. We regret that we did not make this feature of the analysis clear in the original manuscript; in particular, the order in which the synthetic and empirical data embeddings were presented framed this analysis as a stronger hypothesis test than we believe it is. We have therefore revised the ordering of this section of the manuscript, and have added prose to clarify the intent of this analysis.
Additionally, we would like to clarify that we do not use the binary swing/stance variables to estimate limb phases; rather, we estimate limb phases directly from the raw limb positional time series using the discretetime analytic signal method. As described in our response to Essential Revision 3, our investigation of the noise sensitivity of this method using synthetic data shows that measurement error leads to broader peaks in phase distributions but does not change the presence or locations of those peaks.
4) A common feature in the gait transition literature is hysteresis – the transition occurs at a different speed going from slowto fast than fasttoslow. All of the analyses here fundamentally assume that gait is a singlevalued function of speed, however, which may cause for some misinterpretations of the data, since all spread of the distribution would just be assigned to a single mean value. Either adapting the analysis to include this potential effect (or explaining why it is not relevant here) is a necessary addition to the manuscript.
Our primary response is in the comments to Essential Revision 5. Based on our additional analysis, we have not found evidence for hysteresis in our dataset.
Reviewer #2:
In this paper, DeAngelis and colleagues investigate the temporal organization of the fly walking gait. They use a simple tracking algorithm to measure the body orientation and distal leg position of flies walking in a flat, featureless environment. They confirm previous observations that stance duration changes as a function of walking speed. They then extracted swing and stance phases from each leg and used these patterns to classify different walking gaits. They find that flies typically have three feet in stance (tripod gait), but they observe that flies occasionally (and transiently) have 4 or 5 feet in stance, particularly at low walking speeds. Contralateral limbs reliably swing in antiphase, while ipsilateral limb coordination vary a bit more. A lowdimensional embedding visualization of limb coordinates resulted in cyclic (global oscillator phase) and one linear coordinate (stepping frequency). Overall, this suite of analyses suggest that there is no clear differentiation of walking behavior into distinct gaits (e.g., tripod vs. tetrapod). The authors then construct a model that can predict the continuum of walking behavior by tuning of a single parameter. Finally, they find that perturbing walking behavior through activation of the Moonwalker Descending Neuron or visual stimulation alters stance duration, and thus induces a shift along the length of the walking manifold.
Overall, the paper is clearly organized and wellwritten. Although the subject of the fly walking gait is already welltrodden, this is certainly the most exhaustive and definitive treatment to date. I think it is a very strong paper and have only a handful of suggestions to improve it prior to publication.
1) My biggest point of confusion was the coupled phase oscillator model. The paper would be improved if this model were motivated and presented more clearly. The equations are not intuitive at first glance (e.g., the 1/sin terms, the discontinuous indicator functions); is there an easy interpretation of them? And does the specific form of the model matter? How is it related to the other coupled oscillator models in the literature? I would have liked to see a sensitivity analysis that provides some intuition for how the model works and which parameters matter. For example, τ_{stance} is the relevant parameter that the authors tune to get different interlimb coordination patterns, but this is not clearly explained (especially in the main text). Finally, how does this model contribute to our understanding or intuition of the control strategy used by the fly to coordinate walking? The authors should include a deeper discussion of the motivation for and interpretation of this model.
Thank you for this comment. We have addressed this comment in detail in the reply to Essential Revision 4 above.
Reviewer #3:
This paper is notable for collecting a very comprehensive dataset of kinematic gait patterns in Drosophila over a wide range of spontaneous and evoked speed changes. While the optogenetic and visually induced changes are nice as supporting evidence, the core of the paper comes from a detailed and thorough gait analysis. As such the paper has some very good insights when connected to modern dimensionality reduction and some simple models of oscillations, but needs to distinguish itself carefully from the long history of papers that also do gait analysis especially in Drosophila.
1) The structure of the argument is that animals in general are thought to have distinct gaits. The authors cite a number of studies on larger legged animals for this, but in those cases walking is an out of phase oscillation of gravitational potential energy and kinetic energy (the inverted pendulum model) and running is inphase oscillations (the springloaded inverted pendulum). The transition in gaits in many animals from insects to horses corresponds to the transition between these underlying dynamics If this transition is unlikely (although possible?) for Drosophila given their size, it is not clear why we would expect distinct gaits even if larger animals have them.
More critical to the conclusions drawn about Drosophila, previous studies (Wosnitza, et al., 2012, and Mendes, et al., 2013) My reading of these previous studies is that they show that limb coordination patterns in Drosophila are highly variable and that while there were changes with speed that the overall shifts are more continuous than observed in many other animals. These are cited in the paper but the relationship of their conclusions to this paper are only discussed superficially. For example Wosnitza, et al., 2012 conclude that their findings "imply that Drosophila's walking behavior is more flexible than previously thought (Strauss and Heisenberg, 1990): there are no clearly separable gaits and, more specifically, the neural controller producing interleg coordination is not restricted to a fixed tripod pattern." This is very close to the initial conclusions of this paper.
So while some data from those studies do fit tetrapod and tripod gaits (as do some in this study), the conclusion that Drosophila use a variety kinematic patterns that encompass the canonical gaits is already somewhat established. Therefore the hypothesis of distinct gaits in this paper seems like a bit of a strawman. The authors do make some nice analyses of these data and show the point perhaps more convincingly than in previous work. However, prior to the manifold section, the initial presentation of the data did not leave me thinking that I had learned something new. I recognize that the authors may not agree with interpretations in previous works, but the onus is on the authors to clearly articulate what cannot be concluded from the previous work that their work in turn justifies.
More importantly this whole argument of distinct gaits takes away from what for me was the real impact of the paper. The central advance of the paper is the finding that kinematics patterns across a large range of speeds can collapse down onto a manifold with a single parameter governing speed variation. The authors are still making an advance here, but the motivation of the introduction and the conclusion of the first sentence of the Discussion highlight the continuum vs. discrete nature of Drosophila gaits. The emphasis should be about the manifold.
Our primary response to this comment is with Essential Revision 3. We have followed the suggestions of the reviewer and made the manifold structure the primary emphasis of the paper.
2) My concern about the manifold discussion is how this approach is different from the many other dimensionality reductions that have been done on gait. What really sets it apart from the tSNE analyses of behavior that show all walking lying in a low dimensional cluster (Berman, et al., 2014) and the reduced phase and oscillator models of stick insects and cockroaches (e.g. CouzinFuchs, et al., 2015). I think there is novelty here but it must be articulated more in context.
Our main response to this is included with Essential Revision 1. We have updated the prose to expand our discussion of the relative benefits of UMAP over previous dimensionality reduction techniques. In particular, because tSNE is designed to emphasize local rather than global structure, it does not immediately capture the periodic dynamics of legged locomotion, which UMAP highlights (Becht et al., 2019; McInnes et al., 2018). Furthermore, UMAP has been shown to produce embeddings that are more reproduceable across data samples and subsample sizes than those obtained using tSNE (Becht et al., 2019).
To compare the results of applying tSNE with those obtained using UMAP, we used the fast interpolationbased algorithm for tSNE introduced in (Linderman et al., 2017) to generate embeddings of our limb kinematic data. Following standard practices, we used the standard target dimensionality for tSNE (two) and a perplexity of 100, and reduced the dimensionality of the input data to 50 using PCA prior to applying tSNE (Becht et al., 2019; Berman et al., 2014; Linderman et al., 2017; Van Der Maaten, 2014). As shown in Author response image 4, tSNE does capture some of the frequency structure of the manifold produced by UMAP, but, unlike UMAP, it produces an embedding that appears segmented.
We wanted to better understand whether this segmentation arose from features of the data or features of the algorithm. We therefore applied the identical analysis to data generated by our model. By construction, the highdimensional data manifold produced by the model is continuous. However, the tSNE embedding of the model data is also segmented, suggesting that the structure observed in the embedding of experimental data could result from the analysis itself rather than the underlying structure of the data (Author response image 4). This result makes sense given that tSNE is explicitly designed to emphasize local similarities in the highdimensional data at the expense of largescale structure (Van Der Maaten, 2014). This is in some ways opposite to the purpose of UMAP, which is designed to preserve both local and large scale structure in the data. Therefore, unlike UMAP, the application of tSNE to our limb kinematic data does not produce results that are immediately interpretable. Since the goal of our dimensionality reduction analysis is to visualize the global structure of our data, the interpretability of the embedding is key.
Article and author information
Author details
Funding
National Institutes of Health (R01 EY026555)
 Brian DeAngelis
 Damon A Clark
National Science Foundation (GRFP)
 Brian DeAngelis
Richard and Susan Smith Family Foundation
 Brian DeAngelis
 Damon A Clark
Alfred P. Sloan Foundation
 Damon A Clark
Kinship Foundation (Searle Scholar Award)
 Brian DeAngelis
 Damon A Clark
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
Acknowledgements
This work benefitted from many helpful conversations with M Venkadesan, as well as analysis by U Piterbarg. BDD was supported by an NSF GRF and a Gruber Science Fellowship. DAC and this project were supported by the Smith Family Foundation, NIH R01EY026555, a Searle Scholar Award, and a Sloan Fellowship in Neuroscience.
Senior and Reviewing Editor
 Ronald L Calabrese, Emory University, United States
Publication history
 Received: March 2, 2019
 Accepted: June 15, 2019
 Version of Record published: June 28, 2019 (version 1)
 Version of Record updated: December 4, 2020 (version 2)
Copyright
© 2019, DeAngelis et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics

 5,394
 Page views

 596
 Downloads

 48
 Citations
Article citation count generated by polling the highest count across the following sources: Crossref, Scopus, PubMed Central.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading

 Neuroscience
Nicotine intake is likely to result from a balance between the rewarding and aversive properties of the drug, yet the individual differences in neural activity that control aversion to nicotine and their adaptation during the addiction process remain largely unknown. Using a twobottle choice experiment, we observed considerable heterogeneity in nicotinedrinking profiles in isogenic adult male mice, with about half of the mice persisting in nicotine consumption even at high concentrations, whereas the other half stopped consuming. We found that nicotine intake was negatively correlated with nicotineevoked currents in the interpeduncular nucleus (IPN), and that prolonged exposure to nicotine, by weakening this response, decreased aversion to the drug, and hence boosted consumption. Lastly, using knockout mice and local gene reexpression, we identified β4containing nicotinic acetylcholine receptors of IPN neurons as molecular and cellular correlates of nicotine aversion. Collectively, our results identify the IPN as a substrate for individual variabilities and adaptations in nicotine consumption.