iaf_cond_beta

iaf_cond_beta#

class brainpy.state.iaf_cond_beta(in_size, E_L=Quantity(-70., "mV"), C_m=Quantity(250., "pF"), t_ref=Quantity(2., "ms"), V_th=Quantity(-55., "mV"), V_reset=Quantity(-60., "mV"), E_ex=Quantity(0., "mV"), E_in=Quantity(-85., "mV"), g_L=Quantity(16.6667, "nS"), tau_rise_ex=Quantity(0.2, "ms"), tau_decay_ex=Quantity(0.2, "ms"), tau_rise_in=Quantity(2., "ms"), tau_decay_in=Quantity(2., "ms"), I_e=Quantity(0., "pA"), gsl_error_tol=1e-06, V_initializer=Constant(value=-70. mV), g_ex_initializer=Constant(value=0. nS), g_in_initializer=Constant(value=0. nS), spk_fun=ReluGrad(alpha=0.3, width=1.0), spk_reset='hard', ref_var=False, name=None)#

NEST-compatible conductance-based leaky integrate-and-fire neuron with beta-shaped synaptic conductances.

This model implements a conductance-based LIF neuron with beta-function (dual-exponential) synaptic conductances for both excitatory and inhibitory channels. It follows NEST’s iaf_cond_beta implementation, including hard threshold crossing, absolute refractory period, and one-step delayed external current buffering.

1. Mathematical Model

The membrane voltage \(V_m\) evolves according to:

\[C_m \frac{dV_m}{dt} = -I_\mathrm{leak} - I_{\mathrm{syn,ex}} - I_{\mathrm{syn,in}} + I_e + I_\mathrm{stim}\]

where the currents are defined as:

\[\begin{split}I_\mathrm{leak} &= g_L (V_m - E_L) \\ I_{\mathrm{syn,ex}} &= g_\mathrm{ex}(t) (V_m - E_\mathrm{ex}) \\ I_{\mathrm{syn,in}} &= g_\mathrm{in}(t) (V_m - E_\mathrm{in})\end{split}\]

During the refractory period, the membrane voltage is clamped to \(V_\mathrm{reset}\) and \(dV_m/dt = 0\). Outside the refractory period, the effective voltage for synaptic current computation is bounded by \(\min(V_m, V_\mathrm{th})\).

2. Beta-Function Conductance Dynamics

Each synaptic conductance (excitatory and inhibitory) is modeled using two coupled state variables to produce a beta-function (rise-decay) waveform:

\[\begin{split}\frac{d\,dg_\mathrm{ex}}{dt} &= -\frac{dg_\mathrm{ex}}{\tau_{\mathrm{decay,ex}}} \\ \frac{d g_\mathrm{ex}}{dt} &= dg_\mathrm{ex} - \frac{g_\mathrm{ex}}{\tau_{\mathrm{rise,ex}}}\end{split}\]

\[\begin{split}\frac{d\,dg_\mathrm{in}}{dt} &= -\frac{dg_\mathrm{in}}{\tau_{\mathrm{decay,in}}} \\ \frac{d g_\mathrm{in}}{dt} &= dg_\mathrm{in} - \frac{g_\mathrm{in}}{\tau_{\mathrm{rise,in}}}\end{split}\]

Incoming spikes cause instantaneous jumps in \(dg_\mathrm{ex}\) or \(dg_\mathrm{in}\). Positive weights target the excitatory channel; negative weights target the inhibitory channel. Each spike weight (in nS) is multiplied by the beta normalization factor \(\kappa(\tau_\mathrm{rise}, \tau_\mathrm{decay})\) to ensure unit weight produces a 1 nS peak conductance.

The normalization factor is computed as:

\[\kappa = \frac{1/\tau_\mathrm{rise} - 1/\tau_\mathrm{decay}}{\exp(-t_\mathrm{peak}/\tau_\mathrm{decay}) - \exp(-t_\mathrm{peak}/\tau_\mathrm{rise})}\]

