Simultaneous embedding of multiple attractor manifolds in a recurrent neural network using constrained gradient optimization
Authors: Haggai Agmon
Yoram Burak
Venue: NeurIPS 2023 poster (2023)
Keywords: Theoretical Neuroscience
Computational Neuroscience
Recurrent Neural Networks
Attractor models

Abstract: The storage of continuous variables in working memory is hypothesized to be sustained in the brain by the dynamics of recurrent neural networks (RNNs) whose steady states form continuous manifolds. In some cases, it is thought that the synaptic connectivity supports multiple attractor manifolds, each mapped to a different context or task. For example, in hippocampal area CA3, positions in distinct environments are represented by distinct sets of population activity patterns, each forming a continuum. It has been argued that the embedding of multiple continuous attractors in a single RNN inevitably causes detrimental interference: quenched noise in the synaptic connectivity disrupts the continuity of each attractor, replacing it by a discrete set of steady states that can be conceptualized as lying on local minima of an abstract energy landscape. Consequently, population activity patterns exhibit systematic drifts towards one of these discrete minima, thereby degrading the stored memory over time. Here we show that it is possible to dramatically attenuate these detrimental interference effects by adjusting the synaptic weights. Synaptic weight adjustment are derived from a loss function that quantifies the roughness of the energy landscape along each of the embedded attractor manifolds. By minimizing this loss function, the stability of states can be dramatically improved, without compromising the capacity.

Simultaneous embedding of multiple attractor manifolds in a recurrent neural network using constrained gradient optimization
Haggai Agmon
Yoram Burak
The Hebrew University Jerusalem, Israel and Stanford University, USA haggai. agmon@mail hwmpaaill..ji.... ac.
The Hebrew University of Jerusalem, Israel yoram. burak@elsc cusieeisss........... quji.a. ac il
Abstract
The storage of continuous variables in working memory is hypothesized to be sustained in the brain by the dynamics of recurrent neural networks (RNNs) whose steady states form continuous manifolds. In some cases, it is thought that the synaptic connectivity supports multiple attractor manifolds, each mapped to a different context task. For example, in hippocampal area CA3, positions in distinct environments are represented by distinct sets of population activity patterns, each forming continuum. It has been argued that the embedding of multiple continuous attractors in a single RNN inevitably causes detrimental interference: quenched noise in the synaptic connectivity disrupts the continuity of each attrac- tor, replacing it by discrete set of steady states that can be conceptualized as lying on local minima of abstract energy landscape. Consequently, population activity patterns exhibit systematic drifts towards one of these discrete minima, thereby degrading the stored memory time. Here show that is possible to dramatically attenuate these detrimental interference effects by adjusting the synaptic weights. Synaptic weight adjustments are derived from loss function that quantifies the roughness of the energy landscape along each of the embedded attractor manifolds. By minimizing this loss function, the stability of states can be dramatically improved, without compromising the capacity.
Introduction
In the brain, recurrent neural networks (RNNs) involved in working memory tasks are thought to be organized such that their dynamics exhibit multiple attractor states, enabling the maintenance of persistent neural activity even the absence of external stimuli. continuous external variable, these attractors are thought to form a continuous manifold of neural activity patterns. some well studied brain circuits in flies and mammals, neural activity patterns have been shown to robustly and persistently reside along a single, low-dimensional manifold, even when the neural activity dissociated from external inputs to the network 2 21
tasks that require tracking of
[5]. In
other brain regions, however, the same neural circuitry thought to support multiple cow-dimensional manifolds such that activity lies on one of these manifolds at any given time, depending on the context task. The prefrontal cortex, for example, is crucial for multiple working memory and evidence accumulation tasks [14 34 46, 44], and is thus thought that the synaptic connectivity in this brain region supports multiple attractor manifolds. Naively, however, embedding multiple attractor manifolds a single RNN inevitably causes interference that degrades the network performance.
To be specific,
focus here attractor models of spatial representation by hippocampal place cells (28). The synaptic connectivity area CA3 is thought to support auto-associative dynamies, and place cells in this area are often modeled as participating in continuous attractor network
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

