1 Background
============

There is a significant amount of research showing that spike coincidence of neurons encodes information. Singer and colleagues showed not only that networks of neighboring neurons tend to become synchronized, but that tight synchronization occurs at the opposite hemispheres of the brain, in the presence of significant (6 to 8 milliseconds propagation delays) \[[@B1]\]. Zero time-lag synchronization of neural activity has been involved in such important phenomena as ability to recognize objects (by binding different attributes) \[[@B2],[@B3]\], olfactory discrimination \[[@B4]\], and has even been proposed as one of the neural correlates of consciousness \[[@B5]\]. It was found that visuomotor integration is associated with synchronization of signals recorded from the visual and parietal, and motor and parietal areas of the awake cat \[[@B6],[@B7]\].

While much of analysis focuses on network synchrony in the absence of time delays, delays are common in neural networks and there is evidence that under certain parameters \[[@B8]\] or coupling architectures, delays may actually contribute to synchronization. In a recent issue of Science \[[@B9]\], a short article describes a recent laser experiment \[[@B10]\] of mutually coupled lasers in a row. It is noted that when only two lasers are coupled, there is a lag in their phases equal to the amount of time it takes light to pass between them. However, when a third laser is added, the outputs of the outer lasers show zero time-lag synchronization. It has been suggested by Wolf Singer that this phenomena could shed light on synchronization of nerve signals in the brain in the presence of delays \[[@B9]\]. In particular, researchers indicate that this zero-lag synchronization effect in three coupled oscillators could clarify how the hemispheres of the brain synchronize \[[@B9]\].

With this in mind, we propose a three neuron coupling scheme, analogous to the three laser experiment described above. The synchronization phenomena of the three laser system has been recently analyzed \[[@B11]\], elucidating the conditions under which synchronization occurs. It was shown analytically and confirmed with numerics that increasing the delays improves synchronization of the outer lasers \[[@B11]\]. As will be shown in the present paper the dependence of synchronization on parameters is somewhat different in the three coupled neurons case, where it is actually the synaptic time constant, rather than delays in propagation that are important for synchronization of the outer neurons. It is perhaps not surprising that the mechanism behind synchronization in this three neuron model is different, since dynamics at each synapse are modeled by two differential equations, rather than much simpler linear coupling of the three laser scheme.

While the basic coupling architecture we choose to consider is the same as the above described laser experiment, in that we have a symmetric system of three mutually delay coupled oscillators, the projections from the middle to the outer neurons are somewhat stronger (with faster synaptic time constant) than the coupling from the outer to the middle neuron. This difference in coupling strength and synaptic time constant was inspired by the hierarchical structure of many networks in the brain, which to a first approximation is a feed-forward network, modulated by a weaker modulatory feedback \[[@B12]\]. 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 \[[@B3]\]. It has been suggested that exessively strong mutually coupled loops in the brain would promote uncontrollable oscillations, such as in epilepsy \[[@B12],[@B13]\].

