Emergent central pattern generator behavior in chemical coupled two-compartment models with time delay

Physica A 491 (2018) 177–187

Contents lists available at ScienceDirect

Physica A journal homepage: www.elsevier.com/locate/physa

Emergent central pattern generator behavior in chemical coupled two-compartment models with time delay Shanshan Li, Guoshan Zhang, Jiang Wang, Yingyuan Chen, Bin Deng * School of Electrical and Information Engineering, Tianjin University, Tianjin, 300072, PR China

highlights • The modified two-compartment PR model is proved to meet the requirement of CPG and can be used to develop the CPG networks. • The simplest form of CPG is constructed by two inhibitory chemical coupled PR neurons with delay time. • Emergent behaviors of CPG affected by ambient noise, sensory feedback signals, morphology features as well as the coupled delay time are investigated.

article

info

Article history: Received 6 April 2017 Received in revised form 26 June 2017 Available online 28 September 2017 Keywords: Central pattern generator Two-compartment model Ambient noise Sensory feedback signals Morphology Time delay

a b s t r a c t This paper proposes that modified two-compartment Pinsky–Rinzel (PR) neural model can be used to develop the simple form of central pattern generator (CPG). The CPG is called as ‘half-central oscillator’, which constructed by two inhibitory chemical coupled PR neurons with time delay. Some key properties of PR neural model related to CPG are studied and proved to meet the requirements of CPG. Using the simple CPG network, we first study the relationship between rhythmical output and key factors, including ambient noise, sensory feedback signals, morphological character of single neuron as well as the coupling delay time. We demonstrate that, appropriate intensity noise can enhance synchronization between two coupled neurons. Different output rhythm of CPG network can be entrained by sensory feedback signals. We also show that the morphology of single neuron has strong effect on the output rhythm. The phase synchronization indexes decrease with the increase of morphology parameter’s difference. Through adjusting coupled delay time, we can get absolutely phase synchronization and antiphase state of CPG. Those results of simulation show the feasibility of PR neural model as a valid CPG as well as the emergent behaviors of the particularly CPG. © 2017 Elsevier B.V. All rights reserved.

1. Introduction Central pattern generators are neuronal circuits found in both invertebrate and vertebrate animals that can produce rhythmic patterns of neural activity such as walking, breathing, flying and swimming [1–4] in the absence of sensory or descending inputs that carry specific timing information. CPGs play an important role in the formation of repeated oscillatory behaviors and are considered central to their basic survival across much of the animal species. Traditionally, reciprocal synaptic inhibition between two neuronal populations (or two groups of neuronal populations, or two individual neurons [5]) is considered as the standard form of generating CPG behavior in both biological and

*

Corresponding author. E-mail address: [email protected] (B. Deng).

https://doi.org/10.1016/j.physa.2017.08.121 0378-4371/© 2017 Elsevier B.V. All rights reserved.

178

S. Li et al. / Physica A 491 (2018) 177–187

computational systems. The simplest form of CPG is called as a ‘half-center oscillator’, in which two neurons reciprocally inhibit each other. This network was proposed firstly by Brown [6] to explain alternation of extension and flexion phases in cat locomotion, and have subsequently been studied extensively both theoretically [7,8] and experimentally [2]. The half-center CPGs have found in biological models of lamprey [9], in locomotion models such as stick insect locomotion [10]. In Ijspeert’s review [4], several features of CPGs were proposed: (1) typically, CPG models have a few control parameters, which can be used to modulate the rhythmical outputs; (2) CPG models display limit cycles to produce stable rhythmical output and have resistance to noise; (3) CPGs are suitable to integrate sensory feedback signals which provide an opportunity to realize entrainment between CPG and the environment [41]. According to different targets of the research, different neuronal models, such as the Hodgkin–Huxley (HH) model [11], the FitzHugh–Nagumo (FN) model [12,13,42] and Hindmarsh–Rose (HR) neuronal model [14,15], can be used to develop CPG model. Those models are capable of mimicking almost all the behaviors exhibited by real biological neurons, such as spiking, bursting and so on [16]. However, most of neurons forming CPGs are single compartment models. The single compartment model is the simplification of real neurons and do not have physical structures including soma, dendrite, axon et al. In contrast, multi-compartment models have rich firing patterns, and could reflect the spatial location of neural structure. What is more important, the multi-compartment model is built based on physiological data, and is used to mimic physiological experiment and reveal the physiological meaning of simulation results. PR model is a kind of spiking neuronal models with two-compartments—soma and dendrite [17]. The dynamical equations of PR model consist of eight equations, a dozen of functions and tens of parameters, describing membrane potential, activation of sodium current and potassium current, and inactivation of the sodium current and so on. What is more, there are many theoretical studies on the properties of PR model, such as bifurcation [18], analysis of bursting mechanism [19–21] and so on. However, PR model is rarely used to simulate CPG in comparison with other neuronal models, and its application in practice is even less. There is also little systematic analysis about the properties of PR model in the view of CPG network. In this paper, we firstly aim to explore the properties of PR model on the aspect of CPG, and construct a half-center CPG network based on two reciprocally inhibitory PR neurons with chemical synapses. Considering the causal relationship between firing properties and dendrite structures [22], we also investigate the effects of morphological feature (the size ratio of soma and dendrite) of single neural model on the output rhythm of CPG network. The coupling delay time is another important factor influencing the output of CPG [23,24], and it is also considered. The rest of the paper is organized as follows. First, we investigate the spiking properties of modified PR model in the view of CPGs and prove that it can be used to form a CPG network. Then, we construct a half-center CPG network based on modified PR model. Last, we analyze the output rhythm of half-center CPG network affected by ambient noise, sensory feedback signals, and different morphology features of PR neuron as well as the coupled delay time. 2. Model and method 2.1. Mechanism of single neuron The PR model characterizes a typical pyramidal cell as comprising a somatic and a dendritic compartment. Each PR neuron is characterized by eight time-dependent variables, and the equations of membrane potentials are shown in Eq. (1). There are three currents in somatic compartment: transient sodium INa , delayed rectifier potassium IKDR , and leak current. The dendritic compartment contains persistent calcium ICa , calcium activated potassium IKCa , after-hyperpolarisation potassium current IKAHP and leak current. All currents are conductance-based, using the Hodgkin–Huxley formalism [11] of activation and inactivation gates dependent on voltage or intracellular calcium that drive the current. The two compartments are coupled by a coupling current ISD = −IDS = gc (Vd − Vs ). The size of dendrite compartment as a proportion of the entire neuron was given by 1 − p and that of soma compartment as p. Currents and conductances are expressed as densities with units of µA/cm2 and mS/cm2 , respectively. Capacitance Cm is in unit µF/cm2 ; the time unit is ms [17]. Cm Cm

