Emergence of local synchronization in neuronal networks with adaptive couplings

Local synchronization, both prolonged and transient, of oscillatory neuronal behavior in cortical networks plays a fundamental role in many aspects of perception and cognition. Here we study networks of Hindmarsh-Rose neurons with a new type of adaptive coupling, and show that these networks naturally produce both permanent and transient synchronization of local clusters of neurons. These deterministic systems exhibit complex dynamics with 1/fη power spectra, which appears to be a consequence of a novel form of self-organized criticality.


Introduction
Synchronization is ubiquitous in nature [1]. Examples of this phenomenon are: synchronization in arrays of laser [2,3] and microwave oscillators [4], and in superconducting Josephson junctions [5][6][7], the synchronized beating of wings in a flock of birds [8], the chirping of crickets in unison [9], and the synchronized clapping of a crowd [10]. Synchronization is of great importance in neuroscience [11]. There is strong experimental evidence that synchronization of neuronal oscillatory activity is a central mechanism in a variety of cognitive and perceptual processes in the brain, including: the functioning of working memory, the routing of signals across cortical networks, sensory motor awareness, and perceptual grouping [12][13][14]. Synchronization has been observed experimentally between areas of the visual cortex and parietal cortex in awake cats during a visual-motor integration task [14]. The occurrence of synchronization in a cognitive task in humans has also been demonstrated experimentally [15]. Studies have also found that stored memory patterns depend on spontaneously occurring synchrony in neuronal networks [16]. Moreover, abnormally synchronized neural activity has been implicated in a number of clinical disorders, including: schizophrenia, epileptic seizures, and Parkinson's disease [17]. Understanding these phenomena is of fundamental importance in neuroscience, and this motivates the present study on synchronization in neuronal networks.
Previous theoretical studies of neural synchronization have largely focused on the issue of global synchronization, i.e., studying the conditions under which all neurons in the network a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 behave in unison [18]. However, it is known from experimental neuroscience that global synchronization does not occur in normal brains. For example, it has been found experimentally that the visual recognition system of the cat and macaque monkey show local synchronization within brain regions [19]. It has also be found that local synchronization, both prolonged and transient, plays a crucial role in many aspects of memory processes [20]. The purpose of this paper is to study how local synchronization occurs in adaptively coupled networks of neurons. The effect on synchronization of variable couplings in networks of simple oscillators has previously been explored in the literature [21][22][23]. Here we develop this research theme by introducing a new type of adaptive coupling between biologically realistic neurons, and show that it naturally results in the emergence of locally synchronized group of neurons.

