Synchronization of Chemical Synaptic Coupling of the Chay Neuron System under Time Delay

Abstract

This paper studies the chemical synaptic coupling of Chay neurons and the effect of adding time delay on their synchronization behavior. The results indicate that coupling strength stimuli can affect the discharge activity and the synchronization behavior. In the absence of coupling strength, the Chay neurons have chaotic discharge behavior and the system is in a nonsynchronous state. When a certain coupling strength is added, the neurons change from chaotic discharge to ordered periodic discharge, and the system state changes from asynchronous to synchronous. On the other hand, a time lag can alter the coupled system from synchronous to asynchronous.

1. Introduction

Synchronization is a common phenomenon in nonlinear dynamics [1]. In recent years, this phenomenon has received considerable attention due to its huge application prospects in chemistry, biology, physics, and information processing [2,3,4,5,6]. Several theoretical and experimental studies have shown that the synchronization of the nervous system may play an important role both in revealing the brain’s cognitive and sensory processes and in information processing mechanisms. It was discovered that neurological diseases such as Parkinson’s, trembling of the human hands, and epilepsy are caused by some morbid synchronization or non-synchronization [7,8,9]. Therefore, the presence, the absence, or the degree of synchronization can be important parts of the nervous system function or dysfunction. In order to study the nonlinear dynamic behavior and synchronization of single neurons and neural networks, researchers have established neuron mathematical models such as the Hodgkin–Huxley (HH), Morris–Lecar (ML), Hindmarsh–Rose (HR), FitzHugh–Nagumo (FHM), and the Chay representations. The Chay model was established in 1985 as a theoretical formulation with a high degree of uniformity based on different types of excitable cells playing an important role in Ca²⁺-related K⁺ channels, such as neurons, cardiomyocytes, and nerve endings. The model is able to simulate various excitation modes of excitable cells with high biological and physiological characteristics.

The study of synchronous control of biological systems became an important topic having in mind the control and the treatment of some diseases. It is generally believed that there are two ways for synchronizing the system: self-synchronization caused by natural coupling with feedback from various synapses between neurons, and control synchronization, caused by the coupling effect generated by an artificial system applying explicit feedback control [10]. In the first case, synchronization of the coupling neurons of the same type only occurs when a certain critical value of the coupling strength is reached, and this value exceeds the physiological condition of time [11].

Lin Feifei [12] designed an adaptive radial basis function (RBF) neural network and an integer-order parameter adaptive law based on the Lyapunov stability theory. Deng Bin [13] proposed a Proportional and Integral based iterative learning control algorithm for the synchronization control problem of HH neurons. The master-slave neuron synchronization control under four different conditions was simulated. The results showed that after applying the control, HH neurons could quickly track the dynamic behavior of the rest of the neurons.

Xiu Chunbo [14] took hysteretic chaotic neurons and networks as their research object and proposed a control strategy based on filter tracking error to achieve neuron/network stability control. Under the control strategy, the neuron/network gradually stabilized to the extreme point of the optimization function so as to solve the optimization problem. Hettiarachchi [15] adopted auto-disturbance rejection control to study the chaotic synchronization between two neurons coupled through gaps. By adding controllers, it was verified that the neurons were very robust to external disturbances. In biological experiments, the synchronization of two coupled neurons could be achieved by applying an external direct current to the depolarization [16]. With the continuous development of the theory of nonlinear dynamic systems, several methods for achieving synchronization of neurons have been proposed, including feedback control, inversion design, nonlinear control, adaptive control, sliding mode control, and other techniques [17,18,19,20,21].

The time delay has a considerable influence on the dynamic properties of a coupled system. In terms of stability, a time lag often leads to system instability, resulting in various forms of bifurcation, inducing periodic oscillations, multiple steady-state, almost period, and even chaotic motion [22,23,24,25,26,27,28]. Li Hongming [29] studied the Chay neuron model under constant current stimulation based on the stability theory. The results show that with the change of the control parameter I (external constant current), Hopf bifurcation will occur. In terms of chaotic motion, a time lag causes chaos in first-order nonlinear dynamical systems. For systems without a time delay, neither first-order systems, nor second-order autonomous systems, can produce chaos [30]. Zheng Yanhong [31] discussed the influence of time delay on the synchronous dynamics of four ring-coupled chaotic HR neurons. By calculating the synchronization difference, it was found that a proper time lag will induce or improve the complete synchronization between chaotic neurons. The time lag can also cause the system to generate high-dimensional odd, strange attractors with multiple positive Lyapunov exponents, without being limited by the system dimension [32,33]. In terms of synchronization and its transfer dynamics, the time delay can increase or suppress the synchronization of the coupled system. A time delay can induce not only the in-phase and anti-phase synchronization of the coupled system, but also the phase-synchronized transfer of different coupled systems and the transfer of different types of synchronization [34,35,36].

