% ============================================================================
%  Artificial Neural Networks — Methods, Equations and Graphical
%  Representations
%  Author  : Simon-Pierre Boucher — contact@spboucher.ai
%  Chapter 12 : Biologically Inspired Networks (chapters/12-bio-inspired.tex)
% ============================================================================

\chapter{Biologically Inspired Networks: Spiking Neurons, Self-Organizing
         Maps and Reservoir Computing}\label{chap:bio}

Every architecture examined so far shares one commitment: a differentiable
map from inputs to outputs, trained by gradient descent on a global
objective. That commitment is a modelling choice, not a law of nature, and
the brain does not appear to honour it. Cortical neurons communicate with
discrete, asynchronous electrical events rather than real-valued
activations; synapses are modified by locally available information rather
than by a gradient routed backwards from a distant loss; and much of the
structure of sensory cortex emerges from unsupervised competition rather
than from labelled supervision.

This chapter develops three families that each relax one of those
commitments while remaining mathematically precise. \emph{Spiking neural
networks} replace the real-valued activation by a temporal train of events,
and the gradient by a local, timing-dependent plasticity rule.
\emph{Self-organizing maps} replace supervised error correction by
competitive learning, and obtain a topology-preserving projection of the
data as an emergent property. \emph{Reservoir computing} keeps the
recurrent architecture but refuses to train it, learning only a linear
readout — which reduces the entire training problem to convex regression.
For each family we give the state equations, an estimation algorithm in
pseudocode, and a figure of the mechanism.

% ----------------------------------------------------------------------------
\section{Spiking Neural Networks}
% ----------------------------------------------------------------------------

Spiking networks have been called the third generation of neural network
models \cite{maass1997}: the first generation computes with thresholded
binary units, the second with continuous activations, and the third with
sequences of events in continuous time. A spiking neuron does not output a
number; it outputs a set of firing times $\{t^{(1)}, t^{(2)}, \dots\}$,
conventionally written as a sum of Dirac impulses,
\begin{equation}
  S(t) \;=\; \sum_{f} \delta\!\left(t - t^{(f)}\right),
  \label{eq:bio-spiketrain}
\end{equation}
so that the activation of a neuron is a distribution rather than a scalar.
Information may then be carried by how \emph{often} a neuron fires, or by
\emph{when} it fires — a distinction with no counterpart in the
architectures of the preceding chapters.

\subsection{The Leaky Integrate-and-Fire Neuron}

The standard computational unit abstracts the cell membrane as a leaky
capacitor charged by the incoming synaptic current.

\begin{definition}[Leaky integrate-and-fire neuron]\label{def:bio-lif}
A leaky integrate-and-fire (LIF) neuron with membrane time constant
$\tau_m$, resting potential $V_{\mathrm{rest}}$, membrane resistance $R_m$
and threshold $V_{\mathrm{th}}$ evolves according to
\begin{equation}
  \tau_m \, \frac{\mathrm{d}V(t)}{\mathrm{d}t}
  \;=\; -\bigl(V(t) - V_{\mathrm{rest}}\bigr) \;+\; R_m\, I(t),
  \label{eq:bio-lif-ode}
\end{equation}
supplemented by the fire-and-reset rule
\begin{equation}
  \text{if } V(t^-) \ge V_{\mathrm{th}}:
  \quad
  \text{emit a spike at } t,
  \quad
  V(t^+) \leftarrow V_{\mathrm{reset}},
  \quad
  V \equiv V_{\mathrm{reset}} \text{ on } [t,\, t + t_{\mathrm{ref}}],
  \label{eq:bio-lif-reset}
\end{equation}
where $t_{\mathrm{ref}}$ is the absolute refractory period.
\end{definition}

Equation~\eqref{eq:bio-lif-ode} is a first-order linear ODE, so its
behaviour under a constant input current $I(t) \equiv I_0$ is available in
closed form. Starting from $V(0) = V_{\mathrm{reset}}$,
\begin{equation}
  V(t)
  \;=\;
  V_{\mathrm{rest}} + R_m I_0
  \;+\;
  \bigl(V_{\mathrm{reset}} - V_{\mathrm{rest}} - R_m I_0\bigr)
  \, e^{-t/\tau_m},
  \label{eq:bio-lif-sol}
\end{equation}
an exponential relaxation towards the asymptote
$V_\infty = V_{\mathrm{rest}} + R_m I_0$. Two regimes follow immediately.
If $V_\infty < V_{\mathrm{th}}$ the neuron is \emph{subthreshold} and never
fires, however long we wait: the input is filtered out entirely. If
$V_\infty > V_{\mathrm{th}}$, solving \eqref{eq:bio-lif-sol} for the
threshold crossing gives the interspike interval
\begin{equation}
  T_{\mathrm{ISI}}
  \;=\;
  t_{\mathrm{ref}}
  \;+\;
  \tau_m \,
  \ln\!\left(
    \frac{R_m I_0 - (V_{\mathrm{reset}} - V_{\mathrm{rest}})}
         {R_m I_0 - (V_{\mathrm{th}}  - V_{\mathrm{rest}})}
  \right),
  \label{eq:bio-lif-isi}
\end{equation}
and hence the current-to-frequency transfer function, the \emph{f--I curve}
of the neuron,
\begin{equation}
  f(I_0) \;=\; \frac{1}{T_{\mathrm{ISI}}}
  \;=\;
  \left[
    t_{\mathrm{ref}} + \tau_m \ln\!\left(
      \frac{R_m I_0 - V_{\mathrm{reset}} + V_{\mathrm{rest}}}
           {R_m I_0 - V_{\mathrm{th}} + V_{\mathrm{rest}}}
    \right)
  \right]^{-1}.
  \label{eq:bio-lif-fi}
\end{equation}
The refractory term in \eqref{eq:bio-lif-fi} imposes the saturation
$f \le 1/t_{\mathrm{ref}}$: unlike a ReLU, a spiking neuron has a hard
maximum firing rate. Figure~\ref{fig:bio-lif-trace} traces the
corresponding membrane dynamics.

