Abstract

Attractor models are simplified models used to describe the dynamics of firing rate profiles of a pool of neurons. The firing rate profile, or the neuronal activity, is thought to carry information. Continuous attractor neural networks (CANNs) describe the neural processing of continuous information such as object position, object orientation, and direction of object motion. Recently it was found that in one-dimensional CANNs, short-term synaptic depression can destabilize bump-shaped neuronal attractor activity profiles. In this article, we study two-dimensional CANNs with short-term synaptic depression and spike frequency adaptation. We found that the dynamics of CANNs with short-term synaptic depression and CANNs with spike frequency adaptation are qualitatively similar. We also found that in both kinds of CANNs, the perturbative approach can be used to predict phase diagrams, dynamical variables, and speed of spontaneous motion.

1  Introduction

Neurons communicate with each other through the neurotransmitter diffusion initiated by action potentials, or spikes, and the activity of one neuron can excite or inhibit the activity of another neurons. The firing rates of spike trains are thought to carry information, and correlations between the firing rates of neurons depend on the strength of the couplings between those neurons. With some specific settings of couplings, such as Mexican hat couplings (Amari, 1977; Ben-Yishai, Bar-Or, & Sompolinsky, 1995), networks of neurons can support a continuous family of local neuronal activity profiles on a field, which can be used to represent continuous information, such as object position, object orientation, and direction of object motion.

Local neuronal activities associated with continuous information are observed in various brain regions. Typical examples of cells showing such activity are head direction cells (Taube, Muller, & Ranck, 1990; Blair & Sharp, 1995; Zhang, 1996; Taube & Muller, 1998), place cells (O’Keefe & Dostrovsky, 1971; O’Keefe, 1976; O’Keefe & Burgess, 1996; Samsonovich & McNaughton, 1997) and moving-direction cells (Maunsell & van Essen, 1983; Treue, Hol, & Rauber, 2000). These cells have gaussian-like tuning curves as functions of stimulus. Among the numerous models proposed to describe behaviors of these systems are continuous attractor neural networks (CANNs), and a recent study of persistent activity in monkey prefrontal cortices has provided evidence of continuous attractors in the central nervous system (Wimmer, Nykamp, Constantinidis, & Compte, 2014). Persistent neuronal activity in monkeys’ prefrontal cortices was discovered during delayed-response tasks (Funahashi, Bruce, & Goldman-Rakic, 1993). The study by Wimmer et al. (2014) confirmed that pairwise neuronal correlation predicted by theories can be observed in the brain region they investigated (Ben-Yishai et al., 1995; Pouget, Zhang, Deneve, & Latham, 1998; Wu, Hamaguchi, & Amari, 2008).

A one-dimensional (1D) CANN can support a family of gaussian-like neuronal activity profiles. They are attractors of the network. As shown in Figure 1a, attractors have the same shapes and are centered at positions corresponding to different preferred stimuli. These gaussian-like profiles can shift smoothly along the space of attractors. Families of attractor profiles in two-dimensional (2D) CANNs are similar: each attractor profile is a gaussian-like profile and centered at a particular position in the 2D field (see Figure 1b). The dynamics of a network state profile to track a stimulus is widely studied in the literature (Ben-Yishai et al., 1995; Samsonovich & McNaughton, 1997; Wu et al., 2008; Fung, Wong, & Wu, 2010). If couplings between neurons are static over time (quenched), the steady state of neuronal activity profiles will be static because of the homogeneity and translational invariance of CANNs. If they depend on the firing histories of presynaptic neurons, however, the dynamics of CANNs can be different.

Figure 1:

Examples of attractor states of continuous attractor neural networks (CANNs) in (a) a 1D field and (b) a 2D field. (a) Profiles are centered at location z, the location they are representing. (b) x0 is the first coordinate of the preferred stimulus, and x1 is the second coordinate of the preferred stimulus.

Figure 1:

Examples of attractor states of continuous attractor neural networks (CANNs) in (a) a 1D field and (b) a 2D field. (a) Profiles are centered at location z, the location they are representing. (b) x0 is the first coordinate of the preferred stimulus, and x1 is the second coordinate of the preferred stimulus.

Tsodyks and Markram (1997) proposed a model in which the synaptic efficacies between neurons depend on the amount of available neurotransmitters in the presynaptic neurons, and this amount depends on the firing history of the presynaptic neuron (Tsodyks and Markram, 1997; Tsodyks, Pawelzik, & Markram, 1998). This kind of reduction in synaptic efficacies, due to past presynaptic neuronal activity, is called short-term synaptic depression (STD). There are reports that short-term synaptic depression can exhibit rich dynamics in CANNs (York & van Rossum, 2009; Fung, Wong, Wang, & Wu, 2012). Fung, Wong, Wang et al. (2012) reported that in 1D CANNs, short-term synaptic depression can destabilize attractor profiles. Within a broad range of strengths of divisive global inhibition, if we increase the degree of STD to a moderate range, static activity profiles will be translationally destabilized. After some translational perturbations, a bump-shaped profile can move spontaneously along the attractor space. In these scenarios, both bump-shaped static states and moving states can coexist. If we further increase the degree of STD, no static profile can be found. If the degree of STD is too large, the steady state of the system can only be a trivial solution. Similar behavior is reported by York and van Rossum (2009) with a different model. We can see that the instability induced by STD can reshape the intrinsic dynamics of the system, even when no stimulus is presented.

Network response in a CANN with static couplings is always lagging behind a continuously moving stimulus. However, with short-term synaptic depression, due to the translational instability, the underlying dynamics of the network can make the response overtake the actual stimulus. Fung, Wong, and Wu (2012) suggested that this behavior can be used to implement a delay compensation mechanism. Short-term synaptic depression can also induce global (Loebel & Tsodyks, 2002) and local (Fung, Wang, Lam, Wong, & Wu, 2013; Wang, Lam, Fung, Wong, & Wu, 2014) periodic excitements of neuronal activity profiles. These periodic excitements can enhance information processing in the brain. Fung et al. (2013) recently proposed that periodic excitement driven by the intrinsic dynamics can improve the resolution of CANNs. Kilpatrick (2013) also proposed that the neuronal activity pattern may shift between one stimulus and another. These theories suggested that STD may enhance the capability of CANNs.

