aeif_psc_delta

aeif_psc_delta#

class brainpy.state.aeif_psc_delta(in_size, V_peak=Quantity(0., "mV"), V_reset=Quantity(-60., "mV"), t_ref=Quantity(0., "ms"), g_L=Quantity(30., "nS"), C_m=Quantity(281., "pF"), E_L=Quantity(-70.6, "mV"), Delta_T=Quantity(2., "mV"), tau_w=Quantity(144., "ms"), a=Quantity(4., "nS"), b=Quantity(80.5, "pA"), V_th=Quantity(-50.4, "mV"), I_e=Quantity(0., "pA"), gsl_error_tol=1e-06, refractory_input=False, V_initializer=Constant(value=-70.6 mV), w_initializer=Constant(value=0. pA), spk_fun=ReluGrad(alpha=0.3, width=1.0), spk_reset='hard', ref_var=False, name=None)#

NEST-compatible aeif_psc_delta neuron model.

Current-based adaptive exponential integrate-and-fire neuron with delta-shaped synaptic input. Implements NEST models/aeif_psc_delta.{h,cpp} semantics with adaptive RKF45 integration, in-loop spike handling, and optional refractory input buffering.

1. Mathematical Formulation

The model combines exponential spike-initiation current (AdEx), spike-triggered and subthreshold adaptation current \(w\), and delta-function synaptic input that directly jumps the membrane voltage.

Membrane dynamics:

\[C_m \frac{dV}{dt} = -g_L (V - E_L) + g_L \Delta_T \exp\!\left(\frac{V - V_{th}}{\Delta_T}\right) - w + I_e + I_{\mathrm{stim}}.\]

Adaptation dynamics:

\[\tau_w \frac{dw}{dt} = a (V - E_L) - w.\]

Incoming delta spikes are interpreted as instantaneous voltage jumps:

\[V \leftarrow V + J \sum_k \delta(t - t_k).\]

Here \(J\) is the synaptic weight in millivolts.

2. Refractory and Spike Handling

During refractory integration (when refractory_step_count > 0), NEST clamps the effective membrane voltage to V_reset and sets \(dV/dt = 0\). Outside refractory, the RHS uses \(\min(V, V_{\mathrm{peak}})\) as the effective voltage to prevent unbounded exponential growth.

Threshold detection uses:

V_peak if Delta_T > 0 (exponential regime),
V_th if Delta_T == 0 (IAF-like limit).

On each detected spike:

V is reset to V_reset,
adaptation jump w <- w + b is applied immediately,
refractory counter is set to ceil(t_ref / dt) + 1 if t_ref > 0.

Spike handling occurs inside the adaptive RKF45 substep loop. With t_ref = 0, multiple spikes can occur within one simulation step, matching NEST behavior.

3. Refractory Input Buffering

If refractory_input=True (NEST refractory_input flag), delta spikes arriving during refractory are accumulated into refractory_spike_buffer with NEST’s exponential discount factor:

\[\mathrm{buffer} \leftarrow \mathrm{buffer} + J \cdot \exp(-r \cdot \Delta t / \tau_m),\]

where \(r\) is the current refractory step count. The buffered input is applied when the neuron exits refractory.

4. Update Order per Simulation Step

Integrate ODEs on interval \((t, t+\Delta t]\) via adaptive RKF45.
Inside integration loop:
1. Apply arriving delta jump to V (if not refractory).
2. Apply refractory clamp (V <- V_reset if refractory).
3. Apply refractory buffering logic if refractory_input=True.
4. Threshold detection and spike/reset/adaptation handling.
After loop: Decrement refractory counter once.
Apply arriving spike weights directly to V as delta-function pulses.
Store external current input x into one-step delayed I_stim.

5. Numerical Integration Details

The model uses adaptive Runge-Kutta-Fehlberg 4(5) (RKF45) with local error control. The step size integration_step is adjusted per neuron to satisfy the tolerance gsl_error_tol, matching NEST’s GSL solver behavior. Minimum step size is clamped to 1e-8 ms to prevent infinite loops.

The exponential term is computed using the effective voltage \(V_{\mathrm{eff}}\):

\[\begin{split}V_{\mathrm{eff}} = \begin{cases} V_{\mathrm{reset}} & \text{if refractory}, \\ \min(V, V_{\mathrm{peak}}) & \text{otherwise}. \end{cases}\end{split}\]