dVs dt dVd dt

= −ILeak − INa − IKDR −

IDS p

+

Is p

= −ILeak − ICa − IKCa − IKAHP −

ISD 1−p

+

Id 1−p

(1)

.

The activation and inactivation gates of ionic channel evolve as a function of their steady state x∞ and time constant τx , where V represents the membrane potential Vs or Vd , or the intracellular calcium Ca. dx dt

=

x∞ (V ) − x

x∞ (V ) =

τx

;

dCa dt

= −0.13Ca − 0.075Ca

αx (V ) 1 ; τx (V ) = αx (V ) + βx (V ) αx (V ) + βx (V )

(2)

In order to make a detailed mathematical analysis of the PR model from the view of dynamical system, Laura A. Atherton and Krasimira Tsaneva-Atanasova modified the original PR model, and approximated the discontinuous functions by fitting continuous functions directly to the steady state and time activation curves [25]. All activation/inactivation function and parameters of modified PR model can be found in Ref [18]. In the following sections, we will discuss the properties of modified PR model in respect of CPG.

S. Li et al. / Physica A 491 (2018) 177–187

179

Fig. 1. Frequency current (f − I) relation under various electrical and geometric coupling strengths. Frequency is calculated by the number of times Vs increases past 70 mV. The solid triangle represents the spike frequency, the solid circle represents burst frequency and the asterisk represents mixed burst-spike frequency. For standard coupling parameter values (gc = 2.1 mS/cm2 , p = 0.5), there is a transition from resting to periodic very low frequency bursting to somatic spiking to steady depolarization as Is increase. The transition between bursting and spiking involves mixed firing patterns. Increase or decrease the electrical coupling parameter gc , the PR model display spiking. Steady dendritic stimulation can produce periodic bursting and stronger stimulation give rise to mixed burst-spike pattern. Whether stimulating soma or dendrite, increase the soma size, the PR neuron display spiking pattern and decrease the soma size, the PR neuron is in VLF bursting pattern.

2.1.1. Rhythmic behavior of modified PR model As mentioned in previous investigations, the PR model describes a wide range of neural dynamic properties with the variation of Is , Id , gc and p, including stable periodic solutions, undulating spiking pattern, regular spiking activity and very low frequency (VLF) bursting and so on [18]. The salient stimulus–response properties of PR model include the following, just like shown in Fig. 1. For steady somatic injection (Is ), weak Is leads to very low frequency bursting (frequency less than 8 Hz) and the frequency increases modestly with Is . For intermediate Is , there is aperiodic behavior and for strong Is periodic spiking occurs. Steady dendritic stimulation can produce periodic bursting of significant higher frequency (8–20 Hz) than can steady somatic input (<8 Hz). The coupling conductance gc between two compartments is a significant electrotonic parameter. The bursting can occur only for a range of gc intermediate between the extremes of small gc (decoupled compartments) and large gc (an isopotential cell for which a burst is replaced by a composite sodium–calcium spike). When the soma and dendrite are tightly coupled electrically, the PR model reduces to a single compartment and does not burst. When the electronic coupling parameter is infinitesimal, the two compartments are divided into two neurons and do not burst either. Compare with the variation of electrical coupling parameter gc , change the geometric coupling parameter p has analogous impact. When the p infinitely close to 1 (or 0), the two compartments are reduced into individual soma (or dendrite) part, and the neuron is in spike pattern. We note above that PR model can mimic almost all the behavior exhibited in real biological neurons. Specially, we can get several different rhythmic patterns by regulating the external stimulating parameter Is , Id as well as electrical coupling parameter gc and geometrical parameter p. For parameters (gc = 2.1 mS/cm2 , p = 0.5, Id = 0), Is produce bursting in a range of [0.1, 0.6]. For parameters gc = 2.1 mS/cm2 , p = 0.5, Is = 0, Id produce bursting in a range of [0.1, 2.3]. When geometrical parameter p change in the range of [0.1, 0.7], PR model displays bursting pattern, and electrical coupling parameter vary in the range of [1.9, 3], PR model also displays bursting pattern in certain external stimulation.