Studies on CANNs with short-term synaptic depression are mainly on 1D networks. In this article, we discuss the intrinsic dynamics of bump-shaped solutions of two-dimensional CANNs with short-term synaptic depression. Other 2D models possess rich dynamical behaviors such as spiral waves (Kilpatrick & Bressloff, 2010a), breathing pulses (Kilpatrick & Bressloff, 2010a), and collisions of two bump-shaped profiles (Lu, Sato, & Amari, 2011). Here, we focus on the spontaneous motion of a single bump-shaped profile and analyze its stability. We have also studied the influence of STD on the sizes of bump-shaped profiles and their changes in shape.

We study not only CANNs with STD but also CANNs with spike frequency adaptation. Spike frequency adaptation (SFA) is a dynamical feature commonly observed in neurons. Neurons are suppressed after prolonged firing. SFA can be generated by a number of mechanisms (Brown & Adams, 1980; Madison & Nicoll, 1984; Fleidervish, Friedman, & Gutnick, 1996; Benda & Herz, 2003). It can also destabilize the amplitudes and positions of static bumps. SFA-induced destabilization of bump-shaped states in CANNs is reported in the literature (Kilpatrick & Bressloff, 2010b). What we find in our study on CANNs with SFA is similar to what we found in the case with STD. SFA first destabilizes the translational mode and then the amplitudal mode. There is also a parameter region such that both spontaneous-moving-bump solutions and static-bump solutions can coexist.

In this article, in each case, CANNs with STD and CANNs with SFA, we first introduce the model we used to study the problem and then analyze each scenario using the perturbative method proposed by Fung et al. (2010). These sections are followed by a section discussing the comparison between theoretical and simulation results and discussing the limitations of the perturbative method.

2  The Model

In this work we consider a 2D neural field, where neurons are located on a 2D field with positional coordinates . can be interpreted as the preferred stimulus of the neuron sitting at that point in the 2D field. The state of a neural field is specified by the average membrane potential of neurons at at time t. The neuronal activity of a neuron at is given by a nonlinear function of ,
formula
2.1
where is the Heaviside step function and is the divisive global inhibition. The neuronal activity is related to the average firing rate of the neuron. The evolution of the global inhibition is given by
formula
2.2
where k is the parameter controlling the strength of the divisive global inhibition and is the density of neurons over the field. Here, is driven by . This choice of driving term can simplify our calculations. As we will describe in the following, is a weighted sum of . So effectively depends on . Also, at the steady state, the magnitude of is directly proportional to that of .
In the network, neurons are connected by excitatory couplings given by
formula
2.3
where is the norm of the argument. a is the radius of effective excitatory connections and J0 is the intensity of average coupling over the field. The dynamics of is governed by
formula
2.4
is the external stimulus, is the portion of available neurotransmitters of the presynaptic neuron, and is the dynamical variable corresponding to SFA. Since in this work, we are studying intrinsic dynamics in the network, is set to zero throughout the article.
The portion of available neurotransmitters of the presynaptic neuron at , evolves as
formula
2.5
The first term on the right-hand side of this equation is the relaxation of with time constant . The last term is the consumption rate of neurotransmitters. is a parameter proportional to the neurotransmitter consumption due to each spikes. This parameter controls the strength of STD. On the other hand, the dynamics of is given by
formula
2.6
is the timescale of SFA. is the degree of SFA. f is the dependence of SFA variables on average membrane potential, a nondecreasing function. For simplicity, we have chosen
formula
2.7
This choice is convenient for our analytic purpose and should not affect the main conclusion qualitatively as long as the adaptation increases with the average membrane potential, , which is correlated with the average neuronal activity.

In this study, we analyze CANNs with STD and CANNs with SFA separately. For the case of CANNs with STD, we set . For CANNs with SFA, we set .

3  CANNs with STD

3.1  Stationary Solution

In this case, we set and to suppress SFA at the moment. For , , two types of nonzero fixed-point solutions to equation 2.4 exist when :
formula
3.1
where and is an arbitrary location in the field representing the center of local excitation. It is expected that the fixed-point solution with the larger amplitude is stable and the other is unstable. For , consider
formula
3.2
where is the deviation of the profile from the fixed-point solution. Then the differential equation of the deviation is
formula
3.3
Therefore, the smaller fixed-point solutions are unstable, while the larger fixed-point solutions have stable amplitudes. For nonzero , we have to consider the dynamics of . Let , where . The dynamics of the deviation from the fixed-point solution is given by
formula
3.4
Since the trace of the matrix in equation 3.4 is negative for and the sign of the determinant is independent of , the solution with the larger u00 is stable against perturbations in amplitude, while the solution with the smaller amplitude is unstable. The solution of u00 as a function of is shown in Figure 2.
Figure 2:

Solution of as a function of . For each less than 1, the larger solution of is stable and the smaller solution is unstable.

Figure 2:

Solution of as a function of . For each less than 1, the larger solution of is stable and the smaller solution is unstable.

For and , we approximate the attractor profiles of and by nonmoving gaussian distributions:
formula
3.5
formula
3.6
Here, without loss of generality, we consider the case with . They are not the exact solutions of equations 2.2, 2.4, and 2.5. In this ansatz, is assumed to have the same shape as that in the case. The width of is different from because the shape of is similar to at the small limit. The differential equations governing and can be obtained by projecting equation 2.4 onto and projecting equation 2.5 onto .

As in the study by Fung, Wong, Wang et al. (2012), can be replaced by rescaled variables because has a dimension . . And k and can be rescaled by and .