where \(t_\mathrm{peak} = \frac{\tau_\mathrm{rise} \tau_\mathrm{decay}}{\tau_\mathrm{decay} - \tau_\mathrm{rise}} \ln\left(\frac{\tau_\mathrm{decay}}{\tau_\mathrm{rise}}\right)\).

3. Numerical Integration

ODEs are integrated using the Runge-Kutta-Fehlberg (RKF45) adaptive-step method with embedded error control. The integrator maintains a persistent step size estimate (integration_step) across simulation steps, adjusting it based on local truncation error to satisfy a fixed absolute tolerance (gsl_error_tol).

4. Update Order (NEST Semantics)

Each simulation step executes the following operations in order:

Integrate all ODEs on the interval \((t, t+dt]\) using RKF45.
Inside integration loop: apply refractory clamp and spike/reset.
After loop: decrement refractory counter once.
Apply incoming spike weights to \(dg_\mathrm{ex}\) and \(dg_\mathrm{in}\).
Store external current input x into the delayed buffer I_stim (affects next step).

This matches NEST’s ring-buffer semantics: external currents applied at time \(t\) take effect at time \(t + dt\).

5. Design Constraints and Assumptions

Refractory clamping: During refractory period, voltage is fixed at \(V_\mathrm{reset}\) and no integration occurs. NEST uses this approach for consistency with exact spike times.
Beta normalization edge case: When \(\tau_\mathrm{rise} \approx \tau_\mathrm{decay}\), the normalization factor approaches \(e / \tau_\mathrm{decay}\) to avoid division by zero.

Parameters:

in_size (Size) – Population shape, specified as an integer (1D), tuple of integers (multi-dimensional), or brainstate Size object. Determines the shape of all state variables and parameters.
E_L (ArrayLike, optional) – Leak reversal potential. Default: -70 mV. Broadcast to in_size if scalar. Must have units of voltage (mV).
C_m (ArrayLike, optional) – Membrane capacitance. Default: 250 pF. Broadcast to in_size if scalar. Must be strictly positive. Determines voltage response timescale \(\tau_m = C_m / g_L\).
t_ref (ArrayLike, optional) – Absolute refractory period duration. Default: 2 ms. Broadcast to in_size if scalar. Must be non-negative. Converted to discrete grid steps via \(\lceil t_\mathrm{ref} / dt \rceil\).
V_th (ArrayLike, optional) – Spike threshold voltage. Default: -55 mV. Broadcast to in_size if scalar. Must satisfy \(V_\mathrm{reset} < V_\mathrm{th}\).
V_reset (ArrayLike, optional) – Post-spike reset voltage. Default: -60 mV. Broadcast to in_size if scalar. Must be strictly less than V_th. Neuron is clamped to this value during refractory period.
E_ex (ArrayLike, optional) – Excitatory reversal potential. Default: 0 mV. Broadcast to in_size if scalar. Typically positive (depolarizing).
E_in (ArrayLike, optional) – Inhibitory reversal potential. Default: -85 mV. Broadcast to in_size if scalar. Typically more negative than \(E_L\) (hyperpolarizing).
g_L (ArrayLike, optional) – Leak conductance. Default: 16.6667 nS (yields \(\tau_m = 15\) ms with default \(C_m\)). Broadcast to in_size if scalar. Must be strictly positive.
tau_rise_ex (ArrayLike, optional) – Excitatory conductance rise time constant. Default: 0.2 ms. Broadcast to in_size if scalar. Must be strictly positive. Smaller values produce faster rise times.
tau_decay_ex (ArrayLike, optional) – Excitatory conductance decay time constant. Default: 0.2 ms. Broadcast to in_size if scalar. Must be strictly positive. When equal to tau_rise_ex, beta function degenerates to alpha function.
tau_rise_in (ArrayLike, optional) – Inhibitory conductance rise time constant. Default: 2.0 ms. Broadcast to in_size if scalar. Must be strictly positive. Typically slower than excitatory rise for GABA receptors.
tau_decay_in (ArrayLike, optional) – Inhibitory conductance decay time constant. Default: 2.0 ms. Broadcast to in_size if scalar. Must be strictly positive. Determines inhibitory synaptic integration window.
I_e (ArrayLike, optional) – Constant external current. Default: 0 pA. Broadcast to in_size if scalar. Added to membrane current at every time step.
gsl_error_tol (ArrayLike) – Unitless local RKF45 error tolerance, broadcastable and strictly positive.
V_initializer (Callable, optional) – Initialization function for membrane voltage. Default: Constant(-70 mV). Called as V_initializer(varshape) during init_state().
g_ex_initializer (Callable, optional) – Initialization function for excitatory conductance. Default: Constant(0 nS). Called as g_ex_initializer(varshape) during init_state().
g_in_initializer (Callable, optional) – Initialization function for inhibitory conductance. Default: Constant(0 nS). Called as g_in_initializer(varshape) during init_state().
spk_fun (Callable, optional) – Surrogate gradient function for spike generation. Default: ReluGrad(). Applied to scaled voltage \((V - V_\mathrm{th}) / (V_\mathrm{th} - V_\mathrm{reset})\) to produce differentiable spike output for gradient-based learning.
spk_reset (str, optional) – Spike reset mode. Default: 'hard' (stop-gradient reset, matches NEST behavior). Alternative: 'soft' (subtractive reset \(V \leftarrow V - V_\mathrm{th}\)).
ref_var (bool, optional) – If True, create a boolean refractory state variable indicating refractory status. Default: False. Useful for monitoring or conditional computations.
name (str, optional) – Module name for debugging and visualization. Default: None (auto-generated).