Overflow protection: the model validates that \((V_{\mathrm{peak}} - V_{\mathrm{th}}) / \Delta_T\) does not exceed log(DBL_MAX / 1e20) to prevent numerical overflow at spike time.

6. Differences from NEST

Spike output format: NEST emits spike events; brainpy.state returns a binary array (0/1) per simulation step. Internal dynamics (adaptation increments, refractory handling) match NEST exactly.
Surrogate gradients: brainpy.state uses spk_fun (e.g., ReluGrad()) for differentiable spike generation; NEST does not support gradient-based learning.
Spike reset mode: spk_reset='hard' (default) matches NEST; 'soft' is available but non-canonical.

Parameters:

in_size (int, tuple of int) – Shape of the neuron population. Can be an integer (1D) or tuple (multi-dimensional).
V_peak (ArrayLike, optional) – Spike detection threshold in millivolts. Used when Delta_T > 0. Default: 0.0 * u.mV. Can be scalar or array matching in_size.
V_reset (ArrayLike, optional) – Reset potential in millivolts. Default: -60.0 * u.mV. Must satisfy V_reset < V_peak.
t_ref (ArrayLike, optional) – Absolute refractory period in milliseconds. Default: 0.0 * u.ms (no refractoriness). Must be non-negative.
g_L (ArrayLike, optional) – Leak conductance in nanosiemens. Default: 30.0 * u.nS. Must be positive.
C_m (ArrayLike, optional) – Membrane capacitance in picofarads. Default: 281.0 * u.pF. Must be positive.
E_L (ArrayLike, optional) – Leak reversal potential in millivolts. Default: -70.6 * u.mV.
Delta_T (ArrayLike, optional) – Exponential slope factor in millivolts. Default: 2.0 * u.mV. Set to 0.0 for IAF-like limit. Must be non-negative.
tau_w (ArrayLike, optional) – Adaptation time constant in milliseconds. Default: 144.0 * u.ms. Must be positive.
a (ArrayLike, optional) – Subthreshold adaptation coupling in nanosiemens. Default: 4.0 * u.nS.
b (ArrayLike, optional) – Spike-triggered adaptation increment in picoamperes. Default: 80.5 * u.pA.
V_th (ArrayLike, optional) – Spike initiation threshold in millivolts (used in exponential term). Default: -50.4 * u.mV. Must satisfy V_th <= V_peak.
I_e (ArrayLike, optional) – Constant external current in picoamperes. Default: 0.0 * u.pA.
gsl_error_tol (ArrayLike, optional) – RKF45 local error tolerance (unitless). Default: 1e-6. Must be positive.
refractory_input (bool, optional) – If True, accumulate delta spikes arriving during refractory with NEST’s exponential discount factor and apply when refractory ends. Default: False.
V_initializer (Callable, optional) – Membrane potential initializer. Default: braintools.init.Constant(-70.6 * u.mV).
w_initializer (Callable, optional) – Adaptation current initializer. Default: braintools.init.Constant(0.0 * u.pA).
spk_fun (Callable, optional) – Surrogate gradient function for differentiable spike generation. Default: braintools.surrogate.ReluGrad().
spk_reset (str, optional) – Spike reset mode. Options: 'hard' (stop gradient), 'soft' (V -= V_th). Default: 'hard'.
ref_var (bool, optional) – If True, expose self.refractory as a boolean state variable. Default: False.
name (str, optional) – Name of the neuron group.

Parameter Mapping

