Adaptation controls synchrony and cluster states of coupled threshold-model neurons

We analyze zero-lag and cluster synchrony of delay-coupled nonsmooth dynamical systems by extending the master stability approach, and apply this to networks of adaptive threshold-model neurons. For a homogeneous population of excitatory and inhibitory neurons we find (i) that subthreshold adaptation stabilizes or destabilizes synchrony depending on whether the recurrent synaptic excitatory or inhibitory couplings dominate, and (ii) that synchrony is always unstable for networks with balanced recurrent synaptic inputs. If couplings are not too strong, synchronization properties are similar for very different coupling topologies, i.e., random connections or spatial networks with localized connectivity. We generalize our approach for two subpopulations of neurons with nonidentical local dynamics, including bursting, for which activity-based adaptation controls the stability of cluster states, independent of a specific coupling topology.

Figure 1

Stability of synchrony for a homogenous population of coupled aEIF neurons. The maximum Lyapunov exponent (color bar) is plotted as a function of the real ($\alpha$) and imaginary ($\beta$) parts of the eigenvalues of the coupling matrix $\mathbf{C}$. Different panels show the results for strong (top) and weak (bottom) subthreshold adaptation, and for the excitation-dominated (left), balanced (center), and inhibition-dominated (right) regimes. Parameters were selected such that the dynamics of the aEIF model closely matches the regular spiking dynamics of cortical neurons [29, 37]: $\lambda =5$, $\tau =0.3$, ${\tau }_w=10$, ${\tau }_s=0.3$, $V_{th}=5$, $V_r=-5$, $b=2.5$. The external input $I$ was chosen such that the oscillation period was always equal to 5. Circles indicate the unit disk in the eigenvalue plane (insets show blowup).

Figure 2

Eigenvalue spectra of coupling matrices for two biologically relevant connectivity schemes. (a) Schematic diagram of the network: two subpopulations of $N_{\mathcal{E}}$ excitatory (light orange) and $N_{\mathcal{I}}$ inhibitory (dark green) neurons. (b) Circle radius $r$ computed from Eq. (13) for the random connectivity scheme. Lines of equal radius ($r=1$) are plotted in the $(p,c_{\mathcal{E}})$ plane for $\overline{c}=1,0,-1$ (solid green, dashed violet, dotted orange lines). Parameters: $N_{\mathcal{E}}=400$, $N_{\mathcal{I}}=100$. (c) Largest absolute value $r$ of all transverse eigenvalues, averaged over 100 sample matrices each representing a Mexican hat (left) or inverse Mexican hat (right) spatial map connectivity scheme: Excitatory and inhibitory neurons lie uniformly distributed on a two-dimensional sheet with unity edge length and periodic boundary conditions. The probability densities for excitatory and inhibitory connections are normalized Gaussian functions of the spatial distance between neurons with standard deviations ${\sigma }_{\mathcal{E}}$ and ${\sigma }_{\mathcal{I}}$, respectively. Lines of equal absolute value ($r=1$) are plotted with color code as in (b). Parameters: $N_{\mathcal{E}}=400$, $N_{\mathcal{I}}=100$, $p=0.2$; ${\sigma }_{\mathcal{E}}=0.07$, ${\sigma }_{\mathcal{I}}\in [0.07,0.35]$ (left) and ${\sigma }_{\mathcal{E}}\in [0.07,0.35]$, ${\sigma }_{\mathcal{I}}=0.07$ (right). $c_{\mathcal{I}}$ in (b) and (c) via constant row sum condition.

Figure 3

