L ONG

SHORT- TERM MEMORY AND L EARNING - TO LEARN IN NETWORKS OF SPIKING NEURONS

arXiv:1803.09574v1 [] 26 Mar 2018

Guillaume Bellec*, Darjan Salaj*, Anand Subramoney*, Robert Legenstein & Wolfgang Maass
Institute for Theoretical Computer Science
Graz University of Technology
Austria
{bellec,salaj,subramoney,legenstein,maass}@igi.tugraz.at
* = first authors

A BSTRACT
Networks of spiking neurons (SNNs) are frequently studied as models for networks of neurons in the brain, but also as paradigm for novel energy efficient
computing hardware. In principle they are especially suitable for computations
in the temporal domain, such as speech processing, because their computations
are carried out via events in time and space. But so far they have been lacking
the capability to preserve information for longer time spans during a computation,
until it is updated or needed - like a register of a digital computer. This function is
provided to artificial neural networks through Long Short-Term Memory (LSTM)
units. We show here that SNNs attain similar capabilities if one includes adapting
neurons in the network. Adaptation denotes an increase of the firing threshold of
a neuron after preceding firing. A substantial fraction of neurons in the neocortex
of rodents and humans has been found to be adapting. It turns out that if adapting
neurons are integrated in a suitable manner into the architecture of SNNs, the performance of these enhanced SNNs, which we call LSNNs, for computation in the
temporal domain approaches that of artificial neural networks with LSTM-units.
In addition, the computing and learning capabilities of LSNNs can be substantially enhanced through learning-to-learn (L2L) methods from machine learning,
that have so far been applied primarily to LSTM networks and apparently never
to SSNs.
This preliminary report on arXiv will be replaced by a more detailed version in
about a month.

1

I NTRODUCTION

Deep learning in artificial neural networks (ANNs) usually employs feedforward networks. Feedforward networks tend to support computing and learning mechanisms that differ strongly from those
of recurrently connected networks of neurons, which one finds in the brain. Recurrent ANNs have
been applied very successfully for demanding tasks such as speech processing and video recognition,
which require integration of information over time. But virtually all of these applications involve recurrent ANNs that employ specialized modules for preserving information over time, such as Long
Short-Term memory (LSTM) units (Hochreiter & Schmidhuber, 1997). A simpler type of memory
module had been introduced for artificial neural networks by Mikolov et al. (2014). Modules of this
type have been missing in spiking neural network (SNN) models, and as a result, their performance
for tasks that require integration over longer timespans has remained inferior. A first step for alleviating this problem was made through the inclusion of short-term synaptic plasticity in randomly
connected SNNs (Maass et al., 2002). This improved their capabilities for tasks that require integration of information over hundreds of ms, but did not solve the problem for tasks that involve the
behaviourally relevant time scale of seconds and longer.
Neurons in the brain are subject to intracellular processes on many different times scales, from ms
to s, minutes, hours, and days. We show that by including just one of their longer lasting processes,
neuronal adaptation, and by simultaneously considering structured rather than randomly connected
1

network architectures, one arrives at recurrent SNNs with substantially enhanced computing and
learning capabilities. We will refer to these SNNs with longer short-term memory as LSNNs.

2

T HE LSNN

MODEL

We use simple leaky integrate-and-fire neuron models (see Supplementary Information). One population R of regularly firing neurons in the LSNN consists of standard neurons of this type, whereas
a simple adaptation mechanism is included for neurons in a second subpopulation A. Both populations receive spike trains from a population X of external input neurons. Results of computations in
the LSNN are read out by a population Y of external linear readout neurons, see Fig. 1A.
In neuroscience the term neural adaptation refers to a broad range of phenomena and mechanisms. It
was already discovered in the 19th century that the intensity of a constant stimulus tends to cause a
diminishing sensation, see e.g. the Wikipedia article “Neural adaptation”. A closer look shows that
these phenomena take place simultaneously on several time scales, from hundreds of ms to s, hours,
or even days. In the context of neuron models, adaptation refers to changes in the firing threshold or
subthreshold currents of a neuron that reduce its excitability in response to its own firing (Gerstner
et al., 2014). Experimental data on that are documented for example in the form of a distribution of
measured adaptation indices of neurons for large sets of neurons, both from mouse and humans, in
the online Allen Brain Atlas (Allen Institute, 2018). This adaptation index measures the increase of
interspike intervals of a neuron that receives a constant input current. Common ways for modeling
such data are described in (Gerstner et al., 2014; Pozzorini et al., 2015; Gouwens et al., 2018; Teeter
et al., 2018). We are using here the arguably simplest form of integrating adaptation into a neuron
model: We assume that the firing threshold Bj (t) of neuron j increases by some fixed amount β for
each spike of this neuron j, and then decays exponentially back to a baseline value b0j with a time
constant τa,j . The neural dynamics for a discrete time step of δt = 1 reads as follows
Bj (t)
bj (t + 1)

= b0j + βbj (t),
=

ρj bj (t) + (1 − ρj )zj (t),

(1)
(2)

