COBA and PING Dynamics, Discretisation, and Learning

Model family and notation

The COBA and PING classifiers use the same conductance-based leaky integrate-and-fire neuron. Their distinction is architectural. The COBA condition disables the recurrent excitatory-to-inhibitory-to-excitatory pathway, whereas PING enables the loop $𝐸\to 𝐼\to 𝐸$.[1]

Physical time is $𝑡$, the simulation-step index is $𝑘$, and ${𝑡}_{𝑘}=𝑘\mathrm{\Delta }{𝑡}_{\text{sim}}$, where $\mathrm{\Delta }{𝑡}_{\text{sim}}$ is the integration timestep. A presentation contains ${𝑁}_{𝑡}$ timesteps and has duration ${𝑇}_{\text{present}}={𝑁}_{𝑡}\mathrm{\Delta }{𝑡}_{\text{sim}}$. Superscripts on conductances identify AMPA or GABA polarity, while subscripts identify the receiving population.

Continuous membrane dynamics

For neuron $𝑗$ in population $𝑋\in \left\{𝐸,𝐼\right\}$, the membrane equation is

In Equation 1, ${𝑉}_{𝑚}$ is membrane voltage, ${𝐶}_{𝑚}$ is capacitance, ${𝑔}_{𝐿}$ is leak conductance, and ${𝑔}^{𝐸}$ and ${𝑔}^{𝐼}$ are total AMPA and GABA conductances. The reversal potentials are

A conductance $𝑔$ with reversal potential $𝐸$ supplies current $𝐼=𝑔\left(𝐸-{𝑉}_{𝑚}\right)$. Equation 2 therefore makes AMPA normally depolarizing and GABA normally hyperpolarizing. The fixed membrane parameters are

The corresponding passive membrane time constants follow from ${𝜏}_{𝑚}=\frac{{𝐶}_{𝑚}}{{𝑔}_{𝐿}}$:

Continuous synaptic dynamics

The continuous spike train emitted by presynaptic neuron $𝑖$ is

where ${𝑡}_{𝑖,𝑗}$ is event time $𝑗$ and $𝛿$ is the Dirac impulse. A postsynaptic conductance ${𝑔}_{𝑗}$ with decay constant $𝜏$ satisfies

Here ${𝑤}_{\text{event},𝑖,𝑗}$ is the conductance increment produced by one presynaptic event. AMPA and GABA use separate copies of Equation 6 with ${𝜏}_{\text{AMPA}}=2\text{ms}$ and the canonical ${𝜏}_{\text{GABA}}=6\text{ms}$.

Discretising the synaptic equations

Between events, Equation 6 reduces to

Integrating Equation 7 from ${𝑡}_{𝑘}$ to ${𝑡}_{𝑘+1}={𝑡}_{𝑘}+\mathrm{\Delta }{𝑡}_{\text{sim}}$ gives

Exponentiating Equation 8 yields the event-free decay

Integrating Equation 6 across an infinitesimal interval containing one event gives

Combining Equation 9 and Equation 10 gives the timestep recurrence

In Equation 11, $𝒔\left[𝑘\right]$ is a dimensionless binary presynaptic spike vector, $𝑊$ has presynaptic neurons on rows and postsynaptic neurons on columns, and $𝒖\left[𝑘\right]$ is an external conductance increment. The simulator performs exact exponential decay followed by the full, undecayed event increment.[1]

E–I conductance recurrences

Define the dimensionless decay factors

For one E–I layer, Equation 11 expands to

The canonical architecture sets ${𝑊}_{𝐸\to 𝐸}={𝑊}_{𝐼\to 𝐼}=0$. COBA additionally disables ${𝑊}_{𝐸\to 𝐼}$ and ${𝑊}_{𝐼\to 𝐸}$. PING enables those two projections: excitatory spikes recruit inhibitory neurons through Equation 15, and inhibitory spikes subsequently produce GABA conductance through Equation 14. A recurrent spike generated in one iteration is consumed in the next, giving each recurrent projection a one-timestep delay.

Exponential-Euler membrane discretisation

The simulator first updates conductances with Equation 13–Equation 16, then holds those new values constant during the membrane step. Expanding Equation 1 for one neuron gives

Define the total conductance and its conductance-weighted equilibrium voltage by

Substituting Equation 18 into Equation 17 gives

Let $𝑞={𝑉}_{𝑚}-{𝑉}_{\infty }$. Because the conductances and ${𝑉}_{\infty }$ are constant within the timestep,