180

S. Li et al. / Physica A 491 (2018) 177–187

Fig. 2. Frequency and duty change of PR model induced by geometrical parameter p modification. The parameters respectively are (a) p = 0.2 (b) p = 0.3 (c) p = 0.4 (d) p = 0.5 and Is = 1, Id = 0, gc = 2.1 in all the conditions. The black lines represent the reference value of the outputs (zero) and the blue lines are membrane potentials of soma. The red dash-dots are used to calculate the frequency and duty ratio. They denote the transformation of the outputs. The percentages represent the duty rations and the bolded values denote the interbursting interval. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

2.1.2. Frequency and duty cycle of modified PR model From above discussions, we get that the rhythmic pattern can be changed by adjusting external stimulating inputs Is , Id , however, the frequency can only be changed slightly in the bursting range. Besides the frequency change of PR model, we also care about its duration ratio. A desired feature of CPG is that the duty cycle can be regulated for a given frequency. Conveniently, by adjusting geometrical parameter p, we can change the frequency and its duration ratio if the output is bursting, as seen in Fig. 2. In order to calculate the frequency and duty ratio, we calculate the maximal and minimal values of each output (blue line) respectively in the time frames shown in each panel’s abscissa axis, and make the transformation of the outputs. If the output is bigger than zero, we set it to the maximal value and on the contrary. After transformation, the square waves (red lines) are obtained. Through Fig. 2 we get that the period of PR model decreases dramatically (the frequency increases) and the duty ratio increases minimally as p increases. However, the change of frequency can be easily attained but the change of duty ratio is not ideal. In fact, most CPG models cannot realize the diversity of frequency and duration ratio at present except for the one mentioned in Ref. [26]. Zhang et al. proposed an alternative that add some extra components to process the original output of CPG [27]. They set a threshold to the filter output. For the same output, the duration ratio can be changed if we change the threshold. 2.1.3. Sensory feedback underlies entrainment of modified PR model One feature of the CPG is the generation of rhythmic pattern in the absence of sensory inputs, while it should note that the intrinsic frequency of CPG can be entrained by the sensory feedback to the nature frequency of the plant or environment [27]. Although they have similar form, the sensory e is different from external input Is /Id in neuroscience perspective. e denotes the afferent signal from sensory organs, while Is /Id indicates the signals from supraspinal cord or background environment. By adding sensory feedback e to formula (1) we can get the integrated PR model. If the sensory feedback is a series of rhythmic signals, it can entrain the neural oscillator easily, as shown in Fig. 3. In addition, the signal-to-noise ratio for subthreshold sinusoidal signal and noise is quite large, and this implies that the PR neuron can be entrained easily by sinusoidal signal, even in the presence of noise. We can conclude that sensory feedbacks at different frequencies can entrain PR model to follow their own frequency if their amplitude exceed certain values. For the sensory feedback signal with the same period and duration ratio, the number of spikes increases as the amplitude increases and the spacing between spikes in each period decreases. The interburst intervals also shorten. Thus, we can infer that the PR model can easily be entrained by sensory feedback. Besides, we can also tune the parameter of PR model to get relatively reasonable outputs based on the desired gait frequency in practice. 2.1.4. Robustness of PR model In neuronal system, neurons are living in a noisy environment and the noise origins from random transmissions of ions inside and outside the neuron (channel noise) or the stochastic activities of other neurons (synaptic noise). It should be ensured that the PR model is robust to external noises as it can be so easily entrained by sensory feedback. Actually, the PR model can maintain its behavior under influence of disturbances. In consideration of that noise usually interfuses in the afferent signal, such as sensory feedback signal, we add a White Gaussian noise with a power 15 into the sinusoidal sensory feedback signal with a amplitude of 5 and a period of 350 ms. Fig. 4(a) demonstrates the composite signal of white Gaussian noise and sinusoidal sensory feedback signal, and Fig. 4(b) demonstrates the output of PR neuron under noise. We can barely

S. Li et al. / Physica A 491 (2018) 177–187

181

Fig. 3. Entrainment of PR model (Id = 0, gc = 2.1). The blue lines represent the membrane potentials of PR model, and the red lines are the sinusoidal sensory feedbacks. (a) is the membrane potentials without sensory feedback. (b)(c) and (d) are the membrane potentials with a series of sine waves sensory feedback. The sensory feedback of (c) and (d) have the same period of 350 ms and the amplitude is 5 in (c) and 10 in (d). The period of (b) is 200 ms and its amplitude is 5. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 4. Noise-induced oscillations in a single PR neuron for a noise intensity D1 = 15. (a) Sinusoidal sensory feedback signal with white noise; (b) output of the PR neuron.