With these rescaled variables, the differential equations of , and are
formula
3.7
formula
3.8
formula
3.9
As we study the steady-state behavior of the network, we need to solve for the fixed-point solution first. For the first two differential equations, we may let be given. Then we can solve and p00. For a given , and are related by
formula
3.10
The parabolas given by this equation for different s are plotted as dashed lines in Figure 3, where we see that the static fixed-point solution exists only when
formula
3.11
By considering the stability of fixed-point solutions, we obtain the solid line in Figure 3 (also the dotted line in Figure 4). This curve maps the parameter region for the existence of static profiles of . The methodology for studying stability of fixed point solutions can be found in appendix  A.
Figure 3:

Parameter regions of existence of static profiles of . Dashed curves: parabolas given by equations 3.10 for different s. Solid line: boundary separating parameter regions of static profiles of and silent phase.

Figure 3:

Parameter regions of existence of static profiles of . Dashed curves: parabolas given by equations 3.10 for different s. Solid line: boundary separating parameter regions of static profiles of and silent phase.

Figure 4:

Phase diagram of different phases under this theoretical framework. Circles: simulation results on translational stability. Squares: simulation results on amplitude stability. Triangles: simulation results on stability of moving bumps. Curves: corresponding theoretical predictions. Parameter: .

Figure 4:

Phase diagram of different phases under this theoretical framework. Circles: simulation results on translational stability. Squares: simulation results on amplitude stability. Triangles: simulation results on stability of moving bumps. Curves: corresponding theoretical predictions. Parameter: .

3.1.1  Translational Instability

We have simplified equations 2.2, 2.4, and 2.5 by introducing the approximation given by equations 3.5 and 3.6. This simplification, however, is useful for studying only the amplitudal stability of a bump-shaped solution. For the translational stability, we need to consider the stability of the static solution against asymmetric distortions. Displacing originally aligned and profiles of a static solution is a reasonable test, as the dip of is generated by activities of neurons. If the solution is moving, the dip of is always lagging behind. As a result, the asymmetric component of with respect to the center of mass of becomes nonzero. In the calculation, we may drop asymmetric components of for the moment, as we can always choose a frame such that the major asymmetric mode is zero.

Let us assume
formula
3.12
formula
3.13
Due to the symmetry of the preferred stimulus space, here we consider only the distortion along the x0-direction. By substituting these two assumptions into equations 2.2, 2.4, and 2.5, we can study the stability of static solutions against asymmetric distortions. With assumptions equations 3.12 and 3.13, we can derive the stability matrix by calculating the Jacobian matrix at the fixed-point solution. Then the dynamics of distortions of , , , and around the static fixed-point solution is given by
formula
3.14
where the matrix is provided in appendix  A. determines the amplitudal stability of the solution, while M determines translational stability, which is given by
formula
3.15
If , the static solution of two-dimensional CANNs with short-term synaptic depression will obviously be translationally unstable because positional distortion will diverge.
By studying the stability matrix, we found that bump-shaped solutions will be translationally stable only if
formula
3.16
In this case, asymmetric distortions cannot initiate spontaneous motion. When this inequality does not hold, there will be a moving solution such that the bump-shaped profile can move spontaneously with a speed dictated by , , and . This inequality provides a prediction of the boundary separating translationally stable static bumps and translationally unstable static bumps, which is plotted as a solid line in Figure 4.

3.2  Moving Solution

In the one-dimensional situation, when is small and is relatively large, there may be spontaneous moving solutions (Fung, Wong, Wang et al., 2012). To analyze the moving solution, we need to consider higher-order expansions of and . In general, and can be expanded by any basis functions. Here we have chosen eigenstates of quantum harmonic oscillator as basis functions,
formula
3.17
formula
3.18
where
formula
3.19
formula
3.20
In these equations, . ci is the ith component of the velocity of the moving frame. Hn is the nth order physicists’ Hermite polynomial. As in the ansatz for static solutions, equations 3.5 and 3.6, the widths of the basis functions for and are different. This choice makes the perturbative expansion more efficient, although the choice of basis functions can be arbitrary.
After substituting equations 3.17 and 3.18 into equations 2.2, 2.4, and 2.5 and projections on equations 3.19 and 3.20, we obtain
formula
3.21
formula
3.22
formula
3.23
A detailed illustration of the derivation of these equations can be found in appendix  B. Here, s are rescaled dynamical variables: . Coefficients , , , and are defined by
formula
3.24
formula
3.25
formula
3.26
formula
3.27
, , , and can be calculated explicitly. For , , , and with arbitrary k, n, m, and l, the recurrence relations used to generate these coefficients can be found in appendix  B in Fung, Wong, Wang et al. (2012). In practice, we cannot consider infinitely many and . We used finite terms to obtain fairly acceptable results. In Figures 4 and 5, “ Perturbation” means that we consider terms up to .
Figure 5:

(a) Measurements and predictions of . Symbols: measurements from various . Curves: different levels of predictions. (b) Measurements and predictions of . Symbols: measurements from various . Curves: different levels of predictions. (c) Measurements and predictions of . Symbols: measurements from various . Curves: different levels of predictions. (d) Measurements and predictions of the intrinsic moving speed. Symbols: measurements from various . Curves: different levels of predictions. Parameters are the same as Figure 4.

Figure 5:

(a) Measurements and predictions of . Symbols: measurements from various . Curves: different levels of predictions. (b) Measurements and predictions of . Symbols: measurements from various . Curves: different levels of predictions. (c) Measurements and predictions of . Symbols: measurements from various . Curves: different levels of predictions. (d) Measurements and predictions of the intrinsic moving speed. Symbols: measurements from various . Curves: different levels of predictions. Parameters are the same as Figure 4.

