Ben Cao(曹奔), Huaguang Gu(古华光),†, and Yuye Li(李玉叶)

1School of Aerospace Engineering and Applied Mechanics,Tongji University,Shanghai 200092,China

2College of Mathematics and Computer Science,Chifeng University,Chifeng 024000,China

3Institute of Applied Mathematics,Chifeng University,Chifeng 024000,China

Keywords: bifurcation,bursting,excitatory autapse,time delay

1. Introduction

Neuronal electronic activities exhibit complex nonlinear dynamics, for example, the bursting behavior alternating between burst of spikes and quiescent state and bursting pattern transitions or bifurcations,which are involved in information processing,locomotion rhythms,or brain disease.[1–5]In the real nervous system, the neuronal electronic activities are modulated by inhibitory or excitatory synapses including the autapse,which is a special synapse connecting different parts of a same neuron.[6–8]In general, the inhibitory effects play roles opposite to those of excitatory effects: the former one always suppresses the firing activity and the latter one often promotes the firing activity. Recently,the inhibitory and excitatory effects induced-nonlinear phenomena in contrast to the general viewpoint have attract much attention, which are important for both nonlinear physics and neurophysiology.[9–11]For example,the inhibitory stimulus can induce the firing rate increased in biological experiment on auditory system,which is called post-inhibitory facilitation(PIF).[12,13]The PIF phenomenon are very important for the identification of sound location of auditory system and of coincidence detection of the general circuit. In addition,the inhibitory(hyperpolarization)stimulation can induce action potential from steady state,which is called post-inhibitory rebound(PIR)spike,has been observed in real nervous system with hyperpolarization active caution current[14]and simulated in the theoretical model with Hopf bifurcation,[15]which extends the conception of threshold and is related to information process and locomotion rhythm. In theoretical model of neuronal network,the excitatory effect induces suppression of the firing activity[16]and termination of bursting activity,[17]resulting in suppression of firing activity,which is related to the epilepsy.

Recently, normal as well as paradoxical phenomena induced by self-feedback or autapse has attract a lot of attention. The self-feedback mediated by autapse has been identified in various neurons. In biological experiment, inhibitory auatpse can suppress action potential[18]and induce enhanced precision of spike timing,[6]and excitatory autapse can induce persistent firing and enhancement of bursting.[7,19]In the theoretical models, autapse has been identified to play important roles in modulating firing pattern transitions,[20–23]stochastic dynamics,[24–27]and spatiotemporal dynamics of networks.[28–33]In addition to these normal phenomena, the paradoxical phenomenon for both inhibitory and excitatory autapses are acquired. One the one hand, the experimental observation can be explained with theoretical model. For example, in the biological experiment on interneuron which is type-2 excitability and involved in Gamma oscillations,the inhibitory self-feedback can induce spiking from resting state.[34]As suggested by the experimental result, inhibitory autapse with suitable time delay is identified to induce resting state changed to spiking for Hodgkin–Huxley model with type-2 excitability.[35]On the other hand,the prediction of theoretical model can be verified in the experiment. For instance,in theoretical model, inhibitory autapse is identified to promote the synchronization of Gamma oscillation,[36]which is verified in a recent experimental study.[37]Such a result shows that the theoretical study may have guide effect on the experiment. Therefore,the paradoxical effect of inhibitory and excitatory on neural firing patterns and spatiotemporal behaviors of network are investigated in theoretical models. For example, inhibitory autapse is identified to enhance bursting activity[38]and wave propagations in the neuronal network,[39]and the excitatory autapse is identified to suppress spiking[40]and bursting activity.[41]These results extend the contents of nonlinear physics and await experimental demonstration.