Integrating Equation 20 over one timestep gives

Returning to membrane voltage produces the pre-threshold candidate

Equivalently, define

Equation 22 is the update used by the simulator. It is exact for the linear membrane equation under the stated zero-order hold on conductances; it is not an exact solution of the continuously co-evolving neuron and synapses.

Spike, reset, and refractory dynamics

Let $𝑟\left[𝑘\right]$ be the non-negative integer refractory counter. It is decremented before the threshold decision:

A neuron emits the binary spike

where $𝐻$ is the Heaviside step function, ${𝑉}_{\text{th}}=-50\text{mV}$, and $\mathrm{𝟏}$ is an indicator. A neuron that spikes, or remains refractory after Equation 24, is stored at ${𝑉}_{\text{reset}}=-65\text{mV}$. A new spike sets the counter to the population’s refractory step count; otherwise it retains ${𝑟}^{-}\left[𝑘+1\right]$.

The canonical timestep is $\mathrm{\Delta }{𝑡}_{\text{sim}}=0.1\text{ms}$. The excitatory and inhibitory refractory durations, $1.2\text{ms}$ and $0.6\text{ms}$, therefore correspond exactly to 12 and 6 timesteps.

Output readout and classification loss

The canonical classifier uses an output LIF with mean pre-reset voltage. Let $𝒉\left[𝑘\right]$ be the final hidden E-population spike vector and ${𝑊}_{\text{out}}$ the hidden-to-output matrix. With output decay ${𝛽}_{𝑜}=exp\left(-\mathrm{\Delta }\frac{{𝑡}_{\text{sim}}}{{𝜏}_{𝑜}}\right)$, the pre-reset dimensionless output state is

The output spike and subtractive reset are

where ${𝜃}_{𝑜}=1$. Output state and its evidence accumulator begin at zero for each independent image presentation. Hidden E/I voltages, conductances, spikes and refractory counters likewise initialize anew for every independent presentation. The decoder therefore receives presentation boundaries through these per-forward resets. The class score is the mean pre-reset state over the complete readout window:

For minibatch member $𝑏$, softmax-normalized evidence is

It is used for optimization but is not assumed to be a calibrated probability. For true class ${𝑦}_{𝑏}$, the mean cross-entropy training loss is

Some conditions add a one-sided firing-rate ceiling. Let ${𝐾}_{𝑏,𝑙,𝑗}={\sum }_{𝑘=1}^{{𝑁}_{𝑡}}{𝑠}_{𝑏,𝑙,𝑗}\left[𝑘\right]$ be the spike count of hidden E neuron $𝑗$ in layer $𝑙$ for presentation $𝑏$. Its sample-wise population-mean rate is

With tested ceiling ${𝑟}_{𝐸,\text{ceil}}$ and coefficient ${𝜆}_{\text{rate}}$, the penalty is

The canonical reference conditions set ${𝜆}_{\text{rate}}=0$; activity-ceiling conditions use Equation 32 with a positive coefficient.

Backpropagation through time

Collect all differentiable voltages, conductances, spikes, output states and accumulators into the state vector $𝒙\left[𝑘\right]$. One complete simulation step has the abstract form

${𝐹}_{𝑘}$ is a bundle rather than one additional simulator equation. It combines the discrete synapse rule in Equation 11 and the decay factors in Equation 12; the population-specific conductance updates from Equation 13 through Equation 16; the implemented membrane update in Equation 23; the refractory, spike and reset operations in Equation 24 and Equation 25; and the output update in Equation 26 and Equation 27. It also advances the accumulators used to construct Equation 28 and Equation 31. The earlier continuous equations and intermediate derivation steps explain these operations but are not extra updates inside ${𝐹}_{𝑘}$.

Here $𝒖\left[𝑘\right]$ denotes external input, $𝜃$ is the shared trainable parameter vector, and $\ell$ includes the readout and objective defined in Equation 28–Equation 32. At $\mathrm{\Delta }{𝑡}_{\text{sim}}=0.1\text{ms}$, a $200\text{ms}$ presentation applies ${𝐹}_{𝑘}$ 2,000 times. Each copy has a different state but reuses the same parameters.

Define the state sensitivity

Let ${𝑥}_{𝑖\left[𝑘\right]}$ be component $𝑖$ of the current state and ${𝐹}_{𝑘,𝑗}$ component $𝑗$ of ${𝐹}_{𝑘}$. The project evaluates its loss from the terminal readout and rate accumulators. Consequently, for $𝑘<{𝑁}_{𝑡}$, the current state affects the loss through the next state. Applying the multivariable chain rule component by component gives