Stability of cluster states in networks with two groups of excitatory, $\mathcal{E}$, and inhibitory, $\mathcal{I}$, neurons which differ in their local dynamics. (a) Schematic diagram of the network, according to Eqs. (14, 15). (b) Evolution of the membrane voltage $V$ of synchronized excitatory (light orange) and inhibitory neurons (dark green) in the network for three different parameter sets: excitatory neurons: $V_r=-5$, $I=3.725$ (left and center), $V_r=2$, $I=25$ (right); inhibitory neurons: $I=0$ (left and right), $I=1.5$ (center). ${\lambda }^{\mathcal{E}\mathcal{E}}={\lambda }^{\mathcal{I}\mathcal{E}}=5$, ${\lambda }^{\mathcal{E}\mathcal{I}}={\lambda }^{\mathcal{I}\mathcal{I}}=-5$ (left and center); ${\lambda }^{\mathcal{E}\mathcal{E}}={\lambda }^{\mathcal{I}\mathcal{E}}=25$, ${\lambda }^{\mathcal{E}\mathcal{I}}={\lambda }^{\mathcal{I}\mathcal{I}}=-25$ (right). Other parameters: $a=0.5$, $b=2.5$, ${\tau }_w=10$, ${\tau }_s=0.2$, $V_{th}=5$ for the excitatory, and $a=0.1$, $b=0.25$, ${\tau }_w=10$, ${\tau }_s=0.5$, $V_r=-5$, $V_{th}=5$ for the inhibitory neurons; $\tau =0.1$. (c) MSF for each of the three scenarios in (b). The maximum Lyapunov exponent is shown as a function of the real eigenvalues ${\gamma }_1,{\gamma }_2$. White rectangles indicate the unit square in the eigenvalue plane. (d) MSF of the unit square as in (c), but with $a=0.1,0.2,0.3,0.4,0.5$ (top, left to right) and $b=0.5,1,1.5,2,2.5$ (bottom) for the excitatory neurons.

Figure 4

Approximation scheme for the perturbation ${\mathit{\xi }}_i^{+}\left(t_s\right)$ just after the discontinuity of ${\mathbf{x}}_s$, illustrated for scalar-valued functions ${\mathbf{x}}_i$ ($m=1$). Evolution of ${\tilde{\mathbf{x}}}_i\equiv {\mathbf{x}}_s+{\mathit{\xi }}_i$ (dashed dark blue lines), where ${\mathbf{x}}_s$ (solid gray lines) is a component of the synchronous solution of the system Eqs. (3) and (4) and ${\mathit{\xi }}_i$ is a solution component of Eq. (5) with initial condition $\underline{\mathit{\xi }}\equiv {\underline{\mathit{\xi }}}_0$, for $d{t}_i<0$ (left) and $d{t}_i>0$ (right). Note that ${\tilde{\mathbf{x}}}_i$ is not discontinuous at $t_i\ne t_s$, because it is not a solution component of Eqs. (3) and (4). The $i$th solution component of Eqs. (3) and (4) with initial trajectory ${\underline{\mathbf{x}}}_s+{\underline{\mathit{\xi }}}_0$ is indicated (solid light blue lines). Dashed gray lines mark the set $\left\{\mathbf{x}\right|\phi \left(\mathbf{x}\right)=0\}$ which constitutes the threshold hyperplane here. If a solution $\mathbf{x}$ of Eq. (3) reaches this hyperplane from below, a discontinuity occurs: ${\mathbf{x}}^{+}=\mathbf{g}\left(\mathbf{x}\right)$, cf. Eq. (4). Solid red arrows indicate the vectors $\mathbf{F}\left(t_s\right)d{t}_i$ and ${\mathbf{F}}^{+}\left(t_s\right)d{t}_i$ used in the approximation of ${\mathit{\xi }}_i^{+}\left(t_s\right)$. The vectors $\mathbf{F}\left(t_s\right)$ and ${\mathbf{F}}^{+}\left(t_s\right)$ are indicated by dotted red arrows.

Figure 5

Approximation scheme for the perturbation ${\mathit{\xi }}_i^{+}(t_s+\tau )$ just after the kink in ${\mathbf{x}}_s$, illustrated for scalar-valued functions ${\mathbf{x}}_i$ ($m=1$). Evolution of the synchronous solution component ${\mathbf{x}}_s$ (solid gray line), ${\tilde{\mathbf{x}}}_i\equiv {\mathbf{x}}_s+{\mathit{\xi }}_i$ (dashed dark blue lines) and the $i$th solution component of Eqs. (3) and (4) with initial trajectory ${\underline{\mathbf{x}}}_s+{\underline{\mathit{\xi }}}_0$ (solid light blue lines). Solid red arrows indicate vectors that are used in the approximation, as explained in the text. The solid green arrow represents the difference between ${\mathit{\xi }}_i^{+}(t_s+\tau )$ (dashed green arrow) and ${\mathit{\xi }}_i(t_s+\tau )$ (solid dark blue arrow).