where ρj = exp(− τδt
) and zj (t) is the spike train of neuron j. Note that this dynamics of thresha,j
olds of adaptive spiking neurons is similar to the dynamics of the state of context neuron in (Mikolov
et al., 2014).
We optimize the weights, and in some cases also the connections of an LSNN through standard
optimization algorithms for specific ranges of tasks. This optimization algorithm is not claimed to be
biologically realistic. But it might correspond functionally to evolutionary processes that optimized
LSNNs in the brain for tasks that were relevant for survival. For the three experiments that we
report here we used backprop through time (BPTT). BP has already been applied by (Courbariaux
et al., 2016) and (Esser et al., 2016) to feedforward networks of spiking neurons, and by (Huh &
Sejnowski, 2017) to recurrent networks of spiking neurons. The gradient is backpropagated through
spikes by replacing in these approaches the non-existent derivative of the membrane potential at
the time of a spike by a pseudo-derivative that smoothly increases from 0 to 1, and then decays
back to 0. For our applications we needed to modify the pseudo-derivative in order to achieve
stable performance of BPTT for larger time spans, that required backpropagation through several
1000 layers of an unrolled feedforward network of spiking neurons: we reduced (“dampened”) the
amplitude of the pseudo-derivative by a factor < 1.

3

L ONG SHORT- TERM MEMORY IN LSNN S

Store-Recall task:
We first test the performance of LSNNs on a very simple working memory task (see Fig. 1D): The
network has to store a bit – encoded by spikes from input neurons – in working memory if it receives
simultaneously a “store” instruction, which is given in the form of spikes from other input neurons.
The stored bit has to be recalled by the network after a subsequent “recall” instruction, that is also
given in the form of spikes from dedicated input neurons. At the same time other spike inputs should
be ignored by the network.
2

Figure 1: Store-recall task. A An LSNN consists of a regular spiking and an adaptive neural
population. B A standard LSNN architecture that corresponds to the ANN architecture of Mikolov
et al. (2014) C A sample LSNN architecture that resulted from using Deep R, a biologically inspired
technique that prunes and rewires gradually the connectivity matrices during training (Bellec et al.,
2018). D The raster plot shows the activity of a network after convergence. From top to bottom,
the sub-panels show spike rasters from input neurons (X), regular spiking neurons (R), and adaptive
neurons (A) followed by the output of the readout and the dynamics of the adaptive thresholds of
adaptive neurons (this network is based on the standard architecture).

3

More precisely, the input neurons marked in blue and red at the top of Fig. 1D are Poisson neurons
that produce random spike sequences. An input bit to the network is encoded by spiking activity at
50 Hz in the corresponding input channels X0 or X1 for a pattern duration D = 200 ms. “Store”
and “recall” instructions are analogously encoded through firing of populations XS and XR of other
Poisson input neurons at 50 Hz for D = 200 ms. Otherwise input neurons are assumed to be silent
in this simple experiment.
When a “recall” instruction is received, the network should output the bit value received concurrently
with the most recent “store” instruction before that. After a “store”, a “recall” instruction is given
during a subsequent time period of length D with probability π. This results in an expected delay
of D
π between a “store” and “recall” instruction. The next “store” appears similarly in a subsequent
period of length D with probability π.
During a recall cue, the binary output is supposed to be encoded by the sign of the time-averaged
membrane potential of a linear readout neuron Y . The LSNN consisted of 60 regular spiking and
20 adaptive neurons. An example trial after training is shown in Fig. 1D.
A summary of results on performance of trained networks for several versions of the store-recall
tasks is given in Tables 1 and 2 of the Supplement. A SNN without adapting neurons can solve the
store-recall task only if working memory for a few 100 ms suffices. LSNNs work best if the time
constants of their adapting neurons are in the same range as the required durations of working memory in the network. If the adapting neurons of an LSNN have a spread of different time constants,
such LSNN generally performs well for a variety of tasks that require working memory durations
with a similar spread of memory durations. Table 2 shows that the architecture of an LSNN also
influences its performance for these tasks. In general, the standard architecture (Fig. 1 B) works
quite well. But first tests where the network architecture was optimized together with its parameters
suggest that more specific architectures, such as the one shown in Fig. 1 C may perform substantially better for specific tasks (see lowest row of Table 2). We applied there the rewiring method
Deep R (Bellec et al., 2018). Deep R converges theoretically to an optimal network configuration by
continuously updating the set of active connections in a biologically inspired manner (Kappel et al.,
2015; 2017; Bellec et al., 2018).
Sequential MNIST: We next tested the performance of LSNNs on a task that requires continuous
updates of working memory over a longer time span: sequential MNIST. The task is to classify the
handwritten digits of the MNIST datasets when the pixels are presented one after another.
We compare the performance of LSNNs with that of LSTM networks. For fair comparison, the size
of the LSNN is chosen to match the number of parameters of the LSTM network. This leads to
140 regular spiking and 100 adaptive neurons with adaptation time constant 700 ms in the LSNN,
in comparison to 128 LSTM units.
The pixels are presented sequentially in 784 time steps, after which the network is required to output
the class of the presented digit. For LSNNs, the gray level of each pixel is encoded into the spikes
of 80 input neurons. Several different ways of encoding pixel values by spikes can be used. Fig. 2C
shows a version where the gray level of each pixel is encoded by population coding through the firing
probability of the 80 input neurons, 1ms for each pixel. Somewhat better performance is achieved
when each of the 80 input neurons is associated with a particular threshold for the gray value, and
this input neuron fires whenever the gray value of the currently presented pixel crosses its threshold
(two input neurons were used per threshold, one spiked at threshold crossings from below, and one
at the crossings from above; this input convention is chosen for the LSNN results of Fig. 2B). In
either case, an additional input neuron becomes active when the presentation of the 784 pixel values
is finished, in order to prompt an output from the network. The firing of this additional input neuron
is shown at the top right of the top panel of Fig. 2C. The softmax of 10 linear output neurons Y
is trained to produce, during that time segment, the label of the sequentially presented handwritten
digit (see yellow shading around 800 ms of the output neuron for label 6 in the plot of the dynamics
of the output neurons Y in Fig. 2C; this output was correct). In contrast to the store-recall task, we
see that the firing thresholds of adapting neurons undergo complex temporal modifications.
The output of the network is determined by averaging the readout output over the 56 ms following
the presentation of the digit. The network is trained by minimizing the cross entropy between the
softmax of the averaged readout and the label distributions.
4