Unfortunately, the above equations cannot form a complete set of equations sufficient to solve for the fixed-point solution. We also need to consider the self-consistent condition:
formula
3.28
This self-consistent condition helps us to choose the center of mass of to be origin of basis functions so that the perturbative method can be efficient. , , and ci can be solved numerically by regarding as a constant and setting time derivatives equal to zero. For each , and are related by
formula
3.29
Then fixed-point solutions corresponding to each can be solved.
By comparing the fixed-point solution with the simulation results, we found that not only fixed-point solutions but also various measurements of the network state are determined by as shown in Figure 5. In Figures 5b and 5c, we compare the values of the second-order variables between simulations and theoretical predictions by transforming them to the polar coordinates via and . We can convert second-order basis functions for to be combinations of their polar counterparts. Then we can obtain a linear combination of basis functions in polarbr coordinates:
formula
3.30
In order to compare predictions to simulation results without considering the moving direction in the 2D field, we use (in Figure 5b) and (in Figure 5c) for the case of . is the average change in the width of the bump-shaped profile relative to the closest static solution, while is the magnitude of the anisotropic mode. Illustrations of the polar functions in equation 3.30 are shown in Figure 6.
Figure 6:

Distortional modes for .

Figure 6:

Distortional modes for .

Since B is proportional to , an increase of implies a decrease in . In Figure 5a, the slope of is discontinuous at about . For , the slope is significantly larger than that in the region of . This implies that the motion of a bump helps the bump to maintain its magnitude. It also agrees with the tendency of the average membrane potential (or neuronal activity) profile in 1D CANN that the bump tends to move to a region with a higher concentration of neurotransmitters (York & van Rossum, 2009; Fung, Wong, Wang et al., 2012).

Figure 5c suggests that an anisotropic mode happens only when the bump is not static, while Figure 5b shows that the average width of profile increases with the strength of STD. The behavior of the average change in width, , is similar to that in height. The anisotropic mode happens only when the bump is moving. This suggests that the bump get widened unevenly due to the asymmetric profile.

3.3  Phase Diagram

In section 3.2, we showed that the perturbative method can successfully predict different modes of distortions of and the intrinsic speed of spontaneous motion. We also predicted the phase boundary separating translationally stable static bumps and translationally unstable static bumps. To predict different phases in the parameter space, however, we need to study the stability of fixed-point solutions given by equations 3.21 to 3.23. A detailed discussion of the stability issue can be found in appendix  C.

By studying the stability of fixed-point moving solutions, we can obtain a phase boundary separating silent and moving phases, as plotted as a dashed line in Figure 4. In Figure 4, there is a rather complete phase diagram for the model of 2D CANNs with STD we are studying here. For small enough s (below the solid curve in Figure 4), since static bumps are stable in amplitude and translation, this parameter region supports only static bumps. This region is called the static phase.

For parameters between the dotted line and the solid line, static bumps are stable in amplitude but unstable against translational distortions. So in this case, in a rather short time window, both static bumps and moving bumps are able to exist, because once a translationally unstable static bumps get perturbed, it becomes a moving bump. Since both static bumps and moving bumps are observed in this parameter region, we call it the bistable phase. If we further increase , however, even static bumps cannot be observed. And since only moving bumps can be observed in the parameter region, we call it the moving phase. When both and are large, no nontrivial solution can be found in our analysis or simulations.

In this study, we did not find the homogeneous firing patterns reported by York and van Rossum (2009). In the 1D case, however, homogeneous firing can be found in the regime (Wang et al., 2014). Since we are focusing on moderate magnitudes of , the possibility of uniform firing behavior in 2D CANNs with STD near is reserved for future studies.

4  CANNs with SFA

4.1  Stationary Solution

By setting in equation 2.5, the CANN model we are considering becomes a network with SFA only. To study the stationary solution in this case, we assume
formula
4.1
formula
4.2
By substituting equations 4.1 and 4.2 into equations 2.4 and 2.6, we found that the fixed-point solution is given by
formula
4.3
formula
4.4
formula
4.5
is a rescaled variable defined by . Similarly, . We can also rescale and in the same way: and . For , the fixed-point solution with the larger and is stable whenever
formula
4.6
which is labeled by curve L in Figure 7 for . The stability issue can be studied by considering the dynamics of distortions of dynamical variables from their fixed-point solutions. Let
formula
4.7
formula
4.8
formula
4.9
By studying the stability of fixed-point solutions, we obtain a parameter region for static bumps with stable amplitudes. Detailed analysis of the stability issue can be found in appendix  D. Regions for and various s are shown in Figure 7. For , the parameter region over the space becomes smaller as increases. At the limit, the range of to stabilize the static solution is given by
formula
4.10
which agrees with equation 4.6. On the other hand, if , there will be no static solution. In Figure 7, we have shown that if becomes larger, the maximum value of to stabilize the static solution will be smaller.
Figure 7:

Maximum able to stabilize the amplitude of stationary solutions as a function of , given that . L: maximum of the existence of stationary solutions.

Figure 7:

Maximum able to stabilize the amplitude of stationary solutions as a function of , given that . L: maximum of the existence of stationary solutions.

4.2  Translational Stability

To study the translational stability, we consider the lowest-order asymmetric distortions added to the stationary solution. The major concern here is the step function in equation 2.7, but the step function does not have a first-order effect on translational distortions. We let
formula
4.11
formula
4.12
As we are studying the stability issue of the static solution, let us adopt the solution in equations 4.3 to 4.5 and put . Then the dynamics of distortions of dynamical variables near fixed-point static solutions becomes
formula
4.13
where the matrix is provided in appendix  D. determines the amplitudal stability of the static solution. Clearly the variable v10 becomes unstable if . This implies that if the dynamics of SFA is slow enough, the static solution will be translationally unstable.

4.3  Moving Solution

Once is satisfied, spontaneous motion of a bump-shaped solution in 2D CANN with SFA becomes possible. As in the STD case, however, lower-order expansions of and are not sufficient to describe the moving solutions. We have to consider higher-order expansions if we want to obtain good predictions on the behavior of moving solutions. Surprisingly, we found that a limited order of expansion of and can give some fairly good predictions on the moving solutions.