\begin{figure}[htbp]
  \centering
  \begin{tikzpicture}
    \begin{axis}[
        width=13.4cm, height=6.2cm,
        xmin=0, xmax=42, ymin=-0.08, ymax=1.70,
        xlabel={time $t$ (ms)},
        ylabel={membrane potential $V(t)$},
        ylabel style={font=\small}, xlabel style={font=\small},
        tick label style={font=\scriptsize},
        axis lines=left,
        ytick={0,1}, yticklabels={$V_{\mathrm{reset}}$, $V_{\mathrm{th}}$},
        xtick={0,10,20,30,40},
        clip=false,
      ]
      % --- refractory shading (drawn first, behind everything) ---
      \fill[black!8] (axis cs:10.986,0) rectangle (axis cs:12.986,1.40);
      \fill[black!8] (axis cs:23.972,0) rectangle (axis cs:25.972,1.40);
      \fill[black!8] (axis cs:36.958,0) rectangle (axis cs:38.958,1.40);
      % --- threshold and asymptote ---
      \draw[dashed, black!65, thick] (axis cs:0,1)   -- (axis cs:42,1);
      \draw[dotted, black!50, thick] (axis cs:0,1.5) -- (axis cs:42,1.5);
      % --- spikes and resets (unrolled: \foreach variables do not survive
      %     the deferred execution of drawing commands inside an axis) ---
      \draw[coutput, line width=1.1pt]
        (axis cs:10.986,0) -- (axis cs:10.986,1.40);
      \draw[coutput, line width=1.1pt]
        (axis cs:23.972,0) -- (axis cs:23.972,1.40);
      \draw[coutput, line width=1.1pt]
        (axis cs:36.958,0) -- (axis cs:36.958,1.40);
      % --- charging ramps: V = 1.5 (1 - e^{-(t-t_0)/10}) ---
      \addplot[cinput, line width=1.1pt, domain=0:10.986,      samples=60]
        {1.5*(1-exp(-x/10))};
      \addplot[cinput, line width=1.1pt, domain=12.986:23.972, samples=60]
        {1.5*(1-exp(-(x-12.986)/10))};
      \addplot[cinput, line width=1.1pt, domain=25.972:36.958, samples=60]
        {1.5*(1-exp(-(x-25.972)/10))};
      \addplot[cinput, line width=1.1pt, domain=38.958:42,     samples=30]
        {1.5*(1-exp(-(x-38.958)/10))};
      % --- refractory plateaus ---
      \addplot[cinput, line width=1.1pt, domain=10.986:12.986, samples=2]{0};
      \addplot[cinput, line width=1.1pt, domain=23.972:25.972, samples=2]{0};
      \addplot[cinput, line width=1.1pt, domain=36.958:38.958, samples=2]{0};
      % --- annotations ---
      \node[etiquette, anchor=west, text=coutput]
        at (axis cs:11.25,1.20) {spike};
      \node[etiquette, anchor=south] at (axis cs:11.99,0.04)
        {$t_{\mathrm{ref}}$};
      \node[etiquette, anchor=west, fill=white, inner sep=1.5pt]
        at (axis cs:0.8,1.5) {asymptote $V_\infty = R_m I_0$, never reached};
    \end{axis}
  \end{tikzpicture}
  \caption{Membrane potential of an LIF neuron \eqref{eq:bio-lif-ode} driven
    by a constant supra-threshold current, with $\tau_m = 10$~ms,
    $R_m I_0 = 1.5\,V_{\mathrm{th}}$,
    $V_{\mathrm{rest}} = V_{\mathrm{reset}} = 0$ and
    $t_{\mathrm{ref}} = 2$~ms. The potential relaxes exponentially towards
    $V_\infty$, but crosses the dashed threshold first: a spike (red) is
    emitted, the potential is reset, and the shaded refractory window clamps
    it before integration resumes. The period predicted by
    \eqref{eq:bio-lif-isi} is
    $t_{\mathrm{ref}} + \tau_m \ln 3 \approx 12.99$~ms, a firing rate of
    about $77$~Hz.}
  \label{fig:bio-lif-trace}
\end{figure}

For simulation and for training, \eqref{eq:bio-lif-ode} is discretized with
a fixed step $\Delta t$. Writing $\beta = e^{-\Delta t / \tau_m}$ for the
per-step decay, taking $V_{\mathrm{rest}} = 0$, and folding the synaptic
weights $\vect{w}$ into the input current, one obtains the recurrence used
by every modern spiking-network framework:
\begin{align}
  V[t] &= \beta\, V[t-1]
          \;+\; \vect{w}\transp \vect{S}_{\mathrm{in}}[t]
          \;-\; S[t-1]\, V_{\mathrm{th}},
  \label{eq:bio-lif-discrete}\\[2pt]
  S[t] &= \Theta\!\bigl(V[t] - V_{\mathrm{th}}\bigr),
  \qquad
  \Theta(u) = \begin{cases} 1, & u \ge 0,\\ 0, & u < 0,\end{cases}
  \label{eq:bio-lif-heaviside}
\end{align}
where $\vect{S}_{\mathrm{in}}[t] \in \{0,1\}^{n}$ collects the presynaptic
spikes arriving at step $t$. The last term of \eqref{eq:bio-lif-discrete}
implements a \emph{soft reset} by subtraction, which preserves the
supra-threshold excess; the alternative \emph{hard reset} sets
$V[t] \leftarrow 0$ and discards it. Note the structural resemblance
between \eqref{eq:bio-lif-discrete} and the state update of a gated
recurrent unit: an LIF layer is a recurrent network whose forget gate is
frozen at the constant $\beta$ and whose output nonlinearity is a step.

\begin{remark}[Biophysically complete models]\label{rem:bio-hh}
The LIF neuron is a caricature. The Hodgkin--Huxley model instead resolves
the ionic conductances explicitly, through four coupled state variables:
\begin{equation}
  C_m \frac{\mathrm{d}V}{\mathrm{d}t}
  = I_{\mathrm{ext}}
  - \bar{g}_{\mathrm{Na}}\, m^{3} h \,(V - E_{\mathrm{Na}})
  - \bar{g}_{\mathrm{K}}\, n^{4} (V - E_{\mathrm{K}})
  - \bar{g}_{L} (V - E_{L}),
  \label{eq:bio-hh}
\end{equation}
each gating variable $x \in \{m, h, n\}$ obeying first-order kinetics
$\dot{x} = \bigl(x_\infty(V) - x\bigr)/\tau_x(V)$. Here the spike is an
emergent consequence of sodium activation rather than an imposed threshold
rule — at a cost of roughly $10^{3}$ floating-point operations per
millisecond of simulated time.
\end{remark}