This paper introduces a mathematical model consisting of a pair of chemical synaptic-coupled Chay neurons and analyzes the basic dynamic behavior of the coupled neuron system. The paper is organized as follows. Section 1 introduces the synchronization and time delay and reviews the research progress in nonlinear dynamics. Section 2 introduces the Chay neuron model and the modification of a neuron’s electrical activity when some parameters change. Section 3 presents the coupling Chay neuron system, analyzes the electrical activity of the coupled neurons and the synchronization, and evaluates the effect of the time-delay stimulation on the synchronization of the coupled neuron system under a certain coupling strength. Finally, Section 4 summarizes the main results.

2. Chay Neuron Model

The Chay mathematical model describes the characteristics of activation of Ca²⁺-sensitive K⁺ channels when the concentration of Ca²⁺ ions in neurons increases. This model is used to study the abnormal discharge of the nervous system. In the 1980s, Chay proposed six variables and five variable models [37]. In 1985, a three-variable simplified model was developed [38]. This paper adopts the improved three-variable model [39] that was formulated by Chay in 1995. In this model, the nonlinear differential equations are established using the K⁺ channel open probability, Ca²⁺ slow variable, and membrane potential as the dynamic variables, respectively. Six expressions of the open speed and probability variables of K⁺ channel and the speed and probability variables of activation and inactivation of Na⁺ channels are based on electrophysiological properties and relevant experimental data. This description improves the HH model, which involved a considerable number of variables that made it a cumbersome task to handle. More important, the addition of changes in the Ca²⁺-sensitive K⁺ channels explains the nonlinear characteristics of neuronal impulses, spikes, and chaotic discharges. The differential equations are as follows:

$$\frac{\mathrm{dV}}{\mathrm{dt}}={\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{I}}-\mathrm{V}\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}-\mathrm{V}\right){\mathrm{n}}^4{+\mathrm{g}}_{\mathrm{K},\mathrm{C}}\frac{\mathrm{C}}{1+\mathrm{C}}\left({\mathrm{V}}_{\mathrm{K}}-\mathrm{V}\right){+\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_{\mathrm{L}}-\mathrm{V}\right)$$

(1)

$$\frac{\mathrm{dC}}{\mathrm{dt}}=\mathsf{\rho }\left[{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{C}}-\mathrm{V}\right)-{\mathrm{k}}_{\mathrm{C}}\mathrm{C}\right]$$

(2)

$$\frac{\mathrm{dn}}{\mathrm{dt}}=\frac{{\mathrm{n}}_{\infty }-\mathrm{n}}{{\mathsf{\tau }}_{\mathrm{n}}}$$

(3)

Equation (1) describes the variation of neuronal cells. The four terms on the right side of (1) represent the mixed Na⁺–Ca²⁺ channels, K⁺ ion channels with potential-dependent conductance, the K⁺ ion channel current, and the transmembrane leakage current having a conductance that depends on the Ca²⁺ ion concentration rather than the potential. The symbol V_K denotes the reversible potential of the K⁺ ion channels, V_L is the reversible potential of various ions in other channels, V_I is the reversible potential of the Na⁺–Ca²⁺ channels, and V is the equilibrium potential. The parameters g_I, g_K,V, g_K,C, and g_L represent the maximum conductance of each channel, respectively.

Equation (2) describes the variation of the Ca²⁺ ion concentration C in cell membranes. The two terms on the right side of (2) indicate the Ca²⁺ ion channel current in and out of the cell membrane, K_C denotes the ratio constant of intracellular Ca²⁺ ion efflux, ρ is the proportionality constant, and V_C represents the reversible potential of Ca²⁺ ions.

Equation (3) describes the variation of the opening probability for the gate of the K⁺ ion channels with the potential.