In general, we consider
formula
4.14
formula
4.15
Here
formula
4.16
and is the kith-order probabilist’s Hermite polynomial. equations 4.3 and 4.4 become
formula
4.17
formula
4.18
Together with the self-consistent condition, equation 3.28, the fixed-point solution is solvable. Here the definition of is not the same as that in equation 3.24, which is defined by
formula
4.19
can be calculated explicitly. In general, can be obtained by using the recurrence relations given in appendix  E.
If we consider only terms up to and motion along the x0-direction, we can obtain the intrinsic speed of the moving solution (a detailed derivation can be found in appendix  F):
formula
4.20
As the preferred stimulus space is rotationally symmetric, this speed is applicable to motion in any direction. Although this is an approximated solution with relatively few terms, the prediction on the intrinsic speed is fairly good, as shown in Figure 8.
Figure 8:

Speed of spontaneous motions as a function of .

Figure 8:

Speed of spontaneous motions as a function of .

This result suggests that in this case, the terms with are sufficient to give some good predictions on the behavior of bump-shaped solutions. For predictions on measurements of different components of the moving bump, however, we still need to use higher-order perturbation. In Figures 9a to 9c, we have shown that higher-order perturbation is needed to predict if the bump is moving. As in Figures 5b and 5c, in Figures 9b, and 9c, we do not compare , , and directly. Since a moving bump can move in any directions, we compare their projections in polar coordinates. Let us consider and . Then we have
formula
4.21
In Figures 9b and 9c, we compare predictions and simulation results of and for and . is the average change in the width of the bump-shaped profile, while is the magnitude of the anisotropic mode. The graphical illustration of equation 4.21 is omitted because it is similar to the illustrations shown in Figure 6.
Figure 9:

(a) Projection of on . (b) Projection of on rotationally symmetric basis function with . (c) Projection of on anisotropic basis function with . (a) - (c) Parameters: , and .

Figure 9:

(a) Projection of on . (b) Projection of on rotationally symmetric basis function with . (c) Projection of on anisotropic basis function with . (a) - (c) Parameters: , and .

The behavior of moving solutions of CANNs with SFA is similar to that of moving solutions of CANNs with STD. For , its trend is basically the same as in the case with STD shown in Figure 5. The transition happens whenever . Also, anisotropic modes will be available only when the bump is moving, while the width of the profile increases as the strength of SFA increases. The behavior of the average change in width is similar to the behavior of the change in height, . The dependence of anisotropic modes on the strength of SFA, , implies that the widening effect is not uniform when the bump is moving, which is similar to what is seen in a 2D CANN with STD.

4.4  Phase Diagram

A phase diagram similar to Figure 4 for SFA can be obtained in a similar manner used for CANNs with STD. In section 4.3, we showed that perturbative analysis is able to predict dynamical variables of the model. By solving for the moving solutions numerically and testing their stability, we can predict the phase boundary separating parameter regions for moving solutions and trivial solution (silent phase). Using it together with the stability conditions for static bumps given in equations 4.6 and 4.13, we can predict the phase diagram for SFA. We found that the predicted phase diagram matches the phase diagram obtained from computer simulations (see Figure 10). As in Figure 4, there are four phases in the phase diagram: a static phase, a bistable phase, a moving phase, and a silent phase. Their meanings are the same as those of their counterparts in Figure 4. Remarkably, expansions up to can also predict the phase diagram well, which is similar to the prediction of intrinsic speed of moving bumps in CANNs with SFA.

Figure 10:

Phase diagram over the parameter space spanned by . Dotted line: boundary predicted for moving solutions by perturbation. Solid line: boundary predicted for moving solutions by perturbation. Dotted line: boundary predicted for static solutions with stable amplitudes by perturbation. Dot-dashed line: boundary predicted for static solutions with stable amplitudes and translational stability by perturbation. Parameters: same as Figure 9. Symbols: corresponding simulations.

Figure 10:

Phase diagram over the parameter space spanned by . Dotted line: boundary predicted for moving solutions by perturbation. Solid line: boundary predicted for moving solutions by perturbation. Dotted line: boundary predicted for static solutions with stable amplitudes by perturbation. Dot-dashed line: boundary predicted for static solutions with stable amplitudes and translational stability by perturbation. Parameters: same as Figure 9. Symbols: corresponding simulations.

5  Discussion

5.1  Intrinsic Phases and Phase Diagrams

In the parameter spaces of the two models we studied in this article, there are static, bistable, moving, and silent phases. When is large enough, moving bumps can be found in the bistable and moving phases. This behavior can also be seen in the case of SFA. In both the STD and SFA cases, whenever or is not too large (within the bistable phase), static bumps can still exist for a not insignificant period of time. If or is too large, only trivial solutions are stable.

Behaviors of CANNs with STD are very similar to those of CANNs with SFA. This suggests that multiplicative dynamical suppression mechanisms (represented by STD) and subtractive dynamical suppression mechanisms (represented by SFA) can generate similar intrinsic dynamics, especially with regard to spontaneous motion. It also implies that other smooth models having homogeneous couplings and dynamical suppression mechanisms should have a similar phase diagram consisting of static, moving, bistable, and silent phases.

In some studies on CANNs with STD or SFA it is reported that with some parameters, a uniform firing pattern can be found in the network (York & van Rossum, 2009). In this study, uniform firing is not included because the interactive range of the global inhibition in this model is infinite, so uniform firing may be found only in a tight parameter region near . For the 1D STD case, uniform firing can be found in a region and (Wang et al., 2014). Uniform firing should also be found in 2D CANNs with STD or SFA, as should spiral waves and breathing wave fronts. There is a report on spiral waves and breathing wavefronts in 2D CANNs with STD (Kilpatrick & Bressloff, 2010a), but those phenomena are missing from the phase diagram we predicted because spiral waves and breathing wave fronts are not local patterns. Therefore, in our model, they may exist if the magnitude of divisive global inhibition is very small. Richer dynamics of the model for are reserved for future investigations.

5.2  Effectiveness of Perturbative Approach