Figure 2: Sequential MNIST. A The goal of this task is to classify images of handwritten digits
when the pixels are shown sequentially row by row. B The performance of spiking neural networks
is compared for three different setups: without adapting neurons (LIF), an LSNN in the standard
architecture, but with connections from the regular to the adaptive spiking pool (LSNN), and an
LSNN with randomly initialized connectivity that was rewired during training (LSNN rewired). For
comparison, the performance of an artificial fully connected RNN and an LSTM with standard gray
level inputs are shown. C Dynamics of a network classifying a digit as in Figure 1.

5

A performance comparison is given in Figure 2B. LSNNs achieve a 93% classification accuracy on
the test set using the standard architecture with added connections from the regular to the adaptive
spiking pool. A higher performance of 96% is reached if the network is rewired during training
with DEEP R (Bellec et al., 2018). The same performance is achieved by an LSTM receiving
continuous gray levels instead of spikes. The first and fourth bars in Fig. 2B show that this range
of accuracy is out of reach for spiking or nonspiking neural networks without enhanced working
memory capabilities.
Each image presentation with this setup takes 784 ms. To test if the model can also integrate information over a longer time span we expanded this duration to 3.36 s by presenting each pixel for 4
milliseconds instead of 1. An LSNN trained on this expanded presentation duration achieved in a
first experiment 83% classification accuracy (all adaptive neurons had an adaptation time constant
of 3500 ms).

4

L EARNING - TO -L EARN FOR LSNN S

Some of the most intriguing applications of LSTM networks were demonstrations of Learning-toLearn (L2L) or Meta-learning capabilities of this type of neural networks. This research direction
was initiated by (Hochreiter et al., 2001), and has lead in recent years to a host of remarkable
new methods and results, see e.g. (Wang et al., 2016; Gupta et al., 2018). The reason why LSTM
networks have been used for most work on L2L is that they can naturally accommodate two levels of
learning and representation of learned insight: network parameters such as synaptic weights can be
used to encode “innate” or transfer knowledge from preceding learning of more-or-less related tasks.
At the same time, a sufficiently long short-term memory in the neural network can accumulate and
apply knowledge gained during a single learning iteration (generally consisting of many episodes)
for a specific learning task. Of course, many other schemes for encoding, fusing, and applying
results of long-term meta-learning and results from the current learning task are possible and should
be explored.
The results of (Hochreiter et al., 2001) show that besides backprop there are likely to be many other
methods which recurrent neural networks can use to learn nonlinear functions from a teacher. This
supervised learning scenario is also of interest from the perspective of understanding learning in
the brain, since learning to predict future sensory inputs, the consequences of an action, or sensory
input from one modality in terms of sensory inputs from other modalities, are all learning scenarios
where the target output values are given to networks of neurons in the brain – with often some delay
– from the environment or internal processing modules. Another aspect of this L2L approach is
that it produces a host of new ideas and algorithms for implementing learning in recurrent neural
networks. This is of interest from the perspective of neuroscience, because on the basis of the
currently available experimental data it is not clear to what extent synaptic weights are able to
support fast learning by storing knowledge from previous trials for the current learning task. It is
also not clear whether spike-timing-dependent plasticity (STDP) is the primary local learning rule in
neural networks of the brain, and which other local learning rules are there implemented (Lisman &
Spruston, 2010; Titley et al., 2017). In addition, experimental data for STDP question its usefulness
for fast and reliable learning, since STDP is subject to high synapse-to-synapse variability, and
requires dozens of repetitions of pairing of pre- and postsynaptic firing (Froemke et al., 2010).
We show here that the L2L paradigm of Hochreiter et al. (2001) for LSTM networks can be ported to
the biologically more realistic LSNNs. L2L results for SNNs have the advantage that many features
and fingerprints of resulting new learning algorithms for SNNs, such as neural activity (spike trains),
changes in the excitability of neurons and synaptic weights, and changes in the network dynamics
through learning can be compared with corresponding measurements from networks of neurons
in the brain. At the same time, the resulting new methods for learning nonlinear functions are
also of interest for spike-based neuromorphic chips such as Brainscales (Schemmel et al., 2010),
SpiNNaker (Furber et al., 2013), True North (Esser et al., 2016), chips from ETH Zürich such as
(Qiao et al., 2015), and Loihi (Davies et al., 2018). Nonlocal learning rules such as backprop are
real challenges for these neuromorphic devices, and alternative local rules for learning nonlinear
functions are therefore of particular interest.
The standard setup of L2L involves a large, possibly infinitely large, family F of learning tasks
C. Some parameters P of a neural network N , in our case an LSNN, are optimized in an outer
6