The bursting behavior is composed of burst contains multiple spike,[15]which are involved in the physiological functions of various brain neurons.[42]Based on that the bursting behavior appear in the system composed of fast subsystem and slow variable,the bursting behaviors are distinguished into different patterns, according to the bifurcations of the fast subsystem with the slow variable regarded as bifurcation parameter. For example,the“Fold/Homoclinic”bursting pattern,the“Fold/Big Homoclinic” bursting, and the “Circle/Fold cycle”bursting are shown in Figs. 1(a), 1(b), and 1(c), respectively.The bifurcation before the virgule such as the “Fold” corresponds to the beginning phase of the burst,and the bifurcation after the virgule such as the“Homoclinic”corresponds to the ending phase of the burst. Such a method is called fast-slow variable dissection method. With this method, the dynamical mechanism of the suppressed bursting activity induced by excitatory in Ref. [41] is acquired. The bursting belongs to“Fold/Homoclinic”pattern. The autaptic current with suitable time delay and strength take effects before the Fold bifurcation, the quiescent state terminates in advance to transited to burst, resulting in that both durations of quiescent state and burst become shorter and the spikes per burst become less.Because the decrease degree of quiescent state duration is larger than that of the burst duration, the firing frequency increases rather than decreases. Therefore,although the burst activity is suppressed,the firing activity is still enhanced.

For the excitatory autapse with time delay, the“Fold/Homoclinic”bursting exhibits paradoxical phenomenon in burst duration instead of firing rate.[41]An important question arises. Can the excitatory autapse induce not only less spikes per burst but also lower firing frequency for other bursting patterns? Such a question is important for the extension of nonlinear physics, potential functions of excitatory self feedback/autapse,and possible modulation measures to brain neurons. In the present paper,we aim to answer this question. As can be found from experiment,there is time delay between the action potential and autaptic current.[43,44]Therefore,autapse with time delay is used. The lower firing frequency modulated by excitatory autapse cannot be easily simulated in the“Fold/Homoclinic”bursting shown in Fig.1(b),which resembles to Ref.[41]. Due to the results of the“Fold/Homoclinic”bursting resemble those of Rulkov model,[41]in the present paper, we aim to study the influence of excitatory autapse on bursting patterns different from the“Fold/Homoclinic”bursting. The “Circle/Fold cycle” bursting is chosen and both the less spikes per burst and the lower firing frequency can be simulated. Furthermore, reduced bursting activities in both spikes and frequency is well explained with trajectory of bursting combined with the bifurcation points of the fast subsystem. The results present novel examples of paradoxical phenomenon that excitatory autapse induces negative responses of firing frequency, which extends nonlinear dynamics, presents potential functions of auatpse, and possible modulation measure to brain neurons.

The rest of the present paper are organized as follows.The model and method, the results, and conclusions are presented in Sections 2,3,and 4,respectively.

Fig.1. Spike trains:(a)“Fold/Homoclinic”bursting from the Chay model;(b)“Fold/Big Homoclinic”bursting from the modified FitzHugh–Nagumo model(Fig.86 in Ref.[15]);(c)“Circle/Fold cycle”bursting from the modified Morris–Lecar model(Fig.66 in Ref.[15]).

2. Model and method

In the present paper,two theoretical neuron models commonly used to investigate dynamics of different types of bursting patterns are considered, which are modified FitHugh–Nagumo model and modified Morris–Lecar model.[15]Two models exhibit two different types of the bursting patterns.

2.1. The modified FitHugh–Nagumo model

2.1.1. The equations of modified FitHugh–Nagumo(FHN) model

The modified FHN model can simulate the special“Fold/Big Homoclinic”bursting,[15]and the equations are described as follows:

where the variables v and w represent the membrane potential and the membrane recovery variable,and u is a slow variable to modulate dynamics of membrane recovery variable. The first two equations without S(w)=b/(1+exp((c−w)/d))is the origin FHN model, which has been widely used to investigate the dynamics of the spiking and bifurcations of equilibrium points.

The modified FHN model is dimensionless, and the parameter values are as follows: ε =0.15, µ =−0.0005, c=−0.5,b=1.75,d=0.1.

2.1.2. The modified FHN model with autapse

With Eqs.(2) and(3)unchanged, the excitatory autaptic current Iautis introduced to Eq. (1) to form FHN model with autapse described as follows:

where Iaut=−gaut(vpos(t)−vaut)Γ(vpre(t)),[35]and Γ(vpre)=1/(1+exp(−λ(vpre−θaut))). The parameter gautis the autaptic conductance. For an autapse, vpos(t) = v(t) is the postsynaptic membrane potential,and vpre(t)=v(t)is the presynaptic membrane potential. The parameter vaut(=2 in the present paper to ensure the excitatory autapse) is the reversal potential,and θautrepresents the threshold of autaspe. The parameter θaut=0 to ensure that excitatory autaptic current generates when spikes cross the threshold.The parameter λ =30,which indicates a fast release rate of neurotransmitter.

In the present paper, time delay τ of autapse is considered and Iaut=−gaut(v(t)−vaut)Γ(v(t −τ)). The control parameters are the gautand τ,which determine the strength and application phase of the excitatory autaptic current pulse, respectively.

2.1.3. Bifurcations of the fast subsystem

For the modified FHN model, u is the slow variable.Therefore,equations(1)–(2)with u regarded as bifurcation parameters are the fast subsystem,which is shown as follows:

According to bifurcation points to the fast subsystem,[15]the bursting can be classified into different patterns. Firstly, with changing u, the bifurcations of both the fixed points(equilibrium points) and the limit cycles are acquired. Secondly, the bursting trajectory in phase plane(u, v)of the modified FHN model is acquired and plotted with the bifurcations. Last,according to the relation between the bifurcation points and trajectory of bursting, the dynamics or patterns of the bursting are identified.

2.2. The modified Morris–Lecar model

2.2.1. The equations of modified Morris–Lecar (ML)model

The modified ML model can simulate“Circle/Fold cycle”bursting,[15]and the equations are described as follows:

where the first two equations without u is the origin ML model to describe the firing activity of the muscle fiber of the barnacle. The first variable V and the second variable w represent the membrane potential and the gated variable of potassium channel, respectively. The third variable u is a slow dimensionless feedback to modulate membrane current.The reversal potential of three ions, potassium, calcium, and leakage ions are represented by VK, VCa, and VL, respectively. The corresponding maximal conductances are represented by gK, gCa,and gL. The parameter µ is a feedback coefficient. The functions m∞(V)and w∞(V)are steady state of the active probability of Ca2+channel and the gate variable w,respectively,and the function k(V)represents time factor, which are described as follows:

2.2.2. The modified ML model with autapse

With Eqs.(8) and(9)unchanged, the excitatory autaptic current Iautis introduced to Eq. (7) to form ML model with autapse described as follows:

where Iaut=−gaut(Vpos(t)−Vaut)Γ(Vpre(t)),[35]and Γ(Vpre)=1/(1+exp(−λ(Vpre−θaut))). The parameters Vpos(t)=V(t)and Vpre(t) =V(t) are the postsynaptic membrane potential and the presynaptic membrane potential of the ML model,respectively. The meanings and values of other parameters are the same as those in Subsubsection 2.1.2.

2.2.3. Bifurcations of the fast subsystem

For the modified ML model,u is the slow variable.Therefore,equations(7)–(8)with u regarded as bifurcation parameter are the fast subsystem,which is shown as follows:

Similar to Subsubsection 2.1.3,according to bifurcation points to the fast subsystem,[15]the bursting can be classified into different patterns. Firstly,with changing u,the bifurcations of both the fixed points(equilibrium points)and the limit cycles are acquired. Secondly,the bursting trajectory in phase plane(u,V)of the modified ML model is acquired and plotted with the bifurcations. Last, according to the relation between the bifurcation points and trajectory of bursting, the dynamics of bursting are identified.

2.3. Methods

The Euler method is used to integrate the equations of theoretical models,and the time step is 0.01. The bifurcations are acquired with software of XPPAUTO.[45]

3. Results

3.1. The results of FHN model

3.1.1. The “Fold/big Homoclinic” bursting

The “Fold/big Homoclinic” bursting pattern of the modified FHN model without autapse is period-8 bursting, as shown in Fig. 2(a). Different from the “Fold/Homoclinic”bursting shown in Fig.1(a),the“Fold/big Homoclinic”bursting pattern exhibits undershoot characteristic, which means that the minimal value within the burst is small than that of quiescent state. However, the“Fold/Homoclinic”bursting do not manifest undershoot characteristic. The bursting patterns with or without undershoot characteristic widely exist in brain neurons.[42]The period is 500, the quiescent state is 225,and the burst duration is 275. The mean firing frequency is f0=0.0160,which is called intrinsic firing frequency.

