Adaptive Exponential Integrate-and-Fire (aEIF) neuron models¶
-
class
aeif_cond_alpha
: public Archiving_Node¶ - #include <aeif_cond_alpha.h>
Name: aeif_cond_alpha - Conductance based exponential integrate-and-fire neuron model according to Brette and Gerstner (2005). Description:
aeif_cond_alpha is the adaptive exponential integrate and fire neuron according to Brette and Gerstner (2005). Synaptic conductances are modelled as alpha-functions.
This implementation uses the embedded 4th order Runge-Kutta-Fehlberg solver with adaptive step size to integrate the differential equation.
The membrane potential is given by the following differential equation:
\[\begin{split} C_m \frac{dV}{dt} = -g_L(V-E_L)+g_L\Delta_T\exp\left(\frac{V-V_{th}}{\Delta_T}\right) - g_e(t)(V-E_e) \\ -g_i(t)(V-E_i)-w +I_e \end{split}\]and
\[ \tau_w \frac{dw}{dt} = a(V-E_L) - w \]Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
g_ex
nS
Excitatory synaptic conductance
dg_ex
nS/ms
First derivative of g_ex
g_in
nS
Inhibitory synaptic conductance
dg_in
nS/ms
First derivative of g_in
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Synaptic parameters
E_ex
mV
Excitatory reversal potential
tau_syn_ex
ms
Rise time of excitatory synaptic conductance (alpha function)
E_in
mV
Inhibitory reversal potential
tau_syn_in
ms
Rise time of the inhibitory synaptic conductance (alpha function)
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities.
Authors: Marc-Oliver Gewaltig; full revision by Tanguy Fardet on December 2016
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
References:
- 1
Brette R and Gerstner W (2005). Adaptive Exponential Integrate-and-Fire Model as an Effective Description of Neuronal Activity. J Neurophysiol 94:3637-3642 DOI: https://doi.org/10.1152/jn.00686.2005
SeeAlso: iaf_cond_alpha, aeif_cond_exp
-
class
aeif_cond_alpha_multisynapse
: public Archiving_Node¶ - #include <aeif_cond_alpha_multisynapse.h>
Name: aeif_cond_alpha_multisynapse - Conductance based adaptive exponential integrate-and-fire neuron model according to Brette and Gerstner (2005) with multiple synaptic rise time and decay time constants, and synaptic conductance modeled by an alpha function.
Description:
aeif_cond_alpha_multisynapse is a conductance-based adaptive exponential integrate-and-fire neuron model. It allows an arbitrary number of synaptic time constants. Synaptic conductance is modeled by an alpha function, as described by A. Roth and M.C.W. van Rossum in Computational Modeling Methods for Neuroscientists, MIT Press 2013, Chapter 6.
The time constants are supplied by an array, “tau_syn”, and the pertaining synaptic reversal potentials are supplied by the array “E_rev”. Port numbers are automatically assigned in the range from 1 to n_receptors. During connection, the ports are selected with the property “receptor_type”.
The membrane potential is given by the following differential equation:
\[ C dV/dt = -g_L(V-E_L) + g_L*\Delta_T*\exp((V-V_T)/\Delta_T) + I_{syn_{tot}}(V, t)- w + I_e \]where\[ I_{syn_{tot}}(V,t) = \sum_i g_i(t) (V - E_{rev,i}) , \]the synapse i is excitatory or inhibitory depending on the value of \( E_{rev,i}\) and the differential equation for the spike-adaptation current w is:
\[ \tau_w * dw/dt = a(V - E_L) - w \]When the neuron fires a spike, the adaptation current w <- w + b.
Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Delta_T
mV
Slope factor
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
tau_w
ms
Adaptation time constant
Synaptic parameters
E_rev
list of mV
Reversal potential
tau_syn
list of ms
Time constant of synaptic conductance
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities.
Examples:
import nest import numpy as np neuron = nest.Create('aeif_cond_alpha_multisynapse') nest.SetStatus(neuron, {"V_peak": 0.0, "a": 4.0, "b":80.5}) nest.SetStatus(neuron, {'E_rev':[0.0, 0.0, 0.0, -85.0], 'tau_syn':[1.0, 5.0, 10.0, 8.0]}) spike = nest.Create('spike_generator', params = {'spike_times': np.array([10.0])}) voltmeter = nest.Create('voltmeter', 1, {'withgid': True}) delays=[1.0, 300.0, 500.0, 700.0] w=[1.0, 1.0, 1.0, 1.0] for syn in range(4): nest.Connect(spike, neuron, syn_spec={'model': 'static_synapse', 'receptor_type': 1 + syn, 'weight': w[syn], 'delay': delays[syn]}) nest.Connect(voltmeter, neuron) nest.Simulate(1000.0) dmm = nest.GetStatus(voltmeter)[0] Vms = dmm["events"]["V_m"] ts = dmm["events"]["times"] import pylab pylab.figure(2) pylab.plot(ts, Vms) pylab.show()
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
Author: Hans Ekkehard Plesser, based on aeif_cond_beta_multisynapse
SeeAlso: aeif_cond_alpha_multisynapse
-
class
aeif_cond_alpha_RK5
: public Archiving_Node¶ - #include <aeif_cond_alpha_RK5.h>
Name: aeif_cond_alpha_RK5 - Conductance based exponential integrate-and-fire neuron model according to Brette and Gerstner (2005)
Description:
aeif_cond_alpha_RK5 is the adaptive exponential integrate and fire neuron according to Brette and Gerstner (2005). Synaptic conductances are modelled as alpha-functions.
This implementation uses a 5th order Runge-Kutta solver with adaptive stepsize to integrate the differential equation (see Numerical Recipes 3rd Edition, Press et al. 2007, Ch. 17.2).
The membrane potential is given by the following differential equation:
\[\begin{split} C dV/dt= -g_L(V-E_L)+g_L*\Delta_T*\exp((V-V_T)/\Delta_T)-g_e(t)(V-E_e) \\ -g_i(t)(V-E_i)-w +I_e \end{split}\]and\[ \tau_w * dw/dt= a(V-E_L) -w \]Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
g_ex
nS
Excitatory synaptic conductance
dg_ex
nS/ms
First derivative of g_ex
g_in
nS
Inhibitory synaptic conductance
dg_in
nS/ms
First derivative of g_in
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Synaptic parameters
E_ex
mV
Excitatory reversal potential
tau_syn_ex
ms
Rise time of excitatory synaptic conductance (alpha function)
E_in
mV
Inhibitory reversal potential
tau_syn_in
ms
Rise time of the inhibitory synaptic conductance (alpha function)
Numerical integration parameters
HMIN
ms
Minimal stepsize for numerical integration (default 0.001ms)
MAXERR
mV
Error estimate tolerance for adaptive stepsize control (steps accepted if err<=MAXERR). Note that the error refers to the difference between the 4th and 5th order RK terms. Default 1e-10 mV.
Authors: Stefan Bucher, Marc-Oliver Gewaltig.
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
- 1
Brette R and Gerstner W (2005). Adaptive Exponential Integrate-and-Fire Model as an Effective Description of Neuronal Activity. J Neurophysiol 94:3637-3642. DOI: https://doi.org/10.1152/jn.00686.2005
SeeAlso: iaf_cond_alpha, aeif_cond_exp, aeif_cond_alpha
-
class
aeif_cond_beta_multisynapse
: public Archiving_Node¶ - #include <aeif_cond_beta_multisynapse.h>
Name: aeif_cond_beta_multisynapse - Conductance based adaptive exponential integrate-and-fire neuron model according to Brette and Gerstner (2005) with multiple synaptic rise time and decay time constants, and synaptic conductance modeled by a beta function.
Description:
aeif_cond_beta_multisynapse is a conductance-based adaptive exponential integrate-and-fire neuron model. It allows an arbitrary number of synaptic rise time and decay time constants. Synaptic conductance is modeled by a beta function, as described by A. Roth and M.C.W. van Rossum in Computational Modeling Methods for Neuroscientists, MIT Press 2013, Chapter 6.
The time constants are supplied by two arrays, “tau_rise” and “tau_decay” for the synaptic rise time and decay time, respectively. The synaptic reversal potentials are supplied by the array “E_rev”. The port numbers are automatically assigned in the range from 1 to n_receptors. During connection, the ports are selected with the property “receptor_type”.
The membrane potential is given by the following differential equation:
\[ C dV/dt = -g_L(V-E_L) + g_L*\Delta_T*\exp((V-V_T)/\Delta_T) + I_{syn_{tot}}(V, t) - w + I_e \]where:
\[ I_{syn_{tot}}(V,t) = \sum_i g_i(t) (V - E_{rev,i}) , \]the synapse i is excitatory or inhibitory depending on the value of \( E_{rev,i} \) and the differential equation for the spike-adaptation current w is:
\[ \tau_w * dw/dt = a(V - E_L) - w \]When the neuron fires a spike, the adaptation current w <- w + b.
Parameters: The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Delta_T
mV
Slope factor
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
tau_w
ms
Adaptation time constant
Synaptic parameters
E_rev
list of mV
Reversal potential
tau_syn
list of ms
Time constant of synaptic conductance
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities.
Examples:
import nest import numpy as np neuron = nest.Create('aeif_cond_beta_multisynapse') nest.SetStatus(neuron, {"V_peak": 0.0, "a": 4.0, "b":80.5}) nest.SetStatus(neuron, {'E_rev':[0.0,0.0,0.0,-85.0], 'tau_decay':[50.0,20.0,20.0,20.0], 'tau_rise':[10.0,10.0,1.0,1.0]}) spike = nest.Create('spike_generator', params = {'spike_times': np.array([10.0])}) voltmeter = nest.Create('voltmeter', 1, {'withgid': True}) delays=[1.0, 300.0, 500.0, 700.0] w=[1.0, 1.0, 1.0, 1.0] for syn in range(4): nest.Connect(spike, neuron, syn_spec={'model': 'static_synapse', 'receptor_type': 1 + syn, 'weight': w[syn], 'delay': delays[syn]}) nest.Connect(voltmeter, neuron) nest.Simulate(1000.0) dmm = nest.GetStatus(voltmeter)[0] Vms = dmm["events"]["V_m"] ts = dmm["events"]["times"] import pylab pylab.figure(2) pylab.plot(ts, Vms) pylab.show()
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
Author: Bruno Golosio 07/10/2016
SeeAlso: aeif_cond_alpha_multisynapse
-
class
aeif_cond_exp
: public Archiving_Node¶ - #include <aeif_cond_exp.h>
Name: aeif_cond_exp - Conductance based exponential integrate-and-fire neuron model according to Brette and Gerstner (2005).
Description:
aeif_cond_exp is the adaptive exponential integrate and fire neuron according to Brette and Gerstner (2005), with post-synaptic conductances in the form of truncated exponentials.
This implementation uses the embedded 4th order Runge-Kutta-Fehlberg solver with adaptive stepsize to integrate the differential equation.
The membrane potential is given by the following differential equation:
\[\begin{split} C dV/dt= -g_L(V-E_L)+g_L*\Delta_T*\exp((V-V_T)/\Delta_T)-g_e(t)(V-E_e) \\ -g_i(t)(V-E_i)-w +I_e \end{split}\]and
\[ \tau_w * dw/dt= a(V-E_L) -W \]Note that the spike detection threshold V_peak is automatically set to \( V_th+10 mV \) to avoid numerical instabilites that may result from setting V_peak too high.
Parameters: The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
g_ex
nS
Excitatory synaptic conductance
g_in
nS
Inhibitory synaptic conductance
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Spike adaptation parameters
a
nS
Subthreshold adaptation
b
pA
Spike-triggered adaptation
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Synaptic parameters
E_ex
mV
Excitatory reversal potential
tau_syn_ex
ms
Rise time of excitatory synaptic conductance (alpha function)
E_in
mV
Inhibitory reversal potential
tau_syn_in
ms
Rise time of the inhibitory synaptic conductance (alpha function)
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities.
Author: Adapted from aeif_cond_alpha by Lyle Muller; full revision by Tanguy Fardet on December 2016
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
- 1
Brette R and Gerstner W (2005). Adaptive Exponential Integrate-and-Fire Model as an Effective Description of Neuronal Activity. J Neurophysiol 94:3637-3642. DOI: https://doi.org/10.1152/jn.00686.2005
SeeAlso: iaf_cond_exp, aeif_cond_alpha
-
class
aeif_psc_alpha
: public Archiving_Node¶ - #include <aeif_psc_alpha.h>
Name: aeif_psc_alpha - Current-based exponential integrate-and-fire neuron model according to Brette and Gerstner (2005).
Description:
aeif_psc_alpha is the adaptive exponential integrate and fire neuron according to Brette and Gerstner (2005). Synaptic currents are modelled as alpha-functions.
This implementation uses the embedded 4th order Runge-Kutta-Fehlberg solver with adaptive step size to integrate the differential equation.
The membrane potential is given by the following differential equation:
\[\begin{split} C dV/dt= -g_L(V-E_L)+g_L*\Delta_T*\exp((V-V_T)/\Delta_T)-g_e(t)(V-E_e) \\ -g_i(t)(V-E_i)-w +I_e \end{split}\]and
\[ \tau_w * dw/dt= a(V-E_L) -W \]Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
I_ex
pA
Excitatory synaptic current
dI_ex
pA/ms
First derivative of I_ex
I_in
pA
Inhibitory synaptic current
dI_in
pA/ms
First derivative of I_in
w
pA
Spike-adaptation current
g
pa
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Synaptic parameters
tau_syn_ex
ms
Rise time of excitatory synaptic conductance (alpha function)
tau_syn_in
ms
Rise time of the inhibitory synaptic conductance (alpha function)
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities
Author: Tanguy Fardet
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
References:
- 1
Brette R and Gerstner W (2005). Adaptive Exponential Integrate-and-Fire Model as an Effective Description of Neuronal Activity. J Neurophysiol 94:3637-3642. DOI: https://doi.org/10.1152/jn.00686.2005
SeeAlso: iaf_psc_alpha, aeif_cond_exp
-
class
aeif_psc_delta
: public Archiving_Node¶ - #include <aeif_psc_delta.h>
Name: aeif_psc_delta - Current-based adaptive exponential integrate-and-fire neuron model according to Brette and Gerstner (2005) with delta synapse.
Description:
aeif_psc_delta is the adaptive exponential integrate and fire neuron according to Brette and Gerstner (2005), with post-synaptic currents in the form of delta spikes.
This implementation uses the embedded 4th order Runge-Kutta-Fehlberg solver with adaptive stepsize to integrate the differential equation.
The membrane potential is given by the following differential equation:
\[\begin{split} C dV/dt= -g_L(V-E_L)+g_L*\Delta_T*\exp((V-V_T)/\Delta_T)-g_e(t)(V-E_e) \\ -g_i(t)(V-E_i)-w +I_e \end{split}\]and
\[ \tau_w * dw/dt= a(V-E_L) -W \]\[ I(t) = J \sum_k \delta(t - t^k). \]Here delta is the dirac delta function and k indexes incoming spikes. This is implemented such that V_m will be incremented/decremented by the value of J after a spike.
Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables
V_m
mV
Membrane potential
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
tau_w
ms
Adaptation time constant
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_th
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities.
Author: Mikkel Elle Lepperød adapted from aeif_psc_exp and iaf_psc_delta
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
References:
- 1
Brette R and Gerstner W (2005). Adaptive Exponential Integrate-and-Fire Model as an Effective Description of Neuronal Activity. J Neurophysiol 94:3637-3642. DOI: https://doi.org/10.1152/jn.00686.2005
SeeAlso: iaf_psc_delta, aeif_cond_exp, aeif_psc_exp
-
class
aeif_psc_delta_clopath
: public Clopath_Archiving_Node¶ - #include <aeif_psc_delta_clopath.h>
Name: aeif_psc_delta_clopath - Exponential integrate-and-fire neuron model according to Clopath et al. (2010).
Description:
aeif_psc_delta_clopath is an implementation of the neuron model as it is used in [1]. It is an extension of the aeif_psc_delta model and capable of connecting to a Clopath synapse.
Note that there are two points that are not mentioned in the paper but present in a MATLAB implementation by Claudia Clopath [3]. The first one is the clamping of the membrane potential to a fixed value after a spike occured to mimik a real spike and not just the upswing. This is important since the finite duration of the spike influences the evolution of the convolved versions (u_bar_[plus/minus]) of the membrane potential and thus the change of the synaptic weight. Secondly, there is a delay with which u_bar_[plus/minus] are used to compute the change of the synaptic weight.
Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables
V_m
mV
Membrane potential
w
pA
Spike-adaptation current
z
pA
Spike-adaptation current
V_th
mV
Adaptive spike initiation threshold
u_bar_plus
mV
Low-pass filtered Membrane potential
u_bar_minus
mV
Low-pass filtered Membrane potential
u_bar_bar
mV
Low-pass filtered u_bar_minus
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
tau_plus
ms
Time constant of u_bar_plus
tau_minus
ms
Time constant of u_bar_minus
tau_bar_bar
ms
Time constant of u_bar_bar
Spike adaptation parameters
a
nS
Subthreshold adaptation
b
pA
Spike-triggered adaptation
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_peak
mV
Spike detection threshold
V_th_max
mV
Value of V_th afer a spike
V_th_rest
mV
Resting value of V_th
Clopath rule parameters
A_LTD
1/mV
Amplitude of depression
A_LTP
1/mV^2
Amplitude of facilitation
theta_plus
mV
Threshold for u
theta_minus
mV
Threshold for u_bar_[plus/minus]
A_LTD_const
boolean
Flag that indicates whether A_LTD_ should be constant (true, default) or multiplied by u_bar_bar^2 / u_ref_squared (false).
delay_u_bars
real
Delay with which u_bar_[plus/minus] are processed to compute the synaptic weights.
U_ref_squared
real
Reference value for u_bar_bar_^2.
Other parameters
t_clamp
ms
Duration of clamping of Membrane potential after a spike
V_clamp
mV
Value to which the Membrane potential is clamped
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities.
Note:
Neither the clamping nor the delayed processing of u_bar_[plus/minus] are mentioned in [1]. However, they are part of an reference implementation by Claudia Clopath et al. that can be found on ModelDB [3]. The clamping is important to mimic a spike which is otherwise not described by the aeif neuron model.
Author: Jonas Stapmanns, David Dahmen, Jan Hahne
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
References:
- 1
Clopath et al. (2010). Connectivity reflects coding: a model of voltage-based STDP with homeostasis. Nature Neuroscience 13(3):344-352. DOI: https://doi.org/10.1038/nn.2479
- 2
Clopath and Gerstner (2010). Voltage and spike timing interact in STDP – a unified model. Frontiers in Synaptic Neuroscience. 2:25 DOI: https://doi.org/10.3389/fnsyn.2010.00025
- 3
Voltage-based STDP synapse (Clopath et al. 2010) on ModelDB https://senselab.med.yale.edu/ModelDB/showmodel.cshtml?model=144566&file=%2f modeldb_package%2fVoTriCode%2faEIF.m
-
class
aeif_psc_exp
: public Archiving_Node¶ - #include <aeif_psc_exp.h>
Name: aeif_psc_exp - Current-based exponential integrate-and-fire neuron model according to Brette and Gerstner (2005).
Description:
aeif_psc_exp is the adaptive exponential integrate and fire neuron according to Brette and Gerstner (2005), with post-synaptic currents in the form of truncated exponentials.
This implementation uses the embedded 4th order Runge-Kutta-Fehlberg solver with adaptive stepsize to integrate the differential equation.
The membrane potential is given by the following differential equation:
\[\begin{split} C dV/dt= -g_L(V-E_L)+g_L*\Delta_T*\exp((V-V_T)/\Delta_T)-g_e(t)(V-E_e) \\ -g_i(t)(V-E_i)-w +I_e \end{split}\]and
\[ \tau_w * dw/dt= a(V-E_L) -W \]Note that the spike detection threshold V_peak is automatically set to \( V_th+10 \) mV to avoid numerical instabilites that may result from setting V_peak too high.
Parameters:
The following parameters can be set in the status dictionary.
Dynamic state variables:
V_m
mV
Membrane potential
I_ex
pA
Excitatory synaptic current
I_in
pA
Inhibitory synaptic current
w
pA
Spike-adaptation current
Membrane Parameters
C_m
pF
Capacity of the membrane
t_ref
ms
Duration of refractory period
V_reset
mV
Reset value for V_m after a spike
E_L
mV
Leak reversal potential
g_L
nS
Leak conductance
I_e
pA
Constant external input current
Spike adaptation parameters
a
ns
Subthreshold adaptation
b
pA
Spike-triggered adaptation
Delta_T
mV
Slope factor
tau_w
ms
Adaptation time constant
V_t
mV
Spike initiation threshold
V_peak
mV
Spike detection threshold
Synaptic parameters
tau_syn_ex
ms
Rise time of excitatory synaptic conductance (alpha function)
tau_syn_in
ms
Rise time of the inhibitory synaptic conductance (alpha function)
Integration parameters
gsl_error_tol
real
This parameter controls the admissible error of the GSL integrator. Reduce it if NEST complains about numerical instabilities
Author: Tanguy Fardet
Sends: SpikeEvent
Receives: SpikeEvent, CurrentEvent, DataLoggingRequest
References:
- 1
Brette R and Gerstner W (2005). Adaptive Exponential Integrate-and-Fire Model as an Effective Description of Neuronal Activity. J Neurophysiol 94:3637-3642. DOI: https://doi.org/10.1152/jn.00686.2005
SeeAlso: iaf_psc_exp, aeif_cond_exp