\begin{remark}[The Izhikevich compromise]\label{rem:bio-izh}
A two-dimensional reduction of \eqref{eq:bio-hh} retains most of the
qualitative firing repertoire at about $13$ operations per millisecond:
\begin{equation}
  \frac{\mathrm{d}v}{\mathrm{d}t} = 0.04\, v^{2} + 5v + 140 - u + I,
  \qquad
  \frac{\mathrm{d}u}{\mathrm{d}t} = a\,(b v - u),
  \label{eq:bio-izh}
\end{equation}
with the after-spike reset $v \leftarrow c$, $u \leftarrow u + d$ whenever
$v \ge 30$~mV. The four parameters $(a,b,c,d)$ select among regular
spiking, bursting, chattering and fast-spiking regimes.
\end{remark}

\subsection{Neural Coding}

Because a spike train is a point process, the question of what the code
\emph{is} becomes explicit. Under \emph{rate coding} the transmitted
quantity is the count obtained by integrating the spike train
\eqref{eq:bio-spiketrain} over a window of length $T$,
\begin{equation}
  r \;=\; \frac{1}{T}\int_{t}^{t+T} S(u)\,\mathrm{d}u
     \;=\; \frac{\bigl|\{f : t \le t^{(f)} < t+T\}\bigr|}{T},
  \label{eq:bio-ratecode}
\end{equation}
which recovers a real-valued activation and makes spiking networks
approximately equivalent to conventional ones, at the price of a latency
and an energy budget proportional to the spike count. Under \emph{temporal}
or time-to-first-spike coding the quantity is the latency $t^{(1)}$ itself,
so a single spike per neuron suffices — the regime in which neuromorphic
hardware achieves its energy advantage. Phase coding (spike time relative
to a background oscillation) and rank-order coding (the permutation in
which a population fires) are further alternatives.

\subsection{Spike-Timing-Dependent Plasticity}

The canonical learning rule for spiking networks is local in both space and
time: a synapse is strengthened when it participates causally in a
postsynaptic spike, and weakened when it fires just after one. Writing
$\Delta t = t_{\mathrm{post}} - t_{\mathrm{pre}}$ for the signed latency of
a pre/post spike pair, the classical additive exponential window is
\begin{equation}
  \Delta w(\Delta t)
  \;=\;
  \begin{cases}
    \phantom{-}A_{+}\exp\!\left(-\dfrac{\Delta t}{\tau_{+}}\right),
      & \Delta t > 0
      \qquad \text{(pre before post: potentiation)},\\[2.4ex]
    -A_{-}\exp\!\left(\dfrac{\Delta t}{\tau_{-}}\right),
      & \Delta t \le 0
      \qquad \text{(post before pre: depression)},
  \end{cases}
  \label{eq:bio-stdp}
\end{equation}
with time constants $\tau_{\pm}$ of the order of $20$~ms.
Equation~\eqref{eq:bio-stdp} is an unsupervised Hebbian rule that requires
no error signal, no target and no backward pass: it is a pure causality
detector, and every quantity it needs is available at the synapse itself.
Its shape is plotted in Figure~\ref{fig:bio-stdp}.

\begin{figure}[htbp]
  \centering
  \begin{tikzpicture}
    \begin{axis}[
        width=12.4cm, height=6.2cm,
        xmin=-80, xmax=80, ymin=-0.85, ymax=1.22,
        axis lines=middle,
        xlabel={$\Delta t = t_{\mathrm{post}} - t_{\mathrm{pre}}$ (ms)},
        ylabel={$\Delta w$},
        xlabel style={font=\small,
                      at={(axis description cs:0.5,-0.02)}, anchor=north},
        ylabel style={font=\small, at={(ticklabel* cs:1.0)},
                      anchor=south east},
        tick label style={font=\scriptsize},
        xtick={-80,-40,40,80}, ytick={-0.6,-0.3,0.3,0.6,0.9},
        clip=false,
      ]
      \addplot[cgate, line width=1.2pt, domain=0.001:80, samples=120]
        {1.0*exp(-x/17)};
      \addplot[coutput, line width=1.2pt, domain=-80:-0.001, samples=120]
        {-0.6*exp(x/34)};
      \draw[cgate,   dotted, thick] (axis cs:0,0) -- (axis cs:0, 1.0);
      \draw[coutput, dotted, thick] (axis cs:0,0) -- (axis cs:0,-0.6);
      \node[etiquette, anchor=west, text=cgate!55!black, align=left]
        at (axis cs:14,0.88)
        {LTP: pre before post,\\[-1pt]$A_{+}e^{-\Delta t/\tau_{+}}$};
      \node[etiquette, anchor=east, text=coutput!85!black, align=right]
        at (axis cs:-14,-0.56)
        {LTD: post before pre,\\[-1pt]$-A_{-}e^{\Delta t/\tau_{-}}$};
    \end{axis}
  \end{tikzpicture}
  \caption{The pair-based STDP window \eqref{eq:bio-stdp}, with
    $A_{+} = 1.0$, $\tau_{+} = 17$~ms, $A_{-} = 0.6$ and
    $\tau_{-} = 34$~ms. A presynaptic spike arriving shortly \emph{before} a
    postsynaptic one strengthens the synapse (green, right half-plane); the
    reverse order weakens it (red, left half-plane). The rule is
    discontinuous at $\Delta t = 0$ and decays to zero for
    $|\Delta t| \gg \tau_{\pm}$, so only near-coincident pairs are
    consequential. The asymmetry $A_{-}\tau_{-} > A_{+}\tau_{+}$ shown here
    makes the net effect of uncorrelated firing depressing, which keeps
    weights from drifting upwards.}
  \label{fig:bio-stdp}
\end{figure}

Implementing \eqref{eq:bio-stdp} literally would require storing all pairs
of spike times. The standard online form instead keeps one low-pass
\emph{eligibility trace} per terminal, decaying exponentially and
incremented at every spike:
\begin{align}
  \tau_{+}\frac{\mathrm{d}x_{\mathrm{pre}}}{\mathrm{d}t}
    &= -\,x_{\mathrm{pre}}
       + \sum_{f}\delta\!\left(t - t^{(f)}_{\mathrm{pre}}\right),
  \qquad
  \tau_{-}\frac{\mathrm{d}x_{\mathrm{post}}}{\mathrm{d}t}
    = -\,x_{\mathrm{post}}
       + \sum_{f}\delta\!\left(t - t^{(f)}_{\mathrm{post}}\right),
  \label{eq:bio-stdp-traces}\\[2pt]
  \frac{\mathrm{d}w}{\mathrm{d}t}
    &= A_{+}\, x_{\mathrm{pre}}(t)\, S_{\mathrm{post}}(t)
     \;-\; A_{-}\, x_{\mathrm{post}}(t)\, S_{\mathrm{pre}}(t).
  \label{eq:bio-stdp-online}
