Delay-Induced Multistability and Loop Formation in Neuronal Networks with Spike-Timing-Dependent Plasticity

Abstract

Spike-timing-dependent plasticity (STDP) adjusts synaptic strengths according to the precise timing of pre- and postsynaptic spike pairs. Theoretical and computational studies have revealed that STDP may contribute to the emergence of a variety of structural and dynamical states in plastic neuronal populations. In this manuscript, we show that by incorporating dendritic and axonal propagation delays in recurrent networks of oscillatory neurons, the asymptotic connectivity displays multistability, where different structures emerge depending on the initial distribution of the synaptic strengths. In particular, we show that the standard deviation of the initial distribution of synaptic weights, besides its mean, determines the main properties of the emergent structural connectivity such as the mean final synaptic weight, the number of two-neuron loops and the symmetry of the final structure. We also show that the firing rates of the neurons affect the evolution of the network, and a more symmetric configuration of the synapses emerges at higher firing rates. We justify the network results based on a two-neuron framework and show how the results translate to large recurrent networks.

Figure 1

Theoretical prediction of synaptic modification as a function of the initial synaptic strengths in the absence of a frequency mismatch. The colors show the stable fixed point of the time lag Δt^* = χ^*/(2πν) of spiking of the neurons given by Eq. (9) and the vector field shows the direction of the change of the synaptic strengths from Eq. (3). The yellow solid curves denote the simulated synaptic evolution for three different initial values. In each row, the firing frequency ν is constant, but it increases from top to bottom in each column. Based on the initial synaptic coupling, different attractors may be achieved: Bidirectional/unidirectional (left column), unidirectional (middle column), and decoupled/unidirectional (right column) state. The dendritic propagation delay is ${\tau }_{\mathrm{d}}=0.5\mathrm{ms}$ and the axonal delay is denoted above each panel. STDP parameters are A₊ = A₋ = 0.005 and τ₊ = τ₋ = 20 ms.

Figure 2

Time course of simulated synaptic strengths and spiking time lag. (A–C) Correspond to Fig. 1G with ν = 80 Hz and τ_a = 0.3 ms. (D–F) Correspond to Fig. 1I with ν = 80 Hz and τ_a = 1.0 ms. The initially given two-dimensional synaptic strength vector is (g₂₁(0), g₁₂(0)) = (0.6, 0.4), (0.2, 0.7), (0.8, 0.2), (0.7, 0.7), (0.7, 0.3), (0.2, 0.6) for (A–F), respectively. The finally achieved two-dimensional synaptic strength vector (g₂₁, g₁₂) is denoted in the figure.

Figure 3

Theoretical prediction of synaptic modification based on the initial synaptic strengths in the presence of a frequency mismatch. In the presence of a frequency mismatch Δν, the instantaneous fixed point of the time lag is numerically calculated by self-consistently solving Eq. (8) with the Newton-Raphson root-finding method. (A–C) Frequency mismatch Δν = 0.02 Hz with ν₂ > ν₁. (D–F) Frequency mismatch Δν = 0.04 Hz with ν₂ > ν₁. In the presence of frequency mismatch the output synapse of the neuron with higher frequency is more likely to be potentiated, while the reverse synapse is depressed.

Figure 4

Dependence of the two-neuron results on the relative synaptic strength. (Left panels) The interplay between fixed point of the time lag $\mathrm{\Delta }{t}^{⁎}$ and relative synaptic strength Γ determines the final stable connectivity pattern. Colors show the type of the final synaptic connectivity pattern: Bidirectional (red), unidirectional (orange), and decoupled (blue). The blue curve denotes Δt^* calculated from Eq. (9). (A) ν = 80 Hz and τ_a = 0.3 ms. (B) ν = 80 Hz and τ_a = 1.0 ms. The dendritic propagation delay is τ_d = 0.5 ms. (Right panels) The color-coded relative synaptic strength Γ is shown based on the synaptic strengths from Eq. (1). The vector field is the same as in Fig. 1.

Figure 5

Simulation results for a recurrent excitatory neuronal network in the high-frequency regime with heterogeneous initial synaptic strengths and ξ > 0. (A) (Left panel) Time course of the mean synaptic strength $\overline{g}$(t) for different standard deviations σ_g. (Right panel) Initial and final distribution of the synaptic strengths. (B) (Left panel) Time course of the order parameter r(t). Note that the colors indicated in the legend belong to the same σ_g in A and B. (Right panel) Time course of the normalized number of closed loops of length 2 (see Methods), representing the number of bidirectional connections in the network. (C) Final coupling matrices for σ_g = 0.05,0.08,0.10,0.15, respectively. Network and STDP parameters are N = 200, $\nu =80\mathrm{H}\mathrm{z}$, τ_d = 0.5 ms, τ_a = 0.3 ms, A₊ = A₋ = 0.005, and τ₊ = τ₋ = 20 ms.