The second line of Equation 35 sums every route from state component ${𝑥}_{𝑖\left[𝑘\right]}$ to the later loss through next-state component ${𝑥}_{𝑗\left[𝑘+1\right]}$. The third line substitutes Equation 34 and Equation 33. Stacking the component equations produces the transpose-Jacobian vector product in the fourth line; the final line supplies the terminal condition from Equation 33. Because each parameter is reused, its total gradient is

Equation 35 and Equation 36 define backpropagation through time (BPTT). PyTorch autograd evaluates these vector–Jacobian products and accumulates parameter gradients when backward propagation is requested; it need not construct the full Jacobian matrices explicitly.[4]

For example, the local recurrence $𝒈\left[𝑘+1\right]=𝛽𝒈\left[𝑘\right]+𝒔\left[𝑘\right]𝑊$ gives

Equation 37 shows temporal credit assignment directly: every use of the shared weight contributes according to the later loss sensitivity of the conductance it created.

Surrogate spike derivatives

BPTT requires derivatives of the update operations, but the hard spike in Equation 25 has derivative zero away from threshold and is undefined at threshold. Direct use of that derivative would block useful learning signals.

Surrogate-gradient learning retains the binary forward spike while assigning a smooth backward derivative.[2] For $𝜉={𝑉}_{\text{candidate}}-{𝑉}_{\text{th}}$, the chosen proxy is

The approximate equality in Equation 38 denotes an assigned backward rule rather than the true derivative of the Heaviside function. For hidden neurons, $𝜉$ is measured in millivolts and $𝛽$ in reciprocal millivolts. The standard classifiers use $𝛽=1$.

Let $sg\left(𝑧\right)$ denote stop-gradient: it returns $𝑧$ forward but has zero derivative. The simulator constructs

Equation 39 equals $𝐻\left(𝜉\right)$ numerically, while autograd differentiates only $𝑝\left(𝜉\right)$.[5] Boolean reset and refractory decisions remain non-differentiable. Resetting hidden voltage to a constant cuts that voltage-state path, while the floating spike value carries the surrogate derivative into downstream synapses and the readout.

Why voltage gradients can explode

Expanding Equation 35 across many steps produces repeated Jacobian products:

Products in Equation 40 can shrink useful sensitivities or amplify them. This project combines a 2,000-step graph with voltage, conductance, spike and output-state memory. PING adds a closed route by which a voltage perturbation can return to the population from which it began:

The local gain of this route follows from the membrane equations rather than from an extra model quantity. Choose $𝑋=𝐸$ for excitatory conductance or $𝑋=𝐼$ for inhibitory conductance. Differentiating Equation 18 and the decay $𝑎$ in Equation 23 gives

Applying the product rule to $\mathrm{\Delta }𝑉=\left({𝑉}_{\infty }-{𝑉}_{𝑚}\right)\left(1-𝑎\right)$ gives the conductance-to-voltage response

The response ${𝑅}_{𝑋}$ in Equation 43 has units of millivolts per microsiemens. Let ${𝑅}_{𝐸\to 𝐼}$ and ${𝑅}_{𝐼\to 𝐸}$ be diagonal matrices of the appropriate responses at I and E neurons. Let ${𝑆}_{𝐸}$ and ${𝑆}_{𝐼}$ be diagonal matrices containing the local surrogate derivatives from Equation 38. One complete backward traversal of Equation 41 then has return Jacobian

The transposes in Equation 44 convert the simulator’s presynaptic-row, postsynaptic-column weight storage into maps acting on column perturbations.

At a fixed local linearization, a sensitivity returning through the loop $𝑛$ times contains ${𝐽}_{\text{loop}}^{𝑛}𝛿{𝑽}_{𝐸}$. If the spectral radius $𝜌\left({𝐽}_{\text{loop}}\right)>1$, some modes grow geometrically. A singular value above one permits amplification over one traversal, while time-varying, non-normal products can amplify transiently even if each instantaneous eigenvalue has magnitude below one.[3] Inhibitory polarity can reverse a perturbation’s sign without making its magnitude small.

Recurrence therefore does not guarantee exploding gradients. It creates the return operator in Equation 44, which BPTT may traverse many times. Loop-disabled COBA lacks that particular operator, but it still carries gradients through membrane, conductance, spike and readout state across the same long presentation.