In this article, we have shown that in our particular model, the perturbative expansion method is applicable to the study of CANNs with short-term synaptic depression (STD) or spike frequency adaptation (SFA). We found in this study that both STD and SFA can drive similar intrinsic dynamics of the local neural activity profile. In Figures 4 and 10, there are four phases: silent, static, bistable, and moving. Using perturbative expansions on the dynamical variables, we can successfully predict the phase diagrams in both cases with STD and SFA.

As expected, with low-order expansions (i.e., when is small), predictions on s of the static solutions work well in cases with STD and SFA. For moving solutions, higher-order expansions are needed to obtain more accurate solutions. In Figures 5a to 5c and Figure 9, it is shown that low-order perturbative expansions, for STD and for SFA, are able to show the general trend of dynamical variables. Especially, the second-order transitions in all the dynamical variables (i.e., discontinuities of slopes) are observed at this level of perturbative expansion. This suggests that low-order perturbative expansions are good enough for studying general phenomena in different phases.

To predict the behavior of the dynamical variables more accurately, we need to use higher-order perturbative expansions. We have shown in Figures 5a to 5c and Figure 9 that higher-order perturbative expansions, for STD and for SFA, can fit measurements from simulations accurately. Higher-order terms do not, however, significantly improve prediction of the speed of spontaneous motion. This suggests that lower-order perturbative modes are most important to the motion of the local neural activity profile.

5.3  Asymmetric Modes and Moving Solutions

Higher-order perturbative modes are essential for predicting the behavior of dynamical variables of moving solutions accurately because the dynamical variables will become asymmetric about the center of mass of if the neural activity profile is moving. In Figure 11, there are two examples of CANNs with different levels of STD. Figures 11a and 11c are with , which corresponds to the static phase. In this case, both the average membrane potential profile and the fraction of available neurotransmitters are rotationally symmetric about the center of mass of . Therefore, in the static phase, low-order perturbation is good enough for making predictions.

Figure 11:

Snapshots of and in the static phase (a, c) and moving phase (b, d). (a, b) ; (c, d) . Dashed lines: contours for . is the center of mass of . The moving direction of the bump is the direction, shows by the arrows in panels b and d. Parameters: , and for panels a and c, and for panels b and d.

Figure 11:

Snapshots of and in the static phase (a, c) and moving phase (b, d). (a, b) ; (c, d) . Dashed lines: contours for . is the center of mass of . The moving direction of the bump is the direction, shows by the arrows in panels b and d. Parameters: , and for panels a and c, and for panels b and d.

Snapshots of and in the moving phase are shown in Figures 11b and 11d. Here the level of STD is . In this case, seems to be rotationally symmetric about its center of mass. Actually, the shape of the gaussian-like average membrane potential profile is slightly squeezed along the moving direction. For , however, the deformation of the shape is more significant. In Figure 11d, the profile of is strongly biased against the moving direction. This highly skewed profile of makes high-order perturbative mode more important, especially for cases involving stability and predictions of dynamical variables.

5.4  Limitation of Perturbative Approach

In this article, we have shown the perturbative approach is able to successfully predict the phase diagrams of CANNs with STD or SFA, the speed of spontaneous motion, and the dynamical variables (e.g., ). Phenomena richer than spontaneous motion of local neural activity, for example, breathing wave fronts on neural networks, have been considered theoretically (Kilpatrick & Bressloff, 2010a). The perturbative approach is not applicable in that case, however, because the basis functions we used here are local, while the traveling wave fronts are globally spreading. For the same reason, this method cannot be used to analyze spiral waves on a 2D field.

In the case of collisions of two moving bump-shaped profiles of neuronal activity, one may find the dynamics hard to analyze by the perturbative method. This is because there are two centers of mass, one for each bump, and this makes the choice of the origin of basis functions confusing. For example, if one chose one of the centers of mass to be the origin of basis functions, the chosen family of basis functions will not be optimal for the other bump and will make the expansion of that bump less efficient.

6  Conclusion

In this article, we studied models with short-term synaptic depression and spike frequency adaptation based on a two-dimensional CANN model with divisive global inhibition. We found that their intrinsic dynamics were similar. First, there are four phases in each scenario: static, moving, bistable, and silent. In the moving phase, the bump-shaped profile moves spontaneously, while in the static phase, the bump-shaped profile cannot move spontaneously. Interestingly, in the bistable phase, the CANN can support both moving profiles and static profiles.

Second, there are clear phase transitions between static solutions and moving solutions. In Figures 5, 8, and 9 there is a clear discontinuity in slope between the two states. This suggests that spontaneous motion can affect the shapes of the bump-shaped profiles.

The stability of steady states, the shapes of the dynamical variables, and the speed of spontaneous motion can be predicted by perturbative analysis, but the perturbative method cannot be used in some scenarios, including spiral waves, breathing wavefronts and collisions of multiple bumps in 2D fields.

Appendix A:  Amplitudal Stability of Bumps on 2D CANNs with STD

To study the stability of fixed-point solutions of equations 3.7 to 3.9, we assume
formula
A.1
formula
A.2
formula
A.3
Linearizing equations 3.7 to 3.9, we obtain
formula
A.4
where
formula
A.5
and
formula
A.6
formula
A.7
formula
A.8
By calculating the eigenvalues of this matrix for a given fixed-point solution, we can test the stability of each solution. In Figure 3, the parabolas are plotted from to a point such that the real part of one of the eigenvalues turns positive. The solid line in Figure 3 maps the parameter region for the existence of static profiles of .

For the stability matrix in equation 3.14, since for static states, amplitudal stability and translational stability are uncoupled. The in the matrix in equation 3.14 is the same as that given in equation A.5.

Appendix B:  Derivations of Equations 3.21 to 3.23