Parameter	Default	Math equivalent	Description
`in_size`	(required)	—	Population shape
`V_peak`	0 mV	\(V_{\mathrm{peak}}\)	Spike detection threshold (if \(\Delta_T > 0\))
`V_reset`	-60 mV	\(V_{\mathrm{reset}}\)	Reset potential
`t_ref`	0 ms	\(t_{\mathrm{ref}}\)	Absolute refractory duration
`g_L`	30 nS	\(g_{\mathrm{L}}\)	Leak conductance
`C_m`	281 pF	\(C_{\mathrm{m}}\)	Membrane capacitance
`E_L`	-70.6 mV	\(E_{\mathrm{L}}\)	Leak reversal potential
`Delta_T`	2 mV	\(\Delta_T\)	Exponential slope factor
`tau_w`	144 ms	\(\tau_w\)	Adaptation time constant
`a`	4 nS	\(a\)	Subthreshold adaptation coupling
`b`	80.5 pA	\(b\)	Spike-triggered adaptation increment
`V_th`	-50.4 mV	\(V_{\mathrm{th}}\)	Spike initiation threshold (in exponential term)
`I_e`	0 pA	\(I_{\mathrm{e}}\)	Constant external current
`gsl_error_tol`	1e-6	—	RKF45 local error tolerance
`refractory_input`	`False`	—	If True, buffer spikes during refractory with NEST discounting
`V_initializer`	Constant(-70.6 mV)	—	Membrane initializer
`w_initializer`	Constant(0 pA)	—	Adaptation current initializer
`spk_fun`	ReluGrad()	—	Surrogate spike function
`spk_reset`	`'hard'`	—	Reset mode; hard reset matches NEST behavior
`ref_var`	`False`	—	If True, expose boolean refractory indicator

State Variables

Vbrainstate.HiddenState: Membrane potential \(V_m\) in millivolts. Shape: (*in_size,).
wbrainstate.HiddenState: Adaptation current in picoamperes. Shape: (*in_size,).
refractory_step_countbrainstate.ShortTermState: Remaining refractory grid steps (int32). Shape: (*in_size,).
integration_stepbrainstate.ShortTermState: Persistent RKF45 internal step size in milliseconds. Shape: (*in_size,).
I_stimbrainstate.ShortTermState: One-step delayed current buffer in picoamperes. Shape: (*in_size,).
last_spike_timebrainstate.ShortTermState: Last emitted spike time in milliseconds (\(t+\Delta t\) on spike). Shape: (*in_size,).
refractorybrainstate.ShortTermState, optional: Boolean refractory indicator. Only present if ref_var=True. Shape: (*in_size,).

Raises:

ValueError – If V_reset >= V_peak.
ValueError – If Delta_T < 0.
ValueError – If V_peak < V_th.
ValueError – If C_m <= 0.
ValueError – If t_ref < 0.
ValueError – If tau_w <= 0.
ValueError – If gsl_error_tol <= 0.
ValueError – If (V_peak - V_th) / Delta_T exceeds log(DBL_MAX / 1e20) (overflow protection).
ValueError – During integration: if V < -1e3 or |w| > 1e6 (numerical instability).

Examples

Basic usage with delta-function input:

>>> import brainpy.state as bst
>>> import saiunit as u
>>> import brainstate as bs
>>> import jax.numpy as jnp
>>>
>>> # Create 100 AdEx neurons
>>> neu = bst.aeif_psc_delta(100, V_peak=0.*u.mV, V_reset=-60.*u.mV, t_ref=2.*u.ms)
>>>
>>> # Initialize state
>>> with bs.environ.context(dt=0.1*u.ms):
...     neu.init_all_states()
...
...     # Simulate with delta input
...     for i in range(100):
...         # Delta spike input (+1 mV jump)
...         neu.add_delta_input('external', lambda: jnp.ones(100) * 1.0 * u.mV)
...         spk = neu.step_run(i, 0.0*u.pA)
...         print(f"t={i*0.1:.1f}ms: {spk.sum():.0f} spikes")

With refractory input buffering:

>>> # Enable refractory buffering
>>> neu = bst.aeif_psc_delta(
...     100,
...     V_peak=0.*u.mV,
...     V_reset=-60.*u.mV,
...     t_ref=5.*u.ms,
...     refractory_input=True
... )
>>>
>>> with bs.environ.context(dt=0.1*u.ms):
...     neu.init_all_states()
...
...     # Spikes arriving during refractory are buffered and discounted
...     for i in range(100):
...         neu.add_delta_input('external', lambda: jnp.ones(100) * 2.0 * u.mV)
...         spk = neu.step_run(i, 0.0*u.pA)

IAF-like limit (Delta_T = 0):

>>> # Delta_T=0 disables exponential term
>>> neu = bst.aeif_psc_delta(
...     100,
...     Delta_T=0.0*u.mV,
...     V_th=-55.*u.mV,
...     V_peak=-55.*u.mV,  # Must equal V_th when Delta_T=0
...     a=0.0*u.nS,        # No subthreshold adaptation
...     b=0.0*u.pA         # No spike-triggered adaptation
... )

aeif_psc_delta

Contents

aeif_psc_delta#