see the shape of original sinusoidal signal in Fig. 4(a), while the output rhythm is restored and the firing pattern is basically the same as in Fig. 3(c). Thus, we can speculate that the PR neurons have a certain resistance against noise. In addition, for a system subjected to noise alone, there exists noise strength that optimally enhances the coherence in response [28]. When the noise intensity exceeds the threshold, the PR neurons produce chaotic firing pattern.

2.2. Implementation of time-delayed CPG The firing patterns of a single modified PR neuron are well studied in previous sections. However, in order to develop a CPG model, the output rhythm is more important. We need to study more about the mechanism of the coupled PR neurons to obtain proper parameters for more complex networks. In the following section, we construct a CPG network based on two PR neurons coupling with a time-delay chemical synapse.

2.2.1. Chemical synapse In living nervous systems one finds three general types of synaptic connections among neurons [29]: ohmic electrical connections (gap junctions) and two types of chemical connections (excitatory and inhibitory couple). The synapse supplies the current to reshape the membrane potential of multiple electronic neurons. The bursting will be either synchronous (excitatory coupling) with depolarizing current injection or asynchronous (inhibitory coupling) with hyperpolarizing current injection. For our study of the CPG, we implement mutual chemical inhibitory synapses between the two PR neurons [30].

182

S. Li et al. / Physica A 491 (2018) 177–187

Fig. 5. Schematic of two dendrite–dendrite coupled neurons (cells coupled by a chemical synapse).

The inhibitory chemical synapse model is given in [31] as: Isyn = gsyn · r xpre · υ xpost

(

(

dr xpre

)

dt

=

)

(

)

r∞ − r

τr ((

(3)

r∞ = 1 + tanh xpre + 1.2 /0.9 /2 ( ) ( (( ) )) υ xpost = 1 + tanh xpost + 2.81 /0.4 /2

(

)

))