To obtain equations 3.21 to 3.23, we need to deal with differentiations of basis functions and :
formula
B.1
formula
B.2
Combining this result with equations 3.17, and 3.18, we have
formula
B.3
formula
B.4
We now obtain the left-hand sides of equations 2.4 and 2.5. For the right-hand sides, by substituting equations 3.17 and 3.18 into right-hand sides of equations 2.4 and 2.5 and projecting them onto corresponding basis functions, we obtain
formula
B.5
formula
B.6
Here argument t of and is omitted. . Combining expressions of the left-hand and right-hand sides of equations 2.4 and 2.5, we should obtain equations 3.21 and 3.22. Equation 3.23 can be obtained easily by substituting equation 3.17 into equation 2.2 and making use of the orthogonality of s.

Appendix C:  Stability of Moving Bumps on 2D CANNs with STD

To study the stability issue, we define
formula
C.1
formula
C.2
formula
C.3
The stability of fixed-point solutions can be determined by considering the stability matrix given by
formula
C.4
For dynamical variables near their fixed-point solutions, we let
formula
C.5
formula
C.6
formula
C.7
where , , and is the fixed-point solution. Then the dynamics of , , and can be formulated as
formula
C.8
The fixed-point solution is stable only if the maximum of the real parts of eigenvalues of is nonpositive. In general, those eigenvalues can be calculated only by numerical methods. The predicted phase boundary separating the moving and silent phases shown in Figure 4 can be deduced, and the prediction can be verified by simulations.

Appendix D:  Amplitudal Stability of Bumps on 2D CANNs with SFA

For our convenience, we define
formula
D.1
formula
D.2
formula
D.3
Then the stability matrix is given by
formula
D.4
The dynamics of distortions becomes
formula
D.5
If real parts of all eigenvalues of are nonpositive, the static fixed-point solution is stable. For the in equation 4.13, since symmetric terms and asymmetric terms are uncoupled for static bumps, that is the same as the matrix we have shown in this section.

Appendix E:  Recurrence Relation of for CANNs with SFA

We have defined in equation 4.19:
formula
E.1
formula
E.2
where
formula
E.3
To derive the recurrence relations, we first consider
formula
E.4
Here we obtain this identity by using multiple integration by parts. Then we have
formula
E.5
To relate and , we need to use the recursion relations of probabilist’s Hermite polynomials:
formula
E.6
formula
E.7
Then we can derive
formula
E.8
Also,
formula
E.9
formula
E.10
For , we have
formula
E.11
By using equations E.8, E.9, and E.10 and , we can generate any needed for numerical computations.

Appendix F:  Intrinsic Speed of Moving Bumps on 2D CANNs with SFA

Here we consider terms up to . For simplicity we assume . Then,
formula
F.1
formula
F.2
Note that since , terms that are asymmetric along the -direction are dropped. The derived speed should be applicable to any direction because the model is homogeneous.
By projection and orthogonality of the basis functions, equation 2.4 and 2.6 give
formula
F.3
formula
F.4
formula
F.5
formula
F.6
formula
F.7
formula
F.8
formula
F.9
formula
F.10
From these equations we can deduce the ratio between and :
formula
F.11
Combining equations F.4, F.7, and F.10, we have
formula
F.12
formula
F.13
formula
F.14
As shown in Figure 8, this formula can be verified by simulation.

Appendix G:  Stability of Moving Bumps on 2D CANNs with SFA

Using a method similar to that used in appendix  C, we define
formula
G.1
formula
G.2
formula
G.3
Then the matrix concerning the stability issue is
formula
G.4
Similarly, the dynamics of distortions is given by
formula
G.5
Here , , and are defined by
formula
G.6
formula
G.7
formula
G.8
where , and are the fixed-point moving solution to the system. By calculating eigenvalues of , the stability of the fixed-point solution can be determined.

Acknowledgments

This study is partially supported by the Research Grants Council of Hong Kong (grants 604512, 605813 and N_HKUST 606/12).

References