Fig. 2. Dynamics of period-8 bursting in the modified FHN model. (a)Membrane potential; (b) Bifurcations of fast subsystem. H, SH, and SN represent a supercritical Hopf bifurcation point on up branch, a big homoclinic orbit(a saddle-homoclinic)bifurcation,and a saddle-node or fold bifurcation which corresponds to the intersection point between low and middle branch. Stable (solid green) limit cycles with maximal value (upper)and minimal value(lower);(c)panel(b)plotted with trajectory of period-8 bursting(blue solid);(d)the enlargement of panel(c).

The bifurcations of fast subsystem of the modified FHN model with respect to the parameter u are depicted in Fig.2(b).The equilibrium curve of the fast subsystem is an“S”-shaped curve composed of up,middle,and low branches. The middle branch(dotted curve)and low branch(bold solid curve)represent the saddle and stable node,respectively. The intersection point between low and middle branch is saddle-node(SN)or fold bifurcation at u ≈−1.183. There is a super-critical Hopf bifurcation(H)at u ≈−0.748 locating on the up branch. The equilibrium before and after the Hopf point on the up branch is unstable(dashed curve)and stable(thin solid curve)focus,respectively. A stable limit cycle (SLC) appears through the Hopf bifurcation,and evolves from right to left with decreasing parameter u,and terminates when contacts with saddle to form a big homoclinic orbit and a saddle-homoclinic(SH)bifurcation at u ≈−1.215 (The detailed dynamics of the big homoclinic orbit can be found in Ref. [15]). The maximal and minimal values of the SLC are represented with upper and lower green curves, respectively. Different from Chay model shown in Fig. 1 wherein the minimal value V of the SLC is larger than V value of stable node on low branch,the minimal value V of the SLC of the modified FHN model is smaller than V value of stable node.

The trajectory of period-8 bursting (blue solid curve) in(u,V)plane and figure 2(b)are plotted in one figure,as illustrated in Fig.2(c). Figure 2(d)is the enlargement of Fig.2(c).The behavior of the period-8 bursting is transition between the stable node on low branch and the SLC,as can be found from Fig.2(d). The quiescent state and the burst of the bursting correspond to low branch and SLC,respectively.The burst begins from the fold(SN)point,and terminates at the SH point via a big homoclinic orbit. The burst oscillates between maximal and minimal values of the SLC from right to left to form the 8 spikes within a burst. According to Ref. [15], the period-8 bursting belongs to “Fold/Big Homoclinic” bursting. Due to the characteristics of “Fold/Big Homoclinic” bursting behaviors in the modified FHN model,the range of slow variable u is lower than zero.

3.1.2. Less spikes per burst and larger firing frequency

The dependence of the average number of spikes within a burst on the parameters gautand τ is shown in Figs.3(a)and 3(b). The average number larger than 8 appears in the black region,and not larger than 8 in the remained region,as shown in Fig.3(a).To clearly show the decrease or increase of the average number in two colors, the decreased number is labeled with blue and the enhanced number with black, as shown in Fig.3(b),which shows less spikes per burst and resembles the result of Ref.[41].

Fig. 3. The dependence of the indicators of bursting on (τ, gaut). (a) The mean number of spikes per burst (multiple color scales). Number larger than 8 represented by black;(b)the average number of spikes per burst(two color scales). Lower and not lower than 8 in the blue and black regions,respectively;(c)The mean firing frequency f;(d)the ratio f/f0. Larger and not larger than 1 in black and blue regions,respectively.

The dependence of the mean firing frequency f of bursting(Fig.3(c))and the ratio(f/f0)between f and the intrinsic frequency f0(Fig.3(d))on gautand τ are acquired,as depicted in Fig. 3(c) and Fig. 3(d), respectively. The ratio (f/f0) is larger than 1 in whole parameter region,as shown by the black region in Fig. 3(d), which means that the firing frequency is enhanced. Some cases of the reduction of spike number per burst and the increase of mean firing frequency induced by excitatory autapse are introduced as follows.