Parameter Mapping

The following table maps constructor parameters to mathematical notation and NEST equivalents:

Parameter	Default	Math equivalent	Description
`in_size`	(required)	—	Population shape
`E_L`	-70 mV	\(E_\mathrm{L}\)	Leak reversal potential
`C_m`	250 pF	\(C_\mathrm{m}\)	Membrane capacitance
`t_ref`	2 ms	\(t_\mathrm{ref}\)	Absolute refractory duration
`V_th`	-55 mV	\(V_\mathrm{th}\)	Spike threshold
`V_reset`	-60 mV	\(V_\mathrm{reset}\)	Reset potential
`E_ex`	0 mV	\(E_\mathrm{ex}\)	Excitatory reversal potential
`E_in`	-85 mV	\(E_\mathrm{in}\)	Inhibitory reversal potential
`g_L`	16.6667 nS	\(g_\mathrm{L}\)	Leak conductance
`tau_rise_ex`	0.2 ms	\(\tau_{\mathrm{rise,ex}}\)	Excitatory beta rise constant
`tau_decay_ex`	0.2 ms	\(\tau_{\mathrm{decay,ex}}\)	Excitatory beta decay constant
`tau_rise_in`	2.0 ms	\(\tau_{\mathrm{rise,in}}\)	Inhibitory beta rise constant
`tau_decay_in`	2.0 ms	\(\tau_{\mathrm{decay,in}}\)	Inhibitory beta decay constant
`I_e`	0 pA	\(I_\mathrm{e}\)	Constant external current
`gsl_error_tol`	1e-6	—	RKF45 error tolerance
`V_initializer`	Constant(-70 mV)	—	Membrane initializer
`g_ex_initializer`	Constant(0 nS)	—	Excitatory conductance initializer
`g_in_initializer`	Constant(0 nS)	—	Inhibitory conductance initializer
`spk_fun`	ReluGrad()	—	Surrogate spike function
`spk_reset`	`'hard'`	—	Reset mode (`'hard'` matches NEST)
`ref_var`	`False`	—	Expose boolean refractory indicator

Raises:

ValueError – If V_reset >= V_th (reset potential must be below threshold).
ValueError – If C_m <= 0 (capacitance must be strictly positive).
ValueError – If t_ref < 0 (refractory period cannot be negative).
ValueError – If any of tau_rise_ex, tau_decay_ex, tau_rise_in, tau_decay_in are non-positive.

Notes

State Variables