In the sequel Equations (4)–(6) show that the variables ${\mathrm{m}}_{\infty }$, ${\mathrm{n}}_{\infty }$ and ${\mathrm{h}}_{\infty }$ are the stable probability values of the opening of the Ca²⁺, K⁺, and Na⁺ ion channels. In Equation (7), τ_n represents a time scale factor and the smaller its value, the steeper is the spike; the parameter λ_n is the time kinetic constant. The expressions can be written as follows:

$${\mathrm{m}}_{\infty }=\left(0.1\left(25+\mathrm{V}\right)/\left(1-{\mathrm{e}}^{(-0.1\mathrm{V}-2.5)}\right)\right)/\left[\left(0.1\left(25+\mathrm{V}\right)/\left(1-{\mathrm{e}}^{(-0.1\mathrm{V}-2.5)}\right)\right){+4\mathrm{e}}^{(-(\mathrm{V}+50)/18)}\right]$$

(4)

$${\mathrm{n}}_{\infty }=\left(0.01\left(20+\mathrm{V}\right)/\left(1-{\mathrm{e}}^{\left(-0.1\mathrm{V}-2\right)}\right)\right)/\left[\left(0.01\left(20+\mathrm{V}\right)/\left(1-{\mathrm{e}}^{\left(-0.1\mathrm{V}-2\right)}\right)\right)+\left({0.125\mathrm{e}}^{\left(-(\mathrm{V}+30)/80\right)}\right)\right]$$

(5)

$${\mathrm{h}}_{\infty }=\left({0.07\mathrm{e}}^{(-0.05\mathrm{V}-2.5)}\right)/\left[\left({0.07\mathrm{e}}^{(-0.05\mathrm{V}-2.5)}\right)+\left({0.125\mathrm{e}}^{(-(\mathrm{V}+30)/80)}\right)\right]$$

(6)

$${\mathsf{\tau }}_{\mathrm{n}}=1/\left[{\mathsf{\lambda }}_{\mathrm{n}}\left(\left(0.01\left(20+\mathrm{V}\right)/\left(1-{\mathrm{e}}^{\left(-0.1\mathrm{V}-2\right)}\right)\right)+\left({0.125\mathrm{e}}^{\left(-(\mathrm{V}+30)/80\right)}\right)\right)\right]$$

(7)

The parameters adopted herein are depicted in

Table 1. Values of Chay neuron model parameters.

gI = 1800 s−1	gK,V = 1700 s−1	gL = 12 s−1	gK,C = 17 s−1
VK = −68 mV	VL = −40 mV	VC = 200 mV	KC = 3.3/18
λn = 223.8	$\mathsf{\rho }$ = 0.27 mV−1s−1

.

The inter-spike interval bifurcation plot versus the parameter V_I in the Chay neuron model is shown in Figure 1. With the change in V_I in the range [75, 165], the Chay neuron model reveals rich dynamic characteristics. From Figure 1, we find that the value of V_I varies from small to large (marked as region A) and the neuron model enters six cycles and five cycles of periodical discharge patterns. As the V_I value continues to increase (region B), when V_I = 98.2 mV, the neuron model performs four cycles of cluster discharge activity. As the value of V_I increases further, when V_I = 107.8 mV, the neuron model enters a chaotic discharge mode, producing the shaded region marked as C in Figure 1. This is the expression in the inter-spike interval bifurcation diagram. As V_I increases and the value of V_I is greater than 117.6 mV we enter region D. The inter-spike interval sequence is restored to the second-cycle, and the first-cycle discharge mode after an inverse-periodic bifurcation.

3. Synchronization of the Chemical Synaptic-Coupled Chay Neuron System

Synapses have strong adaptability (plasticity) that depend on various external stimuli (such as the action potential of presynaptic neurons). Moreover, chemical synapses can be divided into two types: excitatory and inhibitory [40,41,42,43]. Chemical synapses, which are the most common type of cellular connections in the animal’s nervous system, are highly flexible, making this topic widely studied [44,45,46].

To analyze the synchronous behavior of chemical synaptic-coupled Chay neurons, this paper adopted two chemical synaptic-coupled Chay neuron models as follows:

$$\frac{{\mathrm{dV}}_1}{\mathrm{dt}}={\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{I}}-{\mathrm{V}}_1\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right){\mathrm{n}}_1^4{+\mathrm{g}}_{\mathrm{K},\mathrm{C}}\frac{{\mathrm{C}}_1}{{1+\mathrm{C}}_1}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right){+\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_{\mathrm{L}}-{\mathrm{V}}_1\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_1\right)}{1+\exp (-{\mathsf{\sigma }(\mathrm{V}}_2-\mathsf{\theta }))}$$

(8)

$$\frac{{\mathrm{dC}}_1}{\mathrm{dt}}=\mathsf{\rho }\left[{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{C}}-{\mathrm{V}}_1\right)-{\mathrm{k}}_{\mathrm{C}}{\mathrm{C}}_1\right]$$

(9)

$$\frac{{\mathrm{dn}}_1}{\mathrm{dt}}=\frac{{\mathrm{n}}_{\infty }-{\mathrm{n}}_1}{{\mathsf{\tau }}_{\mathrm{n}}}$$

(10)

$$\frac{{\mathrm{dV}}_2}{\mathrm{dt}}={\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{I}}-{\mathrm{V}}_2\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right){\mathrm{n}}_2^4{+\mathrm{g}}_{\mathrm{K},\mathrm{C}}\frac{{\mathrm{C}}_2}{{1+\mathrm{C}}_2}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right){+\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_{\mathrm{L}}-{\mathrm{V}}_2\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_2\right)}{1+\exp (-{\mathsf{\sigma }(\mathrm{V}}_1-\mathsf{\theta }))}$$

(11)

$$\frac{{\mathrm{dC}}_2}{\mathrm{dt}}=\mathsf{\rho }\left[{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{C}}-{\mathrm{V}}_2\right)-{\mathrm{k}}_{\mathrm{C}}{\mathrm{C}}_2\right]$$

(12)

$$\frac{{\mathrm{dn}}_2}{\mathrm{dt}}=\frac{{\mathrm{n}}_{\infty }-{\mathrm{n}}_2}{{\mathsf{\tau }}_{\mathrm{n}}}$$

(13)

The following numerical values were adopted: synaptic reversible potential V_syn = −65, synaptic threshold θ = −45, and ratio constant σ = 10. Equation groups with subscripts 1 and 2 denote sub-systems 1 and 2, respectively. The variable H_syn represents the coupling strength and V_syn describes the synapse reversible potential of the neuron, which depends on the presynaptic neuron and the recipient. The parameter θ is the synaptic threshold beyond which its presynaptic neurons begin to act on post-synaptic neurons, and σ denotes the rate constant for the onset of excitation or inhibition. The rest of the variables indicate that the significance and value do not change.

In order to better study the synchronization of coupled neurons, we define the synchronization discrepancy between coupled neurons: e₀ = V₂ − V₁, e₁ = n₂ − n₁, e₂ = C₂ − C₁, and e₃ = (|e₀| + |e₁| + |e₂|)/3.

Constructing Lyapunov function as

$$\mathrm{L}=\left({\mathrm{e}}_0^2+{\mathrm{e}}_1^2{+\mathrm{e}}_2^2\right)/2,$$

(14)

then the derivative of L is:

$$\frac{\mathrm{dL}}{\mathrm{dt}}={\mathrm{e}}_0\frac{{\mathrm{de}}_0}{\mathrm{dt}}+{\mathrm{e}}_1\frac{{\mathrm{de}}_1}{\mathrm{dt}}+{\mathrm{e}}_2\frac{{\mathrm{de}}_2}{\mathrm{dt}}$$

(15)

The type of substitution (15) can be obtained as follows:

$$\begin{gathered} \frac{\mathrm{dL}}{\mathrm{dt}}={\mathrm{e}}_0[{\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+{\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}\left({\mathrm{n}}_2^4-{\mathrm{n}}_1^4\right)+{\mathrm{V}}_1{\mathrm{n}}_1^4-{\mathrm{V}}_2{\mathrm{n}}_2^4\right)+{\mathrm{g}}_{\mathrm{K},\mathrm{C}}\left(\frac{{\mathrm{C}}_2}{{1+\mathrm{C}}_2}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right)-\frac{{\mathrm{C}}_1}{{1+\mathrm{C}}_1}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right)\right)+ \\ {\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_2\right)}{1+\exp (-{\mathsf{\sigma }(\mathrm{V}}_1-\mathsf{\theta }))}-\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_1\right)}{1+\exp (-{\mathsf{\sigma }(\mathrm{V}}_2-\mathsf{\theta }))}]+{\mathrm{e}}_1\left(\frac{{\mathrm{n}}_1-{\mathrm{n}}_2}{{\mathsf{\tau }}_{\mathrm{n}}}\right)+{\mathrm{e}}_2\left[\mathsf{\rho }\left({\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+{\mathrm{K}}_{\mathrm{C}}{(\mathrm{C}}_1-{\mathrm{C}}_2)\right)\right] \end{gathered}$$

(16)

If $\frac{dL}{dt}<0$ is established, it needs:

$$\begin{gathered} {\mathrm{e}}_0[{\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_1-{\mathrm{V}}_2\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}\left({\mathrm{n}}_2^4-{\mathrm{n}}_1^4\right)+{\mathrm{V}}_1{\mathrm{n}}_1^4-{\mathrm{V}}_2{\mathrm{n}}_2^4\right)+{\mathrm{g}}_{\mathrm{K},\mathrm{C}}\left(\frac{{\mathrm{C}}_2}{{1+\mathrm{C}}_2}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right)-\frac{{\mathrm{C}}_1}{{1+\mathrm{C}}_1}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right)\right)+ \\ {\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_2\right)}{1+\exp (-{\mathsf{\sigma }(\mathrm{V}}_1-\mathsf{\theta }))}-\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_1\right)}{1+\exp (-{\mathsf{\sigma }(\mathrm{V}}_2-\mathsf{\theta }))}]+{\mathrm{e}}_1\left(\frac{{\mathrm{n}}_1-{\mathrm{n}}_2}{{\mathsf{\tau }}_{\mathrm{n}}}\right)+{\mathrm{e}}_2\left[\mathsf{\rho }\left({\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+{\mathrm{K}}_{\mathrm{C}}{(\mathrm{C}}_1-{\mathrm{C}}_2)\right)\right]<0 \end{gathered}$$

(17)

The substitution of e₀ = V₂ − V₁, e₁ = n₂ − n₁ and e₂ = C₂ − C₁ in (17) gives the following:

$$\begin{gathered} ({\mathrm{V}}_2-{\mathrm{V}}_1)[{\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_1-{\mathrm{V}}_2\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}\left({\mathrm{n}}_2^4-{\mathrm{n}}_1^4\right){+\mathrm{V}}_1{\mathrm{n}}_1^4-{\mathrm{V}}_2{\mathrm{n}}_2^4\right){+\mathrm{g}}_{\mathrm{K},\mathrm{C}}\left(\frac{{\mathrm{C}}_2}{{1+\mathrm{C}}_2}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right)-\frac{{\mathrm{C}}_1}{{1+\mathrm{C}}_1}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right)\right)+ \\ {\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_2\right)}{1+\exp (-\mathsf{\sigma }({\mathrm{V}}_1-\mathsf{\theta }\left)\right)}-\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_1\right)}{1+\exp (-\mathsf{\sigma }({\mathrm{V}}_2-\mathsf{\theta }\left)\right)}]+({\mathrm{n}}_2-{\mathrm{n}}_1)\left(\frac{{\mathrm{n}}_1-{\mathrm{n}}_2}{{\mathsf{\tau }}_{\mathrm{n}}}\right)+ \\ \left({\mathrm{C}}_2-{\mathrm{C}}_1\right)\left[\mathsf{\rho }\left({\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_1-{\mathrm{V}}_2\right)+{\mathrm{K}}_{\mathrm{C}}({\mathrm{C}}_1-{\mathrm{C}}_2)\right)\right]<0 \end{gathered}$$

(18)

When Equation (18) is set up, $\frac{\mathrm{dL}}{\mathrm{dt}}<0$, then we can get the existence and limit of $\begin{array}{c}\lim \\ \mathrm{t}\to \infty \end{array}{{\displaystyle \int }}_0^{\mathrm{t}}{\mathrm{e}}_{\mathrm{i}}\left(\mathsf{\tau }\right)\mathrm{d}\mathsf{\tau },\mathrm{i}=0,1,2$. According to the Lyapunov stability theory, by adjusting the parameters, the $\frac{\mathrm{dL}}{\mathrm{dt}}<0$, the error state e₀, e₁, e₂ and e₃ are asymptotically stable at the origin, namely, $\begin{array}{c}\lim \\ \mathrm{t}\to \infty \end{array}{\mathrm{e}}_{\mathrm{i}}\left(\mathrm{t}\right)=0,\mathrm{i}=0,1,2$, so the system has reached sync.