When τ = 112.05 and gaut= 0.01, a period-2 bursting(black, lower) is induced by the excitatory autaptic current(blue, upper), as shown in Fig. 4(a1). Compared with origin period-8 bursting without autapse, although the number of spikes within a burst is decreased from 8 to 2,the duration of quiescent state is also decreased from 225 to 53.The period of the period-2 bursting pattern is decreased to 121.44,and the mean firing frequency becomes slightly larger,from 0.0160 to 0.0165, which is due to the change of quiescent state plays a relatively dominant role. As can be found from Fig.4(a1),the first autaptic current pulse takes effect at quiescent state before the first spike of the burst and induces the quiescent state changed to period-2 bursting,which is determined by the time delay τ. As can be found from Fig.4(a2),the first pulse of the excitatory autaptic current(arrow)is applied at quiescent state,which is after the sixth spike of the origin period-8 bursting in(u, v) space, therefore, the quiescent state is terminated and changed to spike. After 2 spikes, the burst still terminates at SH point to form the period-2 bursting. Therefore,compared with the origin period-8 bursting, the transition from quiescent state to spikes happens at a phase prior to the SN point of the fast subsystem, which is the cause for the period-2 bursting induced by excitatory inhibitory autaptic current. Such a transition in advance results in reduction of the quiescent state duration,burst duration,and spike number per burst,due to the range of u for both burst and quiescent state becomes narrow.

As τ increases to 171.10 and gaut=0.01,a period-3 bursting appears,as shown in Fig.4(b1),which is due to the application phase of the autaptic current is after the fifth spike of the origin period-8 bursting. The duration of quiescent state is 86, and the period is 179.73. The mean firing frequency is 0.0167,which is slightly larger than the intrinsic frequency f0=0.0160. Compared with the period-2 bursting shown in Figs.4(a1)and 4(a2),the application phase of autaptic current delays,which induces that duration of quiescent state becomes longer and the number of spikes gets larger.

As τ increases to 272.00 with gautfixed to be 0.01, a period-5 bursting is evoked, as shown in Fig. 4(c1). The duration of quiescent state is 130.33, and the period is 281.38.The mean firing frequency is 0.0178, which is larger than the intrinsic frequency f0=0.0160. As can be found from Fig.4(c2),the dynamics for the period-5 bursting is complex.The excitatory autaptic current plays role after the fifth spike of the origin period-8 bursting. Different from period-2 bursting and period-3 bursting wherein the burst terminates at the SH point,the burst for period-5 bursting terminates at a phase left to the SH point. Therefore, 5 spikes instead of 3 spikes appear in a burst,resulting in period-5 bursting.

Fig. 4. The representatives of the reduction of spike number per burst and the increase of mean firing frequency at different τ and gaut values. Left:Membrane potential(black solid)and excitatory autaptic current(blue dash);Right: Bifurcations of fast system,trajectory of bursting corresponding to the left panel(red solid),and trajectory of period-8 bursting(blue dashed). The arrows represent the application phase of the first pulse of the autaptic current. (a1)and(a2)τ =112.05 and gaut=0.01;(b1)and(b2)τ =171.10 and gaut=0.01;(c1)and(c2)τ =272.00 and gaut=0.01.

3.2. The results of ML model

3.2.1. The “Circle/Fold cycle” bursting

The modified ML model without autapse exhibits period-6 bursting with period 365 and firing frequency 0.0164 labeled with f1, as illustrated in Fig. 5(a). There is depolarization block (steady state with high membrane potential following the burst) between burst and quiescent state, which is more complex than both ones without depolarization block depicted in Fig.1. The bursting with depolarization block behavior can be found at the brain neuron at spreading depression.[46]The duration of burst,depolarization block,and quiescent state are about 95,171,and 99,respectively.

Fig.5. Dynamics of period-6 bursting in the modified ML model. (a)Membrane potential;(b)Bifurcations of fast subsystem. H,Flc,and SNIC represent a subcritical Hopf bifurcation point on up branch,a fold point of limit cycles,and a saddle-node bifurcation on an invariant cycle,which corresponds to the intersection point between low and middle branch. Unstable (dotted blue) and stable (solid green) limit cycles with maximal value (upper) and minimal value(lower);(c)panel(b)plotted with trajectory of period-6 bursting(blue solid);(d)The enlargement of panel(c).