Vbrainstate.HiddenState
Membrane potential \(V_m\) with shape (*in_size,) and units mV.
dg_exbrainstate.ShortTermState
Excitatory beta auxiliary state (nS/ms).
g_exbrainstate.HiddenState
Excitatory synaptic conductance with units nS.
dg_inbrainstate.ShortTermState
Inhibitory beta auxiliary state (nS/ms).
g_inbrainstate.HiddenState
Inhibitory synaptic conductance with units nS.
refractory_step_countbrainstate.ShortTermState
Remaining refractory grid steps (integer, dtype int32). Zero when not refractory.
integration_stepbrainstate.ShortTermState
Persistent RKF45 internal step size with units ms. Adapted automatically for numerical stability.
I_stimbrainstate.ShortTermState
One-step delayed external current buffer with units pA. Updated after ODE integration.
last_spike_timebrainstate.ShortTermState
Time of last emitted spike (units ms). Set to \(t + dt\) when spike occurs.
refractorybrainstate.ShortTermState (optional)
Boolean refractory indicator. Only created if ref_var=True.

Performance Considerations:

This model uses per-neuron scalar NumPy integration loops, which are significantly slower than vectorized JAX operations. For large populations, consider using iaf_cond_exp or iaf_cond_alpha with vectorized exponential integrators. The RKF45 method is primarily intended for high-accuracy validation against NEST rather than production simulations.

NEST Compatibility:

This implementation matches NEST 3.9+ iaf_cond_beta semantics, including:

Beta normalization factor computation (exact formula match).
One-step delayed external current handling.
Refractory voltage clamping during integration.
Hard threshold crossing and immediate reset.

Minor differences from NEST:

NEST uses GSL’s RK integrator; this uses a pure-Python RKF45 implementation.
Numerical differences may appear at \(O(10^{-6})\) due to floating-point rounding.

Examples

Basic Usage:

>>> import brainpy.state as bst
>>> import saiunit as u
>>> import brainstate as bstate
>>>
>>> with bstate.environ.context(dt=0.1 * u.ms):
...     neuron = bst.iaf_cond_beta(10, V_th=-50*u.mV, V_reset=-65*u.mV)
...     neuron.init_all_states()
...     # Apply excitatory synaptic input (5 nS conductance jump)
...     neuron.add_delta_input('syn_input', 5.0 * u.nS)
...     spikes = neuron()
...     print(neuron.V.value[:3])  # Membrane voltages of first 3 neurons
[-70. -70. -70.] mV

Comparing Excitatory and Inhibitory Time Constants:

>>> import matplotlib.pyplot as plt
>>> with bstate.environ.context(dt=0.01 * u.ms):
...     fast_ex = bst.iaf_cond_beta(1, tau_rise_ex=0.2*u.ms, tau_decay_ex=2.0*u.ms)
...     slow_in = bst.iaf_cond_beta(1, tau_rise_in=2.0*u.ms, tau_decay_in=10.0*u.ms)
...     fast_ex.init_all_states()
...     slow_in.init_all_states()
...     # Single excitatory spike at t=1ms
...     fast_ex.add_delta_input('spike', 10.0 * u.nS)
...     # Record excitatory conductance
...     g_ex_trace = []
...     for _ in range(500):
...         fast_ex()
...         g_ex_trace.append(fast_ex.g_ex.value[0])
...     plt.plot(g_ex_trace)
...     plt.xlabel('Time (0.01 ms steps)')
...     plt.ylabel('g_ex (nS)')
...     plt.title('Beta-function conductance waveform')

Network with Balanced Excitation and Inhibition:

>>> from brainevent.nn import FixedProb
>>> exc_neurons = bst.iaf_cond_beta(800, E_L=-70*u.mV, V_th=-50*u.mV)
>>> inh_neurons = bst.iaf_cond_beta(200, E_L=-70*u.mV, V_th=-50*u.mV)
>>> exc_neurons.init_all_states()
>>> inh_neurons.init_all_states()
>>> # Create projections (placeholder - requires brainevent)
>>> # exc_proj = FixedProb(exc_neurons, exc_neurons, prob=0.1, weight=2.0*u.nS)
>>> # inh_proj = FixedProb(inh_neurons, exc_neurons, prob=0.2, weight=-5.0*u.nS)

iaf_cond_beta

Contents

iaf_cond_beta#