Synchronized dynamics of cortical neurons with time-delay feedback

The dynamics of three mutually coupled cortical neurons with time delays in the coupling are explored numerically and analytically. The neurons are coupled in a line, with the middle neuron sending a somewhat stronger projection to the outer neurons than the feedback it receives, to model for instance the relay of a signal from primary to higher cortical areas. For a given coupling architecture, the delays introduce correlations in the time series at the time-scale of the delay. It was found that the middle neuron leads the outer ones by the delay time, while the outer neurons are synchronized with zero lag times. Synchronization is found to be highly dependent on the synaptic time constant, with faster synapses increasing both the degree of synchronization and the firing rate. Analysis shows that pre-synaptic input during the inter-spike interval stabilizes the synchronous state, even for arbitrarily weak coupling, and independent of the initial phase. The finding may be of significance to synchronization of large groups of cells in the cortex that are spatially distanced from each other.


I. INTRODUCTION
There is a significant amount of research showing that spike coincidence of neurons encodes information.Synchronization of neural activity has been involved in such important phenomena as ability to recognize objects (by binding different attributes) [9], olfactory discrimination [2], and has even been proposed as one of the neural correlates of consciousness [3].Abnormal neural synchronization has been linked to the onset of epileptic attacks.
Evidence of long range correlation leading to synchronized clusters has been observed in the visual cortex [5] where synchronized oscillatory responses are conjectured to contribute to the binding of distributed features in a visual scene.Further research in [11] on the cat visual cortex supports that oscillatory activity and firing patterns can contribute to long range synchrony among cortical cluster of neurons.Other mechanisms involving gamma oscillations are also conjectured as one underlying cause for long range correlations in the human brain.
On the modeling front, recent models involving self-feedback loops in simple motif networks have shown that coupling architecture may improve long range correlations.Traub et al [16] uses a linear self-feedback loop in the coupling architecture of cortical neurons to inhance synchrony.In [10], theory and experiments on semi-conductor lasers, which exhibit bursting phenomena similar to neurons, show that complete synchronization may occur with self-feedback loops.In addition, lasers coupled in a line without self-feedback have been shown to synchronize the end lasers [6,7], leading to the possibility of long range correlations with delay.
Synchronization of nearby cells is often the result of receiving common input, such as when a pyramidal cell sends projections to a targeted area in the cortex.While pyramidal cells tend to target specific areas, matrix projection cells from the thalamus reach in a diffuse manner into adjacent cortical areas helping to synchronize the activity of large populations of cells.Since many strong synaptic connections between different brain areas tend to be one-way, many of the networks in the brain possess a hierarchical structure, which to a first approximation is a forward network, modulated by a weaker modulatory feedback ( [4]).For example, there is a strong forward projection from the thalamic LGN area to V1 or from V1 to MT, with weaker modulatory feedback projections that modulate the magnitude of the cell's response ( [14]).It has been suggested that exessively strong mutually coupled loops in the brain would promote uncontrollable oscillations, such as in epilepsy ( [13]).
With this in mind, we propose a model of synchronization of cortical cells receiving a common time-delayed input from a different area of the brain, such as the case with projections of pyramidal cells to other cortical areas or the thalamus to the cortex, and receiving a weaker time-delayed feedback.The central dynamical question addressed considers when the synchronous behavior of such a neural network is stable, particularly in regard to the synaptic coupling strengths, and the synaptic time constant.We find that shorter synaptic time constant of the target cells promotes synchronization, and even increases the firing rate, for the same strength of input.
The paper is organized as follows: In Section II, the basic model is set up, (see Fig. 1) and presented along with numerical results that show the dependence of synchronization and firing rate on the strength of coupling and synaptic conductance.Section III, analyzes the synchronization observed for short synaptic times by linearizing about the dynamics of two nearby trajectories.It shows that fast synaptic input tends to synchronize two nearby trajectories, regardless of their phase with respect to the input signal.Section III concludes and summarizes.