Model
Networks of biologically-inspired neurons provide a natural theoretical framework for studying neural synchronization [24]. Here we model individual neurons using the Hindmarsh-Rose equations [25,26]: Eqs (1) and (2) model the fast dynamics associated with a neuron, whereas Eq (3) describes the slow dynamics. The variable x(t) describes the membrane potential; y(t), the fast current, models the transport of potassium and sodium ions across the fast ion channel; and z(t), the slow current, models the transport of other ions through the slow channels. The parameters are as follows: I specifies the membrane input current to the neuron; b allows for switching between spiking and bursting behavior and also controls the frequency of spiking; r controls the rate of change of the slow variable z(t) and, if spiking is present, it determines the spiking frequency, while if bursting is present, it governs the number of spikes per burst; and x 0 determines the resting potential for the neuron. We fix the parameters b, r and x 0 to the standard values b = 3.0, r = 0.006 and x 0 = −1.6 [27]. Furthermore, we set I = 2.8, which corresponds to bursting neurons [27,28], since this behavior seems to be an essential component of information processing in the brain [29,30]. We consider a network Γ of N identical Hindmarsh-Rose neurons coupled through electrical synapses (i.e., via gap junctions, which are known to play a fundamental role in a wide range of neural systems, including the mammalian brain [24,31,32], and to be particularly important in maintaining neural synchronization [33][34][35][36][37][38]). The coupled system of neurons is described by the following equations: Here x i (t), y i (t) and z i (t) denote the x(t), y(t) and z(t) variables of the ith neuron, respectively, k ij (t) (where k ij (t) ! 0, for all i, j and t) represents the coupling between neurons i and j at time t, and A = (A ij ) is the adjacency matrix of the network Γ (i.e., A ij = 1 if neurons i and j are connected by an edge in Γ, and A ij = 0 otherwise). For k ij (t) = k, for all i, j, and Γ the complete network, Eqs (4)- (6) describes N Hindmarsh-Rose neurons, each of which is coupled to all others with common, fixed, coupling strength k [27,[39][40][41].
In this paper we investigate the effects of connection plasticity on the dynamics of the neuronal network defined by Eqs (4)- (6). Previous approaches to modeling connection plasticity between neurons have focused on properties, such as long-term potentiation (LTP) [42,43] and spike-timing-dependent plasticity (STDP) [44], which are particular to chemical synapses. Here, in contrast, we consider connection plasticity in networks of neurons coupled through electrical synapses. We allow the strength of the coupling k ij (t) between adjacent neurons i and j to depend on the states (i.e., the membrane potentials) of the neurons. This assumption is natural since it is known that the electrical coupling between neurons via gap junctions may be modulated by neural activity [45][46][47][48][49]. The form of the dependence that we assume is motivated by Hebb's law that synaptic connections are strengthened between neurons that are active simultaneously [50]. Thus, we assume that the coupling k ij (t) between i and j is adaptive, and will increase in strength if i and j are in (approximately) the same state and will decrease if i and j are in dissimilar states. Precisely, we take the time evolution of the coupling k ij (t) between adjacent neurons i and j to be governed by the equation where α, β and γ are positive parameters. Thus, with k ij determined by Eq (7), the last term in Eq (4) defines an adaptive diffusive coupling in the network of Hindmarsh-Rose neurons. We note some features of the dynamics of Eq (7). (A) If two neurons i and j are unsynchronized (i.e., x i is very different from x j ) then Eq (7) reduces to _ k ij % À gk ij ðk ij þ 1Þ < 0, and the strength of the coupling between i and j decreases, and consequently i and j remain unsynchronized. (B) Note that,k ij ¼ 0 is always an equilibrium of Eq (7). For 0 < k ij ( 1, Eq (7) becomes _ k ij % k ij ðae À bðx i À x j Þ 2 À gÞ. For such values of k ij , if neurons i and j are synchronized (i.e., x i = x j ) then Eq (7) further reduces to _ k ij % k ij ða À gÞ; thus we require that α > γ to ensure that _ k ij > 0 (and thusk ij is unstable) for synchronized neurons, ensuring that the coupling between i and j increases, maintaining their synchronization. (C) We also note that, for synchronized neurons i and j (i.e., x i = x j ), Eq (7) has a second equilibrium, k Ã Eqs (4)-(7) model the adaptively coupled neuronal network that we study here. We observe that the form of the variation in coupling strength between neurons that we have postulated has the effect of increasing the synchronization between neurons that are already in similar states and of further decreasing the synchronization between neurons in dissimilar states. It would be very interesting to consider the full stability analysis of Eqs (4)-(7), however, this appears to be a highly non-trivial enterprise which is beyond the scope of the current work. While it would be most illuminating to understand the fixed points of this system and the possible bifurcation that can occur as the parameters α, β and γ vary, achieving such an understanding would seem to be quite ambitious as even the fixed-point equations constitute a high dimensional system of coupled transcendental equations. As a possibly more manageable first step for future work in this direction it may be worth considering the stability and bifurcation analysis of a single pair of adaptive coupled neurons satisfying Eqs (4)-(7).

Results
We now present results obtained by numerically integrating this neuronal network. We consider the time evolution of the network for initial coupling strengths k ij chosen uniformly randomly from (0, 1). Here we fix α = 1, γ = 0.5 (so k ij 2 [0, 1], for all i, j) and treat β as a control parameter. The results described here are for β = 12, although we obtain qualitatively similar results for the range of values β 2 [10, 12.5].
We first consider the case in which Γ is the complete network (i.e., A ij = 1, for all i 6 ¼ j). Thus, every neuron is potentially coupled to every other neuron. Fig 1 shows the time evolution of the coupling strengths for a network of N = 100 neurons. We observe that over time the coupling k ij (t) between any pair of neurons i, j evolves into one of the following three classes: undergoes sustained oscillations between 0 and 1. This coupling dynamics induces a corresponding synchronization dynamics for the neurons: (1) pairs of neurons i, j for which k ij (t) ! 1 as t ! 1 become completely synchronized (i.e., |x i (t) − x j (t)| ! 0 as t ! 1); (2) i, j for which k ij (t) ! 0 as t ! 1 become completely unsynchronized; and (3) i, j for which k ij (t) continuously varies between 0 and 1 become synchronized when k ij (t) is high (greater than approximately 0.8) and then lose synchronization as k ij (t) falls. The variations in the couplings result in the emergence over time of clusters of both permanently and transiently synchronized neurons. Within each cluster of permanently synchronized neurons, the couplings between neurons evolve dynamically to 1. For transiently synchronized neurons, the couplings between such neurons continuously change between 0 and 1 under the dynamics, resulting in the emergence of clusters of neurons that are synchronized for a period of time when the couplings are high and then become desynchronized as the couplings decrease. In addition, the couplings between some neurons evolve dynamically to 0, resulting in these neurons being permanently unsynchronized. The nature of the local synchronization process may be elucidated by considering the longterm dynamics of the couplings between neurons in the network. We define the time-average coupling hk ij i between neurons i and j over the interval [T 0 , T] to be hk ij i ¼ 1 The perpetual variation in the total coupling strength K(t), apparent in Fig 5(a), reflects the dynamic in which groups of neurons transiently synchronize when the couplings between them are strong, subsequently desynchronize as the coupling strengths fall, and later resynchronize as the coupling strengths rise again. We note, from Fig 5(b), that the power P(f ) associated with frequency f satisfies, to a good approximation, the power-law P(f ) ∝ 1/f η , where η = 2.628 ± 0.002. The fact that the power spectrum satisfies a 1/f η power law implies that variations in K(t) occur on all time scales. We find the same 1/f η power-law for all values of β 2 [10, 12.5], the range of βvalues for which local synchronization occurs.
We have also investigated the time series and power spectrum of the total membrane potential XðtÞ ¼ P N i¼1 x i ðtÞ, which is shown in Fig 5(c) and 5(d), respectively. The power spectrum of X(t) has greater deviations from a power law than that of K(t), however, X(t) also exhibits variations on a wide range of time scales. The dynamics of the emergence of permanently and transiently synchronized groups of neurons in the adaptively coupled network is considerably elucidated by considering an appropriate order parameter. Here we define the order parameter wðtÞ ¼ 1 N j P N j¼1 e 2pix j ðtÞ j, wherê x j ðtÞ ¼  Fig 6(a). The continual variation in χ(t), apparent in Fig 6(a), vividly illustrates the dynamic of repeated synchronization and desynchronization of clusters of transiently synchronized neurons in the adaptively coupled network. The power spectrum of χ(t), shown in Fig 6(b), shows greater deviations from a power law than that for K(t), but still exhibits variations on a large range of time scales, indicative of the complex dynamics of χ(t). The distributions of the complex numbers fe 2pix j ðtÞ g    The system subsequently evolves in time to less synchronized states, followed by more synchronized states, and perpetually repeats this pattern of variation.
We have also studied the coupled neuronal network introduced here for a large variety of model and empirical network topologies Γ. Here we shall only briefly describe certain of these results deferring a more detailed account of these cases to a subsequent publication.
We find that the coupled neuronal network defined on a wide diversity of network topologies Γ, including small-world [51] and scale-free [52] topologies, exhibit exactly analogous behavior to that found above for the complete network topology. We again see that the time evolution of the coupling k ij (t) between every pair of adjacent neurons i and j in Γ falls into one of three classes: (1) k ij (t) ! 1 as t ! 1, resulting in the neurons i and j becoming completely synchronized; (2) k ij (t) ! 0 as t ! 1, results in i and j becoming completely unsynchronized; or (3) k ij (t) maintains continued oscillations between 0 and 1, results in i and j being transiently synchronized. The variations in the coupling between adjacent neurons results in the formation over time of clusters of both the permanently and transiently synchronized neurons in Γ.
The transiently synchronized neurons result in a complex dynamics in which the power spectrum of the total coupling K(t) satisfies to an excellent approximation a power-law P(f ) ∝ 1/f η . The power spectra for networks Γ of N = 100 neurons with small-world and scale-free topologies  are shown in Fig 7. For small-world topologies with rewiring probability 0.3 and mean degree 10 and 20, the power spectra, shown in Fig 7(a), satisfy the power-laws P(f ) ∝ 1/f η , with η = 2.013 ± 0.002 and η = 2.070 ± 0.002, respectively. The power spectra for scale-free topologies with mean degree 10 and 20, shown in Fig 7(c), satisfy the power-laws P(f ) ∝ 1/f η , with η = 2.094 ± 0.002 and η = 2.028 ± 0.002, respectively. For both small-world and scale-free topologies the power spectra of the total membrane potential X(t), shown in Fig 7(b) and 7(d), display greater deviations from power-law behavior than K(t), but nevertheless exhibit variations on a  It is an intriguing result, worthy of further study, that the power spectrum of K(t) conforms to a power-law much more closely on complex networks than on the complete network. It is also potentially significant that the power-law exponents are almost equal for quite different complex network topologies, suggesting a possible universality underlying the dynamics of these coupled complex networks of neurons.
complex numbers that define the order parameter on the unit circle in the complex plane at different times. At time t = 13 (shown in (c)) the system is largely unsynchronized, and subsequently evolves in time to successively more synchronized states at t = 3825 (shown in (d)), t = 3910 (shown in (e)), and t = 4000 (shown in (f)). https://doi.org/10.1371/journal.pone.0178975.g006

Discussion
Most theoretical studies of synchronization in neuronal networks have focused on understanding the conditions under which global synchronization of neurons occur [18]. However, it is clear empirically that many important aspects of neural behavior depend on local, and often transient, synchronization of groups of neurons in the brain [19]. Here we have studied adaptively coupled neuronal networks, in which the coupling between two neurons is determined dynamically by the states of the neurons, and which has the key feature that the coupling between pairs of adjacent neurons that are in a similar state is strengthened, while that between pairs of adjacent neurons in dissimilar states is weakened. The dynamics of these networks result in the spontaneous emergence of locally synchronized groups of neurons, some of which are permanently synchronized while other show transient synchronization. This pattern of synchronization is found for a significant parameter range in the model. The existence of transient local synchronization (in addition to permanent local synchronization) results in the neuronal couplings having a complex dynamics, which is characterized by the power spectrum P(f ) of the sum of the couplings obeying a power-law P(f ) ∝ 1/f η . For the complete network topology, we find the same power law for all values of the parameter β that yields local synchronization. We have also studied the model for other values of α and γ, and find that for a wide variety of these parameters there exists a significant range of β values for which local synchronization occurs and for which the total coupling exhibits a power law spectrum. The exponent of the power law varies with α and γ, but for fixed α and γ is independent of β. For complex network topologies, including small-world and scale-free topologies, we find no significant variation in the power-law exponents for a wide range of the parameters α, β and γ, which suggests an intriguing universality underlying the dynamics of these adaptively coupled neuronal networks. The 1/f η power-law in this deterministic model does not depend on any fine tuning of parameters and appears to be a consequence of an interesting new type of selforganized criticality.
The mechanism responsible for driving continuous oscillations in the couplings between pairs of transiently synchronized neurons appears to be a type of dynamical frustration. In this frustration process distinct groups of completely synchronized neurons form in the network, and since the oscillations between different synchronized groups are typically out of phase this results in pairs of neurons spanning different groups having a constantly changing difference in their membrane potentials, which in turn leads to a continuously changing coupling strength between such pairs of neurons. The operation of this process between all pairs of neurons spanning different completely synchronized groups results in many pairs of transiently synchronized neurons with couplings varying on a wide range of time scales. Thus, this deterministic dynamical system apparently self-organizes into a state in which variations in the couplings between transiently synchronized neurons occur on all time scales and, therefore, the emergence of a power-law spectrum in this deterministic model, which does not depend on any fine tuning of parameters, seems to be a consequence of a novel form of self-organized criticality [53]. Concepts of self-organized criticalty, while having been applied to a variety of problems in neuroscience [54][55][56][57], including certain aspects of synchronization [58], appears not to have been previously connected to synchronization in networks of neurons adaptively coupled through electric synapses.
The frustration process that we have suggested here to underlie the emergence of a powerlaw spectrum in an adaptively coupled network of Hindmarsh-Rose neurons is novel in that it depends upon frustration occurring between dynamically formed clusters of neurons, where the neurons within a given cluster are completely synchronized, but neurons lying in different clusters are unsynchronized. We note here that a different, and more conventional type of geometric frustration process plays an important role in other types of neuronal networks [59][60][61][62]. Obtaining a deeper understanding of the dynamical frustration process we have suggested for adaptively coupled neuronal networks, and elucidating the connection to power-law spectra seems to be an important task for future research.
In this paper we have considered adaptively coupled networks of Hindmarsh-Rose neurons, where the membrane input current I is set to produce bursting behavior for the neurons, as this is known to be a key aspect of information processing in the brain [30,63]. We have also explored the behavior of our adaptively coupled networks when the membrane input current I is set to give spiking behavior for the neurons. We find that the outcome for spiking neurons is very similar to the case of bursting neurons, with the spontaneous emergence of groups of completely synchronized neurons and also of clusters of transiently synchronized neurons, the latter of which result in a complex dynamical process analogous to that found for bursting neurons.
We should also like to mention that the type of adaptively coupled neuronal networks that we have considered here may be relevant to the problem of pattern selection. It is known that complex pattern formation, including spiral waves, can occur in some neuronal networks and it would be an interesting topic for future research to investigate whether such pattern formation occurs in the adaptively coupled Hindmarsh-Rose networks that we have considered here [64][65][66][67].
Finally, we remark that while we have focused here on neuronal synchronization it is possible that the type of adaptive coupling we have postulated may have broader applications to the study of local synchronization in other networked systems.

Author Contributions
Conceptualization: TK PI BS.