The fast subsystem of the modified ML model exhibits complex bifurcations with respect to u, as illustrated in Fig. 5(b). The curve of fixed or equilibrium points is composed of low branch (stable node, black bold curve), middle branch corresponding to saddle (black dotted curve), and up branch. The low branch insects with middle branch to form a saddle-node bifurcation on an invariant cycle(SNIC)appearing at u ≈−0.071. At u ≈−0.040 on the up branch, a Hopf bifurcation(subcritical,H)appears. Stable focus locates prior to the H point and unstable focus after the H point. An unstable(dotted blue)and stable(solid green)limit cycle with maximal value(upper)and minimal value(lower)contact with each other to form a fold or saddle-node bifurcation of limit cycles(Flc)at u ≈−0.091.

After plotting(u,V)trajectory(blue solid)of bursting to Fig.5(b)containing bifurcations of fast subsystem,figure 5(c)is acquired, and figure 5(d) is local enlargement around the burst of Fig. 5(c). The burst starts from the SNIC point, oscillates between the maximal value and minimal value of the stable limit cycle,evolves from right to left,and terminates at the Flc to form burst with 6 spikes. After the burst, the trajectory runs around first and then along the stable focus from left to right,runs across the H point,further runs along the unstable focus to form the depolarization block,runs to a phase wherein the repelling effect of unstable focus(middle branch)and attraction of stable node(low branch)can induce transition to the low branch,and runs along the low branch from right to left to form the quiescent state. Such a bursting is called“Circle/Fold cycle”bursting in Ref.[15].Different from“Fold/Big Homoclinic” bursting, for the “Circle/Fold cycle” bursting in the modified ML model, there is depolarization block behavior, which induces the maximal u value increased to a large extent to be higher than zero, and the range of u value corresponding to bursting is from values less than zero to larger zero.

3.2.2. Less spikes per burst and lower firing frequency

Fig.6. The dependence of the indicators of bursting on(τ,gaut). (a)The mean number of spikes per burst(multiple color scales). Number larger than 6 represented by black;(b)The average number of spikes per burst(two color scales). Lower and not lower than 6 in the blue and black regions,respectively;(c)The mean firing frequency f;(d)The ratio f/f1. Smaller and not smaller than 1 in blue and black regions,respectively.

The dependence of the average number of spikes within a burst on the parameters gautand τ is shown in Figs.6(a)and 6(b). The average number larger than 6 appears in the black region,and not larger than 6 in the remained region,as shown in Fig.6(a).To clearly show the decrease or increase of the average number in two colors, the decreased number is labeled with blue and the enhanced number with black, as shown in Fig.6(b),which shows less spikes per burst and resembles the result of Ref.[41].

The dependence of the mean firing frequency f of bursting (Fig. 6(c)) and the ratio (f/f1) between f and the intrinsic frequency f1(Fig. 6(d)) on gautand τ are acquired.The mean firing frequency f larger than the intrinsic firing frequency f1= 0.0164 appears in the black region, which shows enhancement of the bursting activity, as depicted in Fig. 6(d). However, in the blue region, f is smaller than f1=0.0164,which shows that the excitatory autapse induces the reduction of firing frequency.The reduced firing frequency presents novel example of paradoxical phenomenon different from Ref.[41].

3.2.3. Four cases of bursting behaviors

According to both indicators shown in Fig. 6, the spike number and the firing frequency,the bursting activities can be divided into four cases shown in Fig.7.

Fig.7. The distributions of Case 1(gray),Case 2(orange),Case 3(green),and Case 4(blue)of bursting in the parameter plane(τ,gaut).

Case 1 (gray): Both the mean of spike number and firing frequency decrease,which represents the paradoxical phenomenon and the novel bursting activity. Case 1 mainly appears when gautis relatively small.

Case 2 (orange): Both indicators increase, which is in accord with the common view. Case 2 mainly appears when gautand τ are relatively small(down-left region)or gautand τ are relatively large(up-right region).