II. BASIC MODEL AND NUMERICS
Neurons can be largely divided into two classes, dependent on their spiking properties.
Class I neurons can be stimulated to fire at an arbitrarily low frequency, due to a saddle-node bifurcation, with increasing frequency as the magnitude of the stimulating current increases.
Class II neurons, on the other hand, only begin to fire at relatively high frequency, with their limit cycle resulting from a subcritical Hopf bifurcation.Class II neurons are well represented by the squid axon, while a large majority of the mammalian neurons are of the Class I type.Dynamics of a human neo-cortical neuron in the absence of synaptic connections are well approximated by the following Class I neuron equations: [12] where V and R are voltage and recovery variables, respectively, and τ R = 5.6 ms.The dR/dt equation is written as a sum of a linear term, for the normal Na + and K + currents, and a quadratic terms in V to approximate the the transient potassium current contributions [17].The above model has been optimized to provide an accurate quantitative fit to the shape of a regular spiking neuron potentials obtained from human neocortical neurons [19].
Figure 2 shows the limit cycle of a cortical neuron in Eq. ( 1) when the applied current, I, is above the bifurcation value, resulting in a saddle-node bifurcation.Due to a quadratic term in the recovery variable, the spike rate of a cortical neuron can be arbitrarily low, for low currents, and increases as I is increased.The above neurons can be coupled by adding an additional term to dV /dt in Eq. ( 1) proportional to g (V − E syn ), where g is the synaptic conductance variable, E syn = 0 for excitatory synapses and −0.92 for inhibitory.
The synaptic conductance, g, is obtained from the following equations, commonly used for synaptic coupling [12], where H step (x) = 1, if x > 0 and zero if x < 0, V pre is the voltage of the presynaptic neuron and τ syn is the synaptic conductance.In numerical simulation, Ω = −0.20 mV was chosen [12].The reason for using two synaptic equations is that, depending on τ syn , the conductance will peak after V pre , continuing to depolarize the membrane after the end of the presynaptic spike (see Figure 3).This type of response is consistent with physiological data.For brief stimulus spike at t = 0, the two synaptic equations produce a response proportional to [13] Figure 3 shows g for different values of τ syn .In each case, the area under the curve is the same (equal to one), with τ syn determining the width of the post-synaptic spike.
Using Eqs. ( 1) and ( 2) we can now set up the basic model of three mutually coupled neurons.The basic set-up is shown in Figure 1, where 3 neurons are coupled in a line.The middle neuron sends an excitatory time-delayed signal to the two outer neurons which are coupled to the middle one via weaker delayed coupling with a longer synaptic time constant.
This model reflects the typically observed hierarchy in the brain (described in the introduction), where the neurons in one brain area, such as the thalamus, send strong excitatory connections to cortical neurons, receiving weak modulatory feedback.Since the equations of the outer neurons are identical, receiving the same synaptic input from the middle cell, their synchronization would indicate that a group of identical neurons fire synchronously when subjected to a particular synaptic input.Using Eq. ( 1), the equations for the coupled neurons are given by: where i is the index of each neuron, with outer neurons receiving the same current, I 1,3 ≡ I s , and the same synaptic strength, δ 1 = δ 3 ≡ δ.The current input to the middle neuron is higher than to the outer ones: I 2 > I s , leading to the lower uncoupled spiking rate for the outer cells.This was done so that the higher activity of the inner cell drives the outer ones, as might be the case in a typical hierchical network.The synaptic input to the middle neuron is weaker, δ 2 < δ.The synaptic time-constant in Eq. ( 2) is the same for the outer neurons, τ syn1 = τ syn3 ≡ τ syn , but longer for the middle neuron, to model the affect of a slowly-varying modulatory feedback.The conductances, g i , are obtained from Eq. ( 2), with presynaptic voltage to the middle cell given by: and to the outer cells by The delay in the propagation of the signal is given by τ d .