\end{align}
Equations \eqref{eq:bio-stdp-traces}--\eqref{eq:bio-stdp-online} reproduce
\eqref{eq:bio-stdp} exactly for isolated spike pairs while requiring $O(1)$
memory per synapse, which is why they, and not the pairwise window, are
what neuromorphic hardware implements.
Algorithm~\ref{alg:bio-snn-stdp} assembles the discretized simulation loop;
substituting a richer neuron model such as \eqref{eq:bio-izh} changes only
the membrane update on line~4, leaving the plasticity block untouched.

\begin{algorithm}[htbp]
\caption{Clocked simulation of an LIF layer with online STDP}
\label{alg:bio-snn-stdp}
\begin{algorithmic}[1]
\Require step $\Delta t$; horizon $T$; weights
  $\mat{W} \in \R^{n_{\mathrm{post}} \times n_{\mathrm{pre}}}$; decays
  $\beta = e^{-\Delta t/\tau_m}$, $\beta_{\pm} = e^{-\Delta t/\tau_{\pm}}$;
  rates $A_{\pm}$; bound $w_{\max}$
\State $\vect{V} \leftarrow \vect{0}$;\;
       $\vect{S}_{\mathrm{post}} \leftarrow \vect{0}$;\;
       $\vect{x}_{\mathrm{pre}} \leftarrow \vect{0}$;\;
       $\vect{x}_{\mathrm{post}} \leftarrow \vect{0}$;\;
       $\vect{\rho} \leftarrow \vect{0}$
       \Comment{$\vect{\rho}$: refractory counters}
\For{$t = 1, \dots, T/\Delta t$}
  \State $\vect{S}_{\mathrm{pre}} \leftarrow$ presynaptic spikes at step $t$
  \State $\vect{V} \leftarrow \beta\,\vect{V}
          + \mat{W}\,\vect{S}_{\mathrm{pre}}
          - V_{\mathrm{th}}\,\vect{S}_{\mathrm{post}}$
          \Comment{leak, integrate, soft reset \eqref{eq:bio-lif-discrete}}
  \State $\vect{V} \leftarrow \vect{V} \odot \mathbf{1}[\vect{\rho} = 0]$
          \Comment{clamp refractory units \eqref{eq:bio-lif-reset}}
  \State $\vect{S}_{\mathrm{post}} \leftarrow
          \Theta(\vect{V} - V_{\mathrm{th}}) \odot \mathbf{1}[\vect{\rho} = 0]$
          \Comment{\eqref{eq:bio-lif-heaviside}}
  \State $\vect{\rho} \leftarrow \max\!\bigl(\vect{\rho} - 1,\, 0\bigr)
          + \lceil t_{\mathrm{ref}}/\Delta t \rceil\,
            \vect{S}_{\mathrm{post}}$
  \Statex \hspace{\algorithmicindent}\textit{plasticity: traces first, then
          the two-term update}
  \State $\vect{x}_{\mathrm{pre}} \leftarrow
          \beta_{+}\,\vect{x}_{\mathrm{pre}} + \vect{S}_{\mathrm{pre}}$;\quad
         $\vect{x}_{\mathrm{post}} \leftarrow
          \beta_{-}\,\vect{x}_{\mathrm{post}} + \vect{S}_{\mathrm{post}}$
          \Comment{\eqref{eq:bio-stdp-traces}}
  \State $\mat{W} \leftarrow \mat{W}
          + A_{+}\,\vect{S}_{\mathrm{post}}\,\vect{x}_{\mathrm{pre}}\transp
          - A_{-}\,\vect{x}_{\mathrm{post}}\,\vect{S}_{\mathrm{pre}}\transp$
          \Comment{rule \eqref{eq:bio-stdp-online}}
  \State $\mat{W} \leftarrow \mathrm{clip}\bigl(\mat{W},\, 0,\, w_{\max}\bigr)$
          \Comment{additive STDP is otherwise unbounded}
\EndFor
\State \Return $\mat{W}$
\end{algorithmic}
\end{algorithm}

\subsection{Surrogate Gradients}

STDP is unsupervised and cannot by itself solve a supervised task. To train
a deep spiking network by gradient descent one must differentiate
\eqref{eq:bio-lif-heaviside}, and the derivative of the Heaviside step is a
Dirac impulse: zero almost everywhere, unbounded at the threshold. Gradient
descent through it is impossible. The \emph{surrogate gradient} method
keeps the step function in the forward pass but substitutes a smooth,
bounded pseudo-derivative in the backward pass, for instance the
fast-sigmoid derivative
\begin{equation}
  \frac{\partial S}{\partial V}
  \;\approx\;
  \sigma'_{\mathrm{surr}}(V)
  \;=\;
  \frac{1}{\bigl(1 + \gamma\,|V - V_{\mathrm{th}}|\bigr)^{2}},
  \qquad \gamma > 0,
  \label{eq:bio-surrogate}
\end{equation}
which peaks at the threshold and decays away from it, so that only neurons
close to firing receive appreciable credit. With \eqref{eq:bio-surrogate}
in place the discretized network \eqref{eq:bio-lif-discrete} becomes an
ordinary recurrent network over the time axis, trainable by
backpropagation through time (Algorithm~\ref{alg:bio-surrogate}). The
steepness $\gamma$ controls a bias--variance trade-off: large $\gamma$
approaches the true, useless derivative, while small $\gamma$ leaks
gradient to neurons far from threshold.

\begin{algorithm}[htbp]
\caption{Surrogate-gradient training of a deep spiking network}
\label{alg:bio-surrogate}
\begin{algorithmic}[1]
\Require encoded input spike trains $\vect{S}^{(0)}[1{:}T]$; targets
  $\vect{y}$; layers $\ell = 1, \dots, L$; steepness $\gamma$; rate $\eta$