Typically, neurons are first arranged on an abstract neural sheet which mapped their preferred firing location in the environment. The connectivity between any two neurons is then assumed to depend on their distance on the neural sheet, with excitatory synapses between nearby neurons and effectively inhibitory synaptic connections between far away neurons. This architecture constrains the population activity dynamics to express localized self-sustaining activity pattern, or 'bump', which can be centered anywhere along the neural sheet. Thus, the place cell network possesses two-dimensional continuum of possible steady states. As the animal traverses its environment, this continuum of steady states can be mapped in a one-to-one manner to the animal's 's position.
(CAN) 7 41 40 31
49
36
II 13
17
The same place cells, however, participate in representations of multiple distinct environments, collectively exhibiting global remapping [27]: the relation between activity patterns of a pair of place cells in one environment cannot predict the relation of their activity patterns in a different environment and thus each environment has its unique neural population code. To account for this phenomenon in the attractor framework, multiple, distinet continuous attractors which rely on the same place cells are embedded in the connectivity, each representing distinct single environment [36 5 25 II.
Conceptually, this is similar to the Hopfield model [18] but instead embedding discrete memory patterns, each embedded memory pattern is continuous manifold that corresponds to different spatial map. Despite the network's 's ability to represent multiple environments, it has been argued that the embedding of multiple continuous attractors inevitably produces frozeen (or quenched) noise that
eliminates the continuity of steady states in each of the discrete attractors [36 32
19 25 26), leading to detrimental drifts of the neural activity. Unlike stochastic dynamical noise which is expressed in instantaneous neural firing rates, time-independent quenched noise is hard-wired in the connectivity it breaks the symmetry between the continuum of neural representations.
From a qualitative perspective, the bump of activity can be conceptual- ized as residing on a minimum of an abstract energy landscape in the N dimensional space of neural activ- ity (where N is the number of neu- rons). This energy landscape is pre- cisely flat in the single map case, in- dependently of the bump's ' position (Fig. [a). However, contributions the synaptic connectivity included to support additional maps distort this flat energy landscape into a sim- ilar, yet wrinkled energy landscape, and the continuous attractors are re- placed by localized discrete attrac- for states (Fig. [p). Consequently,
b
Energy
Energy
Figure 1: Energy landscape. Schematic illustration of energy surfaces along a 1-dimensional attractor when a single (a). and multiple (b). maps are embedded. Cyan traces show energy landscapes along bottom of surfaces.
the system can no longer stably represent continuous manifold of bump states for any of its embed- ded attractors. Instead, activity patterns will systematically drift into one of the energy landscape's s discrete minima, thereby degrading the network's 's ability to sustain persistent memory. This highly detrimental effect has been explored in previous works and was recognized as an inherent property of
such networks (25 26 [. External sensory inputs can pin and stabilize the representation [26] D. but it remained unknown whether internal mechanisms could potentially attenuate the quenched noise in the system and stabilize the representations, independently of external sensory inputs. These insights raise the question. whether it is possible to embed multiple continuous attractors in single RNN while eliminating interference effects and retaining a precisely flat energy landscape for each attractor, but without simply reducing the network capacity.
Here, we first formalize the concept of an energy landscape in the multiple attractor case. We then minimize an appropriately defined loss function to flatten the energy landscape for all embedded attractors. We show that the energy landscape can be made nearly flat, which reflected by an attenuation of the drifts, and dramatic improvement in the stability of attractor states across the multiple maps. These results provide proof of principle that internal brain mechanisms can support maintenance of persistent representations in neural populations which participate in multiple continuous attractors simultancously, with much higher stability than previously and naively expected.

---

Results
We consider the coding of spatial position in multiple environments, or maps, by hippocampal place cells. Within each map. the connectivity between place cells produces a continuous attractor. We adopt a simple CAN synaptic architecture with spatially periodic boundary conditions. In the one-dimensional analogue of this architecture, this connectivity maps into the ring attractor model
[7 41 31 49| neurons, functionally arranged on a ring, excite nearby neighbors, while global inhibitory connections elicit mutual suppression of activity between distant neurons. This synaptic architecture leads to a bump of activity, which can be positioned anywhere along the ring.
The overall synaptic connectivity J is expressed as a sum over contributions from different maps, J= J', where L. is the number of embedded maps (see Supplementary Material, SM). To mimic the features of global remapping. distinct spatial map is generated independently for each environment, by choosing a random permutation that assigns all place cells to a set of preferred firing locations that uniformly tile this environment. Consequently, the network possesses a discrete set of continuous ring attractors, each mapping a distinet environment to a low-dimensional manifold of but adapted here to the
population activity patterns. This is similar to previous models 36 5 25 formalism of a dynamical rate model [
Since the connectivity is symmetric, it is guaranteed that the network dynamics will converge to a fixed point [TO]. Qualitatively, these fixed points can be viewed as the minima of an abstract high-dimensional energy landscape in which the population activity will eventually lie at steady state: each population activity pattern is associated with an energy value, and the population dynamics are constrained to follow only trajectories that cannot increase the energy and must ultimately settle in a local minimum.
When only a single map is embedded in the connectivity, the representations span true continuum of steady states and can be accurately read out. From the energy perspective, this continuity corresponds to completely flat energy landscape where a continuum of steady states share an identical energy value. Initialization of the network from any arbitrary state in this case will always converge to steady state which is a localized idealized bump of activity. Idealized bumps have a smooth symmetric structure, and are invariant as they can evolve at any position along the attractor.
is sufficient, however, to embed only one additional map to distort the continuity of true steady states. Embedding multiple maps introduces quenched noise in the system which eliminates the ability to represent a true continuum of steady states in each of the attractors. From the energy perspective, the quenched noise is reflected in distortions of the energy landscape, it becomes wrinkled with multiple minima (Fig. Tb). Thus, the discrete attractors lose their ability to represent a true continuum of steady states, and the memory patterns consequently accrue detrimental drifts (Fig. 2).
The idealized bumps from the single map case multiple maps are embedded. Instead, activity patterns converge on distorted versions of idealized bumps. Nevertheless, as long the network is below its capacity, this activity is still localized. When initializing the network from an idealized bump at a random position along any of the discrete attractors, it will instantaneously become slightly distorted and will then systematically drift the energy minimum point in its basin of attraction (Fig. (F).. These effects are enhanced as the number of embedded maps (and thus the magnitude of quenched noise) is increased and until capacity is breached when the exhibited activity patterns are not expected to be localized in any of the attractors.
longer steady states of the dynamics when
Flattening the energy landscape
To attenuate the systematic drifts that emerge upon embedding of multiple maps, we focused on flattening of the wrinkled energy landscape. The methodology throughout this work relied on the following two-step procedure: first, the wrinkled energy landscape was evaluated. Then, based on this evaluation, small modifications M were added to the original synaptic weights to flatten the energy landscape (SM). Importantly, the goal was to flatten the energy landscape for all maps simultaneously, which rules out the trivial solution of flattening subset of maps by unembedding the others.
Generally, the energy value of each arbitrary state for any RNN with symmetric connectivity given by the following Lyapunov function 110:
3

---

a
map
map
map #3
mss
Firing rate [H22]
10L
10
1500 ms
10
ma
Bump position
500
8.
400
300
200
100
200
500
400
500
500
100
200
300
400
600
500
0.5
1.5
2.5
Bump position
Bump position
Time (sec]
maps
16
10 maps
map #5
Firing rate [Hz]
12
1500
12
Bump position
500
400
300
200
100
100
200
-900
400
500
500
100
200
300
400
500
500
0.5
1.5
8.5
Bump position
Bump position
Time [sec]
Figure 2: Embedding multiple maps distorts the idealized bump and induces systematic drifts. a, Two examples (left and right) showing snapshots of population activity at three time points (0 ms, 10 ms and 1500 ms). Activity was initialized example. Population activity remains stationary when a single map is embedded (top), but distorts and drifts when maps are embedded (bottom). Superimposed bump positions versus time of ten independent idealized bump initialization (red), in two different embedded maps (top map #3, bottom- map #5), out ten total embedded maps. The stable bump position obtained when single map is embedded is plotted for reference (dashed traces). Systematic drift is evident when ten maps embedded.
random position along the attractor in each
N (%/² dzzzz (zs) hh) B(I, W Σ
h6 (T:)
Σ j=1
d (Is) W Wajd (1i)
(I)
where E is the energy value, I is the population synaptic activity. W is the connectivity, h is an external constant current, and is the neural transfer function. To flatten the wrinkled energy landscape, it was first necessary to evaluate it along each of the approximated attractor manifolds. Since h, o. and number of neurons (N) are assumed to be fixed, the energy value (Eq, [ of each state depends only on the population activity (ī) and the synaptic connectivity (W) which, in case, is decomposed into the embedded maps and the weight modifications, namely, W = E + M.
We first took a perturbative approach in which when only a single map is embedded in the connectivity, while treating all the contributions to the connectivity that arise from the other maps inducing small perturbation to the energy. In Eq. ( corrections to the energy of stationary bump states in map 1 arise from two sources: First, a contribution arising directly from the synaptic weights associated with the other maps (J" where I' * I), as well as from the modification weights M (third term in Eq. I.. This term, leading order, is linear in the synaptic weights. Second, corrections arising from deformation of the bump state. Because deformed bump drifts slowly when the perturbation
examined the idcalized bumps that emerge
the synaptic weights is weak, it is possible to conceptually define bump states along the continuum of positions, which are nearly stationary (a precise way to do will be introduced later on). This contribution to the energy modification is quadratic in the deformations, because the idealized bump states are minima of the unperturbed energy functional, and therefore the energy functional is locally quadratic near these minima. Hence, to leading order in the perturbation, we neglect this modification in this section. Overall, the energy Ei of idcalized bump (ib) state centered around position & in map ( is,

---

a
3117"
map #3
2307
10 mape; map #3
600
map
Energy (a0)
Bump Doallion
Ineuinn Untts)
400
Energy
200
100
200
900
400
600
2000
500
400
500
600
Bump position
Bump position
Time [sed]
110"
A10"
maps: maps
R5
200
Energy [aui
Energy [au]]
Buno posilion
(mouurnn urtti units)
000
100
100
200
300
100
6000
600
1200
1800
2400
3000
5600
4500
6800
5400 0000
0.5
7.5
Bump position
Concatenated bump position
Time [sec]
Figure 3: Energy landscape is distorted when multiple maps are embedded but can be flattened to reduce drifts. a, Evaluated energy landscape using idealized bumps, along two representative embedded maps (top map #3, bottom- map #5), out of the total I. = 10 embedded maps (orange) used Fig. 2 The energy landscape obtained when L = 1 is superimposed in blue. For visibility purposes, energy landscapes are shown after subtraction of the mean across all positions and maps (SM). b, Re-eevaluated energy landscape using idcalized bumps after weight modifications were added (orange). The landscape is precisely flat. Note that this is only approximation to the actual energy landscape, due to the use of idealized bumps (see text). Top: energy along map #3 (blue trace is identical orange trace in panel top). Bottom: same top, for all embedded maps concatenated. e, Bump trajectories as in Fig. 2p, but with added weight modifications. Qualitatively, the magnitude of the drifts is reduced (compare with Fig. ))
E/by Eu
Σ 1'A
pk k.1" J"
1.1 2"0
Mr
(2)
where Eg the energy value of an idealized bump state in the case of a single embedded map (which idealized bump around position in map I.
independent of le and I), and ro.. ④ (irk.))
As first step, we evaluated the first two terms in Eq. 2 which represent the Lyapunov energy of idealized bump states in the absence of the weight modifications M. In each one of the maps, the energy landscape was evaluated using these idcalized bumps at N uniformly distributed locations, centered around the N preferred firing locations of the neurons. Since this resolution is much smaller than the width of the bump, the continuous energy landscape was densely sampled. As expected, precisely flat energy landscape blue traces), as the attractor is truly continuous in this case. However, when multiple maps were embedded. a wrinkled energy landscape observed (Fig. 33, orange traces), as consequence of the quenched noise in the system.
observed when only single map was embedded (Fig. 3p,
We next sought small modifications in the synaptic weights which could potentially flatten the energy landscape with only mild effects on the population activity patterns. Under the approximations considered above, flattening the energy landscape requires that the right hand side of Eq. 2 is constant, independent
and Because the energy depends on the connectivity a linear fashion to leading order, obtain set of linear equations. Using the discretization described above and when N > 2L + 1, such a system is underdetermined with more unknown weight modifications than sampled energy evaluations along the attractors. Out of the infinite space of solutions, the least squares solution, which has the overall minimal L2 norm of connectivity modifications was chosen.
As expected, re-evaluating the energy using idealized bumps with the addition of the weight mod- ifications yielded a precisely flat energy landscape (Fig. 3b, orange traces). This does not imply, however, that drifts will vanish since evaluated the true energy landscape only up first order
5

---

a
Baseline
Modified
2.5
Bump score [au]
85
LE (au)
-2
[ees] Tinl
80-
25
20
0.5
10
C
Drift
[nouron
Baseline Modited
Baseel
Mociffdd
10
20
Figure 4: Effects of weight modifications obtained using idealized bumps scheme on the stability of bump states. a, Energy difference between consecutive network states without (left) and with (right) the weight modifications a function of L, the total number of embedded maps. In both cases, the changes in energy of network states are negligible after 3 seconds, indicating that steady states were reached (see also Supplementary Fig. 1). b, The bump score (SM) without (blue) and with (orange) the weight modifications as a function of L. c, Mean squared change of all neuron's 's firing rates [Hz] between time points 0 and 3 sec, without (blue) and with (orange) the weight modifications as a function of L. Note that the scale is logarithmic. d, Mcasured drifts without (blue) and with (orange) the weight modifications as function of L. Error bars are #1.96 SEM in panels (SM).
in the deviations from idealized bump state. To test whether the weight modifications led to decreased drifts, we independently initialized the network using an idealized bump, centered each realization at one of ten uniformly distributed positions in each attractor. First, changes in the energy of the population activity through time were monitored, to verify that steady states were reached. As expected for Lyapunov energy function, this quantity decreased monotonically, and stabilized almost completely for all initial conditions within three seconds, indicating that the convergence to steady state was nearly complete (Fig. 4.). To further validate that activity was nearly settled within three seconds, we measured the drift rates and the rate of mean squared change in the firing rate of all the neurons as defined below. Only subtle changes were still observed beyond this time (Supplementary Fig. 1), justifying the consideration of the states achieved after three seconds as approximated steady states. We next evaluated bump score that quantifies the similarity between the population activity with any one of the idealized bumps states, placed at all possible positions along all the attractors (SM). In the large N limit, sharp transition expected when capacity breached, but this effect is smoothed for smaller and finite number of neurons, as used in this study (SM). Nevertheless, adding the weight modifications did not decrease the bump score (Fig. 4 indicating that the network can still reliably represent its multiple embedded maps. Our focus on networks below the capacity, which we roughly estimate to be ~0 maps for N 600.
To quantify the stability of the network, we measured the mean squared change in the firing rate of all neurons, from initialization and until approximated steady states were reached (at 3s). This measure of stability improved dramatically due to the synaptic weight modifications (Fig. 4e). A second measure of stability was the drift of the bumps, obtained by examining the mean distance bumps traveled. The bump location was defined as the position that maximized the bump score. We found that the distance bumps traveled over 38 from initialization decreased significantly when the weight modifications were added to the connectivity (Fig- ft). Since the spatial resolution of the energy landscape evaluation identical to that in which neurons tile the environment, drifts lower than single unit practically correspond to perfectly stable network. These results (Fig. 4b-d) demonstrate

---

that flattening the approximated energy landscape across all maps by adding the synaptic weight modifications indeed increased the stability of the multiple embedded attractors.
Despite the approximation of the energy function to first order in the weight modifications, our approach achieved dramatic improvement in the network stability. Two limitations of this approach can be clearly identified: first, the energy landscape was evaluated for idealized bump states, even though these activity patterns are not the true steady states when multiple maps are embedded in the connectivity, Second, weight modifications were designed to correct the energy only up to the first order. Once introduced, these weight modifications produce slight changes in the structure of the steady states, which in turn generate higher-order corrections to the energy that were neglected in the approximation. Next, we describe how these limitations can be overcome by defining a more precise loss function for the roughness of the energy landscape.
Iterative constrained gradient descent optimization
As discussed above, when L 1 the unmodified system does not possess true steady states at all positions: following initialization, an idcalized bump state is immediately distorted, and then systematically drifts to local minimum of the energy functional (Fig. 2 order to define an energy landscape over continuum of positions, it is thus necessary to first formalize the concept of states that have reached an energy minimum but also span a continuum of positions. To do so, we considered the minima of the energy functional under a constraint on the center of mass of the bump. In practice, such minima are found by idealized bump initialization, followed by gradient descent on the energy functional under the constraint. A distorted bump then dynamically evolves to minimize the Lyapunov function, while maintaining a fixed center of mass at its initial position (Supplementary Fig, 2). See SM for details on the constrained optimization procedure and its validation.
Since the goal of the constraint is to parametrize all states along the approximated continuous attractor. its exact definition is not crucial as long as the optimization is applied along a dense representation of positions. By finding the minima of the energy functional along such dense sample of positions, it is possible to precisely measure the energy landscape of the system along each one of the attractors. As expected, the evaluated energy landscape obtained using gradient optimization achieved lower energy values compared to the energy landscape obtained using idealized bumps (Fig. 5a. orange traces).
Using this formal definition of the position dependent energy, it is possible to formulate a learning objective for the modification weights M, designed to flatten the energy landscape. We denote by E the constrained minimum of the energy at position k in map 4, evaluated using the constrained gradient minimization scheme described above. We approached our goal by attempting to equate the values of E for all i and I, using an iterative gradient- based scheme (SM). It is straightforward to evaluate the gradient of Edd with respect to M, since up to first order in AM,
1 R AT (AM) R
E (M+ - E (M) = -
where li is the firing rate population vector at the constrained minimum corresponding to E We note that the constrained minima themselves are modified due to the change in M, but these modifications contribute only terms of order (AM)² to the change in the energy (see SM).
In each iteration of the weight modification scheme, E&d and R k,l were first evaluated numerically for all k and I using constrained gradient optimization and the values of M from the previous iteration (starting from vanishing modifications in the first iteration). Next, in iteration i, we sought adjustments AM, the weight modifications M, that equate the energy across all l: and I, to leading order. As in the previous section, Eq. 3 yields an underdetermined set of linear equations (assuming that N > 2L + 1), where the least square solution was chosen in each iteration.
Fig. 5p-c shows results obtained over several gradient optimization iterations. The mean absolute value and standard deviation of weight modifications decreased with the iteration number (Fig. 5p), and the flatness of the energy landscape improved systematically in successive iterations (Fig. 15
top, and Supplementary Fig. 3). The standard deviation of the energy across locations and maps decreased by almost four orders of magnitude after four iterations, and the maximal deviation of the energy from its mean decreased dramatically as well (Fig. 5. bottom).
We next examined the consequences of improved flatness of the energy landscape on the network's stability. Qualitatively, individual states exhibited much less drift after five iterations (Fig. (b)
7

---

b
x10"
maps; map #3
Mean absolute
value of wrights
10
Standard deviation of landscape
Energy (au)
10
Idealited bumps
Grudient oprimcatioon
AM.
AM
10
Ixroice
400
10
10
20
50
60
Bump position
X10'
10 maps: map #5
10
10"
Standard deviation
10
Maximal energy
weights
difference from the
10
10"
Energy
10
Orginsl
-10
1003lood bumps
opranication
450
500
10
20
30
AMM
AM.
40
50
AM,
AM.
AM.
60
10
10
20
30
Famulion
Tanulion-4
Twendinn
50
60
Bump position
Figure 5: Energy landscape evaluation and modifications using constrained gradient optimiza- tion. a, Energy landscape evaluated using constrained gradient optimization (orange) for maps #3 (top) and #5 (bottom), out of L 10 embedded maps as used in Figs. 2 and For reference, energy landscapes evaluated using idealized bumps are plotted as well (blue traces, identical to the orange traces from Fig. Bh. As expected, the energy after gradient optimization (orange trace) is lower. Re-evaluating the energy landscape to leading order in the weight modifications yielded precisely flat energy landscape when evaluated using the pre-modified network steady states (yellow). b, Top: mean absolute value of the weights in the original connectivity matrix (blue) and in modification adjustments for each gradient optimization iteration, as a function of L, the total number embedded maps. Bottom: corresponding standard deviation. Error bars are #1 11.6 SEM (SM). Top: standard deviation of the energy landscape without (blue) and with the weight modification for each gradient optimization iteration function of L. Bottom: same as top, but showing the peak absolute difference between all energy values and their mean. Error bars are +1.96 SEM.
compared to the unmodified scenario (Fig. 2p), and also in comparison with the scheme based on idealized bumps (Fig. Bp.. To systematically quantify the improved stability, first validated that approximated steady states were reached three seconds after initialization (Fig. (p), and that the bump score not significantly affected by these modifications (Fig. (b) in the idealized bump
approach (Fig. 4)..). Next, we examined two measures of stability (as in Fig. --- The mean squared change in the firing rate of all the neurons, across 3s from initialization was dramatically improved, by almost three orders of magnitude, compared to the pre-modified networks for the smaller values of I. (Fig. 6t... Measured drifts were attenuated significantly well (Fig- (Ge) Note that, as described above, drifts which are helow or equal to the discretization precision of a single unit correspond to perfect stability. For completeness, drifts were also measured using the phase of the population activity vector and demonstrated similar results (Supplementary Fig. 4). Taken together, these results show that with few constrained gradient optimization iterations, dramatic increase can be achieved in the inherent stability of network forming a discrete set of continuous attractors.
Discussion
Interference between distinct manifolds highly detrimental for the function of RNNs designed to represent multiple continuous attractor manifolds. Naively, it is sufficient to embed only second manifold network destroy the continuity of the energy landscapes in each attractor. This leads to the loss of the continuity of steady states, which are replaced by a limited number of discrete attractors. In tasks that require storage a continuous parameter working memory, these effects are manifested by systematic drifts that degrade the memory over time. Sensory inputs can, hypothetically, pin and stabilize the representation (26] I but in the absence of correcting external inputs the representation is destined to systematically drift to the nearest local energy minima, thereby
8

---

a
map #3
mup
Baseline
Bump score [au]
idcalized bumps
Bump position
Gradient optimizaioon
Ineuron untte
400
ME (au)
[c0d] SUIL
1.5
100
30
40
00
60
60
L
Time [eBc]
map #5
d
e
500-
10
Bump poeition
Ineuron unit]]
400-
Incuron units]
Drift
-Bassiin
-dediizedd bumps
GGadien optimizaiion
Baseline - Ideealized bumpo Gradient aptimizaiion
0.5
25
15
Time [98c]
L
Figure 6: Effects of weight modifications obtained using constrained gradient optimization on the stability bump states. a, Bump trajectories as in Figs. [ and and but with added weight modifications obtained after five constrained gradient optimization iterations. Qualitatively, the drifts almost vanish (compare with Figs. 2 3nd and note cyclic boundary conditions). b, Same as Fig. T. but after five iterations of the constrained gradient optimization scheme. c-e, Same as Fig. 4b-d, but with superimposed results obtained after five constrained gradient optimization iterations (yellow traces). Adding weight modifications did not reduce the bump score (panel c). The network stability improved significantly (panels d-e) compared to the unmodified network (blue traces) and to the modified network using the idealized bumps scheme (orange traces). Error bars are + +1.96 SEM (SM).
impairing network functionality. This has been thought to be a fundamental limitation of such networks.
Here we reexamined this question by adding small weight modifications, tailored to flatten the high dimensional energy landscape along the approximate attractor, to the naive form of the connectivity matrix. We found that appropriately chosen weight modifications can dramatically improve the stability of states along the embedded attractors. Furthermore, the objective of obtaining flat energy landscape appears to be largely decoupled from the question of capacity, which is also limited by interference between multiple attractors. Indeed, introducing the weight modifications did not qualitatively affect the dependence of the bump score on L (Fig. (b). In particular, the kink in this function occurs at approximately L 20 (for N 600) both in the pre- and post-modified networks. Theoretically, the position of the kink is expected to roughly match the location of a sharp transition in the large N limit, while keeping L/N fixed. A recent work (6) suggests that the capacity may depend logarithmically, and thus weakly, on the prescribed density of steady states. It will be interesting to explore the relation of this result to framework
was obtained for RNNs composed of binary neurons, with linear kernel support vector machine based learning rule.
We explored our schemes in 1D to reduce the computational cost, but it is straightforward to extend our approach to 2D. One notable difference between ID and 2D environments is that the number of neurons required to achieve a good approximation to a continuous attractor, even for single map, scales in proportion to the area in 2D, as opposed to length in ID. However, for a given number of neurons, there is substantial difference between the two cases in terms of the complexity of the problem: the number of equations scales as NL, and the number of parameters (synaptic weights) scales as 2 Since the random permutations completely unrelated to the spatial organization of the firing fields, quenched (frozen) noise expected to behave similarly in the two cases.
For simplicity, we assumed that each cell is active in each environment, but it is straightforward to adapt the architecture to one which the participation ratio, p, defined as the average fraction of maps in which each cell participates, is smaller than unity. Measurements in CA3, performed

---

in the same cells in multiple environments, indicate that CA3 cells are active only in a subset of environments, with p ~ 15% a rough estimate [. Clearly, the quenched noise in the system increases in proportion to the number of embedded maps T., and it will be interesting in future work to assess how the roughness of the energy landscape depends on N. L, and p. We note that even though small p implies low interference, it also implies a reduction in the number of neurons participating in each map, and in each bump state. This is expected to reduce the resilience of each attractor to the quenched noise. Hence, the overall effect of p the roughness of the energy landscape (when keeping N fixed) is non-trivial.
Previous works have suggested that synaptic scaling 32 or short term facilitation 19 38 mecha- nisms can stabilize working memory in networks with heterogeneous connectivity. These mechanisms were achieved, however, in networks which are equivalent to the single map case but with random added heterogeneity.
this respect, the source of heterogeneity in these networks is different from the one in this work, where heterogeneity arises from the embedding of multiple coexisting attractors. would be interesting to investigate whether a synaptic facilitation mechanism [9]] can be implemented in the multiple map case and further improve the network's stability alongside the energy-based approach proposed here. However, since in a discrete set of continuous attractors each neuron is participating in multiple representations, a synaptic scaling mechanism 32] is unlikely to achieve such stabilization. In addition, it may also be interesting to explore other potential methods for the design of synaptic connectivity that supports multiple, highly continuous attractors, other than the one explored here. These might include generalization of approaches that were based on the escudo-inverse learning rule for embedding of correlated memories in binary networks (29), or training of rate networks to possess a near-continuum of persistent states by requiring stability of a low-dimensional readout variable [12
Our approach is based on the minimization of a loss function that quantifies the energy landscape's s flatness. This raises several important questions from biological standpoint. First, the connectivity structure obtained after training is fine-tuned, raising the question of whether neural networks the brain can achieve a similar degree of fine tuning. A similar question applies very broadly to CAN models, yet, there is highly compelling evidence for the existence of CAN networks in the brain 2 48 46 39 15). We also note that out of an infinite space of potential weight modifications only a specific set (least-squares) was chosen, leaving many unexplored solutions to the energy flattening problem which may relax the fine-tuning requirement. Second, our gradient-based learning rule for the minimization of the loss function was not biologically plausible. In recent years, however, many important insights were obtained computation in biological neural networks by training RNN models using gradient based learning 43 33 44 24
37 16 22 42). The (often implicit) assumption that biological plasticity in the brain can reach similar connectivity structures as those obtained from gradient based learning, even if the learning rules are not yet fully understood. Thus, our results should be viewed as a proof of principle that multiple embedded manifolds, previously assumed to be inevitably unstable, can be restructured through learning to achieve highly stable, nearly continuous attractors. It will be of great interest to seek biologically plausible learning rules that could shape the neural connectivity into similar structures. Such rules might be derived from the goal of stabilizing the neural representation during memory maintenance, either based solely on the neural dynamics, or on corrective signals arising from a drift of sensory inputs relative to the internally represented memory 23
30 20
Acknowledgmentts
The study was supported by the European Research Council Synergy Grant no. 951319 ("KILONEU- RONS"), and by grant nos. 10s.1978/11 1978/13, and 1745/18 from the Israel Science Foundation. We further acknowledge support from the Gatsby Charitable Foundation. Y.B. is the incumbent of the William N. Skirball Chair in Neurophysics. This work is dedicated to the memory of Mrs. Lily Safra, great supporter of brain research.
10

---

References
[1] Haggai Agmon and Yoram Burak. A theory of joint attractor dynamics in the hippocampus and the entorhinal cortex accounts for artificial remapping and grid cell field-to-field variability. eLife, 9ce56894, 2020.
[2] Emre Aksay, Itsaso Olasagasti, Brett D Mensh, Robert Baker, Mark S Goldman, and David W Tank. Functional dissection of circuitry neural integrator. Nature neuroscience, 10(4))494 504, 2007.
[3] Charlotte B Alme, Chenglin Miao, Karel Jezek, Alessandro Treves, Edvard Moser, and May- Britt Moser. Place cells in the hippocampus: eleven maps for eleven rooms. Proceedings of the National Academy of Sciences, 111(52)18428818444 2014.
[4] Omri Barak. Recurrent neural networks as versatile tools of neuroscience research. Current opinion in neurobiology, 46:1-6, 2017.
[5] F. P. Battaglia and A. Treves. Attractor neural networks storing multiple space representations: model for hippocampal place fields. Physical Review E, 586))77338-7777 dee 1998.
[6] Aldo Battista and Rémi Monasson. Capacity-resolution trade-off in the optimal learning of multiple low-dimensional manifolds by attractor neural networks. Physical review letters, 124(4):48302, 2020.
[7] R. Ben- Yishai, R. L. Bar-Or, and H. Sompolinsky. Theory of orientation tuning in visual cortex. Proceedings of the National Academy of Sciences, 92(9):3844-3844 1995.
[8] Y. Burak and I. R. Fiete. Accurate Path Integration in Continuous Attractor Network Models of Grid Cells. PLoS Computational Biology, 5(2)xe1000291 feb 2009.
[9] Rishidev Chaudhuri, Berk Gerçek, Biraj Pandey, Adrien Peyrache, and Ila Fiete. The intrinsic attractor manifold and population dynamics of canonical cognitive circuit across waking and sleep. Nature Neuroscience, 22(9):1512-1500 2019.
[10] Michael A. Cohen and Stephen. Grossberg. Absolute Stability of Global Pattern Formation and Parallel Memory Storage by Competitive Neural Networks. IEEE Transactions Systems, Man and Cybernetics, 13:815-826 826, 1983.
[11] Albert Compte, Nicolas Brunel, Patricia S Goldman-Rakic, and Xiao-Jing Wang. Synaptic mechanisms and network dynamics underlying spatial working memory in a cortical network model. Cerebral cortex, 10(9):910-995 923, 2000.
[12] Ran Darshan and Alexander Rivkind. neous neural networks. Cell Reports, 39(1): 10612, 2022.
Learning to represent continuous variables in heteroge-
[13] M. C. Fuhs and D. S. Touretzky. A spin glass model of path integration in rat medial entorhinal cortex. Journal of Neuroscience, 26 (6(16):4266-4 4276, 2006.
[14] Shintaro Funahashi, Matthew V Chafee, and Patricia S Goldman- Rakic. Prefrontal neuronal activity in rhesus monkeys performing delayed anti-saccade task. Nature, 365(6448):753-756 1993.
[15] Richard J Gardner, Erik Hermansen, Marius Pachitariu, Yoram Burak, Nils A Baas, Benjamin A Dunn, May-Britt Moser, and Edvard I Moser. Toroidal topology of population activity in grid cells. Nature, pages 1 1-6, 2022.
[16] Mikhail Genkin and Tatiana A Engel. Moving beyond generalization to accurate interpretation of flexible models. Nature machine intelligence, 2 2(11):674-63, 2020.
[17] A. Guanella, D. Kiper, and P. Verschure. A model of grid cells based on a twisted torus topology. International Journal of Neural Systems, 17(04)2311-240, aug 2007.
[18] J. Hopfield. Neural networks and physical systems with emergent collective computational abilities. Proceedings of the national academy of sciences, 79(8):2554-255 1982.

---

[19] Vladimir Itskov, David Hansel, and Misha Tsodyks. Short-term facilitation may stabilize parametric working memory trace. Frontiers in Computational Neuroscience, 5:40. 2011.
[20] Ingmar Kanitscheider and lla Fiete. Training recurrent networks to generate hypotheses about how the brain solves hard navigation problems. Advances Neural Information Processing Systems, 30, 2017.
[21] Sung Soo Kim, Hervé Rouault, Shaul Druckmann, and Vivek Jayaraman. Ring attractor dynamics in the Drosophila central brain. Science, 356(6340):849-883, 2017.
[22] Thomas Langlois, Haicheng Zhao, Erin Grant, Ishita Dasgupta, Tom Griffiths, and Nori Jacoby. Passive attention in artificial neural networks predicts human visual selectivity. Advances in Neural Information Processing Systems, 14:27094-27100 2021.
[23] David MacNeil and Chris Eliasmith. Fine-tuning and the stability of recurrent neural networks. PLoS ONE, 6())x22285,, 2011.
[24] Valerio Mante, David Sussillo, Krishna Shenoy, and William Newsome, Context-dependent computation by recurrent dynamics in prefrontal cortex. nature, 503(7474):78- 84, 2013.
[25] R. Monasson and S. Rosay. Crosstalk and transitions between multiple spatial maps in an attractor neural network model of the hippocampus: Phase diagram. Physical Review E Statistical, Nonlinear, and Soft Matter Physics, 87(6):062813, 2013.
[26] R. Monasson and S. Rosay. Crosstalk and transitions between multiple spatial maps in an attractor neural network model of the hippocampus: Collective motion of the activity. Physical Review E, 89())033803, mar 2014.
[27] R. U. Muller and J. L. Kubie. The effects of changes in the environment on the spatial firing of hippocampal complex-spike cells. Journal of Neuroscience, 777):1951-1988, 1987.
[28] J. O' Keefe and J. Dostrovsky. The hippocampus as spatial map: preliminary evidence from unit activity in the freely-moving rat. Brain Research, 34:171-175, 1971.
[29] L Personnaz, Guyon, and G Dreyfus. Collective computational properties of neural networks: New learning mechanisms. Physical Review A, 34(5):4217, 1986.
[30] Kanaka Rajan, Christopher D Harvey, and David W Tank. Recurrent network models of sequence generation and memory. Neuron, 90(1))122- 142, 2016.
[31] A. David Redish, Adam N. Elga, and David S. Touretzky. A coupled attractor model of the rodent Head Direction system. Network: Computation in Neural Systems, 7(4):671-685, 1996.
[32] Alfonso Renart, Pengcheng Song, and Xiao Jing Wang. Robust spatial working memory through homeostatic synaptic scaling in heterogeneous cortical networks. Neuron, 38(3):473 485, 2003.
[33] Blake A Richards, Timothy P Lillicrap, Philippe Beaudoin, Yoshua Bengio, Rafal Bogacz, Amelia Christensen, Claudia Clopath, Rui Ponte Costa, Archy de Berker, and Surya Ganguli. A deep learning framework for neuroscience. Nature neuroscience, 22(11):1761-1770) 2019.
[34] Ranulfo Romo, Carlos D Brody, Adrián Hernández, and Luis Lemus Neuronal correlates of parametric working memory in the prefrontal cortex. Nature, 399(6735):470-473 1999.
[35] Erik Rybakken, Nils Baas, and Benjamin Dunn. Decoding of neural data using cohomological feature extraction, 2019.
[36] A.
Samsonovich and B. L. McNaughton. Path integration and cogmitive mapping in a continuous attractor neural network model. Journal of Neuroscience, 17(15):5900-5920, 1997.
[37] Friedrich Schuessler, Francesca Mastrogiuseppe, Alexis Dubreuil, Srdjan Ostojic, and Omri Barak. The interplay between randomness and structure during learning in RNNs. Advances in Neural Information Processing Systems, 33:13352-11366. 2020.
12

---

[38] Alexander Secholzer, Moritz Deger, and Wulfram Gerstner. Stability of working memory in con- tinuous attractor networks under the control of short-term plasticity. Stability of working memory in continuous attractor networks under the control of short-term plasticity, 2018.
[39] Johannes D. Seelig and Vivek Jayaraman. Neural dynamics for landmark orientation and angular path integration. Nature, 521(7551): 2015.
[40] H. Sebastian Seung. How the brain keeps the eyes still. Proceedings of the National Academy of Sciences, 93(23):13339-1334 1996.
[41] W. E. Skaggs, J. J. Knierim, H. S. Kudrimoti, and B. L. McNaughton. A model of the neural basis of the rat's sense of direction. Advances in neural information processing systems, pages 173-180, 1995.
[42] Ben Sorscher, Gabriel C Mel, Samuel A Ocko, Lisa M Giocomo, and Surya Ganguli. A unified theory for the computational and mechanistic origins of grid cells. Neuron, 11(112121-137, 2023.
[43] David Sussillo. Neural circuits computational dynamical systems. Current opinion in neurobiology, 25:156-163, 2014.
[44] David Sussillo and Larry F Abbott. Generating coherent patterns of activity from chaotic neural networks. Neuron, 63(4):544-557 2009.
[45] J. Wills, C. Lever, Cacucci, N. Burgess, and Oeefe.. Attractor dynamics in the hippocampal representation of the local environment. Science (New York, N.Y.J, 308(5723):873- 876, may 2005.
[46] Klaus Wimmer, Duane Q Nykamp, Christos Constantinidis, and Albert Compte. Bump attractor dynamics in prefrontal cortex explains behavioral precision in spatial working memory. Nature neuroscience, 17(3):31-439, 2014.
[47] Guangyu Robert Yang, Madhura Joglekar, H Francis Song, William T Newsome, and Xiao- Jing Wang. Task representations in neural networks trained to perform many cognitive tasks. Nature neuroscience, 22(2):297-306, 2019.
[48] K. Yoon, M.A. Buice, C. Barry, R. Hayman, N. Burgess, and L.R. Fiete. Specific evidence of low-dimensional continuous attractor dynamics in grid cells. Nature neuroscience, 16(8):1077, 2013.
[49] K. Zhang. Representation of spatial orientation by the intrinsic dynamics of the head-direction cell ensemble: A theory. Journal of Neuroscience, 16(6):2112-211 1996.
13


Reviews:
[{"review_id": "PTshLyt4Jt", "reviewer": "Reviewer_m6Fs", "summary": "The authors present a new technique for embedding multiple attractor manifolds into RNNs. To do so, they first randomly choose a number of attractor manifolds, embed these into an RNN, and then make weight adjustments to smooth out the interference created by multiple manifolds. The authors propose two strategies for weight adjustment, which consider first- and then second-order interference effects. These adjustments are shown to iteratively improve several intuitive metrics.", "questions": "the phrase \"quenched noise\" is used a lot - can this be defined more clearly upon introduction? My idea of \"quenched\" is equivalent to \"reduced\", which does not seem to align with how this phrase is being used.\n\nL15: typo: adjustment -> adjustments\nL151: typo: \"emfbedded\"\n\nFig 1b top: only 9 red lines - need to expand axis limits?\n\nFig 2b top: it is interesting that modification - while flattening the energy landscape - also raises the overall energy of the system. It would be instructive if the authors could offer some intuition for why this is the case. Also, in Fig 4a it appears that the energy after second-order modification is actually less than the single map case - will this always be the case? How should I think about this result?\n\nFigs 3-5: what are error bars? The captions say mean+/-SEM, but what are these statistics computed over? Are there multiple network instantiations for each value of L?", "limitations": "One limitation that is already pointed out by the authors is that their approach is not a biologically plausible learning algorithm; however, I agree that this proof-of-principle is a solid first step and that investigating other mechanisms for weight updates is an interesting direction for future work.", "rating": 8, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "This paper addresses and interesting and well-defined problem. This is not my area of expertise, but the results seem quite general to understanding RNN function, and thus significant. \n\nThe paper is well written and the results clearly presented. The stability metrics are intuitive.\n\nThe presented solutions are strikingly effective!", "weaknesses": "While the results are very strong, the paper only explores a single task (embedding 1D ring attractors). Do the results hold when moving to, say, two dimensions? How does the number of neurons in the network interact with task complexity?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": ">While the results are very strong ... task complexity?\n\nIn order to reduce the computational cost, we explored our schemes in 1D, but it is straightforward to extend our approach and implement it in 2D. Conceptually, we do not expect a qualitative difference when doing so. One notable difference between 1D and 2D environments is that the number of neurons required to achieve a good approximation to a continuous attractor, even for a single map, scales in proportion to the area in 2D, as opposed to length in 1D. But for a given number of neurons, there is no substantial difference between the two cases in terms of the complexity of the problem: the number of equations scales as N*L, and the number of parameters (synaptic weights) scales as N^2. Since the random permutations are completely unrelated to the spatial organization of the firing fields, quenched (frozen) noise is expected to behave similarly in the two cases.\n\nPlease see also our last response to reviewer fsQk.\n\n>the phrase \"quenched noise\" ... being used.\n\nThe term 'quenched noise' is used in statistical physics to describe frozen disorder in many-particle systems (as opposed to noise and disorder arising from dynamic fluctuations). We defined this term in our context in the introduction (line 54): it describes a hard-wired (and thus also referred to as frozen) noise in the connectivity which is independent of time. This is in contrast to dynamical noise that affects instantaneous neural firing rates (see also Refs. [23,5]).\n\nWe will add the following clarification in the final version:\n“Unlike stochastic dynamical noise which is expressed in instantaneous neural firing rates, time-independent quenched noise is hard-wired in the connectivity – it breaks the symmetry between the continuum of neural representations.”\n\n>L15: typo: adjustment -> adjustments L151: typo: \"emfbedded\"\n\nBoth fixed.\n\n>Fig 1b top: only 9 red lines - need to expand axis limits?\n\nYes! Fixed, thanks.\n\n>Fig 2b top: it is interesting ... why this is the case.\n\nThe (flat) modified energy landscapes in Fig. 2b (top and bottom, orange traces) are not the true energy landscapes, because the energy is evaluated using idealized bumps and these are not true steady states of the multi-map network. These precisely flat energy traces are shown as sanity check, to demonstrate that the added weight modifications indeed produce the expected outcome of a flat energy landscape, when re-evaluating the energy landscape using idealized bumps. For this reason, the comparison between the absolute energy values of the pre-modified (blue) and modified (orange) landscapes in Fig. 2b is not fully justified.\n\nSubsequently in the manuscript, we performed a precise evaluation of the energy landscape, by replacing the idealized bumps with steady states obtained through the constrained gradient optimization scheme. The appropriate figure to look at is therefore Supplementary Fig. 3b, where we show the precise energy landscape across multiple iterations. Here it is evident that the overall (mean) energy does increase, as suggested by the reviewer.\n\nIn our scheme, we do not seek to reach a prescribed value of the energy in each iteration, but only to equalize the energy across states at a value which is determined implicitly in the optimization scheme. This introduces a bit more freedom in the choice of the correction weights than if targeting for a prescribed energy. We suspect that the increase of the energy in our scheme is specific to the choice of the single-map connectivity that we worked with, and that other choices of this connectivity could elicit a decrease of the energy during optimization. We did not check this systematically, however, since uniform shifts of the energy landscape are inconsequential for the stability and dynamics of the bump states.\n\n>Also, in Fig 4a ... about this result?\n\nWe did not show how the energy depends on the number of embedded maps. We realize now that Fig. 2a might be confusing in this respect: throughout the manuscript, we uniformly shifted the energy by a constant in order to allow for the landscapes to be shown on the same plot for a single map and for 10 maps (Fig. 2a). This constant was selected such that the mean energy of idealized bump states across all states and maps (without any weight modification) is zero: thus, the mean of the blue trace in the bottom panel of Fig. 2b is zero (as well as single map blue traces, Fig. 2a). Therefore, comparison of absolute mean energies in Fig. 2a between single map and 10 maps traces is not meaningful. This is now explicitly explained in the final version.\n\nWe can, however, compare the absolute mean energy value for landscapes with a varying number of maps before shifting and centering them around 0: when doing so, we find that as more maps are embedded in the connectivity, the absolute mean energy of the landscape decreases. To understand why this happens one should look at the expression for the energy (Eq. 1). The network parameters (Supplementary Material) were chosen such that the sum of rows/columns of the basic connectivity matrix is a constant value. Therefore, each additional map embedded in the connectivity will contribute, when averaging over a random permutation, a constant shift to the energy through the third term of Eq. 1. This constant is negative for our choice of the single-map connectivity matrix.\n\nIn summary, there is a systematic shift of the mean energy with addition of new maps that depends on the specific choice of the single-map connectivity matrix, and occurs even without any weight modifications. The weight modifications introduce additional shifts in the energy, but these are fairly subtle compared to the pre-modified dependence of the mean energy on the number of maps.\n\n>Figs 3-5: ... value of L?\n\nYes, please see Supplementary Material lines 138-141. We added in the final version a reference to the Supplementary Material in the captions of Figs 3-5 (after mentioning SEM).\n"}, {"review_id": "6EhBMpjKJ2", "reviewer": "Reviewer_vC81", "summary": "This work tackles the problem of interference between continuous attractors when they are held in a single RNN.  The authors adopted a Lyapunov function as a depiction of energy of the network and tried to flatten the energy landscape of attractors by adding a modulation term to the original connection matrix.  The modulation term was derived by two methods respectively: a first-order approximation of the energy function and a constrained iterative gradient-based method.  They showed that by constraining bumps at their initial position, the gradient descent achieved a better result.", "questions": "1. Since the learning method is used to determine M, why the authors use W=\\sum_i J_i+M and learn the modification term, why not learn W from the first space? ", "limitations": "The limitations are discussed.", "rating": 4, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "The work achieved the goal of encoding multiple continuous attractors into a single recurrent network. ", "weaknesses": "1. The recent experimental data actually showed that in the remapping of cognitive maps in hippocampus, place cells encoding different maps actually have little overlap, i.e., the hippocampus recruits different groups of neurons to form different continuous attractors. In other words, the interference between multiple continuous attractors is only a mathematical problem, not a biological problem. This limits the contribution of this study to neuroscience. \n2. There are several issues in this study whose biological plausibility are not justified, the energy function (note the real neuronal connections are not symmetry), the modification of synapse strengths based on the attractors the network has stored, the gradient-based learning method. Overall, the insight of this study to neuroscience is rather limited. \n3. Some important references are missed, such as the work of Misha Tsodyks et al. on stroring multiple continuous attractors (PLoS Computational Biology?)\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": ">The recent experimental data actually showed that in the remapping of cognitive maps in hippocampus, place cells encoding different maps actually have little overlap, i.e., the hippocampus recruits different groups of neurons to form different continuous attractors. In other words, the interference between multiple continuous attractors is only a mathematical problem, not a biological problem. This limits the contribution of this study to neuroscience.\n\nRemapping of CA3 (and CA1) place fields in distinct environments (“global remapping”) has been documented in a large number of studies (see, for example, the following two review articles: \nhttps://doi.org/10.1016/j.tins.2008.06.008, https://doi.org/10.3389/fnbeh.2017.00253).\n\nTypically, in a given small environment, only a subset of CA3 place cells are active. Therefore, some cells will be active only in one environment and not in another. However, in any two environments, there is a significant fraction of cells that express firing fields under both conditions. The hallmark of global hippocampal remapping is that the spatial relationship between these fields is unrelated in the different environments. \n\nAs an example, see:\nhttps://doi.org/10.1038/nature05601\nFigure 1 shows firing rate maps of a few CA3 cells under global remapping, and Supporting Figure 2 (panels A and B) shows rate maps of many additional cells.\n\nWe are therefore not sure why the reviewer comments that the hippocampus recruits distinct, non-overlapping groups of neurons to represent different environments. We may have misunderstood the comment of the reviewer, and therefore we kindly ask that the reviewer will reexamine this claim.\n\n>There are several issues in this study whose biological plausibility are not justified, the energy function (note the real neuronal connections are not symmetry), the modification of synapse strengths based on the attractors the network has stored, the gradient-based learning method. Overall, the insight of this study to neuroscience is rather limited.\n\nProgress in the theoretical understanding of the brain has benefited from various approaches and styles of research: among these approaches, those based on mathematical abstraction have often proved to be useful, because this facilitates the identification of general principles. For example, the Hopfield model of associative memory was formulated in terms of binary units, with symmetric connectivity, and using many other gross idealizations of the biological reality. Despite these choices – and, in fact, perhaps because of them, the model had tremendous influence on the thinking about short-term and long-term memory, and the relation between neural network organization and function. Numerous later studies expanded on this model. Some examined how the same principles can be realized with more biological realism, while others used the model in unexpected contexts not foreseen originally. \n\nLikewise, most models of head-direction, grid-cell, and hippocampal networks to date assume symmetric connectivity, as in our study. This is not meant to suggest that the question of non-symmetric connectivity is unimportant: however, by providing a proof or principle that multiple continuous attractors can be embedded without compromising the stability (as presumed previously) we provide a new insight – even if we do so at this stage only for symmetric connectivity. It is likely that our work will motivate future studies, in which some of the limitations mentioned by the reviewer (and discussed in the section of our manuscript on limitations) will be confronted. \n\nFurthermore, the value of rigorous mathematical theory is often realized in unexpected contexts. The ring attractor model was initially proposed for orientation selectivity on V1, yet it has motivated the formulation of theoretical models of head direction cells, and has been recently extremely valuable in understanding the physiology and structure of the fly’s central complex, where surprisingly tight correspondence has been observed between network organization (and dynamics) and the seemingly idealized mathematical theories that were developed in previous studies for head-direction cells.\n\nWe therefore believe that our study provides a significant development in the theory of attractor networks in the brain, and that it fits well into the scope of NeurIPS.\n\n>Some important references are missed, such as the work of Misha Tsodyks et al. on stroring multiple continuous attractors (PLoS Computational Biology?)\n\nWe were not sure which paper was meant here, and wonder whether it is this one:\nhttps://doi.org/10.1371/journal.pcbi.1000869.\nThis article proposes a network architecture in which two or more maps are embedded in a single network, in the case where the patterns of activity associated with the maps are correlated. The paper does not address the problem of achieving precise stability, yet it is thematically related to our work in a broad sense because of the study of simultaneously embedded maps. We will be happy to cite it, and will be grateful if the reviewer could provide us with other concrete suggestions for additions to the reference list.\n\n>Since the learning method is used to determine M, why the authors use W=\\sum_i J_i+M and learn the modification term, why not learn W from the first space?\n\nWe first started from the established naive form of connectivity J, which was identified and used in previous studies for embedding multiple attractors. This initial connectivity can be viewed as arising from a Hebbian form of learning, as discussed in previous studies. We implemented a *perturbative* approach in which we started from this known approximate solution to the problem, where the energy landscape (as we show) is not precisely flat. The perturbative approach is key to our methodology (please see lines 137-152), and the weight corrections that are required to rescue the continuity of the attractor are indeed small (Fig. 4b).\n"}, {"review_id": "7yyPWl9rxM", "reviewer": "Reviewer_NH89", "summary": "A new method is proposed to allow the simultaneous embedding of multiple attractors in an RNN through minimization of the energy function corresponding to the dynamics.\nTwo different methods to achieve this are considered, the first one based on the linearized energy function and the second on  constrained optimization. The second method optimizes the weights for a flat energy landscape with the constraint that there is no change in the center of mass of the bumps.\n", "questions": "174: To analyse the modified network, one can use the linearized network or Lyapunov exponents to identify stability. Try to find the equilibria of the modified network and assess the stability of these points. Do the equilibrium states have a Lyapunov exponent close to 0? (This is expected from the ring attractor structure.)\n\n\nCANNs can implement path integration. Is it possible to do that reliably in the proposed network for the different embedded maps?\n\n\nHowever, the proposed idealized bump method to do this does not guarantee that the acquired energy landscape is truly flat. Can you plot the energy landscape for values in between bump positions?\n\nSI 111: Do you mean that the function $f$ is defined \\emph{in terms of} the weighted averaged position of the rates? The definition that follows is different.\n\n132-133: Could you further explain how these steps follow from each other? The gradient of $I^{k,l}$ in S24 is not dependent on the constraint, why couldn't the derivative of it with respect to $M_{ij}$  change in another  direction? \n\n190: Could you further assess the stability of the network by introducing small displacements in the activity? \nWhat are the Lyapunov exponents of the network? Do they correspond to a continuous attractor?\n\n194: Could a different definition of the location lead to a different conclusion from the analysis? \nThe actual attractor states are not the idealized attractors any more necessarily.\nTherefore, the mean square change defined like this might be biased.", "limitations": "\nThey are adequately discussed.\n", "rating": 5, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "The aim to flatten the energy landscape of an attractor network with multiple attractors embedded in it is novel.\nThe main aim and the methods are clearly described.\nThe discussed problem of detrimental inference is definitely very important for theoretical neuroscience.", "weaknesses": "\nThere is a lack of assurance that the modified network actually has the equalized energy function maintains the bumps as minima of the energy landscape.\n\n\nThe procedure also doesn't necessarily contribute to a flat energy landscape, it just ensures that all the evaluated states have the same energy value.\n\n\nA more thorough way to enforce flattening energy functions for RNNs is described in: Noorman, Marcella, et al. \"Accurate angular integration with only a handful of neurons.\" bioRxiv (2022): 2022-05.\n\n\nThe time to implement the algorithm is also really long 72 hours for $L=60$ and does not seems to be practical.\n\nThe reason to use the bump score is not fully justified. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": ">There is ... landscape.\n\nExistence of the Lyapunov function rigorously guarantees that network dynamics will converge to stationary states which are local minima of the energy (lines 96-97). We have precisely characterized these minima, and carefully verified our numerical scheme by checking that the energy function is locally quadratic near these minima. Nevertheless, we agree that there is value also in checking independently that convergence takes place, as additional verification of the numerical scheme. Please see Supp Figs. 1 and 2b.\n\n>The procedure ... value.\n\nThis is true only for the first scheme in which we evaluate the energy of idealized bump states. Next, we find the true minima of the energy, constrained on the center of mass of the bump thereby precisely evaluating the energy landscape as a function of bump position. Our optimization scheme then aims to flatten this landscape, and we succeed in doing so as shown in Fig. 4c and Supp Fig. 3b.\n\n>A more ... (2022):\n\nNoorman et al aimed to generate a *single* continuous ring attractor with a small number of neurons. Their approach, like ours, is based on the flattening of a Lyapunov energy landscape. Their work treats a very specific form of connectivity with cosine weights that enables the derivation of some results analytically, but is difficult to generalize to other forms of connectivity. In this sense, our approach is more general. In fact, we successfully used our approach to address the problem of Noorman et al using other forms of initial connectivity. In the *multi-map* problem addressed in our submission, the random permutations yield highly irregular weight profiles, and the problem of flattening the landscape is considerably more difficult.\n\n>The time ... practical.\n\nOur focus is on asking a basic theoretical question: whether it is possible to embed multiple flat continuous attractors in a single network? We did not put effort into optimizing our code or running it on state-of-the-art hardware and it is likely that there are more efficient ways to do so. For example, much of the ~72 hours mentioned above is due to our lack of an automated pipeline for gathering each iteration results and generating the next one. The key value of our work is in showing that the flattening is possible, and not in the efficiency of the algorithm used to obtain this result.\n\n>The reason ... justified.\n\nIn the literature on Hopfield based models, the retrieval of a memory is obtained by computing the intuitive overlap measure between the memory pattern and the network state. This is equivalent to the computation of position using our bump score. Previous works that studied multiple embedded maps have used this measure as well [32,23,1]. For completeness, we also used the population vector to measure the position of localized activity. Both measures have shown a dramatic improvement in the network drift (Supp Fig. 4). We now clarify this choice in the SM.\n\n>174: ... structure.)\n\nSince rate networks with symmetric weights are guaranteed to settle on stationary steady states, it’s sufficient to examine the Jacobian of the dynamics to assess stability. If the attractor is continuous, one expects to observe an eigenvalue close to zero. We didn’t check this explicitly, but the slow dynamics over a continuum of near-steady states and the near-flatness of the energy function (Figs. 4c and 5a,d,e, and Supp Figs. 1 and 3b) implies so.\n\n>CANNs ... maps?\n\nIt has been argued that path integration within the hippocampus would be difficult, because the synaptic connectivity required to do so would need to be tailored for each embedded map. Thus, it was suggested that path integration occurs elsewhere - perhaps in the entorhinal cortex. This idea was recently examined in Ref. [1], demonstrating that path integration can reliably be implemented despite the interference between the maps. The flattening procedure discussed in our work is expected to substantially improve the accuracy of this computation.\n\n>However, ... positions?\n\nAs the reviewer notes, Fig. 2a,b for 10 maps only shows an approximation to the energy landscape, because the idealized bump states are not true minima of the energy. We later evaluate the precise energy landscape by finding states that minimize the Lyapunov function under a constraint on their center of mass. This is done at a dense set of positions (Fig. 4a and Supp Fig. 3b), with single neuron resolution. Even the classical ring attractor composed of N neurons is not precisely continuous, but has N precisely steady states. This is an issue of practical significance only when the number of neurons is small (Noorman et al).\n\n>SI 111: ... different.\n\nYes, fixed.\n\n>132-133: ... direction?\n\nI^kl in Eq. S24 is defined as the state that minimizes the energy under the constraint. This is where the dependence on the constraint is coming from. Any modification of I^kl due to changes in M must remain in the subspace in which the constraint is obeyed. We modified the SM to clarify this point. In addition, note that we numerically verified that the second term in Eq. S24 vanishes (Supp Fig. 2b).\n\n>190: ... attractor?\n\nThe steady states that we identify are minima of the energy under the constraint on the center of mass, and therefore the manifold of states that we identify is stable. The stability of the approximate attractor is evident in Fig. 5a,d,e. See above response regarding the Lyapunov exponent.\n\n>194: ... biased.\n\nThis is true for our idealized bump scheme (Fig. 3c), where the mean squared change (MSC) is compared between the converged and the initialized idealized bump states. However, in Fig. 5 we quantify the stability of states that are true minima of the energy under a constraint on the position of the bump: we start from these deformed bump states and follow their dynamics (SM lines 84-85). Since the MSC is compared between deformed converged states of consecutive iterations there is no bias in relation to idealized bumps."}, {"review_id": "8uccfZ6QgD", "reviewer": "Reviewer_fsQk", "summary": "The paper studies the storage of multiple continuous attractors in a recurrent neural network. Specifically, the authors tackle the interference between attractors and its effect on activity bump drift. By using a perturbative approach, they compute a correction to the connectivity that reduces the drift dramatically.\nContinuous attractors (e.g., ring model) are important models in neuroscience, and understanding them is an important task. Multiple attractors are relevant, for example, in the CA3 region of the hippocampus, where remapping of place cells in different environments is common. Nevertheless, a naïve connectivity that is a superposition of several ring-connectivities results in only approximate continuous attractors. The result is a few stable points in each attractor, to which dynamics converge.\nThe authors use the perspective of an energy function (Lyapunov), and examine how the interference renders this function non-flat. They then calculate the perturbation to leading order, and solve for a change in connectivity that will flatten the energy. Furthermore, using gradient descent, they are able to achieve even greater precision.\n", "questions": "L146 energy IS quadratic\nL178 “vast majority” the sentence isn’t very clear. Was the intention “close to a fixed point” or something similar?\nL146 – I didn’t fully understand the argument on why the first order vanishes.\nCan you say something about higher dimensions? Place cells are relevant in 2D. The results of [5,23] suggest differences when dimension increases.\n", "limitations": "yes", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "Continuous attractors are a fundamental building block in the study of recurrent neural networks in neuroscience contexts. There are relatively few studies tackling multiple such attractors. The method to reduce interference is novel.", "weaknesses": "First, as the authors note, the resulting connectivity is extremely fine-tuned. This is a known problem with continuous attractors that is not addressed here.\nSecond, the problem and solution are highly related to similar problems in discrete attractors. Interference in Hopfield networks was tackled using pseudo inverse rules, either approximated online or as a global formula. There is no discussion of the relation to these works. The SVM approach of Battista and Monasson (Ref 5) is perhaps a similar example in continuous attractors.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": ">First, as the authors note, the resulting connectivity is extremely fine-tuned. This is a known problem with continuous attractors that is not addressed here.\n\nWe agree (please see also our 7'th response to reviewer MLg2, regarding fine-tuning).\n\n>Second, the problem and solution are highly related to similar problems in discrete attractors. Interference in Hopfield networks was tackled using pseudo inverse rules, either approximated online or as a global formula. There is no discussion of the relation to these works. The SVM approach of Battista and Monasson (Ref 5) is perhaps a similar example in continuous attractors.\n\nThe potential relation to the works on pseudo-inverse based variants of the Hopfield network architecture is interesting, and we thank the reviewer for bringing this up. While it will be interesting to consider this relation in depth, we would like to point out two potential differences from the scheme explored in our work: First, the goal of the pseudo-inverse rule is to precisely embed a set of prescribed patterns as attractor states of the dynamics, whereas we are willing to accept distortions in the bump states and only care about (i) keeping them localized and (ii) precisely equalizing their energy. We suspect that an attempt to embed the idealized bumps as perfect steady states (using perhaps a generalization of the pseudo-inverse rule) would prove more difficult than our goal, yet this is an interesting question and our intuition on this is only an initial thought. The second difference is that the key benefit of the pseudo-inverse method is in increasing the capacity of the network. There are some conceptual issues in relating a definition of capacity in our problem to the definition in the discrete and binary case, but if we judge the capacity by the number of maps at which delocalized states start to appear it seems that our approach does not strongly affect the capacity (Fig. 5c), whereas its main outcome is in the flattening of the energy landscape.\n\nWe now briefly mention in the Discussion the potential interest in exploring a generalization of the pseudo-inverse rules, as a means of enhancing the stability of states along simultaneously embedded continuous attractors.\n\n>L146 energy IS quadratic\n\nFixed, thanks.\n\n>L178 “vast majority” the sentence isn’t very clear. Was the intention “close to a fixed point” or something similar?\n\nYes, this was the intention, and we thank the reviewer for pointing out that the phrasing was unclear.\n\nWe now revised it by: “...indicating that the convergence to a steady state was nearly complete (Fig. 3a).”\n\n>L146 – I didn’t fully understand the argument on why the first order vanishes.\n\nWe now revised the sentence in lines 146-148 as follows:\n\n“The precise structure of these states affects the second and third terms on the right hand side of Eq. (2) via the deviation of the state from an idealized bump. This contribution to the energy modification is quadratic in the deformations, because the idealized bump states are minima of the unperturbed energy functional, and therefore the energy functional is locally quadratic near these minima.”\n\nWe hope that this helps clarify the argument. Please let us know if further clarification is required.\n\n>Can you say something about higher dimensions? Place cells are relevant in 2D. The results of [5,23] suggest differences when dimension increases.\n\nRef. [23] has shown that unlike in 1D, activity bumps can bypass an energy barrier in 2D when a constant external input (or force) is applied. Intuitively, this happens since in 1D there is only a single possible direction for the bump’s motion along the applied force while in 2D there may be multiple directions (each with an overall smaller force), in which the bump could bypass the energy barrier. The question of bypassing energy barriers given external inputs is, however, decoupled from the question of systematic drifts that we address here, which are independent of external inputs. Comparing, perhaps, the typical magnitude of these systematic drifts in a naive 1D vs 2D network is indeed interesting and will be explored in future work. Nevertheless, the proof of principle provided here in 1D should be easily extended to 2D, and so it is expected that emerging systematic drifts in a naive 2D network will be attenuated using the methodology we demonstrated in 1D.\n\nThere are N true steady states in a continuous attractor formed by N neurons, independently from a 1D or 2D organization. In order to produce an approximate continuous attractor, these states should densely tile the represented space, and this implies a different scaling with the size of the environment in 2D (quadratic) vs 1D (linear). Once the number of neurons is sufficiently large to allow for the single-map attractor to be nearly continuous, we do not see any conceptual difference in the problem of flattening of the energy landscapes for simultaneously embedded attractors between 1D and 2D. In both cases, the number of equations that need to be solved scales as N*L, whereas the number of weights scales as N^2. While the flattening of the energy landscape is expected to behave very similarly, we do not know how the number of true minima of the energy landscape (which is never precisely flat, even in our scheme) will behave. We did not quantify this quantity, as our main focus was on stability, or in other words, on reducing the systematic drifts to minimum. This latter question is more directly related to the quantifications in Ref. [5], where the key interest was in the number of true attractors, and not on the speed of convergence when starting at intermediate states. Relating our results, in a rate model network, to those of Ref. [5] (in a binary network) is interesting, as mentioned in the Discussion.\n\nPlease see also our first response to reviewer m6Fs."}, {"review_id": "mFAdhGMOt0", "reviewer": "Reviewer_MLg2", "summary": "The paper investigates the challenge of embedding multiple continuous attractor manifolds within a single RNN, with a focus on hippocampal place cells. The issues arise due to the presence of discrete steady states, visualized as minima on an abstract energy landscape, which disrupt the continuity of network activity patterns. This disruption prompts systematic drift of population activity patterns towards these discrete states, resulting in degraded memory over time. Past studies have considered the stabilizing influence of external stimuli; however, solutions in their absence remain unclear. The authors address this issue by modifying the synaptic weight to flatten the energy landscape, showcasing through simulations how this significantly stabilizes the activity pattern.", "questions": "o\tHow applicable is your approach if connection weights are not symmetrical?\n\no\tHow robust is the pattern stability if the modified weights are slightly perturbed?\n\no\tCan your approach provide biological insights into the structure of weight modification when new maps are added?\n\no\tDoes the resulting weight (after the modification) violate Dale’s law?\n\no\tWhile it's evident that an energy landscape embedded with minima leads to activity patterns settling into one of them, thereby degrading the network’s ability to sustain persistent memory, does a flat energy landscape guarantee pattern stability?\n\no\tHere are some additional related works that could be included in the discussion [1-2].\n\no\tReferences:\n\n\t[1] Genkin and Engel, Nature Machine Intelligence, 2020: https://www.nature.com/articles/s42256-020-00242-6 \n\n\t[2] Whittington et al., Cell, 2020:\nhttps://www.sciencedirect.com/science/article/pii/S009286742031388X \n", "limitations": "The authors have acknowledged certain limitations, such as the unclear biological plausibility of their approach. ", "rating": 7, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "o\tThe methodology and problem formulation are clearly articulated.\n\no\tSimulations indicate that weight modification significantly improves the stability of activity patterns.\n\no\tThe authors provide the code in the appendix. \n", "weaknesses": "o\tWhile the manuscript is generally well-written, some areas could benefit from improvement. A diagrammatic illustration could help elucidate the issue of discretized states in the context of the energy landscape. Further details of the simulations (e.g., duration and specific inputs to the network) need to be provided. The appendix, which helps clarify the methodology, should be referred to more often in the Results section.\n\no\tStrengthening the paper could involve additional simulations and/or discussions (see questions below).\n\no\tThere is a minor typo: \"emfbedded map\" should be \"embedded map.\" \n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": ">While the manuscript is generally well-written, some areas could benefit from improvement. A diagrammatic illustration could help elucidate the issue of discretized states in the context of the energy landscape.\n\nWe agree that it will be helpful to include a schematic illustration to clarify the concept of discrete minima in the energy landscape. Since the final version can include one additional page, we plan to include such illustration if accepted. We thank the reviewer for suggesting this.\n\n>Further details of the simulations (e.g., duration and specific inputs to the network) need to be provided.\n\nThe Supplementary text includes full details of the network dynamics, parameters, the numerical integration scheme used to simulate the dynamics, the optimization scheme, and quantifications. We will be grateful if the reviewer would let us know in case something specific seems to be missing. Note that there are no training inputs in our study.\n\n>The appendix, which helps clarify the methodology, should be referred to more often in the Results section.\n\nWe thank the reviewer for this comment and now refer to the Supplementary Material more frequently in the Results section.\n\n>Strengthening the paper could involve additional simulations and/or discussions (see questions below).\n\nSee responses below.\n\n>There is a minor typo: \"emfbedded map\" should be \"embedded map.\"\n\nFixed, thanks!\n\n>How applicable is your approach if connection weights are not symmetrical?\n\nOur approach is based on the existence of a Lyapunov function [9], and this relies on symmetric connections. All existing models for embedding multiple attractors in CA3 adhere to this rule. Furthermore, one of the key values of our work, in our view, is that it provides a proof of principle for the existence of a more stable version for a discrete set of continuous attractors than previously expected, which was assumed to be inevitably unstable. In future work, it will be highly interesting to seek biologically plausible learning rules, and to relax the assumption of symmetric connectivity.\n\n>How robust is the pattern stability if the modified weights are slightly perturbed?\n\nOne of the main criticisms of continuous attractor networks is that they require fine tuning of the weights: small perturbations to the weights destroy the continuity of the attractor. This is also why the embedding of multiple maps, using the standard prescription, leads to a wrinkled energy landscape. Considering this, it is a priori expected that the recovery of a flat energy landscape would rely on weights that require fine tuning. One of the interesting features of our results is that the required correction weights are small (see Fig. 4b). A hypothetical learning rule that is based on stability of the attractor states would thus need to explore only small local perturbations to the initial weights, whereas the initial weights could be learned based on a simple Hebbian rule. Even though the correction weights require fine tuning along some dimensions, note that out of an infinite space of potential weight modifications only a specific set (least-squares) was chosen, leaving many unexplored solutions to the energy flattening problem which may relax the fine-tuning requirement to some extent.\n\n>Can your approach provide biological insights into the structure of weight modification when new maps are added?\n\nSee above response. The main correction to the weights is the one arising from the naive embedding scheme. The modifications that we identify in our study are rather subtle. Therefore, attempts to relate plasticity under exposure to new maps to theoretical constructs should first focus on the predictions of the naive scheme.\n\n>Does the resulting weight (after the modification) violate Dale’s law?\n\nLike many theoretical models of attractor dynamics in the hippocampus (and elsewhere in the brain) we abstract away the distinction between excitatory and inhibitory connections. One imagines that inhibitory connections are realized via interneurons.\n\nThe unmodified synaptic connections in our model can be thought of as consisting of all-to-all inhibitory connections (mediated by interneurons), and excitatory connections which are specific to neuron pairs that have similar tuning in any one of the maps. Therefore, it is easy to reformulate the model in accordance with Dale’s law.\n\nSince in our model each neuron forms ~100 excitatory connections associated with each map, almost all neurons end up being connected to each other with an excitatory connection when the number of embedded maps exceed 6 (since N=600). The small modifications to the weights that are required to flatten the energy landscape could thus be implemented as a subtle increase (or decrease) in the strength of the excitatory weights. Alternatively, the weight modifications could be applied to both excitatory and inhibitory synapses. It will be interesting to flesh out these ideas in a follow-up work.\n\n>While it's evident that an energy landscape embedded with minima leads to activity patterns settling into one of them, thereby degrading the network’s ability to sustain persistent memory, does a flat energy landscape guarantee pattern stability?\n\nYes. The existence of a Lyapunov function guarantees that dynamics will settle on local minima of the function. When minima form a manifold (a continuous set of states, all sharing the same energy), stability of each state along this manifold is guaranteed. The property of this Lyapunov function is that whenever dI/dt≠0 → dE/dt<0, and thus a flat energy landscape guarantees pattern stability.\n\n>Here are some additional related works that could be included in the discussion [1-2].\nReferences:\n[1] Genkin and Engel, Nature Machine Intelligence, 2020: https://www.nature.com/articles/s42256-020-00242-6\n[2] Whittington et al., Cell, 2020: https://www.sciencedirect.com/science/article/pii/S009286742031388X\n\nThanks, we will incorporate these references."}]

Author rebuttal:
We thank the reviewers for their careful reading of our submission and for their insightful comments. Please see our point-by-point responses to each review. We will highly appreciate any additional comments or requests for clarification that may arise in response to our answers to the questions. 


Meta-review: This paper presents a new method for embedding multiple attractor manifolds into recurrent neural networks. Naive embedding of such manifolds lead to interference between them. The method described here fine tunes the synaptic weights in a way to reduce this interference. This is an important contribution that will be helpful in building mechanistic models of brain function.
Decision: Accept (poster)