The effect of delays on correlation and synchronization
While the length of the synaptic delay, τ d , does not seem to affect the degree of synchronization, it has a substantial affect on correlations and phase relations between neurons.
Figure 4 shows correlations in spiking output between the middle and outer neuron, and between outer neurons for τ d = 10.The x-axis indicates the time-shift at which the correlations function was computed.The outer neurons are synchronized, since C 13 = 1 at t = 0.It can be seen that the greatest correlations between the inner and outer neurons occur when t is shifted by the delay time, τ d .The delay also creates spikes in correlation at intervals of twice the delay time, as can be seen in Fig. 4. Since, in this case, the outer neurons are synchronized, the spikes indicate that there are self-correlations in the time series of a single neuron at intervals of 2τ d .This is the round-trip, or feed-back, time, since its the minimum time that it would take for a signal to travel from one of the neurons, affect the target, and get back to that same neuron.It follows that time delays lead to self-correlations in the spike trains that are not observed when delays are absent, and which may lead to more regular patterns in the time-series data.
For reasons explained in the introduction, δ > δ 2 and τ syn < τ syn2 were used to model the affect of a stronger forward and a weaker modulatory feedback.For this type of coupling, the middle neuron leads the outer ones by the delay time, τ d .This type of phase-locking behavior has been observed in other time-delay systems, such as lasers [18,20].Thus the time-shift between two correlated spike trains should be directly related to the delays in transmission, with the input leading the output by the delay time.
The effect of the coupling strength on firing rate and synchronization Increasing the synaptic coupling in general increases the firing rate of the neurons for excitatory synapses.Above a certain value of the coupling strength, there is phase-locking between the inner and the outer neurons, whereby all firing rates are equal.The bifurcation value of synaptic coupling at which 1 : 1 frequency locking occurs depends on the current injected into each cell.In the absence of coupling, a significant difference in injected current for the outer and inner neurons, (I s = 0.22 and I 2 = 0.5 in simulations) leads to a big difference in firing rate.Figure 5 shows a typical voltage trace for the uncoupled, δ = 0, and coupled neurons, δ = 4, where coupling is sufficiently strong to cause phase-locking.
Phase-locking or 1 : 1 frequency locking occurs for δ > 3.4.This type of behavior has been observed for mutually coupled neurons in the absence of delays [12].For lesser values of the coupling strength, different frequency locked behaviors are observed, with frequency ratio between 1 and f i /f o , where f i and f o are the frequencies of inner and outer neurons, respectively, in the absence of synaptic coupling.Figure 6 shows 1 : 4 frequency locking that occurs for weak synaptic coupling: δ = 0.5, δ 2 = 0.3.When the coupling is too weak, the outer neurons become desynchronized, as shown in Figure 7.The bifurcation value of δ, however depends on the synaptic time constant, τ syn , of the outer-neurons.Thus Figs. 6 and   8, which have a very short time constant of τ syn = 0.03, shows synchronization at a lower coupling strength of δ = 0.5, compared to a minimum δ = 1.03 needed for synchronization of neurons in Fig. 7, where τ syn = 0.5 (a more realistic value for fast synapses).This sensitive dependence of synchronization on the coupling strength and synaptic time constant of the outer-neurons is explored analytically in Section III.
The effect of synaptic time-constant on synchronization and firing rate Figure 9 shows an increase in firing rate of the outer neurons as their synaptic time constant, τ syn , is decreased.This increase in firing rate may be surprising, since the time constant only controls the width of the conductance spike and not the area under the curve, as shown in Fig. 3. Thus the contribution of a pre-synaptic spike to a change in postsynaptic voltage during an inter-spike interval is largely independent of the synaptic time constant.This can be seen by integrating the g (V i − E syn ) term in Eq. ( 4) over the interval of conductance change.Figure 10 shows a fluctuation in conductance for a system given by Eqs. ( 2) and ( 4), with τ syn = 1 and τ syn2 = 3.A lower firing rate that occurs for slower synapses may be the result of the decay of any increase in voltage during the inter-spike interval back to the limit cycle trajectory (see the bottom of Fig. 8).A pre-synaptic spike from the inner neuron can trigger a spike from the outer one, when it is delivered toward the end of the inter-spike interval.Thus a more narrow jump in conductance and the resultant jump in voltage, V 1,3 , may mean that there is less decay before the critical threshold is reached, thereby increasing the likelihood of a spike.
The most noticeable affect of the decrease in τ is greater synchronization.This is shown in Fig. 11, where the outer neurons become progressively synchronized as τ syn is decreased from 0.5 to 0.2.It can be seen that for relatively weak coupling of δ = 1, δ 2 = 0, and τ syn = 0.2, the outer neurons are completely synchronized, after the transients die out.The next section analyzes the affect of a presynaptic spike on a fast synapse (small τ syn ) in synchronizing two nearby trajectories of the outer neurons.