The present paper investigates a simple model for 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 sending 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 2, the basic model is set up, (see Fig. [1](#F1){ref-type="fig"}) and presented along with numerical results that show the dependence of synchronization and firing rate on the strength of coupling and synaptic conductance. Section 3, 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 4 concludes and summarizes.

![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 *τ*~*syn*2~for the outer and the middle neurons, respectively.](1753-4631-1-2-1){#F1}

2 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: \[[@B14]\]

$$\begin{matrix}
{\frac{dV}{dt} = - \left\{ 17.81 + 47.58V + 33.8V^{2} \right\}(V - 0.48) - 26R(V + 0.95) + I = F(V,R) + I} \\
{\frac{dR}{dt} = \frac{1}{\tau_{R}}\left( {- R + 1.29V + 0.79 + 3.3{(V + 0.38)}^{2}} \right) = G(V,R),} \\
\end{matrix}$$

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 *N a*^+^and *K*^+^currents, and a quadratic terms in *V*to approximate the transient potassium current contributions \[[@B15]\]. 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 \[[@B16]\]. Figure [2](#F2){ref-type="fig"} 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 \[[@B14]\],

![Limit cycle of an uncoupled cortical neuron, given by Eq. (1). *I*= 0.5 nA.](1753-4631-1-2-2){#F2}

$$\begin{matrix}
{\frac{df}{dt} = \frac{1}{\tau_{syn}}( - f + H_{step}(V_{pre} - \Omega))} \\
{\frac{dg}{dt} = \frac{1}{\tau_{syn}}( - g + f)} \\
\end{matrix}$$

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 time constant. In numerical simulation, Ω = -0.20 mV was chosen \[[@B14]\]. 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](#F3){ref-type="fig"}). 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 \[[@B13]\]

![Conductance, *g*, after a single pre-synaptic spike. *τ*~*syn*~= 1, 3, and 5 ms.](1753-4631-1-2-3){#F3}

$$g = \left( {t/\tau_{syn}^{2}} \right)exp\left( {- t/\tau_{syn}} \right)$$

Figure [3](#F3){ref-type="fig"} 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](#F1){ref-type="fig"}, 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 reflect 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 weaker modulatory feedback. Since the equations of the outer neurons are identical (resulting from a fit with experimental data) and they receive the same synaptic input from the middle cell, their synchronization would indicate that a group of cortical neurons with the same type of synaptic coupling will fire synchronously when subjected to a particular synaptic input. Using Eq. (1), the equations for the coupled neurons are given by:

$$\begin{matrix}
{\frac{dV_{i}}{dt} = F(V_{i},R_{i}) + I_{i} - \delta_{i}g_{i}(V_{i} - E_{syn})} \\
{\frac{dR_{i}}{dt} = g(V_{i},R_{i})} \\
\end{matrix}$$

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 hierarchical network \[[@B3]\]. The synaptic input to the middle neuron is weaker, *δ*~2~\<*δ*. The synaptic time-constant in Eq. (2) is the same for the outer neurons, *τ*~*syn*1~= *τ*~*syn*3~≡ *τ*~*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: *V*~*pre*2~= *V*~1~(*t*- *τ*~*d*~) + *V*~3~(*t*- *τ*~*d*~), and to the outer cells by *V*~*pre*1,3~= *V*~2~(*t*- *τ*~*d*~). The delay in the propagation of the signal is given by *τ*~*d*~.

2.1 The effect of delays on correlation and synchronization
-----------------------------------------------------------

While the length of the synaptic delay, *τ*~*d*~, does not seem to effect the degree of synchronization, it has a substantial effect on correlations and phase relations between neurons. Figure [4](#F4){ref-type="fig"} shows correlations in spiking output between the middle and outer neuron, and between outer neurons for *τ*~*d*~= 10. This delay time is similar, for instance, to the propagation delay between the two hemispheres of the brain \[[@B1]\]. 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](#F4){ref-type="fig"}. 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.

![Correlations of spiking output. Top: between the outer and the middle neuron. Bottom: Between outer neurons. *τ*~*d*~= 10 ms, *δ*= 4, *δ*~2~= 2, *τ*~*syn*~= 1 ms, *τ*~*syn*2~= 2 ms, *I*~1,3~= *I*~*s*~= 0.22 nA, *I*~2~= 0.5 nA.](1753-4631-1-2-4){#F4}

For reasons explained in the introduction, *δ*\> *δ*~2~and *τ*~*syn*~\<*τ*~*syn*2~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 \[[@B17],[@B18]\]. 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.

2.2 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](#F5){ref-type="fig"} 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 \[[@B14]\]. 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](#F6){ref-type="fig"} 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](#F7){ref-type="fig"}. The bifurcation value of *δ*, however depends on the synaptic time constant, *τ*~*syn*~, of the outer-neurons. Thus Figs. [6](#F6){ref-type="fig"} and [8](#F8){ref-type="fig"}, which have a very short time constant of *τ*~*syn*~= 0.03, show synchronization at a lower coupling strength of *δ*= 0.5, compared to a minimum *δ*= 1.03 needed for synchronization of neurons in Fig. [7](#F7){ref-type="fig"}, 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.

![Left: spiking in the absence of coupling. Right: coupling, *δ*= 4, *δ*~2~= 2, *τ*~*syn*~= 1, *τ*~*syn*2~= 2. In both cases, *I*~*s*~= 0.22, *I*~2~= 0.5.](1753-4631-1-2-5){#F5}

![Spiking voltage for weak coupling, *δ*= 0.5, *τ*~*syn*~= 0.03, *δ*~2~= 0.3. There is 1 : 4 frequency locking between the outer and the middle neurons.](1753-4631-1-2-6){#F6}

![Sensitive dependence of synchronization on coupling strength. Top: *δ*= 1.02, Bottom: *δ*= 1.03, *τ*~*syn*~= 0.5, *δ*~2~= 0.3.](1753-4631-1-2-7){#F7}

![Synchronization for weakly coupled, fast synapses. *δ*= 0.5, *τ*~*syn*~= 0.03, *δ*~2~= 0.3. Top: synchronization after the transients die out. Bottom: Limit cycle of one of the outer neurons. Same parameters as in Fig. 6.](1753-4631-1-2-8){#F8}

2.3 The effect of synaptic time-constant on synchronization and firing rate
---------------------------------------------------------------------------

Figure [9](#F9){ref-type="fig"} 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](#F3){ref-type="fig"}. Thus the contribution of a pre-synaptic spike to a change in post-synaptic 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](#F10){ref-type="fig"} shows a fluctuation in conductance for a system given by Eqs. (2) and (4), with *τ*~*syn*~= 1 and *τ*~*syn*2~= 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](#F8){ref-type="fig"}). 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.

![Increase in firing rate as *τ*~*syn*~decreases. Left: *τ*~*syn*~= 2. Right: *τ*~*syn*~= 0.1. In both cases, *δ*= 1, *δ*~2~= 0.2, *τ*~*syn*2~= 4, *I*~*s*~= 0.22, *I*~2~= 0.5.](1753-4631-1-2-9){#F9}

![Conductance for *τ*~*syn*~= 1, *τ*~*syn*2~= 3. All else as in Fig. 4.](1753-4631-1-2-10){#F10}

The most noticeable affect of the decrease in *τ*is greater synchronization. This is shown in Fig. [11](#F11){ref-type="fig"}, 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.

![Dependence of synchronization on the synaptic time constant, *τ*~*syn*~. Top: *τ*~*syn*~= 0.5, Middle: *τ*~*syn*~= 0.3, Bottom: *τ*~*syn*~= 0.2, In all three cases, *δ*= 1, *δ*~2~= 0.](1753-4631-1-2-11){#F11}

3 Analysis of synchronization for weakly coupled, fast 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: {*V*~1~= *V*~3~= *V*~0~, *R*~1~= *R*~3~= *R*~0~}, is a solution \[[@B19]\]. 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: $\widetilde{V}$ = *V*~1~- *V*~3~and $\widetilde{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:

$$\frac{d\widetilde{V}}{dt} = - n(V_{0})\widetilde{R} - h(V_{0},R_{0})\widetilde{V} - \delta\left( {t_{s}/\tau_{syn}^{2}} \right)exp\left( {- t_{s}/\tau_{syn}} \right)\widetilde{V}$$

$$\frac{d\widetilde{R}}{dt} = \frac{1}{5.6}\left( {- \widetilde{R} - m(V_{0})\widetilde{V}} \right)$$

where *n*(*V*~0~) = 26 (*V*~0~+ 0.95), *h*(*V*~0~, *R*~0~) = (-5.03 + 101.4$V_{0}^{2}$ - 32.45*V*~0~+ 26*R*~0~), and *m*(*V*~0~) = -(3.8 + 6.6*V*~0~). *t*~*s*~denotes the time since the arrival of the last spike from the center to outer neurons. Equation (5) is valid for a short synaptic time constant, when only the effect of the last presynaptic spike is significant (see Fig. [3](#F3){ref-type="fig"}). 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](#F2){ref-type="fig"}, 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 {$\widetilde{V}$, $\widetilde{R}$} denote the difference between the two nearby trajectories along the limit cycle, there is a relationship between the two variables given by

$$\widetilde{R}(t) = - l(V_{0})\widetilde{V}(t)$$

where *l*is the negative of the slope of the limit cycle at {*V*~0~(*t*), *R*~0~(*t*)}. From Fig. [2](#F2){ref-type="fig"}, *l*(*V*~0~) ≈ 2 during the inter-spike interval. The above equation is valid as long as two nearby trajectories remain on the limit cycle.

The dynamics of $\widetilde{V}$ along a limit cycle can now be calculated by substituting Eq. (7) into Eq. (5), dividing by $\widetilde{V}$, bringing *dt*over to the right-hand side of the equation and integrating. For a small synaptic time constant, *τ*~*syn*~≤ 1, the width of a duration of a post-synaptic spike, Δ*t*, is short compared to the time-scale of neuronal dynamics of outer cells during the inter-spike interval. Here, Δ*t*measures the duration of a conductance spike, (*t*~*s*~/$\tau_{syn}^{2}$) *exp*(-*t*~*s*~/*τ*~*syn*~), which is quite narrow for sufficiently fast synapses (see Fig. [3](#F3){ref-type="fig"}). Thus most of the change in $\widetilde{V}$ during the narrow post-synaptic spike is due to the spike itself. This can be seen by considering a voltage difference of two nearby trajectories at time *t*+ Δ*t*in the **absence**of synaptic input: ${\widetilde{V}}_{u}(t + \Delta t) \approx \widetilde{V}(t)\exp(\Delta t(nl - h))$, with *n*and *h*defined after Eq. (6) and *l*given in Eq. (7). This expression was obtained by substituting Eq. (7) into Eq. (5) and integrating over Δ*t*while assuming that *V*~0~stays almost constant over Δ*t*. Since Δ*t*is small, we have \|Δ*t*(*nl*- *h*)\| ≪ {1, *δ*}. So that only the conductance term in Eq. (5) makes a significant contribution to $\widetilde{V}$ during the duration, Δ*t*, of a post-synaptic spike. Dividing Eq. (5) by $\widetilde{V}$ and integrating, we get

$$\widetilde{V}(t + \Delta t) \approx \widetilde{V}(t)e^{- \delta}$$

The above equation gives a change in $\widetilde{V}$ following a synaptic change in conductance. Using Eqs. (6)-(8), the dynamics of $\widetilde{R}$ immediately following a narrow post-synaptic spike can be approximated as

$$\frac{d\widetilde{R}}{dt} = \frac{1}{5.6}\left( {- 1 + \left( \frac{m(V_{0})}{l(V_{0})} \right)e^{- \delta}} \right)\widetilde{R}$$

From Eq. (8), a single presynaptic spike from a center neuron acts to decrease the voltage difference, $\widetilde{V}$, of the outer neurons by a factor of *exp*(-*δ*). This presynaptic spike also decreases $\widetilde{R}$ by decreasing the positive contribution from the *m*/*l*\> 0 term in Eq. (9). Thus the immediate effect of a presynaptic input is to decrease the perturbation from synchronized state of the outer neurons, pushing their trajectories closer in phase-space.

Eqs. (8) and (9) 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, $\widetilde{V}$, $\widetilde{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 trajectory collapses onto a limit cycle. 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](#F8){ref-type="fig"} 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](#F6){ref-type="fig"}, 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](#F8){ref-type="fig"}. 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 \[[@B14]\].

Eqs. (8) and (9) show a sensitive dependence of synchronization on the coupling strength from the middle to the outer neurons. This is confirmed by numerics. Figure [7](#F7){ref-type="fig"} 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.

4 Conclusion
============

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 \[[@B3]\]. 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 \[[@B3],[@B20]\]. The three neuron scheme investigated in this paper is a simple model for studying these types of hierarchical networks, since it incorporates this phenomena of synchronization of certain areas of the brain (outer neurons) due to common input from a different area (middle neuron), and since synchronization of the two outer neurons would indicate synchronization of many neurons, if coupled to the middle neuron in the same way as the two outer neurons in a three neuron model.

This 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 firing 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 of outer neurons, *τ*~*syn*~was also highly significant. It was found that shorter synaptic constant substantially improves correlations, leading to zero-lag synchronization of end neurons even when the coupling strength is very weak. A short synaptic constant was also able to significantly increase firing rate, to the point of inducing 1 : 1 frequency locking with the middle neuron, at a much weaker mutual coupling than would otherwise occur for the input currents used. Analysis of dynamics for fast synapses showed that fast 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 suficiently 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.

Acknowledgements
================

This work was supported by a grant from the Office of Naval Research. ASL is currently a post doctoral fellow with the National Research Council.