Figure 6

Simulation results for a recurrent excitatory neuronal network in the high-frequency regime with heterogeneous initial synaptic strengths and ξ > 0. (A) (Left panel) Time course of the mean synaptic strength $\overline{g}$(t) for different standard deviations σ_g. (Right panel) Initial and final distribution of the synaptic strengths. (B) (Left panel) Time course of the order parameter r(t). (Right panel) Time course of the normalized number of closed loops of length 2, representing the number of bidirectional connections in the network. (C) Final coupling matrices for σ_g = 0.05,0.08,0.10,0.15, respectively. Network and STDP parameters are N = 200, ν = 80 Hz, ${\tau }_{\mathrm{d}}=0.5\mathrm{ms}$, τ_a = 1.0 ms, A₊ = A₋ = 0.005, and τ₊ = τ₋ = 20 ms.

Figure 7

Illustration of the accordance between two-neuron theoretical framework and network simulation results. (A–D) The probability P_asym as a measure of asymmetry in the two-neuron motif (solid) and the asymmetry index $C_{\mathrm{net}}$ for the neuronal network (dotted curves) are calculated in the presence of different axonal propagation delays with τ_a = 0.3 ms (red) and τ_a = 1.0 ms (blue curves). (A) The effect of the firing frequency on P_asym (C_net) in the two-neuron motif (network simulation) with inhomogeneity in the initial distribution of the synaptic strengths represented by σ_g. In the figure σ_g = 0.05. (B) Same as A, but in the presence of inhomogeneity in the firing frequencies indicated by Δν (${\sigma }_{\nu }$) in the two-neuron motif (network simulation). In the figure Δν = σ_ν = 0.01 Hz. (C) The effect of inhomogeneity in the frequencies Δν (σ_ν) on P_asym (C_net). In the figure ν = 80 Hz and σ_g = 0.05. (D) Same as C, but for inhomogeneity in the initial distribution of the synaptic strengths σ_g. In the figure ν = 80 Hz and Δν = σ_ν = 0.01 Hz. (a₁–d₁) Samples of final coupling matrix indexed by the number of pre- (j) and postsynaptic (i) neurons correspond to a₁–d₁ markers in A–D, representing the value of the asymmetry index C_net = 0.31,0.00,0.75,0.64,0.62,0.80 in the simulated network, respectively.

Figure 8

Simulation results for a recurrent excitatory neuronal network in the high-frequency regime with heterogeneous firing frequencies and ξ > 0. (A) (Left panel) Time course of the mean synaptic strength $\overline{g}$(t) for different standard deviations σ_ν. (Right panel) Distribution of the firing frequencies and the final synaptic strengths. (B) (Left panel) Time course of the order parameter $r(t)$. (Right panel) Time course of the normalized number of closed loops of length 2, representing the number of bidirectional connections in the network. (C) Final coupling matrices for σ_ν = 0.01, 0.02, 0.03, 0.04 Hz, respectively. Network and STDP parameters are N = 200, $\overline{\nu }$ = 80 Hz, τ_d = 0.5 ms, ${\tau }_{\mathrm{a}}=0.3\mathrm{ms}$, A₊ = A₋ = 0.005, and τ₊ = τ₋ = 20 ms.

Figure 9

Simulation results for a recurrent excitatory neuronal network in the high-frequency regime with heterogeneous firing frequencies and ξ > 0. (A) (Left panel) Time course of the mean synaptic strength $\overline{g}$(t) for different standard deviations ${\sigma }_{\nu }$. (Right panel) Distribution of the firing frequencies and the final synaptic strengths. (B) (Left panel) Time course of the order parameter r(t). (Right panel) Time course of the normalized number of closed loops of length 2, representing the number of bidirectional connections in the network. (C) Final coupling matrices for σ_ν = 0.01, 0.02, 0.03, 0.04 Hz, respectively. Network and STDP parameters are N = 200, $\overline{\nu }$ = 80 Hz, τ_d = 0.5 ms, ${\tau }_{\mathrm{a}}=1.0\mathrm{ms}$, A₊ = A₋ = 0.005, and τ₊ = τ₋ = 20 ms.