\For{each mini-batch}
  \State $\vect{V}^{(\ell)} \leftarrow \vect{0}$,\;
         $\vect{S}^{(\ell)} \leftarrow \vect{0}$ for all $\ell$
  \For{$t = 1, \dots, T$}
        \Comment{forward pass: the exact step function is retained}
    \For{$\ell = 1, \dots, L$}
      \State $\vect{V}^{(\ell)}[t] \leftarrow
              \beta\, \vect{V}^{(\ell)}[t-1]
              + \mat{W}^{(\ell)} \vect{S}^{(\ell-1)}[t]
              - V_{\mathrm{th}}\, \vect{S}^{(\ell)}[t-1]$
      \State $\vect{S}^{(\ell)}[t] \leftarrow
              \Theta\bigl(\vect{V}^{(\ell)}[t] - V_{\mathrm{th}}\bigr)$
    \EndFor
  \EndFor
  \State $\hat{\vect{y}} \leftarrow
          \softmax\bigl(\textstyle\sum_{t} \vect{V}^{(L)}[t]\bigr)$
          \Comment{or a spike-count readout \eqref{eq:bio-ratecode}}
  \State $\Loss \leftarrow \mathrm{CE}(\hat{\vect{y}}, \vect{y})$
  \State backpropagate through time, replacing every factor
         $\partial \vect{S}^{(\ell)}[t] / \partial \vect{V}^{(\ell)}[t]$
         by $\sigma'_{\mathrm{surr}}$ of \eqref{eq:bio-surrogate}
  \State $\mat{W}^{(\ell)} \leftarrow \mat{W}^{(\ell)}
          - \eta\, \nabla_{\mat{W}^{(\ell)}} \Loss$ for all $\ell$
\EndFor
\end{algorithmic}
\end{algorithm}

% ----------------------------------------------------------------------------
\section{Self-Organizing Maps}
% ----------------------------------------------------------------------------

The self-organizing map \cite{kohonen1990} is an unsupervised, competitive
algorithm that projects a high-dimensional data space onto a
low-dimensional discrete lattice — usually a two-dimensional grid — while
preserving topology: inputs close in $\R^{n}$ are mapped to units close on
the lattice. Unlike the projection produced by an autoencoder, the output
space is prescribed and discrete, and the mapping is learned without any
gradient.

Each unit $i$ of the lattice carries two objects: a fixed position
$\vect{r}_i$ \emph{on the lattice}, and a trainable codebook vector
$\vect{m}_i \in \R^{n}$ living \emph{in the data space}. Learning
alternates two steps for each presented sample $\vect{x}(t)$.

\paragraph{Competition.} The best matching unit (BMU) is the nearest
codebook vector,
\begin{equation}
  c(t) \;=\; \argmin_{i} \bigl\| \vect{x}(t) - \vect{m}_i(t) \bigr\|,
  \label{eq:bio-som-bmu}
\end{equation}
a nearest-prototype rule identical to the assignment step of $k$-means.

\paragraph{Cooperation.} What distinguishes the SOM from $k$-means is that
the winner does not adapt alone. A neighbourhood kernel, measured in
\emph{lattice} coordinates rather than data coordinates, spreads the update
to nearby units:
\begin{equation}
  h_{ci}(t)
  \;=\;
  \exp\!\left(
    -\,\frac{\bigl\|\vect{r}_c - \vect{r}_i\bigr\|^{2}}{2\,\sigma^{2}(t)}
  \right)
  \;\in\; (0, 1],
  \label{eq:bio-som-neigh}
\end{equation}
so that $h_{cc} = 1$ and $h_{ci}$ decays with lattice distance. The
codebook update is then a shrinkage of every unit towards the current
sample, weighted by that kernel:
\begin{equation}
  \vect{m}_i(t+1)
  \;=\;
  \vect{m}_i(t)
  \;+\;
  \alpha(t)\, h_{ci}(t)\,
  \bigl[\vect{x}(t) - \vect{m}_i(t)\bigr].
  \label{eq:bio-som-update}
\end{equation}
Because \eqref{eq:bio-som-update} couples lattice neighbours to a common
data point, units adjacent on the grid are dragged towards adjacent regions
of the data space. This coupling, and not any explicit constraint, is what
produces topology preservation.

Both the learning rate and the neighbourhood radius are annealed,
\begin{equation}
  \alpha(t) = \alpha_0\, e^{-t/\lambda_{\alpha}},
  \qquad
  \sigma(t) = \sigma_0\, e^{-t/\lambda_{\sigma}},
  \label{eq:bio-som-schedule}
\end{equation}
which splits training into two phases with distinct roles. While $\sigma$
is large the kernel \eqref{eq:bio-som-neigh} is nearly flat and the whole
map moves almost rigidly, unfolding into the data cloud and fixing the
global ordering. As $\sigma$ shrinks towards one lattice unit the updates
become local and the phase turns into fine quantization. Ordering, once
achieved, is never undone — which is the practical reason the annealing
must be slow. Figure~\ref{fig:bio-som} shows one such update.