optimization loop to support fast learning of a randomly drawn task C from F. Like in (Hochreiter
et al., 2001) we let all synaptic weights W of N belong to the set of parameters P that are optimized
through the outer loop, i.e., through meta-learning. Hence the network is forced to encode all results
from learning the current task C in its internal state, in particular in its short-term memory. Thus the
synaptic weights W of the neural network N encode an efficient algorithm for learning arbitrary
tasks C from F, rather than results of learning a particular task C. The teacher input for learning a
particular function C from F is given in (Hochreiter et al., 2001) in a delayed manner: The target
output value of C is given after N has provided its guessed output value for the preceding input
values. This avoids that N learns to cheat by simply using the teacher input as its guessed output
value. We are using the same convention in our experiment.
The learning methods that can possibly be encoded by the synaptic weights W range over a substantially larger range of local update rules than commonly considered for network learning. Note
for example that most commonly considered local update rules are formulated as functions where
a parameter update (e.g., a weight update) depends in a particular but smooth manner on previous
parameter values and specific aspects of the recent network activity, such as spike times or rates of
particular neurons. As long as this function is smooth, it can be approximated by a neural network.
But importantly, ANY smooth update function – even a function that was previously never considered for neural network learning can be approximated by a neural network. The L2L approach of
(Hochreiter et al., 2001) allows in principle in the outer optimization loop a search over the best
performing ones among all these potential local update rules. They are only restricted to those that
can be implemented by the neural network N . In particular, they are constrained to local learning
rules if the neural network N only supports local rules for computing and changing the content of
short-term memories, as is the case for the LSNNs N that we are considering here. For example,
the network N is free to use the delayed teacher input for computing the error that it made on the
preceding step, and to transmit this error through its internal spiking activity to other neurons in N ,
in particular also to adapting neurons.
One of the more challenging family F of tasks considered by Hochreiter et al. (2001) for LSTM
networks was the class of all quadratic functions over two real-valued variables. We found that
L2L for this family F works also very well for LSNNs. We present here instead a potentially even
more interesting variation of this L2L task, where F consists of all functions of two real-valued
variables that can be realized by a minimal class of artificial neural networks: small two-layer target
neural networks TN, i.e., feedforward networks with one hidden layer, consisting of sigmoidal (nonspiking) neurons. We found that LSNNs can learn to learn quickly to reproduce the input-output
behavior of any such ANN TN, by just using their network dynamics and their internal states for
learning and emulating this ANN. Note that BPTT was only used for the outer learning loop, in
order to determine the synaptic weights W of the LSNN. But these could not be optimized for a
particular target network TN. Constraints of compute time did not yet allow us to emulate learning
for larger target networks TN, and we also did not try to optimize the architecture of the LSNN for
this task. Hence the results below provide just a proof of concept.
The target network TN represents a function C(x1 , x2 ), with 2 inputs x1 and x2 and 2 hidden
neurons, with a sigmoid activation at both the hidden and output layers (shown in Fig. 3A). The
outputs of this function, which serve as the targets for the output Y of the LSNN, are scaled to be
between [0, 1].

7

Figure 3: (Caption next page.)

8

Figure 3: Learning to Learn for SNNs. A Illustration of the two-layer feed-forward target network
(TN) used to generate the targets for the LSNN to learn. B The architecture of the LSNN used for
the learning to learn task. It received as inputs at each step x1 , x2 and the target from the previous
step C(x01 , x02 ). Its prediction Ĉ(x1 , x2 ) was output through a linear readout. C Performance during
training and testing. A new function (target network TN) was randomly chosen for each episode.
The LSNN is enabled by its weights W to learn approximations to the input-output behavior of
target networks TN that was significantly better than the best linear approximation to its inputoutput behavior. D The first 100 steps in a sample test episode where the network initially had high
error, but after few steps learned to predict the new targets successfully. E The correlation of the
output of the LSNN vs target outputs at the start, middle and end of training. F Sample spike raster
during testing, of the input neurons (every third neuron, top panel), and the regular (middle panel)
and adaptive neuron populations (bottom panel) for the same episode as in (D). The point where the
error goes to 0 is marked by a red cross corresponding to the one in (D).