Case 3 (green): The average number of spikes decreases while the mean firing frequency increases, which is paradoxical phenomenon in spike number and resembles the result of Ref.[15]. Case 3 mainly appears when gautis relatively strong and τ is not long.

Case 4 (blue): The average number of spikes remains unchnaged while the mean firing frequency decreases, which presents novel result of paradoxical phenomenon in firing frequency. Case 4 appears in several small regions.

Along the horizontal dashed line for gaut=0.015,each of Cases 1, 2, 3, and 4 appears more than one time in different discrete narrow τ regions when 0 ＜τ ＜155, and Case 1 appears for large τ (155 ＜τ ＜365). Along the vertical dashed line with τ =125, Case 2, Case 1, and Case 3 appear in sequence with increasing gaut.

3.2.4. The dynamical mechanism for the four cases

The four cases of bursting behaviors are explained with representatives as follows:

When τ =40 and gaut=0.015,the period-3 bursting appears and belongs to Case 1. The autaptic current (upper,blue)and spike trains of bursting(lower, black)are shown in Fig. 8(a1). The autaptic current with τ =40 (from the hollow cycle to the peak of the arrow)takes effect near the third trough, which results in the period-3 bursting. From period-6 to period-3, the spikes within a burst becomes less, correspondingly,the burst duration becomes short from 95 to 47.72,while the depolarization block duration slightly increases from 171 for to 177.4,and the duration of quiescent state increases from 99 to 116. Therefore, although the spike number per burst and burst duration decreases to a large extent (about 50%), while the bursting period decreases to a small extent(from 365 to 341.12, about 6.5%), due to the increase of depolarization block and quiescent state, which shows that the decrease of burst duration and spikes per burst are dominant.Therefore,the firing frequency decreases from f1=0.0164 to 0.0088.

The dynamics of the period-3 (red solid) and period-6 bursting (blue dashed) and the bifurcations of fast subsystem are shown in Fig. 8(a2). For the period-3 bursting, the burst terminates near the 3rd through,and the trajectory transits stable focus. Therefore, the burst duration and spike number within a burst decrease. Compared with period-6 bursting,the initial phase and the termination phase of the depolarization block of period-3 bursting move to right,resulting in increase of the depolarization block duration.The initial phase(i.e.,the termination phase of the depolarization block)of the quiescent state moves to right,and the termination phase of the quiescent state remains unchanged(i.e.,the SNIC point). Therefore,the quiescent state duration becomes long.

For τ = 125 and gaut= 0.005, Case 2 behavior corresponding to period-8 bursting appears, as depicted in Figs. 8(b1) and 8(b2). During the quiescent state, the positive autaptic current generates,resulting in the increase of the membrane potential and larger than the stable node on the low branch. Therefore, the quiescent state terminates at a phase right to the SNIC point and changes to spikes, which leads more spikes within a burst to form period-8 bursting with a longer duration 138.05 (95 for period-6 bursting). During the depolarization block, the excitatory autaptic current exhibits multiple oscillations, which induces larger oscillations of the membrane potential to result in left shift of termination phase of the depolarization block. Therefore,the depolarization block duration decreased to 154.96 (171 for period-6 bursting). The initial and termination phases of quiescent state move to left and right, respectively, to form a shortened quiescent state duration 82.3 (99 for period-6 bursting). From period-6 bursting to period-8 bursting,the spikes per burst increase 33.3%, which is dominant, and period of bursting increases 2.8%(from 365 to 375.31). Therefore,the mean firing frequency become large from f1=0.0164 to 0.0213.

Fig. 8. The representatives of four cases at different τ and gaut values. Left: Membrane potential (black solid) and excitatory autaptic current (blue dash); Right: Bifurcations of fast system, trajectory of bursting shown in left figure(red solid), and trajectory of period-6 bursting(blue dashed). (a1)and(a2)τ=40 and gaut=0.015: Case 1;(b1)and(b2)τ=125 and gaut=0.005: Case 2;(c1)and(c2)τ=125 and gaut=0.07: Case 3;(d1)and(d2)τ =70 and gaut=0.015: Case 4.