\begin{figure}[htbp]
  \centering
  \begin{tikzpicture}[font=\small]
    % ---- input vector -------------------------------------------------
    \node[ninput, minimum size=9mm] (x) at (-1.1, 4.5) {$\vect{x}$};
    \foreach \r in {1,...,5}{
      \draw[black!35, -{Stealth[length=1.6mm]}] (x) -- (1.02, 1.5*\r);
    }
    \node[etiquette, anchor=north, align=center, text width=2.6cm]
      at (-1.1, 3.9) {compared with\\all units,
                      \eqref{eq:bio-som-bmu}};
    % ---- neighbourhood radius (drawn first, so the opaque unit fills
    %      hide the arcs that would otherwise cross their labels) -------
    \draw[dashed, cgate!70!black, thick] (3.0,3.0) circle (2.25cm);
    \node[etiquette, fill=white, inner sep=1.5pt, text=cgate!50!black]
      at (3.58, 0.83) {$\sigma(t)$};
    % ---- lattice ------------------------------------------------------
    \foreach \c/\r/\h/\f in {%
      1/1/0.64/53, 1/2/0.80/64, 1/3/0.64/53, 1/4/0.33/31, 1/5/0.11/16,
      2/1/0.80/64, 2/3/0.80/64, 2/4/0.41/37, 2/5/0.14/17,
      3/1/0.64/53, 3/2/0.80/64, 3/3/0.64/53, 3/4/0.33/31, 3/5/0.11/16,
      4/1/0.33/31, 4/2/0.41/37, 4/3/0.33/31, 4/4/0.17/20, 4/5/0.06/12,
      5/1/0.11/16, 5/2/0.14/17, 5/3/0.11/16, 5/4/0.06/12, 5/5/0.02/9}{
      \node[neuron, minimum size=9mm, font=\tiny,
            fill=chidden!\f, draw=chidden!80!black]
        (u\c\r) at (1.5*\c, 1.5*\r) {\h};
    }
    % ---- best matching unit -------------------------------------------
    \node[neuron, minimum size=9mm, font=\tiny, fill=chidden!78,
          draw=coutput, line width=1.1pt] (bmu) at (3.0, 3.0) {1.00};
    % ---- profile of the update along the winner's row -----------------
    \foreach \c/\l in {1/0.80, 2/1.00, 3/0.80, 4/0.41, 5/0.14}{
      \draw[cinput!85!black, line width=2.4pt]
        (1.5*\c - \l*0.5, 0.15) -- (1.5*\c + \l*0.5, 0.15);
    }
    \node[etiquette, anchor=east, align=right, text width=2.9cm]
      at (0.5, 0.15) {update size\\$\alpha(t)\,h_{ci}(t)$};
    % ---- legend -------------------------------------------------------
    \node[etiquette, anchor=west, align=left] at (8.2, 3.15)
      {\textcolor{coutput}{thick ring}:\\BMU $c$};
    \node[etiquette, anchor=west, align=left] at (8.2, 2.15)
      {numbers: $h_{ci}(t)$};
  \end{tikzpicture}
  \caption{One SOM update on a $5 \times 5$ lattice. The sample $\vect{x}$
    is compared with every codebook vector, and the winner $c$ (red ring) is
    the nearest one, \eqref{eq:bio-som-bmu}. The number inside each unit is
    its neighbourhood weight $h_{ci}$ from \eqref{eq:bio-som-neigh} with
    $\sigma = 1.5$ lattice units, and the fill intensity is proportional to
    it; the dashed circle marks the radius $\sigma(t)$. The bars at the
    bottom give the profile of the update magnitude
    \eqref{eq:bio-som-update} along the winner's row: the whole
    neighbourhood moves towards $\vect{x}$, by an amount decaying with
    lattice distance. Note that distances in \eqref{eq:bio-som-neigh} are
    measured on the lattice, whereas the competition
    \eqref{eq:bio-som-bmu} is decided in the data space.}
  \label{fig:bio-som}
\end{figure}

A batch variant removes the dependence on presentation order by recomputing
every codebook vector as a kernel-weighted mean of the assigned samples,
\begin{equation}
  \vect{m}_i
  \;=\;
  \frac{\sum_{t} h_{c(t)\,i}\; \vect{x}(t)}
       {\sum_{t} h_{c(t)\,i}},
  \label{eq:bio-som-batch}
\end{equation}
which is to \eqref{eq:bio-som-update} what Lloyd's algorithm is to online
$k$-means. Map quality is assessed by two complementary scalars, the
quantization error and the topographic error,
\begin{equation}
  \mathrm{QE} = \frac{1}{N}\sum_{t=1}^{N}
      \bigl\|\vect{x}(t) - \vect{m}_{c(t)}\bigr\|,
  \qquad
  \mathrm{TE} = \frac{1}{N}\sum_{t=1}^{N}
      \mathbf{1}\!\left[
        \vect{r}_{c_1(t)} \text{ and } \vect{r}_{c_2(t)}
        \text{ non-adjacent}
      \right],
  \label{eq:bio-som-errors}
\end{equation}
where $c_1$ and $c_2$ denote the first and second best matching units. The
pair matters: QE alone is minimized by abandoning topology altogether, and
TE alone by collapsing every codebook vector onto the same point.
Algorithm~\ref{alg:bio-som} states the complete online procedure; replacing
its inner loop by \eqref{eq:bio-som-batch} yields the batch variant.

\begin{algorithm}[htbp]
\caption{Self-organizing map — online training}
\label{alg:bio-som}
\begin{algorithmic}[1]
\Require data $\{\vect{x}(t)\}_{t=1}^{N}$; lattice positions
  $\{\vect{r}_i\}_{i=1}^{M}$; $\alpha_0$, $\sigma_0$, $\lambda_\alpha$,
  $\lambda_\sigma$; epochs $E$
\State initialize $\vect{m}_i$ by random sampling of the data, or on the
       plane of the two leading principal components
\For{$e = 1, \dots, E$}
  \For{each sample $\vect{x}(t)$, in random order}
    \State $c \leftarrow \argmin_{i}\|\vect{x}(t) - \vect{m}_i\|$
           \Comment{competition \eqref{eq:bio-som-bmu}}
    \State $\alpha \leftarrow \alpha_0 e^{-t/\lambda_\alpha}$;\quad
           $\sigma \leftarrow \sigma_0 e^{-t/\lambda_\sigma}$
           \Comment{annealing \eqref{eq:bio-som-schedule}}
    \For{each unit $i = 1, \dots, M$}
      \State $h \leftarrow \exp\bigl(-\|\vect{r}_c - \vect{r}_i\|^{2} /
             (2\sigma^{2})\bigr)$
             \Comment{lattice distance \eqref{eq:bio-som-neigh}}
      \State $\vect{m}_i \leftarrow \vect{m}_i
             + \alpha\, h\, \bigl(\vect{x}(t) - \vect{m}_i\bigr)$
             \Comment{cooperation \eqref{eq:bio-som-update}}
    \EndFor
  \EndFor
  \State monitor QE and TE of \eqref{eq:bio-som-errors}
\EndFor
\State \Return $\{\vect{m}_i\}_{i=1}^{M}$
\end{algorithmic}
\end{algorithm}

% ----------------------------------------------------------------------------
\section{Reservoir Computing and Echo State Networks}
% ----------------------------------------------------------------------------

Training a recurrent network by backpropagation through time is expensive
and numerically delicate. Reservoir computing takes a radical shortcut:
keep a large, sparse, randomly initialized and permanently frozen recurrent
network — the reservoir — and train only a linear readout on top of its
states \cite{jaeger2001}. The reservoir serves as a fixed,
high-dimensional, nonlinear expansion of the input history, and the
learning problem that remains is linear regression.

The echo state network state update, with leak rate $\alpha \in (0,1]$, is
\begin{equation}
  \vect{h}_t
  \;=\;
  (1 - \alpha)\,\vect{h}_{t-1}
  \;+\;
  \alpha\,\tanh\!\bigl(
    \mat{W}_{\mathrm{in}}\,\vect{u}_t
    + \mat{W}\,\vect{h}_{t-1}
    + \vect{b}
  \bigr),
  \label{eq:bio-esn-state}