Amari
,
S.
(
1977
).
Dynamics of pattern formation in lateral-inhibition type neural fields
.
Biological Cybernetics
,
27
,
77
87
.
Benda
,
J.
, &
Herz
,
A.V.M.
(
2003
).
A universal model for spike-frequency adaptation
.
Neural Computation
,
15
,
2523
2564
.
Ben-Yishai
,
B.
,
Bar-Or
,
R. L.
, &
Sompolinsky
,
H.
(
1995
).
Theory of orientation tuning in visual cortex
.
Proceedings of the National Academy of Sciences of the United States of America
,
92
(
9
),
3844
3848
.
Blair
,
H. T.
, &
Sharp
,
P. E.
(
1995
).
Anticipatory head direction signals in anterior thalamus: Evidence for a thalamocortical circuit that integrates angular head motion to compute head direction
.
Journal of Neuroscience
,
15
,
6260
6270
.
Brown
,
D. A.
, &
Adams
,
P. R.
(
1980
).
Muscarinic supression of a novel voltage-sensitive K+ current in a vertebrate neuron
.
Nature
,
183
,
673
676
.
Fleidervish
,
I. A.
,
Friedman
,
A.
, &
Gutnick
,
M. J.
(
1996
).
Slow inactivation of Na+ current and slow cumulative spike adaptation in mouse and guineapig neocortical neurones in slices
.
Journal of Physiology
,
493
,
83
97
.
Funahashi
,
S.
,
Bruce
,
C. J.
, &
Goldman-Rakic
,
P. S.
(
1993
).
Dorsolateral prefrontal lesions and oculomotor delayed-response performance: Evidence for mnemonic “scotomas.”
Journal of Neuroscience
,
13
(
4
),
1479
1497
.
Fung
,
C. C. A.
,
Wang
,
H.
,
Lam
,
K.
,
Wong
,
K. Y. M.
, &
Wu
,
S.
(
2013
).
Resolution enhancement in neural networks with dynamical synapses
.
Frontiers in Computational Neuroscience
,
7
,
73
.
Fung
,
C. C. A.
,
Wong
,
K. Y. M.
,
Wang
,
H.
, &
Wu
,
S.
(
2012
).
Dynamical synapses enhance neural information processing: Gracefulness, accuracy and Mobility
.
Neural Computation
,
24
,
1147
1185
.
Fung
,
C. C. A.
,
Wong
,
K. Y. M.
, &
Wu
,
S.
(
2010
).
A moving bump in a continuous manifold: A comprehensive study of the tracking dynamics of continuous attractor neural networks.
Neural Computation
,
22
,
752
792
.
Fung
,
C. C. A.
,
Wong
,
K. Y. M.
, &
Wu
,
S.
(
2012
).
Delay compensation with dynamical synapses
. In
P.
Bartlett
,
F. C. N.
Pereira
,
C. J. C.
Burges
,
L.
Bottou
, &
K. Q.
Weinberger
. (Eds),
Advances in Neural Information Processing Systems
,
25
(pp.
1097
1105
).
Red Hook, NY
:
curran
.
Kilpatrick
,
Z. P.
(
2013
).
Short term synaptic depression improves information transfer in perceptual multistability
.
Frontiers in Computational Neuroscience
,
7
,
85
.
Kilpatrick
Z. P.
, &
Bressloff
P. C.
, (
2010a
).
Spatially structured oscillations in a two-dimensional excitatory neuronal network with synaptic depression
.
Journal of Computational Neuroscience
,
28
(
2
),
193
209
.
Kilpatrick
,
Z. P.
, &
Bressloff
,
P. C.
(
2010b
).
Effects of synaptic depression and adaptation on spatiotemporal dynamics of an excitatory neuronal network
.
Physica D
,
239
,
547
560
.
Loebel
,
A.
, &
Tsodyks
,
M.
(
2002
).
Computation by ensemble synchronization in recurrent networks with synaptic depression
.
Journal of Computational Neuroscience
,
13
,
111
124
.
Lu
,
Y.
,
Sato
,
Y.
, &
Amari
,
S.-I.
(
2011
).
Traveling bumps and their collisions in a two-dimensional neural field
.
Neural Computation
,
23
(
5
)
1248
1260
.
Madison
,
D. V.
, &
Nicoll
,
R. A.
(
1984
).
Control of the repetitive discharge of rat CA1 pyramidal neurones in vitro
.
Journal of Physiology
,
354
,
319
331
.
Maunsell
,
J. H.
, &
van Essen
,
D. C.
(
1983
).
Functional properties of neurons in middle temporal visual area of the macaque monkey. I. Selectivity for stimulus direction, speed, and orientation
.
Journal of Neurophysiology
,
49
,
1127
1147
.
O’Keefe
,
J.
(
1976
).
Place units in the hippocampus of the freely moving rat
.
Experimental Neurology
,
51
(
1
),
78
109
.
O’Keefe
,
J.
, &
Burgess
,
N.
(
1996
).
Geometric determinants of the place fields of hippocampal neurons
.
Nature
,
381
(
6581
),
425
428
.
O’Keefe
,
J.
, &
Dostrovsky
,
J.
(
1971
).
The hippocampus as a spatial map: Preliminary evidence from unit activity in the freely-moving rat
.
Brain Research
,
34
(
1
),
171
175
.
Pouget
,
A.
,
Zhang
,
K.
,
Deneve
,
S.
, &
Latham
,
P.E.
(
1998
).
Statistically efficient estimation using population coding
.
Neural Computation
,
10
,
373
401
.
Samsonovich
,
A.
, &
McNaughton
,
B.
(
1997
).
Path integration and cognitive mapping in a continuous attractor neural network model
.
Journal of Neuroscience
,
17
,
5900
5920
.
Taube
,
J. S.
, &
Muller
,
R. U.
(
1998
).
Comparisons of head direction cell activity in the postsubiculum and anterior thalamus of freely moving rats
.
Hippocampus
,
8
(
2
),
87
108
.
Taube
,
J. S.
,
Muller
,
R. U.
, &
Ranck
,
J. B. Jr
, (
1990
).
Head-direction cells recorded from the postsubiculum in freely moving rats. II. Effects of environmental manipulations
.
Journal of Neuroscience
,
10
(
2
),
436
447
.
Treue
,
S.
,
Hol
,
K.
, &
Rauber
,
H.-J.
(
2000
).
Seeing multiple directions of motion—physiology and psychophysics
.
Nature Neuroscience
,
3
(
3
),
270
276
.
Tsodyks
,
M. V.
, &
Markram
,
H.
, (
1997
).
The neural code between neocortical pyramidal neurons depends on neurotransmitter release probability
.
Proceedings of the National Academy of Sciences of the United States of America
,
94
(
2
),
719
723
.
Tsodyks
,
M. S.
,
Pawelzik
,
K.
, &
Markram
,
H.
(
1998
).
Neural networks with dynamic synapses
.
Neural Computation
,
10
,
821
835
.
Wang
,
H.
,
Lam
,
K.
,
Fung
,
C. C. A.
,
Wong
,
K. Y. M.
, &
Wu
,
S.
A rich spectrum of neural field dynamics in the presence of short-term synaptic depression
. (
2014
).
Unpublished manuscript
.
Wimmer
,
K.
,
Nykamp
,
D. Q.
,
Constantinidis
,
C.
, &
Compte
,
A.
(
2014
).
Bump attractor dynamics in prefrontal cortex explains behavioral precision in spatial working memory
.
Nature
,
17
(
3
),
431
439
.
Wu
,
S.
,
Hamaguchi
,
K.
, &
Amari
,
S.
(
2008
).
Dynamics and computation of continuous attractors
.
Neural Computation
,
20
(
4
),
994
1025
.
York
,
L. C.
, &
van Rossum
,
M. C. W.
(
2009
).
Recurrent networks with short term synaptic depression
.
Journal of Computational Neuroscience
,
27
,
607
620
.
Zhang
,
K.-C.
(
1996
).
Representation of spatial orientation by the intrinsic dynamics of the head-direction cell ensemble: A theory
.
Journal of Neuroscience
,
16
,
2112
2126
.