The training procedure was as follows: Network training was divided into training episodes. At
the start of each training episode, a new set of input pairs x1 , x2 ∈ [−1, 1] and new weights
(between [−1, 1]) for the target network TN were randomly chosen and used to generate targets
C(x1 , x2 ) ∈ [0, 1]. These input pairs and targets were used as training data and presented to the
LSNN during the episode. The LSNN parameters were updated using BPTT to minimize the mean
squared error between the LSNN output and the target in the training set, using gradients computed
over batches of 10 such episodes. In other words, each weight update included gradients calculated
on the input/target pairs from 10 different TNs, which prevented the LSNN from specializing on
predicting the output of one single TN.
In each episode, for the given TN, 500 different input/target pairs were presented, one per step.
Each step lasted for 20 timesteps, which corresponds to 20 ms. The LSNN received at each step
three input values as input: x1 , x2 and C(x01 , x02 ), where x01 , x02 were the inputs for the preceding
step (i.e., C(x01 , x02 ) was the target of the previous step; it was set to 0 in the first step) during both
the training and testing process, as illustrated in Fig. 3B. While testing, the weights of the LSNN
remained fixed, and it was required to predict targets from functions Ctest ∈ F that it had never
seen before.
The input population consisted of three subpopulations of 100 neurons each, one population for each
of the inputs x1 , x2 and C(x01 , x02 ). These input populations produced Poisson spike trains based on
a Gaussian population-rate encoding of the input values, between 0 and 200 Hz. The range of each
input was split evenly among all the neurons responsible for that input, and for each given analog
input value, the firing rate of the neurons had a Gaussian activity profile with mean at that analog
value, and a variance of 0.1.
We used an LSNN architecture that contained, in addition to the connections of the standard architecture, connections from the regular-spiking to the adaptive population (see Fig. 3B), since performance of the standard architecture was much worse. The network contained 100 in each of the
regular spiking and adaptive populations. A sample spike raster is shown in Fig. 3F (which is for
the same episode as shown in Fig. 3D). The top panel shows the input values.
During training, there was an initial steep drop in error (see Fig. 3C) when the network learned to
output values in the same range as the target functions, see Fig. 3E left and middle. The performance
of the network was compared to a linear baseline: for each TN, an optimal linear regressor was
computed using the same input/target pairs presented to the LSNN. The mean squared error for the
predictions of this linear regressor on these input/target pairs was the linear baseline.
During testing, even though the network was shown target functions Ctest (x1 , x2 ) it had not seen
before, the test error remained low. The temporal evolution of test errors for the first 100 steps in a
sample test episode for a new target function C is shown in Fig. 3D. It can be seen that the network
had a high error in the first few steps of the episode, but it quickly learned to output the correct
targets, even though its weights W remained fixed. Thus the LSNN had learnt to learn in other
ways than through changes of its synaptic weights.
9

5

D ISCUSSION

We have shown that the inclusion of adapting neurons significantly enhances computing and learning
capabilities of recurrent networks of spiking neurons. In particular we have shown that the resulting
LSNN is able to integrate and maintain information over several seconds in a working memory.
In fact, it is reasonable to assume that any duration for working memory can be achieved with
corresponding large time constants of adapting neurons. Furthermore, our tests so far suggest that
LSNNs are able to approach the computational power and learning capability of the generally best
performing ANN model for temporal processing tasks: LSTM networks. BPTT tends to work well
for LSNNs because the error gradient can be propagated through the slowly changing firing threshold
of an adapting neuron without causing vanishing or exploding gradients. Hence adapting neurons
appear to provide “highways” for propagating errors backwards even over 100s and 1000s of layers
of an unrolled LSNN.
A remarkable feature of a realization of working memory in SNNs through increased firing thresholds is that this method is very energy-efficient: When a memory is stored through an increased firing
threshold, this tends to reduce firing activity in the SNN. In contrast, working memory was previously modelled in SNNs via persistent firing of a subset of neurons, i.e., through increased firing
activity. In neuroscience there is an ongoing debate to what extent working memory is implemented
in the brain through persistent firing or through activity-silent mechanisms (Stokes, 2015). Interestingly, experimental data show that also one form of long-term memory, familiarity, is expressed in
some brain regions (perirhinal and entorhinal cortex) through reduced activity during the presentation of a familiar stimulus (Eichenbaum et al., 2007). These results are consistent with a memory
model based on longer lasting increases of the firing threshold of specific subsets of neurons, see
also (Titley et al., 2017).
The implementation of short-term memory through neural adaptation is somewhat analogous to an
interesting variant of long short-term memory in artificial neural networks (Mikolov et al., 2014).
They introduced a particular population of neurons, termed context neurons, whose state changes on
a slower time scale than that of the other neurons. If one interprets the firing threshold of adapting
neurons as their “state”, then the dynamics of this state is analogous to the dynamics of context
neurons in their approach. However, in the non-spiking setting of (Mikolov et al., 2014) one cannot
see an advantage of this model with regard to energy-efficiency of a physical realization: the state
of context neurons has to be read out at every discrete time step, like for all other neurons. In other
words, a clocked network cannot make use of the possibility to convey information by not sending
a message, as one can do in SNNs. This option, to maintain and use information by not firing, is
extensively used by LSNNs.
In (Costa et al., 2017) an LSTM module was modelled through a recurrent network of non-spiking
neurons, i.e., of sigmoidal neurons whose analog outputs are interpreted as firing rates. Their model
is consistent with models for working memory based on persistent firing. It achieved remarkable
results for several datasets, including the same sequential MNIST data that we have considered for
the case of spike-based network inputs and spike-based network computations.
We have shown in Fig. 3 that the inclusion of adapting neurons in SNNs opens the door to porting
methods and results on L2L from LSTM networks to SNNs. Our result suggests new ways of storing
the “program” of a computation in SNNs: In the states of adapting neurons, i.e., in the working
memory of the LSNN. Similarly as previously shown for artificial neural networks in (Hochreiter
et al., 2001), this approach simultaneously opens the door to the exploration of new classes of
learning algorithms for SNNs that could previously not be considered. In particular, our result
shows that LSNNs can learn the input/output behaviour of feedforward artificial neural network
with a hidden layer by just using readily available local mechanisms. This shows that LSNNs are
able to emulate the function of backprop learning in quite different ways than considered so far.
We will show in (Bellec et al., 2018) that LSNNs also pave the way for porting new ideas for learning
to learn from rewards (Wang et al., 2016) to networks of spiking neurons. In particular, one can see
that LSNNs are able to learn nontrivial sequential behaviors, and that they can make use of control
transfer learning, as shown there for LSTM networks.
10