Where gsyn is the maximal synaptic (conductance, xpost and xpre are the membrane potential of the postsynaptic neuron ) determined by the presynaptic activity, and and presynaptic neuron, respectively. r xpre is the synaptic activation ) ( variable τr is the characteristic time constant of the synapse (τr ≈ 100). υ xpost is a nonlinear function of the membrane potential of the postsynaptic neuron that reduces the strong effect of the large spikes in postsynaptic neuron. 2.2.2. Architecture of CPG network with time-delay The CPG network studied in this article is a network of two coupled modified-PR models. The coupling form belongs to inhibitory chemical couple, and the coupled site is dendrite-soma between neurons. The schematic of this CPG is shown in Fig. 5. The mathematical model of the CPG is described by the following equations: Cm Cm

dVs1 (t ) dt dVd1 (t ) dt dVs2 (t )

= −ILeak − INa − IKDR −

IDS1 p

+

= −ILeak − ICa − IKCa − IKAHP −

Is1

+ gsyn · r (Vd2 (t − τ21 )) · υ (Vs1 (t ))

p ISD1

+

Id1

1−p 1−p Is2 Cm = −ILeak − INa − IKDR − + + gsyn · r (Vd1 (t − τ12 )) · υ (Vs2 (t )) dt p p dVd1 (t ) ISD2 Id2 Cm = −ILeak − ICa − IKCa − IKAHP − + dt 1−p 1−p IDS2

(4)

gsyn is the synaptic coupling strength, and it summarizes how information is distributed between neurons. τji is the net time delay-the time for the action potential propagate along the axon connecting the pre-synaptic neuron j to post-synaptic neuron i. The mutual delay in the coupling is motivated by the propagation delay of action potentials between the two neurons [32,33]. 3. Results The dynamics of all networks depend on the ongoing interplay between the intrinsic properties of the neurons that make up networks and the strength, time course, and time-dependent properties of the synapses among them. CPGs can exhibit robust oscillatory behavior. In the following part, we investigate the influence of ambient noise, sensory feedback signals, morphology parameter (the size ratio of soma and dendrite), as well as the coupling delay time on the output rhythm of CPG. 3.1. Stochastic synchronization of CPG Noise can induce oscillations even though the fixed point is stable [34]. The noise sources then play the role of stimulating the excitable subsystems. Even if only one subsystem is driven by noise, it induces oscillations of the whole system through

S. Li et al. / Physica A 491 (2018) 177–187

183

Fig. 6. (a) The phase plane of noise-induced oscillations in a PR neuron (gsyn = 0) for a noise intensity for D1 = 2. (b) Time series of two inhibitory chemical coupled CPG system for different noise intensities D1 : (i) D1 = 0, (ii) D2 = 10, (iii) D2 = 40. The red, blue and black curves represent to the membrane potentials of V1 , V2 and their sum V1 + V2 , respectively. The other parameters are p1 = p2 = 0.5, gsyn = 2.5, D2 = 0.2. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

the chemical coupling. We shall first consider two coupled PR neuron as in Eq. (1) without delay in the coupling (τ = 0). Since the neurons can be coupled to different parts of the neural network or of the environment, the noise intensities in the two neurons can be quite different. We set that each neuron is driven by Gaussian white noise ξi (t ) (i = 1, 2) with zero mean and unit variance. The noise densities are denoted by the parameters D1 and D2 . In our study, we set the noise intensity D2 in the second subsystem equals to a small value D2 = 0.09 in order to model some background noise level. Depending on the coupling strength gsyn and the noise intensity D1 in the first subsystem, the two neurons show cooperative dynamics. Fig. 6(a) depicts the temporal dynamics of a single PR neuron (gsyn = 0) due to stochastic input. One can see that the system is excitable because it performs large excursions in the phase plane. Fig. 6(b) shows the temporal dynamics of two coupled CPG system for increasing noise intensities D1 , where the red, blue and black curves correspond to V1 , V2 and their sum V1 + V2 , respectively. For D1 = 0, the first neuron is enslaved and emits a burst every time the second neuron does. The synchronization decreases with the increase of D1 . For strong noise intensities D1 the dynamics of the first neuron is independently dominated by its own stochastic input. 3.2. Entrainment of CPG network In previous section, we get that the sinusoidal sensory feedback signal can entrain the rhythmic output of PR neuron. In fact, the sinusoidal sensory feedback signals also entrain the output of CPG network. Even if only one subsystem is driven by the sensory feedback, it produces the synchrony of the whole system through the coupling. In the following, we add a rhythmic sensory feedback signal to the CPG consisting of two non-identical neurons without time delay in the coupling (τ = 0). The different of two neurons are indicted by geometrical parameters p1 = 0.3 and p2 = 0.5. The neuron 1 is stimulated by sinusoidal sensory feedback signal with different periods and the rhythmic outputs of the CPG are demonstrated in Fig. 7. The amplitude of sinusoidal signals is 5 µA. With the increase of period, the time needed to achieve synchronization increases, and phase-lock ratio decreases. We also observe that increase the amplitude to 10 µA in Fig. 7(e), the regularity of synchronization change; the couple strength between neurons increase results in the prolongation of time needed to achieve stable synchronization. The results above illustrate that the CPG network can be entrained by sensory feedback signal, and we can obtain the required biological rhythms through adjusting the amplitude and period of sensory feedback signals. 3.3. Effects of Morphology on outputs of CPG The basic structure of the single neuron is typically consisted of cell body (i.e., soma), axon and dendrites. The soma usually connects with more than one dendrite, but only one axon. The dendrites mainly receive and conduct synaptic inputs from other cells to the soma. The PR model was originally formulated over 20 years as a two-compartment reduction of the 19-compartment Traub CA3 cell model developed earlier [35]. The reduced model segregates the fast current for sodium spiking into a proximal, soma-like, compartment and the slower calcium and calcium-mediated currents into a dendritelike compartment. The parameter p denotes the size ratio of soma and dendrite. p = 1 represents that the neuron consist of soma part, while p = 0 represents that the neuron consist of dendrite part. In this part, we investigate the output of CPG corresponding to varying geometrical parameter p in a range of [0.1, 0.9]. Here we consider the phase of each neuron and propose phase synchronization to describe the output of CPG network. Firstly, it requires the definition of phase ϕ1 (t ) and ϕ2 (t ) for each neuron and comparison between them. In this paper, we introduce the phase for each system as

184

S. Li et al. / Physica A 491 (2018) 177–187

Fig. 7. Output rhythm of the CPG under sinusoidal sensory feedback signal. The red line represents the output of neuron 1 stimulated by sinusoidal sensory feedback with amplitude 5, and the blue line represents the output of neuron 2 stimulated by constant somatic input. The synchronization and phase-lock between two PR neurons vary with the increase of feedback signal periods. (a) the period of sensory feedback is 100 ms, and the two PR neurons synchrony with phase-lock 2:1 (b) and (c) have same phase-lock ratio 1:1 while different synchronization index ; (d) the phase-lock ratio is 2:3. Comparing with (a), the amplitude of sinusoidal signal is 10 in (e) and the phase-ratio change a little; the phase-ratio in (f) is still 2:1, while the time needed to achieve stable synchronization increase because of strong chemical couple between two neurons. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 8. (a) Phase synchronization indexes of CPG network under standard electrical coupling parameter gsyn = 2.5. (b) Phase synchronization index for different p1 and p2 = 0.1 : 0.05 : 0.9. The first row represents the phase synchronization index for p2 = 0.8; the second row represents the phase synchronization index for p2 = 0.45; the third row represents the phase synchronization index for p2 = 0.25.

t −t

ϕj (t ) = 2π t −ti−1 + 2π (i − 1) , j = 1, 2, . . .. Where ti is the time at which we observe a spike in the respective system’s i i−1 realization. We define n : m phase synchronization to occur if the phase difference ∆ϕn,m (t ) = ϕ1 (t ) − m ϕ (t ) exhibits n 2 horizontal plateaus of sufficient duration. Usually, if n : m synchronization takes place, ∆ϕn,m (t ) demonstrates plateaus occasionally interrupted by 2π jumps. On the plateaus, ∆ϕn,m (t ) usually oscillates around some local average level. As time grows, ∆ϕn,m (t ) drifts to plus or minus infinity [34]. In Ref. [36] several measures phase synchronization √⟨to characterize ⟩2 ⟨ ⟩2 were introduced. Here, we choose to estimate the synchronization index γn,m = cos ∆ϕn,m (t ) + sin ∆ϕn,m (t ) , γn,m can vary between 0 (no synchronization) and 1 (perfect n:m phase synchronization). Fig. 8(a) displays the phase synchronization index under standard coupling strength. Fig. 8(b) displays the synchronization indexes for three different geometrical parameter for neuron 1. The maximal value of phase synchronization occurs for two identical coupling neurons, and the synchronization index decreases with the increase of difference between two coupling neurons. When the geometrical parameters of two coupling neurons equal, the two neurons have same firing frequency and

S. Li et al. / Physica A 491 (2018) 177–187

185

Fig. 9. (a) Phase synchronization indexes for different coupling delay times. τ12 represents the coupling delay time for neuron 1 acting upon neuron 2, and τ21 represents the coupling delay time for neuron 2 acting upon neuron 1. (b) two different synchronous firing patterns. i denotes the complete synchronization for τ12 = τ21 ; ii denotes the antiphase oscillations for τ12 = 0, τ21 = 100, gsyn = 10.

phase. Therefore, it is much easier to achieve synchronization. The bigger of difference between two neuron’s structures, the bigger of difference between the intrinsic frequency, and it is more difficult to achieve synchronization. 3.4. Effects of coupled delay time on outputs of CPG For networks, both propagation delays of the electrical signals connecting different neurons and local neurovascular couplings lead to time delays [37–39]. Signal transmission time delays in a network of nonlinear oscillators are known to be responsible for a variety of interesting dynamic behaviors including phase-flip transitions leading to synchrony or out of synchrony [40]. While the chemical synaptic time delays are small (∼2 ms), the axonal conduction delays, which depend on the distance between neurons in the brain, can reach up to tens of milliseconds. Tuning the time delays may provide us with additional ‘degree of freedom’ for suppressing or changing the rhythmic behavior in biological systems. The time-delays can vary the strength of the control actions and modulate the timing of the actions as well. By appropriate choice of the delay time, synchronization can be either enhanced or suppressed. Fig. 9(a) displays phase synchronization indexes for different coupling delay times. The two neurons are identical completely, and the coupling strength gsyn equals to 10. It is obviously that phase synchronization index equals to 1 when delay times between neurons are similar. The phase synchronization indexes decrease with the increase of difference between delay times. Fig. 9(b) denote two different synchronous firing patterns: (i) absolutely phase synchronization for τ12 = τ21 and any coupling strength; (ii) antiphase oscillations for τ12 = 0, τ21 = 100, gsyn = 10. So, we can know that antiphase oscillations in delay-coupled neurons can be induced for sufficiently large delay and coupling strength. 4. Conclusions and discussions Infer from the study conducted above, we can get the conclusion that the modified PR model is capable of being a CPG. Firstly, PR model can generate rhythmic signals, which conforms to the basic property of CPG in appearance. Secondly, the PR model has several parameters controlled, including external input Is , electrical coupling parameter gc , geometrical parameter p and so on. Tuning the Is and p can easily regulate the period and duty ratio of output rhythm. Thirdly, the PR model can maintain its behavior under influence of disturbances below certain strength. Lastly, we have shown that the sensory feedback signals can be easily incorporated into PR model, and the entrainment is explored in single neuron and coupled neural network separately. In a word, PR model is suitable to be used as a CPG. Of course, the PR model is not perfect. It has an irregularing bursting area and chaotic patterns that are undesired when generating rhythmic patterns in practice. We also investigate the behavior of half-center CPG consisting of two PR model coupled with inhibitory chemical synapses, and consider the stochastic synchronization, the entrainment of CPG by sensory feedback signals, as well as the effects of different geometrical parameters and coupling delay time on the rhythmical output of CPGs. Even only one neuron is driven by noise, the two neurons of half-center CPG are synchronous oscillation, and the synchronization index decreases with the increasing of noise intensities. In addition, the CPG is easily entrainment by sensory feedback signals. We can obtain the required biological rhythms through adjusting the amplitude and period of sensory feedback signals. The morphological features of single neuron play a crucial role in the rhythmical output of CPG. When two neurons have identical structure, it is much easier to achieve synchronization for their similar frequencies; with the increasing of difference between two neuron’s structures, the synchronization index decrease and it is more difficult to achieve synchronization. The coupling delay time is another factor affecting the output of CPG. Proper coupling delay time can maintain the output of CPG as the desired pattern.

186

S. Li et al. / Physica A 491 (2018) 177–187

The formation and modulation of rhythmical signals generated by CPGs for essential behavior remains an open field. Thus, it is of great importance to explore the theoretical possibilities, so that we can provide biological experimentalists and roboticists ample paradigms from which to choose. While we have demonstrated a novel type of CPG, and found several important features affecting the output rhythmic, there may still exists other features playing a crucial role for the generation of rhythmical signals that are not yet addressed by the CPG literature or the current study. Thus, it is important to keep an open mind about the mechanism of how a CPG producing rhythmical output. Acknowledgments This work was supported by the National Natural Science Foundation of China (Grant No. 61471265, 61372010) and the Youth Fund of the National Natural Science Foundation of China (Grant No. 61401312). References [1] A.A. Hill, S.D.V. Hooser, R.L. Calabrese, Half-center oscillators underlying rhythmic movements, Nature 261 (2002) 146–148. [2] E. Marder, D. Bucher, Central pattern generators and the control of rhythmic movements, Curr. Biol. 11 (2001) (2001) 986–996. http://dx.doi.org/10. 1016/S0960-9822(01)00581-4. [3] J.L. Feldman, C.A. Del Negro, Looking for inspiration: new perspectives on respiratory rhythm, Nat. Rev. Neurosci. 7 (2006) 232–242. http://dx.doi.org/ 10.1038/nrn1871. [4] A.J. Ijspeert, Central pattern generators for locomotion control in animals and robots: a review, Neural Netw. 21 (2008) 642–653. http://dx.doi.org/10. 1016/j.neunet.2008.03.014. [5] A.A. Sharp, F.K. Skinner, E. Marder, Mechanisms of oscillation in dynamic clamp constructed two-cell half-center circuits, J. Neurophysiol. 76 (1996) 867–883. [6] T.G. Brown, on the nature of the fundamental activity of the nervous centers; together with an analysis of the conditioning of rhythmic activity in progression, and a theory of the evolution of function in the nervous system, J. Physiol. 48 (1914) 18–46. http://dx.doi.org/10.1113/jphysiol.1914. sp001646. [7] X.J. Wang, J. Rinzel, Alternating and synchronous rhythms in reciprocally inhibitory model neurons, Neural Comput. 4 (1992) 84–97. http://dx.doi.org/ 10.1162/neco.1992.4.1.84. [8] X.J. Wang, J. Rinzel, Spindle rhythmicity in the reticularis thalami nucleus: synchronization among mutually inhibitory neurons, Neuroscience 53 (1993) 899–904. http://dx.doi.org/10.1016/0306-4522(93)90474-T. [9] S. Grillner, P. Wallen, Central pattern generators for locomotion, with special reference to vertebrates, Annu. Rev. Neurosci. 8 (1985) 233–261. http://dx.doi.org/10.1146/annurev.ne.08.030185.001313. [10] U. Bassler, A. Buschges, Pattern generation for stick insect walking movements—multisensory control of a locomotor program, Brain. Res. Rev. 27 (1998) 65–88. http://dx.doi.org/10.1016/j.cub.2015.09.044. [11] A.L. Hodgkin, A.F. Huxley, A quantitative description of membrane current and its application to conduction and excitation in nerve, J. Physiol.-London 117 (1952) 500–544. http://dx.doi.org/10.1016/S0092-8240(05)80004-7. [12] A.D. Sánche, G.G. Izús, M.G. Dell’Erba, R.R. Deza, Hub-enhanced noise-sustained synchronization of an externally forced FitzHugh-Nagumo ring, Physica A 468 (2017) 289–298. http://dx.doi.org/10.1016/j.physa.2016.10.096. [13] A.D. Sánche, G. Cascallares, A.D. Sánche, Noise-sustained synchronization between electrically coupled fitzhugh–nagumo networks, Physica A 433 (2015) 356–366. http://dx.doi.org/10.1016/j.physleta.2016.04.017. [14] N. Dahasert, Ismail Öztürk, Recai Kiliç, Experimental realizations of the hr neuron model with programmable hardware and synchronization applications, Nonlinear Dyn. 70 (2012) 2343–2358. http://dx.doi.org/10.1007/s11071-012-0618-5. [15] N. Korkmaz, Ismail. Öztürk, R. Kiliç, The investigation of chemical coupling in a hr neuron model with reconfigurable implementations, Nonlinear Dyn. (2016) 1–14. http://dx.doi.org/10.1007/s11071-016-2996-6. [16] E.M. Izhikevich, Which model to use for cortical spiking neurons?, IEEE Trans. Neural. Networks. 15 (2004) 1063–1070. http://dx.doi.org/10.1109/TNN. 2004.832719. [17] P.F. Pinsky, J. Rinzel, Intrinsic and network rhythmogenesis in a reduced traub model for CA3 neurons, J. Comput. Neurosci. 1 (1994) 39–60. http://dx.doi.org/10.1007/BF00962717. [18] A.A. Laura, Y.P. Luke, T.A. Krasimira, Bifurcation analysis of a two-compartment hippocampal pyramidal cell model, J. Comput. Neurosci. 41 (2016) 91–106. http://dx.doi.org/10.1007/s10827-016-0606-8. [19] K.C. Bittner, C. Grienberger, S.P. Vaidya, A.D. Milstein, J.J. Macklin, J. Suh, S. Tonegawa, J.C. Magee, Conjunctive input processing drives feature selectivity in hippocampal CA1 neurons, Nat. Neurosci. 18 (2015) 1133–1142. http://dx.doi.org/10.1038/nn.4062. [20] A. Bose, V. Booth, Bursting in 2-compartment neurons: a case study of the Pinsky-Rinzel model, Bursting: The Genesis of Rhythm in the Nervous System 0 (2005) 123–144. http://dx.doi.org/10.1142/9789812703231_0005. [21] J.T. Brown, A.D. Randall, Activity-dependent depression of the spike after-depolarization generates long-lasting intrinsic plasticity in hippocampal CA3 pyramidal neurons, J. Physiol.-London. 587 (2009) 1265–1281. http://dx.doi.org/10.3410/f.1158243.618370. [22] Z.F. Mainen, T.J. Sejnowski, Influence of dendritic structure on firing pattern in model neocortical neurons, Nature 382 (1996) 363–366. http: //dx.doi.org/10.1038/382363a0. [23] J.P. Tseng, A novel approach to synchronization of nonlinearly coupled network systems with delays, Physica A 452 (2016) 266–280. http://dx.doi.org/ 10.1016/j.physa.2016.02.025. [24] J. Feng, F. Yu, Y. Zhao, Exponential synchronization of nonlinearly coupled complex networks with hybrid time-varying delays via impulsive control, Nonlinear. Dyn. 85 (2016) 621–632. http://dx.doi.org/10.1007/s11071-016-2711-7. [25] J. Hellgren, S. Grillner, A. Lansner, Computer simulation of the segmental neural network generating locomotion in lamprey by using populations of network interneurons, Biol. Cybern 68 (1992) 1–13. http://dx.doi.org/10.1007/BF00203132. [26] L. Righetti, A.J. Ijspeert, Design methodologies for central pattern generators: an application to crawling humanoids, P. Robotics: Syst. Sci. (2006) 191–198. http://dx.doi.org/10.15607/RSS.2006.II.025. [27] D.G. Zhang, X.Y. Zhu, L. Lan, K.Y. zhu, Mathematical study on ionic mechanism of lamprey central pattern generator model, Int. J. Neural Syst. 19 (2009) 409–424. http://dx.doi.org/10.1142/S0129065709002117. [28] J. Baltanas, J. Casado, Noise-induced resonances in the hindmarsh-rose neuronal model, Phys. Rev. E. 65 (2002) 041915. http://dx.doi.org/10.1103/ PhysRevE.65.041915. [29] M.S. Gordon, The Synaptic Organization of the Brain, in: M.S. Gordon (Ed.), Neurobiology, Oxford Univ. Press, New York, 1994.

S. Li et al. / Physica A 491 (2018) 177–187

187

[30] R.D. Pinto, P. Varona, A.R. Volkovskii, Synchronous behavior of two coupled electronic neurons, Phys. Rev. E. 62 (2000). http://dx.doi.org/10.1103/ PhysRevE.62.2644. [31] R. Huerta, P. Varona, M.I. Rabinovich, H.D.I. Abarbanel, Topology selection by chaotic neurons of a pyloric central pattern generator, Biol. Cybern 84 (2001) 1–8. http://dx.doi.org/10.1007/PL00007976. [32] E. Schöll, G. Hiller, P. Hövel, M.A. Dahlem, Time-Delayed feedback in neurosystems, Philos. Trans. Math. Phys. Eng. Sci. 367 (2008) 1079–1096. http://dx.doi.org/10.1098/rsta.2008.0258. [33] N. Buric, D. Todorovic, Dynamics of FitzHugh–Nagumo excitable systems with delayed coupling, Phys. Rev. E. 67 (2003) 066222. http://dx.doi.org/10. 1103/PhysRevE.67.066222. [34] B. Hauschildt, N.B. Janson, A. Balanov, E. Schöll, Noise-induced cooperative dynamics and its control in coupled neuron models, Phys. Rev. E. 74 (2006) 051906. http://dx.doi.org/10.1103/PhysRevE.74.051906. [35] R.D. Traub, R.K.S. Wong, R. Miles, H. Michelson, A model of a CA3 hippocampal pyramidal neuron incorporating voltage-clamp data on intrinsic conductances, J. Neurophysiol. 66 (1991) 635–650. [36] M.G. Rosenblum, A.S. Pikovsky, J. Kurths, C. Schäfer, P. Tass, Dectecting unstable periodic orbits in biological systems, in: F. Moss, S. Gielen (Eds.), The Handbook of Biological Physics, Elsevier Science, Amsterdam, 2001, pp. 279–321. [37] H.R. Wilson, Spikes, Decisions and Actions: The Dynamical Foundations of Neuroscience, Oxford University Press, New York, 1999. [38] W. Gerstner, W. Kistler, Spiking Neuron Models, Cambridge University Press, UK, 2002, pp. 277–280. [39] Haken H., Brain Dynamics: Synchronization and Activity Patterns in Pulse-Coupled Neural Nets with Delays and Noise, Springer, Berlin, Germany, 2006. [40] B.M. Adhikari, A. .Prasad, M. Dhamala, Time-delay-induced phase-transition to synchrony in coupled bursting neurons, Chaos 21 (2011) 023116. http://dx.doi.org/10.1063/1.3584822. [41] D.G. Zhang, Q. Zhang, X.Y. Zhu, Exploring a type of central pattern generator based on hindmarsh-rose model: from theory to application, Int. J. Neural Syst. 25 (2015). http://dx.doi.org/10.1142/S0129065714500282. [42] R. FitzHugh, Impulses and physiological states in theoretical models of nerve membrane, Biophys. J. 1 (1961) 445–466. http://dx.doi.org/10.1016/ S0006-3495(61)86902-6.