As illustrated in Figs. 8(c1) and 8(c2), a period-4 bursting appears when τ =125 and gaut=0.07, which belongs to Case 3 appears. The burst duration decreases to 50.62. The burst with 4 spikes terminates near the 4th trough, which is induced by excitatory autaptic current pulse. The autaptic current within the depolarization block of the period-4 bursting is about 0.158, as illustrated in Fig. 8(c1), which is stronger than 0.029 for the period-3 bursting (Fig. 8(a1)), because gaut= 0.07 is much stronger than gaut= 0.015. Therefore,the membrane potential of depolarization block of period-4 bursting is higher than that of the stable focus, as shown in Figs.8(c2)and 8(a2),respectively.The large distance between the depolarization block of the period-4 bursting and the focus induces the termination phase of the depolarization block moved to left to a large extent,resulting in a shorted duration of the depolarization block 54.43. For the quiescent state,the beginning phase moves to left while the ending phase remains at the SNIC point. Therefore,the quiescent state duration decreases to 21.93. Although spikes per burst decrease 33.3%,the period of bursting decreases to a larger extent(65.2%,from 365 to 126.98),because the large decrease of both depolarization block and quiescent state duration is dominant. Therefore, the mean firing frequency increases to 0.0315 (0.0164 for period-6 bursting).

When τ=70 and gaut=0.015,the bursting is still period-6, which belongs to Case 4, as illustrated in Figs. 8(d1) and 8(d2). The autaptic current is positive at last half phase of the burst, resulting in a slightly decreased burst duration 92.81.Positive autaptic current plays role within the duration of depolarization block, resulting in large oscillations near the termination phase of depolarization block, right shift of the termination phase of depolarization block,and slightly increased depolarization block duration 182.94. The initial phase of the quiescent state, i.e., the termination phase of depolarization block, moves to right, resulting in the increased quiescent state duration 122.92. The period of bursting increases to 398.67. Therefore,the mean firing frequency becomes lower,from 0.0164 to 0.015.

4. Conclusion

In the present paper, the less spikes per burst and the higher firing frequency of bursting induced by excitatory autapse for“Fold/Big Homoclinic”and the less spikes per burst and the lower firing frequency of bursting induced by excitatory autapse for “Circle/Fold cycle” are acquired. The former resemble to Ref. [41], while the latter is different from Ref. [41], which presents novel example of the paradoxical phenomenon that excitatory autapse can reduce firing frequency. Such a results has significance in three aspects. First,it extends the nonlinear science. The excitatory self-feedback can induce negative response for bursting behavior, which is different from the negative responses of the spiking activity[40]or networks[16,17]induced by excitatory effects. Such a result presents an extension to the common viewpoint that excitatory self-feedback always induces positive response. Second,for real brain neurons with excitatory autapse and bursting behaviors, which have been widely observed in the experiment, the result presents potential function of the excitatory autapse,awaited experimental demonstration. Last,the result presents modulation measure to bursting activity,which can be implemented in circuit design, due to that circuit implement has been a feasible manner to realize some novel modulation measures.[9,10,47]

The cause for reduction of both spike number and firing frequency of “Circle/Fold cycle” bursting are acquired.Compared with “Fold/Homoclinic” bursting in Ref. [41] and“Fold/Big Homoclinic”bursting in the present paper,the“Circle/Fold cycle” bursting is more complex due to it has depolarization block. As excitatory autapse can induce less spikes per burst, the reduction of spike number is not dominant for the “Fold/Homoclinic” and “Fold/Big Homoclinic” bursting but dominant for“Circle/Fold cycle”bursting because the existence of depolarization block. Therefore, reduction of firing frequency appears for“Circle/Fold cycle”bursting,which can be explained with the relation between the trajectory of bursting and bifurcation points of fast subsystem. In future,the influences of noise[48–51]on the paradoxical phenomenon for various bursting patterns should be investigated and extended to networks models.[16]As suggested in Refs.[52–54],when the coupled neurons with time delay exhibit complete synchronization and nearly complete synchronization,the synchronous behavior resembles the behavior of uncoupled single neuron to a large extent. Therefore, it can be speculated that the single cells in coupling state have the phenomenon similar to that of single neuron in the present paper,which will be studied in future.