ACKNOWLEDGMENTS
Research leading to these results has in parts been carried out on the Human Brain Project PCP
Pilot Systems at the Juelich Supercomputing Centre, which received co-funding from the European
Union (Grant Agreement #604102). We gratefully acknowledge Sandra Diaz, Alexander Payser
and Wouter Klijn from the Simulation Laboratory Neuroscience of the Jülich Supercomputing Centre for their support. This research was carried out under partial support by the Human Brain Project
of the European Union #720270 and #785907.

R EFERENCES
Allen Institute. c 2018 Allen Institute for Brain Science. Allen Cell Types Database, cell feature search.
Available from: celltypes.brain-map.org/data. 2018.
Guillaume Bellec, Anand Subramoney, Robert Legenstein, and Wolfgang Maass. Networks of spiking neurons
learn to learn. In Preparation, 2018.
Rui Costa, Ioannis Alexandros Assael, Brendan Shillingford, Nando de Freitas, and Tim Vogels. Cortical
microcircuits as gated-recurrent neural networks. In Advances in Neural Information Processing Systems,
pp. 272–283, 2017.
Matthieu Courbariaux, Itay Hubara, Daniel Soudry, Ran El-Yaniv, and Yoshua Bengio. Binarized neural networks: Training deep neural networks with weights and activations constrained to+ 1 or-1. arXiv preprint
arXiv:1602.02830, 2016.
Mike Davies, Narayan Srinivasa, Tsung-Han Lin, Gautham Chinya, Yongqiang Cao, Sri Harsha Choday, Georgios Dimou, Prasad Joshi, Nabil Imam, Shweta Jain, et al. Loihi: A neuromorphic manycore processor with
on-chip learning. IEEE Micro, 38(1):82–99, 2018.
Howard Eichenbaum, Andrew P Yonelinas, and Charan Ranganath. The medial temporal lobe and recognition
memory. Annu. Rev. Neurosci., 30:123–152, 2007.
Steven K. Esser, Paul A. Merolla, John V. Arthur, Andrew S. Cassidy, Rathinakumar Appuswamy, Alexander
Andreopoulos, David J. Berg, Jeffrey L. McKinstry, Timothy Melano, Davis R. Barch, Carmelo di Nolfo,
Pallab Datta, Arnon Amir, Brian Taba, Myron D. Flickner, and Dharmendra S. Modha. Convolutional networks for fast, energy-efficient neuromorphic computing. Proceedings of the National Academy of Sciences,
113(41):11441–11446, November 2016. ISSN 0027-8424, 1091-6490. doi: 10.1073/pnas.1604850113.
Robert C Froemke, Dominique Debanne, and Guo-Qiang Bi. Temporal modulation of spike-timing-dependent
plasticity. Frontiers in synaptic neuroscience, 2:19, 2010.
Steve B Furber, David R Lester, Luis A Plana, Jim D Garside, Eustace Painkras, Steve Temple, and Andrew D
Brown. Overview of the spinnaker system architecture. IEEE Transactions on Computers, 62(12):2454–
2467, 2013.
Wulfram Gerstner, Werner M. Kistler, Richard Naud, and Liam Paninski. Neuronal dynamics: From single
neurons to networks and models of cognition. Cambridge University Press, 2014.
Nathan W Gouwens, Jim Berg, David Feng, Staci A Sorensen, Hongkui Zeng, Michael J Hawrylycz, Christof
Koch, and Anton Arkhipov. Systematic generation of biophysically detailed models for diverse cortical
neuron types. Nature communications, 9(1), 2018.
Abhishek Gupta, Russell Mendonca, YuXuan Liu, Pieter Abbeel, and Sergey Levine. Meta-reinforcement
learning of structured exploration strategies. arXiv preprint arXiv:1802.07245, 2018.
Geoffrey Hinton, Nitish Srivastava, and Swersky Kevin. Lecture 6d - a separate, adaptive learning rate for each
connection. Slides of Lecture Neural Networks for Machine Learning, 2012.
Sepp Hochreiter and Jürgen Schmidhuber. Long short-term memory. Neural computation, 9(8):1735–1780,
1997.
Sepp Hochreiter, A Steven Younger, and Peter R Conwell. Learning to learn using gradient descent. In
International Conference on Artificial Neural Networks, pp. 87–94. Springer, 2001.
Dongsung Huh and Terrence J Sejnowski. Gradient descent for spiking neural networks. arXiv preprint
arXiv:1706.04698, 2017.

11

David Kappel, Stefan Habenschuss, Robert Legenstein, and Wolfgang Maass. Network Plasticity as Bayesian
Inference. PLOS Computational Biology, 11(11):e1004485, 2015.
David Kappel, Robert Legenstein, Stefan Habenschuss, Michael Hsieh, and Wolfgang Maass. Reward-based
stochastic self-configuration of neural circuits. arXiv preprint arXiv:1704.04238, 2017.
Diederik P Kingma and Jimmy Ba.
arXiv:1412.6980, 2014.

Adam: A method for stochastic optimization.

arXiv preprint

John Lisman and Nelson Spruston. Questions about STDP as a general model of synaptic plasticity. Frontiers
in synaptic neuroscience, 2:140, 2010.
Wolfgang Maass, Thomas Natschläger, and Henry Markram. Real-time computing without stable states: A new
framework for neural computation based on perturbations. Neural computation, 14(11):2531–2560, 2002.
Tomas Mikolov, Armand Joulin, Sumit Chopra, Michael Mathieu, and Marc’Aurelio Ranzato. Learning longer
memory in recurrent neural networks. arXiv preprint arXiv:1412.7753, 2014.
Christian Pozzorini, Skander Mensi, Olivier Hagens, Richard Naud, Christof Koch, and Wulfram Gerstner.
Automated high-throughput characterization of single neurons by means of simplified spiking models. PLoS
computational biology, 11(6):e1004275, 2015.
Ning Qiao, Hesham Mostafa, Federico Corradi, Marc Osswald, Fabio Stefanini, Dora Sumislawska, and Giacomo Indiveri. A reconfigurable on-line learning spiking neuromorphic processor comprising 256 neurons
and 128k synapses. Frontiers in neuroscience, 9:141, 2015.
Johannes Schemmel, Daniel Brüderle, Andreas Grübl, Matthias Hock, Karlheinz Meier, and Sebastian Millner.
A wafer-scale neuromorphic hardware system for large-scale neural modeling. In Circuits and systems
(ISCAS), proceedings of 2010 IEEE international symposium on, pp. 1947–1950. IEEE, 2010.
Mark G. Stokes. ‘Activity-silent’ working memory in prefrontal cortex: a dynamic coding framework. Trends
in Cognitive Sciences, 19(7):394–405, 2015.
Corinne Teeter, Ramakrishnan Iyer, Vilas Menon, Nathan Gouwens, David Feng, Jim Berg, Aaron Szafer,
Nicholas Cain, Hongkui Zeng, Michael Hawrylycz, et al. Generalized leaky integrate-and-fire models classify multiple neuron types. Nature communications, 1(1):1–15, 2018.
Heather K Titley, Nicolas Brunel, and Christian Hansel. Toward a neurocentric view of learning. Neuron, 95
(1):19–32, 2017.
Jane X Wang, Zeb Kurth-Nelson, Dhruva Tirumala, Hubert Soyer, Joel Z Leibo, Remi Munos, Charles
Blundell, Dharshan Kumaran, and Matt Botvinick. Learning to reinforcement learn. arXiv preprint
arXiv:1611.05763, 2016.

12

S UPPLEMENTARY I NFORMATION
Neuron model: In continuous time the spike trains xi (t) and zj (t) are formalized as sums of
Dirac pulses. Neurons are modeled according to a standard adaptive leaky integrate-and-fire model.
A neuron j spikes as soon at its membrane potential Vj (t) is above its threshold Bj (t). At each
spike time t, the membrane potential Vj (t) is reset by subtracting the current threshold value Bj (t).
Importantly at each spike the threshold Bj (t) of an adaptive neuron is increased by a constant β.
Then the threshold decays back to a baseline value b0j . Between spikes the membrane voltage Vj (t)
and the threshold Bj (t) are following the dynamics
τm V̇j (t)
τa,j Ḃj (t)

= −Vj (t) + Rm Ij (t)
=

b0j

− Bj (t),

(3)
(4)

where τm is the membrane time constant, τa,j is the adaptation time constant and Rm is the membrane resistance. The input current Ij (t) is defined as the weighted sum of spikes from external
inputs and other neurons in the network:
X
X
in
rec
Ij (t) =
Wji
xi (t − din
Wji
zi (t − drec
(5)
ji ) +
ji ),
i

i

in
rec
where Wji
and Wji
denote respectively the input and the recurrent synaptic weights and din
ji
rec
and dji the corresponding synaptic delays. All network neurons are connected to a population of
out
. When network neuron j spikes, the output synaptic strength
readout neurons with weights Wkj
out
Wkj is added to the membrane voltage yk (t) of all readout neurons k. yk (t) also follows the

dynamics of a leaky integrator τm ẏk (t) = −yk (t). The physical units of voltages, currents, delays,
Dirac pulses and synaptic weights are respectively Volt, Ampere, Second, Hertz and Coulomb.
Implementation in discrete time: Our simulations were performed in discrete time with a time
step δt = 1 ms. In discrete time, the spike trains are modeled as binary sequences xi (t), zj (t) ∈
{0, 1}. Neuron j emits a spike at time t if it is currently not in a refractory period, and its membrane
potential Vj (t) is above its threshold Bj (t). During the refractory period following a spike, zj (t) is
fixed to 0. The dynamics of the threshold is defined by Bj (t) = b0j + βbj (t) where β is a constant
which scales the deviation bj (t) from the baseline b0j . The neural dynamics in discrete time reads as
follows
Vj (t + 1) =
bj (t + 1) =

αVj (t) + (1 − α)Rm Ij (t) − Bj (t)zj (t)
ρj bj (t) + (1 − ρj )zj (t),

(6)
(7)

where α = exp(− τδtm ) and ρj = exp(− τδt
). The term Bj (t)zj (t) implements the reset of the
a,j
membrane voltage after each spike. The current Ij (t) is the weighted sum of the incoming spikes.
The definition of the input current in equation (5) holds also for discrete time, with the difference
that spike trains now assume values in {0, 1}.
Parameter values: For adaptive neurons, we used βj = 1.7, and for regular spiking neurons we
used βj = 0 (i.e. Bj is constant). The baseline threshold voltage was b0j = 0.01 and the membrane
time constant τm = 20 ms. Initial network weights were drawn from a Gaussian distribution Wji ∼
√w0 N (0, 1), where nin is the number of afferent neurons in the considered weight matrix (i.e.,
nin
the number of columns of the matrix), N (0, 1) is the zero-mean unit-variance Gaussian distribution,
and w0 is a weight scaling factor chosen to be w0 = R1Vm δt . For this choice, a spike transmitted
with synaptic weight w0 resulted in an increase of post-synaptic potential of 1 − α. In this way,
the standard initialization of synaptic weights initialized the networks in a stable regime producing
realistic firing rates.
In store-recall the test performance is computed over sequences of 32 consecutive store-recall pairs
and averaged over 32 random input sequences. During training with BPTT the depth of the unrolled
network is adapted to the distribution of store-recall delays. For expected delays of 2, 4, and 6 s
it resulted in 3200, 5200 and 8000 unrolled time steps respectively. All networks were trained for
500 iterations. We used the Adam optimizer (Kingma & Ba, 2014) with default parameters and a
13

learning rate of 0.01, with sequences presented in batches of sizes indicated in Table 1 and 2 (as
“batch size”). To avoid unrealistically high firing rates, the loss function contains a regularization
term that minimizes the distance of average firing rate between individual neurons and a target firing
rate of 10 Hz.
In sequential MNIST all networks were trained for 36, 000 iterations with approximately 66, 000
parameters, with a batch size of 256, using the default Adam optimizer, and a learning rate initialized
at 0.01 and decayed by a factor 0.8 every 2500 iterations. The baseline artificial RNN contains 128
hidden units and uses hyperbolic tangent as the activation function. The LIF network is formed by a
fully connected population of 220 regular spiking neurons.
In the learning to learn task, the network was trained for 5, 000 iterations with a batch size of 10
episodes, using the RMSProp (Hinton et al., 2012) optimizer, and a learning rate of 0.001.
Propagation of gradients in recurrent networks of LIF neurons: In artificial recurrent neural
networks such as LSTMs, gradients can be computed with back propagation through time (BPTT).
For BPTT in spiking neural networks, complications arise from the non-differentiability of the output of spiking neurons, and from the fact that gradients need to be propagated either through continuous time or through many time steps if time is discretized. Therefore, in (Courbariaux et al., 2016;
Esser et al., 2016) it was proposed to use a pseudo-derivative.
dzj (t)
:= max{0, 1 − |vj (t)|},
dvj (t)
V (t)−B (t)

(8)

where vj (t) denotes the normalized membrane potential vj (t) = j Bj (t)j . This made it possible
to train deep feed-forward networks of deterministic binary neurons (Courbariaux et al., 2016; Esser
et al., 2016). We observed that this convention tends to be unstable for very deep (unrolled) recurrent
networks of spiking neurons (see Fig. 1E). To achieve stable performance we dampened the increase
of back propagated errors through spikes by using a pseudo-derivative of amplitude γ < 1 (typically
γ = 0.3):
dzj (t)
:= γ max{0, 1 − |vj (t)|}.
(9)
dvj (t)
Note that in adaptive neurons, gradients can propagate through many time steps in the dynamic
threshold. This propagation is not affected by the dampening.

14

expected delay between “store”
and “recall” instructions
D (ms)
batch size
LSNN τa = 200ms
LSNN τa = 2s
LSNN τa = 4s
LSNN τa = 6s
LSNN τa = 8s
LSNN τa spread over [0, 8]s
LIF

200ms

2s

4s

6s

uniform in [0-8]s

50
24
0.11
0.24
0.30
0.58
0.40
0.41
0.08

200
24
51
1.59
7.62
7.24
16.3
3.08
49

200
24
48
1.97
3.77
9.17
45.1
3.25
49

200
24
48
15.7
8.27
5.06
34.0
8.12
51

200
256
49
16.9
6.96
3.44
8.33
4.36
49

Table 1: Test error (%) on store-recall tasks with different delays, and different time constants
of adapting neurons in LSNNs with the standard architecture. Mean test errors (over 3 random
network initializations on test batch of size 512) reached by different networks after 500 training
iterations. Best and second best results per task (column) are marked with boldface (D: pattern
presentation time). For comparison, the performance of a fully connected recurrent SNN without
adapting neurons is shown in the last row (LIF).

expected delay
D (ms)
batch size
Standard
Fully connected
Deep R

4s
200
24
3.25
6.39
0.41

uniform in [0-8]s
200
256
4.36
31.8
5.17

Table 2: Test error (%) of different LSNN architectures on two store-recall tasks. Mean test
errors (over 3 random network initializations on test batch of size 512) reached by different networks
after 500 training iterations. Rows show results for SSNs with different architectures. The adaptation
time constants were spread for all these networks uniformly over [0, 8] s, like in the second to last
row of Table 1.

15