Global gradient clipping cannot by itself repair this internal propagation. It acts only after BPTT has formed the complete parameter gradient. Internal sensitivities may already be enormous, non-finite, or dominated by the amplifying paths in Equation 40 before clipping occurs.

Voltage-gradient damping

The project therefore intervenes inside the backward recurrence while preserving the physical forward trajectory. For dimensionless divisor ${𝑑}_{\text{grad}}>0$, define

Forward evaluation of Equation 45 gives ${\mathrm{𝒟︀}}_{𝛼\left(𝑧\right)}=𝑧$, but its backward derivative is

The hidden E and I membrane updates replace Equation 23 by

Equation 45 and Equation 47 leave voltages and spikes unchanged forward. Backward, each selected conductance-to-voltage response is multiplied by ${𝛼}_{\text{grad}}$. The PING return term in Equation 44 consequently becomes

The standard paired classifiers use ${𝑑}_{\text{grad}}=1000$, so ${𝛼}_{\text{grad}}={10}^{-3}$ and this particular complete loop term is scaled by ${10}^{-6}$. This does not scale the entire state Jacobian because Equation 45 is inserted only on the $\mathrm{\Delta }𝑉$ branch of Equation 47. For one hidden neuron, the three local derivatives entering its candidate voltage are therefore

Thus the conductance-to-voltage blocks are multiplied by ${𝛼}_{\text{grad}}$, but the direct identity route from ${𝑉}_{𝑚\left[𝑘\right]}$ to ${𝑉}_{\text{candidate}\left[𝑘+1\right]}$ retains coefficient one. The conductance decay and spike-injection blocks, the refractory-counter update, and the local output recurrence contain no damping operator and are not directly rescaled. Downstream composite derivatives can nevertheless inherit attenuation when their path passes through a damped hidden-voltage increment. There is therefore no scalar ${𝛼}_{\text{grad}}$ for which the complete step Jacobian is simply ${𝛼}_{\text{grad}}{𝐽}_{\text{step}}$.

Holding conductances fixed before bounds and reset, Equation 49 and Equation 23 give

Without damping, Equation 50 reduces to $𝑎$. Strong damping moves the direct old-voltage derivative toward one while attenuating sensitivities entering through the membrane increment. It can therefore suppress recurrent-drive amplification without simply erasing membrane memory.

Voltage-gradient damping is an altered backward rule, not a biophysical process or the exact derivative of the forward simulator. It deliberately biases temporal credit assignment and may suppress useful sensitivities along with unstable ones. The divisor is therefore an empirical training choice. Standard classifiers apply it to hidden E and I membrane increments, before voltage bounds, spike detection and reset; output-neuron updates are undamped.[1]

Gradient clipping, AdamW, and Dale projection

After BPTT has accumulated the batch gradient $𝒈={\nabla }_{𝜃}{𝐿}_{\text{total}}$, the trainer computes its global Euclidean norm

With clipping threshold $𝑐=1$, it replaces the gradient by

A non-finite norm causes the optimizer step to be skipped. Otherwise AdamW forms first- and second-moment estimates

The square in Equation 53 is elementwise. Bias correction gives

For learning rate $𝜂$, decoupled weight decay ${𝜆}_{\text{wd}}$, and numerical stabilizer $𝜀$, AdamW applies

The standard training recipe used $𝜂=0.0004$, ${𝜆}_{\text{wd}}=0$, minibatches of 256 presentations and 50 epochs. With zero weight decay, Equation 55 is numerically the Adam update for the same hyperparameters.

Dale-constrained trainable connection weights are finally projected onto the non-negative orthant:

The canonical input and readout matrices are trainable, while recurrent E–I matrices are fixed unless an experiment explicitly makes them trainable. Thus one optimization cycle evaluates the forward dynamics, computes Equation 32, applies surrogate-gradient BPTT with the internal rule in Equation 47, clips by Equation 52, updates by Equation 55, and finishes with Equation 56.

Scope and interpretation

The dynamical discretisation is a numerical approximation to the continuous hybrid system, although Equation 9 and Equation 22 exactly solve their respective subproblems under the stated within-step assumptions. Surrogate spike derivatives and voltage-gradient damping introduce additional backward rules: they make learning possible and tractable but do not turn BPTT into the exact derivative of a smooth physical neuron.

The loop derivation isolates one interpretable amplification route rather than the full state Jacobian, and it does not imply that gamma oscillation itself causes gradient explosion. A paired empirical comparison can vary loop engagement and damping independently.