SYNAPSES
The three neuron model described by Eqs. ( 2) and (4) possesses internal symmetry, e.g. its equations of motion do not change if the variables {V 1 , R 1 } and {V 3 , R 3 } are interchanged.
It follows that the synchronized regime, where the symmetric variables are exactly equal: [1].For an uncoupled system, this solution would not be a stable one, since any perturbation along the limit cycle would result in a phase-difference.A spiking input, V 2 , from the middle neuron stabilized the synchronous state when the synaptic time-constant, τ syn is sufficiently short.To study the affect of a single pre-synaptic spike on two identical neurons, with nearby trajectories along the limit cycle, introduce new variables: Ṽ = V 1 − V 3 and R = R 1 − R 3 .These new variables correspond to a perturbation transverse to the synchronized state: {V 0 , R 0 }.Using Eqs.
(1), ( 3) and ( 4), the linearized dynamics are: where n (V 0 ) = 26 (V 0 + 0.95), h (V 0 , R 0 ) = (−5.03+ 101.4V 2 0 − 32.45V 0 + 26R 0 ), and m (V 0 , R 0 ) = − (3.8 + 6.6V 0 ).t s denotes the time since the last spike of the center neuron.Equation ( 5) is valid for a short synaptic time constant, when only the effect of the last presynaptic spike is significant (see Fig. 3).During an inter-spike interval of the outer cells, n, h and m are all positive.This can be easily confirmed by using Fig. 2, where −0.8 < V < −0.6 and .18< R < 0.4 during the inter-spike interval, and calculating the lowest possible values of n, h and m for a given range of V 0 and R 0 .Since { Ṽ , R} denote the difference between the two nearby trajectories along the limit cycle, there is a relationship between the two variables given by R(t) = −l(V 0 ) Ṽ (t) (7) where l is the negative of the slope of the limit cycle at {V 0 (t), R 0 (t)}.From Fig. 2, l(V 0 ) ≈ 2 during the inter-spike interval.
A change in Ṽ following a post-synaptic spike can now be calculated by substituting Eq. ( 7) into Eq.( 5), dividing by Ṽ and integrating over △t, the duration of the conductance spike, t s /τ 2 syn exp (−t s /τ syn ).For a short synaptic time constant, τ syn ≤ 1, the width of △t is small compared to the time-scale of neuronal dynamics during the inter-spike interval of the outer cells.A change in Ṽ following a post-synaptic spike is now given by, where Ṽu (t + △t) is the voltage difference of two nearby trajectories at time t + △t in the absence of synaptic input: Ṽu (t + △t) ≈ Ṽ (t) exp (△t (nl − h)), with n and h defined after Eq. ( 6) and l given in Eq. ( 7).For small △t, so that |△t (nl − h) | ≪ {1, δ}, most of the change in voltage during a post-synaptic spike is caused by the spike itself so that Eq. ( 8) can be approximated as Using Eqs. ( 6)-( 8), the dynamics of R immediately following a narrow post-synaptic spike can be approximated as From Eqs. ( 8) and ( 9), a single presynaptic spike from a center neuron acts to decrease the voltage difference, Ṽ , of the outer neurons by a factor of exp (−δ).This presynaptic spike also decreases R by decreasing the positive contribution from the m/l > 0 term in Eq. (10).
Thus the effect of the presynaptic input is to decrease the perturbation from synchronized state of the outer neurons, pushing their trajectories closer in phase-space.
Eqs. ( 8)- (10) show that a spiking input from the middle neuron has a stabilizing affect on the synchronized state when the synaptic time constant is short.Since the difference in trajectories, Ṽ , R is taken along a limit cycle, the maximum Lyapunov exponent in the absence of synaptic coupling would be zero (corresponding to the displacement along a trajectory in phase-space), and the transverse exponents must be negative since the systems is dissipative.It follows that, for sufficiently small τ syn , a common synaptic input acting during the inter-spike interval should eventually synchronize the neurons, even for weak synaptic coupling.
Figure 8 shows the synchronization of the outer neurons for a very short synaptic time constant, τ syn = 0.03 and weak coupling, δ = 0.5.After the transients die out, the outer neurons become synchronized, thereby falling on a straight line in the V 1 vs V 3 plot.As can be seen in Figure 6, for weak coupling, there is a big difference in firing rate between the outer and the middle neuron, due to differences in injected current.It follows that the synchronization is not due to phase-locking between the middle and the outer cells.The effect of the spiking input on the limit cycle trajectory can be seen at the bottom of Figure 8).There is an integer ratio between the inner and the outer frequencies, whereby the onset of a spike in the outer neuron is triggered by spiking input from the middle one, delivered toward the end of the inter-spike interval.The phenomena is similar to the subharmonic resonance where the limit cycle responds at a subharmonic of the stimulus frequency [12].
Eqs. ( 8)- (10) show a strong dependence of synchronization on the coupling strength.This is confirmed by numerics.Figure 7 shows V 1 vs V 3 for two slightly different values of δ, δ = 1.02 and δ = 1.03.A slight change in δ leads to an onset of synchronization between the outer neurons.