In a first phase, let us consider the discharge activities of two Chay neuron models without the coupling strength. When the coupling strength is H_syn = 0, solving the two sets of equations leads to the (V₁, V₂) locus of the two synaptic-coupled Chay neurons shown in Figure 2. The irregular movements indicate that the system is in an asynchronous state. The time history of V₁ and V₂ are represented in Figure 3. We verify that the values of V₁ and V₂ are different and that the discharge of the neurons is in a state of chaos; in other words, the system is in a chaotic state. In addition, the synchronization difference map of the coupled neurons is depicted in Figure 3. We observe that the synchronization discrepancy for all time instants indicates that the two neuron system is not in synchronization. In summary, it can be concluded that when H_syn = 0, the two neuron systems are in a state of chaotic discharge and are in an asynchronous state.

In a second phase, let us consider the coupling strength H_syn = 0.8 that leads to the (V₁, V₂) diagram of the chemical synaptic-coupled Chay neuron system shown in Figure 4. The plot is a straight line with a slope of 1, indicating that the system is in a synchronized state. The time evolution of V₁ and V₂ depicted in Figure 5 shows also that the curves of V₁ and V₂ coincide. Similarly, it can also be seen from the synchronization discrepancy diagram of Figure 5 that the difference of the synchronization discrepancy with time is always zero, indicating that the coupled neuron system is in a synchronized state. In summary, when the coupling strength is H_syn = 0.8, the system is in a synchronized state.

It can be concluded after the above analysis that appropriate changes in the coupling strength can induce the synchronization of the chemical synaptic-coupled Chay neuron system.

4. Time-Delay Inhibition of Synchronization of the Chemical Synaptic-Coupled Chay Neuron System

The coupling modes between neurons are divided into two types: electrical and chemical synaptic coupling. The chemical synapses depend on the release of neurotransmitters and are connected by means of dendrites and axons. Furthermore, chemical synapses dominate throughout the neuronal system. Due to the presence of synaptic gaps in the nervous system, there is a significant delay in the transmission of neural signals through the synapse. This effect results in a time delay in the transmission of signals, thereby causing a time lag in the nervous system [47,48,49]. Due to the presence of a time lag, the coupled neurons become an infinite dimensional dynamic system that exhibit a richer nonlinear behavior.

We consider two chemical synaptic-coupled Chay neuron models with time delays as follows:

$$\frac{{\mathrm{dV}}_1}{\mathrm{dt}}={\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{I}}-{\mathrm{V}}_1\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right){\mathrm{n}}_1^4{+\mathrm{g}}_{\mathrm{K},\mathrm{C}}\frac{{\mathrm{C}}_1}{{1+\mathrm{C}}_1}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_1\right){+\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_{\mathrm{L}}-{\mathrm{V}}_1\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_1\right)}{1+\exp (-\mathsf{\sigma }({\mathrm{V}}_2(\mathrm{t}-\mathsf{\tau })-\mathsf{\theta }\left)\right)}$$

(19)

$$\frac{{\mathrm{dC}}_1}{\mathrm{dt}}=\mathsf{\rho }\left[{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{C}}-{\mathrm{V}}_1\right)-{\mathrm{k}}_{\mathrm{C}}{\mathrm{C}}_1\right]$$

(20)

$$\frac{{\mathrm{dn}}_1}{\mathrm{dt}}=\frac{{\mathrm{n}}_{\infty }-{\mathrm{n}}_1}{{\mathsf{\tau }}_{\mathrm{n}}}$$

(21)