\end{equation}
and the trained readout is the affine map
\begin{equation}
  \vect{y}_t \;=\; \mat{W}_{\mathrm{out}}\,
  \begin{bmatrix} 1 \\ \vect{u}_t \\ \vect{h}_t \end{bmatrix}.
  \label{eq:bio-esn-readout}
\end{equation}
Only $\mat{W}_{\mathrm{out}}$ is learned: $\mat{W}_{\mathrm{in}}$,
$\mat{W}$ and $\vect{b}$ are drawn once at random and never touched. For
this to be sensible the reservoir must forget, so that its state is
asymptotically a function of the input history alone and not of an
arbitrary initial condition.

\begin{property}[Echo state property]\label{prop:bio-esp}
The reservoir \eqref{eq:bio-esn-state} has the echo state property when,
for any two initial states $\vect{h}_0$ and $\vect{h}'_0$ driven by the
same input sequence, $\|\vect{h}_t - \vect{h}'_t\| \to 0$ as
$t \to \infty$. A sufficient condition is $\sigma_{\max}(\mat{W}) < 1$; in
practice one uses the far less conservative spectral-radius heuristic
\begin{equation}
  \rho(\mat{W}) \;=\; \max_{i} \bigl|\lambda_i(\mat{W})\bigr| \;<\; 1,
  \label{eq:bio-esn-spectral}
\end{equation}
tuning $\rho$ just below unity. The reservoir is then a \emph{fading
memory} filter: larger $\rho$ lengthens the memory, smaller $\rho$ shortens
it, and $\rho$ slightly above one often still works on strongly driven
tasks.
\end{property}

Condition \eqref{eq:bio-esn-spectral} is enforced by construction rather
than by optimization. One draws a sparse random matrix
$\mat{W}_{\mathrm{raw}}$ — typically one to five per cent of entries
nonzero — computes its spectral radius numerically, and rescales:
\begin{equation}
  \mat{W}
  \;\leftarrow\;
  \frac{\rho^{\star}}{\rho(\mat{W}_{\mathrm{raw}})}\; \mat{W}_{\mathrm{raw}},
  \qquad \rho^{\star} \lesssim 1 .
  \label{eq:bio-esn-rescale}
\end{equation}

The decisive consequence of freezing the recurrence is that the training
objective becomes convex. Collecting the extended states of all retained
time steps into $\mat{H} \in \R^{(1 + n_u + n_h) \times T}$ and the targets
into $\mat{Y} \in \R^{n_y \times T}$, the ridge-regularized least-squares
solution is available in closed form:
\begin{equation}
  \mat{W}_{\mathrm{out}}
  \;=\;
  \argmin_{\mat{W}}
  \left\{
    \bigl\|\mat{W}\mat{H} - \mat{Y}\bigr\|_{F}^{2}
    + \lambda \bigl\|\mat{W}\bigr\|_{F}^{2}
  \right\}
  \;=\;
  \mat{Y}\mat{H}\transp
  \bigl(\mat{H}\mat{H}\transp + \lambda \mat{I}\bigr)^{-1}.
  \label{eq:bio-esn-ridge}
\end{equation}
There is no local minimum, no learning rate, no vanishing gradient and no
backward pass through time: one matrix inversion of size
$(1 + n_u + n_h)$ replaces the entire training procedure. The ridge term is
not optional — reservoir states are strongly correlated, so
$\mat{H}\mat{H}\transp$ is badly conditioned.
Algorithm~\ref{alg:bio-esn} gives the complete recipe, and
Figure~\ref{fig:bio-esn} the architecture.

\begin{figure}[htbp]
  \centering
  \begin{tikzpicture}[font=\small]
    % ---- reservoir boundary -------------------------------------------
    \draw[dashed, black!55, thick]
      (4,0) ellipse [x radius=2.4cm, y radius=1.8cm];
    % ---- reservoir units ----------------------------------------------
    \node[nhidden, minimum size=6mm] (n1) at (2.55, 0.55) {};
    \node[nhidden, minimum size=6mm] (n2) at (2.55,-0.65) {};
    \node[nhidden, minimum size=6mm] (n3) at (3.60, 1.15) {};
    \node[nhidden, minimum size=6mm] (n4) at (3.75,-1.15) {};
    \node[nhidden, minimum size=6mm] (n5) at (4.10, 0.05) {};
    \node[nhidden, minimum size=6mm] (n6) at (4.95, 1.05) {};
    \node[nhidden, minimum size=6mm] (n7) at (5.15,-0.85) {};
    \node[nhidden, minimum size=6mm] (n8) at (5.55, 0.10) {};
    % ---- fixed sparse recurrent connections ---------------------------
    \begin{scope}[black!40, -{Stealth[length=1.6mm]}, semithick]
      \draw (n1) -- (n5);
      \draw (n5) to[bend left=15]  (n3);
      \draw (n3) -- (n6);
      \draw (n6) to[bend left=12]  (n8);
      \draw (n8) -- (n7);
      \draw (n7) -- (n4);
      \draw (n4) to[bend left=12]  (n2);
      \draw (n2) -- (n1);
      \draw (n5) to[bend right=15] (n7);
      \draw (n6) to[bend left=18]  (n5);
      \draw (n2) to[bend right=12] (n5);
      \draw (n5) to[out=-115, in=-65, looseness=6] (n5);
    \end{scope}
    % ---- input --------------------------------------------------------
    \node[ninput, minimum size=9mm] (u) at (0,0) {$\vect{u}_t$};
    \draw[fleche, black!55] (u) -- (n1);
    \draw[fleche, black!55] (u) -- (n2);
    \node[etiquette, anchor=east] at (1.68, 0.80)
      {$\mat{W}_{\mathrm{in}}$ (fixed)};
    % ---- readout ------------------------------------------------------
    \node[noutput, minimum size=9mm] (y) at (8.7,0) {$\vect{y}_t$};
    \draw[fleche, coutput, line width=1.1pt] (n6) -- (y);
    \draw[fleche, coutput, line width=1.1pt] (n8) -- (y);
    \draw[fleche, coutput, line width=1.1pt] (n7) -- (y);
    \node[etiquette, anchor=south, text=coutput!85!black] at (7.4, 0.80)
      {$\mat{W}_{\mathrm{out}}$ (\textbf{trained})};
    % ---- titles -------------------------------------------------------
    \node[anchor=south, font=\small] at (4, 2.05)
      {Reservoir: $n_h$ fixed random units};
    \node[etiquette, anchor=north, align=center] at (4, -2.05)
      {$\mat{W}$ sparse, random, rescaled to
       $\rho(\mat{W}) = \rho^{\star} < 1$, \eqref{eq:bio-esn-rescale}};
  \end{tikzpicture}
  \caption{An echo state network. The input weights
    $\mat{W}_{\mathrm{in}}$ and the sparse recurrent matrix $\mat{W}$ (grey)
    are drawn at random and \emph{never trained}; they merely expand the
    input history into a high-dimensional dynamical state
    \eqref{eq:bio-esn-state}, whose fading memory is guaranteed by the
    spectral-radius condition \eqref{eq:bio-esn-spectral}. Only the readout
    $\mat{W}_{\mathrm{out}}$ (red) is learned, and it is learned in closed
    form by ridge regression \eqref{eq:bio-esn-ridge}. For legibility only
    three readout connections are drawn; in fact every reservoir state, and
    the input itself, feeds the readout \eqref{eq:bio-esn-readout}.}
  \label{fig:bio-esn}