IV. CONCLUSION
A model of three mutually coupled cortical neurons with delays was studied using analysis and numerical simulation.The outer neurons were stimulated with smaller current and had a much lower frequency in the uncoupled case, with their frequency significantly increasing depending on the strength of synaptic coupling with the middle neuron.At higher values of the synaptic coupling constant, typical phase-locked behavior and 1 : 1 frequency locking was found between the middle and the outer neurons, with different frequency locking ratios as the synaptic strength was lowered.It was found that delays affected the time-series data by introducing correlations at the time-scale of the delay.While the spiking behavior in the synchronized case was fairly regular, this effect would be interesting to explore for a more complicated, chaotic spike train that can be achieved by incorporating slow adaptation currents into the neuron model.In the case of phase or frequency locking, the middle neuron leads the outer by the delay time, τ d .
While synchronization of outer neurons was sensitive to the synaptic strength, the synaptic time constant for outer neurons, τ syn was also highly significant.It was found that shorter synaptic constant substantially improves synchronization, able to synchronize the neurons even when the coupling strength is very weak.A short synaptic constant was also able to induce 1 : 1 frequency locking at a much weaker mutual coupling than this type of dynamics might be expected, for the input currents used.Analysis of dynamics for fast synapses showed that spiking synaptic input during the inter-spike interval stabilized the synchronization manifold, even for arbitrarily weak coupling, and independent of the phase relationship between the inner and outer cells.This indicates that even a very weak synaptic input can synchronize cells, as long as the synaptic time constant is sufficiently short.The finding may have significance in synchronizing large groups of cells in the cortex via weak synaptic input from other areas, such as the thalamus, or other areas in the cortex proper.

FIG. 1 :FIG. 2 :FIG. 3 :
FIG. 1: The basic model of three synaptically coupled neurons.The strength of the coupling from the middle to the outer neurons is given by δ and from the outer to the middle by δ 2 .The injected currents are shown.All synaptic coupling has propagation delay of τ d and a synaptic time constant of τ syn and τ syn2 for the outer and the middle neurons, respectively.