Spike patterning of a stochastic phase model neuron given periodic inhibition

We present a phase model of a repetitively firing neuron possessing a phase-dependent stochastic response to periodic inhibition. We analyze the dynamics in terms of a stochastic phase map and determine the invariant phase distribution. We use the latter to compute both the distribution of interspike intervals (ISIs) and the stochastic winding number (mean firing rate) as a function of the input frequency. We show that only low-order phase locking persists in the presence of weak phase dependence, and is characterized statistically by a multimodal ISI distribution and a nonmonotonic variation in the stochastic winding number as a function of input frequency.

Figure 1

(a) $p_{AB}\left(\tau \right)$ for the equations ${\varphi }_{n+1}=T_B+{\varphi }_n+a_0+z$ with $z\sim N(0,\sigma )$, $a_0=-0.2$, $\sigma =0.05$, and $T_B=0.9$. Notice that even with $R\left(\theta \right)=a_0$ (no phase-dependent coupling) the distribution of relative timings still has a single mode of more likely input times. (b) Schematic diagram showing the construction of $p_{AB}\left(\tau \right)$ from the terms $p_{AB}^j\left(\tau \right)$: The initial phase $\psi$ at time $t=0$ occurs with probability $q^{*}\left(\psi \right)d\psi$. Some of these phase trajectories, comprising the density $G_1\left(q^{*}\right)\left(\tau \right)=p_{AB}^1\left(\tau \right)$, spike at time $T_B-\tau$ before the next input at $t=T_B$ [top row of panel (b)]. Others starting at the initial phase $\psi$ spike after $t=T_B$ but before $t=2T_B$. These phase trajectories comprise the density $G_1\left[G_0\right(q^{*}\left)\right]\left(\tau \right)=p_{AB}^2\left(\tau \right)$, shown in the bottom row of panel (b). Phase trajectories that spike in the time interval $[2T_B,3T_B)$ (not shown) comprise the density $G_1\left\{G_0\right[G_0\left(q^{*}\right)\left]\right\}\left(\tau \right)=p_{AB}^3\left(\tau \right)$, etc.

Figure 2

(a) $p_{AA}(T\mid \tau )$ for the equations ${\varphi }_{n+1}=T_B+{\varphi }_n+a_0+\xi$ with $\xi \sim N(0,\sigma )$, $a_0=-0.1$, $\sigma =0.025$, $\tau =0.25$, and $T_B=0.5$ $({\omega }_B=2)$. The discontinuities occur at $T=k{T}_B+\tau$, $k=1,2,\dots$, and reflect the separation of oscillators that fire after $k$ inputs. (b) Below the graph is an idealized representation of likely spike trains (output times) in response to inputs (input times), given the first input was at $\tau =0.25$. For this particular case the most likely interspike intervals are found in the intervals $[2T_B+\tau ,3T_B+\tau )$ and $[3T_B+\tau ,4T_B+\tau )$, associated with the respective probability densities $p_{AA}^2(T\mid \tau )$ and $p_{AA}^3(T\mid \tau )$.

Figure 3

(a) $q^{*}\left(\psi \right)$ calculated as the dominant eigenvector of $G$ (solid line) and the perturbative solution (dashed line). Parameters are $f(\psi ,\xi )=\psi +T_B+a_0+ϵ\sin \left(2\pi \psi \right)+\xi$ with $T_B=1.25$, $ϵ=0.01$, $\xi \sim N(0,\sigma )$, $\sigma =0.025$, and $a_0=-0.2$. (b) With the same parameters except $ϵ=0.05$ $q^{*}\left(\psi \right)$ can only be calculated as the dominant eigenvector of $G$ as the perturbative solution is not valid for $ϵ$ this large.

Figure 4

Calculations of the stochastic winding number ${\Omega }_{A,\mathrm{stoch}}$ for the sine circle map through a frequency sweep from ${\Theta }^{-1}=0.9$ to ${\Theta }^{-1}=1.1$ with model parameters $a_0=0$, $ϵ=0.01$, and $\sigma =0.05$, 0.1, and 0.15 represented in panels (a), (b), and (c), respectively. The thick solid line is the actual stochastic winding number; the thick dashed line is the perturbative solution. For reference, the thin dashed line is the deterministic winding number that was numerically calculated as in Eq. (16) for large $n$. In all three plots 1:1 statistical phase locking produces the paradoxical region. Increased noise squashes the magnitude and broadens the width of the paradoxical region. Notice that the amplitude of the singularly perturbed solution is not a good match to the actual solution but the minima and maxima provide a good approximation to the size of the paradoxical region.

Figure 5

Calculations of $p_{AA}\left(T\right)$ for the noisy circle map with no phase dependent coupling $(ϵ=0)$ through a frequency sweep from ${\omega }_B=0.7$ to ${\omega }_B=1.6$. Model parameters are $a_0=-0.2$ and $\sigma =0.025$.

Figure 6

Calculations of $p_{AA}\left(T\right)$ for the noisy circle map through a frequency sweep from ${\omega }_B=0.7$ to ${\omega }_B=1.6$. Model parameters are $a_0=-0.2$, $\sigma =0.025$, and $ϵ=0.1$.

Figure 7

Calculations of ${\Omega }_{A,\mathrm{stoch}}$ for the noisy circle map through a frequency sweep from ${\omega }_B=0.7$ to ${\omega }_B=1.6$. Model parameters are $a_0=-0.2$ and $ϵ=0.1$. Solid line represents $\sigma =0.025$, dashed line $\sigma =0.1$, and dot-dashed line $\sigma =0.2$

Figure 8

Calculations of $p_{AA}\left(T\right)$ (left column) and the associated $q^{*}\left(\psi \right)$ (right column) for the low noise level $(\sigma =0.025)$ sine circle map through a frequency sweep from ${\omega }_B=0.7$ through ${\omega }_B=1$. Model parameters are $a_0=-0.2$ and $ϵ=0.1$.

Figure 9

Calculations of $p_{AA}\left(T\right)$ (left column) and the associated $q^{*}\left(\psi \right)$ (right column) for the higher noise level $(\sigma =0.1)$ sine circle map through a frequency sweep from ${\omega }_B=0.7$ through ${\omega }_B=1$. Model parameters are $a_0=-0.2$ and $ϵ=0.1$.

Figure 10

Calculations of $p_{AA}\left(T\right)$ (left column) and the associated $q^{*}\left(\psi \right)$ (right column) for the low noise level $(\sigma =0.025)$ sine circle map through a frequency sweep from ${\omega }_B=1.35$ through ${\omega }_B=1.5$. Model parameters are $a_0=-0.2$ and $ϵ=0.1$.