\end{figure}

\begin{algorithm}[htbp]
\caption{Echo state network — construction and closed-form training}
\label{alg:bio-esn}
\begin{algorithmic}[1]
\Require sequence $\{(\vect{u}_t, \vect{y}^{\star}_t)\}_{t=1}^{T}$;
  reservoir size $n_h$; sparsity $s$; target radius $\rho^{\star}$; leak
  $\alpha$; ridge $\lambda$; washout $T_0$
\Statex \textit{1. build the reservoir, once, at random}
\State draw $\mat{W}_{\mathrm{in}} \in \R^{n_h \times n_u}$ with i.i.d.\
       entries uniform on $[-a, a]$
\State draw $\mat{W}_{\mathrm{raw}} \in \R^{n_h \times n_h}$ with a fraction
       $s$ of nonzero i.i.d.\ entries
\State $\mat{W} \leftarrow
       \bigl(\rho^{\star}/\rho(\mat{W}_{\mathrm{raw}})\bigr)\,
       \mat{W}_{\mathrm{raw}}$
       \Comment{enforce \eqref{eq:bio-esn-spectral} by
                \eqref{eq:bio-esn-rescale}}
\Statex \textit{2. drive the reservoir and harvest its states}
\State $\vect{h}_0 \leftarrow \vect{0}$
\For{$t = 1, \dots, T$}
  \State $\vect{h}_t \leftarrow (1-\alpha)\,\vect{h}_{t-1}
          + \alpha \tanh\bigl(\mat{W}_{\mathrm{in}}\vect{u}_t
          + \mat{W}\vect{h}_{t-1} + \vect{b}\bigr)$
          \Comment{\eqref{eq:bio-esn-state}}
  \If{$t > T_0$}
    \State append $[\,1;\, \vect{u}_t;\, \vect{h}_t\,]$ as a column of
           $\mat{H}$, and $\vect{y}^{\star}_t$ as a column of $\mat{Y}$
  \EndIf
\EndFor
\Statex \textit{3. solve one convex problem}
\State $\mat{W}_{\mathrm{out}} \leftarrow \mat{Y}\mat{H}\transp
        \bigl(\mat{H}\mat{H}\transp + \lambda\mat{I}\bigr)^{-1}$
        \Comment{ridge regression \eqref{eq:bio-esn-ridge}}
\State \Return $\mat{W}_{\mathrm{in}}$, $\mat{W}$, $\mat{W}_{\mathrm{out}}$
\end{algorithmic}
\end{algorithm}

The discarded prefix of length $T_0$ in Algorithm~\ref{alg:bio-esn} is not
a detail: by Property~\ref{prop:bio-esp} the influence of $\vect{h}_0$
decays only asymptotically, so the earliest states are contaminated by the
arbitrary initial condition and would bias the regression. Echo state
networks are competitive on chaotic time-series prediction, at a training
cost orders of magnitude below that of a fully trained recurrent network.

\begin{remark}[Liquid state machines]\label{rem:bio-lsm}
The liquid state machine is the spiking counterpart of the echo state
network: the reservoir is a recurrent population of leaky
integrate-and-fire neurons in the sense of Definition~\ref{def:bio-lif} —
the liquid — and the readout is trained on its instantaneous state.
Real-time universal computing power follows from two conditions: a
\emph{separation property} for the liquid, meaning that different input
histories drive measurably different states, and an \emph{approximation
property} for the readout. The idea has outlived its original setting in
\emph{physical} reservoir computing, where the fixed nonlinear dynamics are
supplied by photonic, memristive or even mechanical substrates rather than
by simulated neurons.
\end{remark}

% ----------------------------------------------------------------------------
\section{What Is Learned, What Is Fixed, What Emerges}
% ----------------------------------------------------------------------------

The three families of this chapter occupy three distinct positions on the
question of where the adaptivity of a network should reside;
Table~\ref{tab:bio-comparison} summarizes them.

\begin{table}[htbp]
  \centering
  \small
  \caption{The three biologically inspired paradigms compared. None of them
    uses a global gradient of a supervised loss in its native form.}
  \label{tab:bio-comparison}
  \begin{tabular}{@{}lllll@{}}
    \toprule
    & \textbf{Fixed} & \textbf{Learned} & \textbf{Learning signal}
    & \textbf{Emergent property} \\
    \midrule
    SNN with STDP
      & neuron dynamics
      & all synapses
      & local spike timing
      & causality detection \\
    SOM
      & lattice topology
      & codebook $\vect{m}_i$
      & competition, no error
      & topology preservation \\
    ESN
      & $\mat{W}_{\mathrm{in}}$, $\mat{W}$
      & readout only
      & convex least squares
      & fading memory \\
    \bottomrule
  \end{tabular}
\end{table}

\begin{remark}[The common thread]\label{rem:bio-thread}
Each paradigm buys something by giving up end-to-end differentiability.
Spiking networks buy energy: under temporal coding a neuron may transmit
its message with a single event, which is why these models map so naturally
onto event-driven neuromorphic silicon. Self-organizing maps buy
interpretability: the lattice is a fixed, inspectable output space, so the
learned representation can be read off directly rather than probed.
Reservoir computing buys convexity, and with it a training procedure that
is exact, fast and free of the pathologies of backpropagation through time.
The price, uniformly, is expressive power on large supervised problems —
precisely where the gradient-trained architectures of the preceding
chapters remain unmatched. The surrogate-gradient method of
\eqref{eq:bio-surrogate} is the most successful attempt so far to have
both, restoring differentiability to the spiking forward pass without
abandoning its event-driven semantics.
\end{remark}