$$\frac{{\mathrm{dV}}_2}{\mathrm{dt}}={\mathrm{g}}_{\mathrm{I}}{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{I}}-{\mathrm{V}}_2\right){+\mathrm{g}}_{\mathrm{K},\mathrm{V}}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right){\mathrm{n}}_2^4{+\mathrm{g}}_{\mathrm{K},\mathrm{C}}\frac{{\mathrm{C}}_2}{{1+\mathrm{C}}_2}\left({\mathrm{V}}_{\mathrm{K}}-{\mathrm{V}}_2\right)+{\mathrm{g}}_{\mathrm{L}}\left({\mathrm{V}}_{\mathrm{L}}-{\mathrm{V}}_2\right)+\frac{{\mathrm{H}}_{\mathrm{syn}}\left({\mathrm{V}}_{\mathrm{syn}}-{\mathrm{V}}_2\right)}{1+\exp (-\mathsf{\sigma }({\mathrm{V}}_1(\mathrm{t}-\mathsf{\tau })-\mathsf{\theta }\left)\right)}$$

(22)

$$\frac{{\mathrm{dC}}_2}{\mathrm{dt}}=\mathsf{\rho }\left[{\mathrm{m}}_{\infty }^3{\mathrm{h}}_{\infty }\left({\mathrm{V}}_{\mathrm{C}}-{\mathrm{V}}_2\right)-{\mathrm{k}}_{\mathrm{C}}{\mathrm{C}}_2\right]$$

(23)

$$\frac{{\mathrm{dn}}_2}{\mathrm{dt}}=\frac{{\mathrm{n}}_{\infty }-{\mathrm{n}}_2}{{\mathsf{\tau }}_{\mathrm{n}}}$$

(24)

where t and τ denote time and time lag, respectively. In Equations (19)–(24), the significance and values of other parameters and variables are consistent with Section 2.

Figure 6 shows the (V₁, V₂) locus for a coupling strength H_syn = 1.55 and time lag τ = 0 in the case of the two coupled neuron system (19)–(24). The values of V₁ and V₂ are identical and the system is in synchronization state. Similarly, the time history diagram of V₁ and V₂ reveals that the two curves of V₁ and V₂ overlap. The synchronization difference diagram is also shown in Figure 7, indicating that the two neuron system is in a synchronized state. In summary, when the coupling strength H_syn = 1.55 and the time lag τ = 0, the two neuron coupling system is in a synchronized state.

In a third phase, we discuss the synchronization of two chemical synapse-coupled neurons with the coupling strength H_syn = 1.55 and τ = 3.4. Figure 8 reveals that the locus (V₁, V₂) has an irregular behavior, indicating that the system is in an asynchronous state. The time history of V₁ and V₂, illustrated in Figure 9, shows that the two curves do not overlap. Similarly, Figure 9 demonstrates that the synchronization difference between V₁ and V₂ is always present, revealing that the system does not reach the synchronized state. In summary, the two chemical synaptic coupling Chay neuron system is in an asynchronous state when the coupling strength H_syn = 1.55 and time lag τ = 3.4.

In order to further verify the effect of coupling strength and time delay on the neuron system, this paper fixed the time delay and changed the coupling strength: τ =2.5 and H_syn = 2.7. By observing Figure 10, the curve of random motion of the system in phase plane shows that the neuron has chaotic discharge activity and the neuron system is in the unsynchronized state. In Figure 11, the time history diagram does not overlap the motion curve of the neurons and the difference of the synchronous difference map is not zero, which indicates that the system is in an asynchronous state.

Finally, when the coupling strength increases to H_syn = 5.3 with τ = 2.5, the phase plane diagram of the coupled neuron system (Figure 12) shows that the trajectory of the system is still an irregular random curve, indicating that the system does not reach the synchronous state. The time history diagram and the synchronous difference diagram in Figure 13 also verify that the system is in an asynchronous state.

From the above discussion, we conclude that the addition of time lags in the two chemical synaptic-coupled Chay neuron system yields a system change from periodic to chaotic discharge activity, leading to the modification from synchronous to non-synchronized state and destroying the original synchronization state.

5. Conclusions

In this paper, the chemical synaptic coupling of a Chay neuron system was studied. The discharge characteristics of the Chay neuron model were analyzed through the state locus, synchronization discrepancy plot, and inter-spike interval bifurcation diagram. The effect of time delay and coupling strength on the coupled Chay neuron system was discussed. The model simulations showed that in a chemical synaptic-coupled Chay neuron system, coupling strength can promote synchronization behavior. Moreover, adding a time delay can inhibit the occurrence of its synchronous behavior, changing from a synchronous to a non-synchronized state. The results have given a deeper insight into the dynamic characteristics of neuron models and provide a basis for further research of neuronal system synchronization.