LETTER

Communicated by Walter Senn

Optimality Model of Unsupervised Spike-Timing-Dependent
Plasticity: Synaptic Memory and Weight Distribution
Taro Toyoizumi
taro.toyoizumi@brain.riken.jp
School of Computer and Communication Sciences and Brain-Mind Institute,
Ecole Polytechnique Fédérale de Lausanne, CH-1015 Lausanne EPFL, Switzerland,and Department of Complexity Science and Engineering, Graduate School of Frontier
Sciences, University of Tokyo, Tokyo 153-8505, Japan

Jean-Pascal Pfisterjean-pascal.pfister@epfl.ch
School of Computer and Communication Sciences and Brain-Mind Institute,
Ecole Polytechnique Fédérale de Lausanne, CH–1015 Lausanne EPFL, Switzerland

Kazuyuki Aihara
aihara@sat.t.u.-tokyo.ac.jp
Department of Information and Systems, Institute of Industrial Science, University ofTokyo, 153–8505, Tokyo, Japan, and Aihara Complexity Modelling Project, ERATO,
JST, Shibuya-ku, Tokyo 151–0064, Japan

Wulfram Gerstner
wulfram.gerstner@epfl.ch
School of Computer and Communication Sciences and Brain-Mind Institute,
Ecole Polytechnique Fédérale de Lausanne, CH-1015 Lausanne EPFL, Switzerland

We studied the hypothesis that synaptic dynamics is controlled by threebasic principles: (1) synapses adapt their weights so that neurons caneffectively transmit information, (2) homeostatic processes stabilize themean firing rate of the postsynaptic neuron, and (3) weak synapses adapt
more slowly than strong ones, while maintenance of strong synapses iscostly. Our results show that a synaptic update rule derived from theseprinciples shares features, with spike-timing-dependent plasticity, is sensitive to correlations in the input and is useful for synaptic memory.
Moreover, input selectivity (sharply tuned receptive fields) of postsynaptic neurons develops only if stimuli with strong features are presented.
Sharply tuned neurons can coexist with unselective ones, and the distribution of synaptic weights can be unimodal or bimodal. The formulationof synaptic dynamics through an optimality criterion provides a simple
graphical argument for the stability of synapses, necessary for synaptic
memory.
Neural Computation 19, 639–671 (2007)


C 2007 Massachusetts Institute of Technology

640

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

1 Introduction
Synaptic changes are thought to be involved in learning, memory, and cortical plasticity, but the exact relation between microscopic synaptic properties
and macroscopic functional consequences remains highly controversial. In
experimental preparations, synaptic changes can be induced by specificstimulation conditions defined through pre- and postsynaptic firing rates(Bliss & Lomo, 1973; Dudek & Bear, 1992), postsynaptic membrane potential(Kelso, Ganong, and Brown, 1986), calcium entry (Malenka, Kauer, Zucker,& Nicoll, 1988; Lisman, 1989), or spike timing (Markram, Lübke, Frotscher,& Sakmann, 1997; Bi & Poo, 2001. In the theoretical community, conditions
for synaptic changes are formulated as synaptic update rules or learning
rules (von der Malsburg, 1973; Bienenstock, Cooper, & Munroe, 1982; Miller,
Keller, & Stryker, 1989) (for reviews, see Gerstner & Kistler, 2002; Dayan &Abbott, 2001; Cooper, Intrator, Blais, & Shouval, 2004), but the exact featuresthat make a synaptic update rule a suitable candidate for cortical plasticity
and memory are unclear.
From a theoretical point of view, a synaptic learning rule should be(1) sensitive to correlations between pre- and postsynaptic neurons (Hebb,
1949) in order to respond to correlations in the input (Oja, 1982); theyshould (2) allow neurons to develop input selectivity (e.g., receptive fields)(Bienenstock et al., 1982; Miller et al. 1989), in the presence of strong input
features, but (3) distribution of synaptic strength should remain unimodal
otherwise (Gütig, Aharonov, Rotter, & Sompolinsky, 2003). Furthermore (4)synaptic memories should show a high degree of stability (Fusi, Drew, &Abbott, 2005) and nevertheless remain plastic (Grossberg, 1987). Moreover,experiments suggest that plasticity rules are (5) sensitive to the presynaptic
firing rate (Dudek & Bear, 1992), but (6) depend also on the exact timing of
the pre- and postsynaptic spikes (Markram et al., 1997, Bi & Poo, 2001).
Many other experimental features could be added to this list, for example, the role of intracellular calcium and of NMDA receptors, but we will
not do so (see Bliss & Collingridge, 1993, and Malenka & Nicoll, 1993, forreviews).
The items in the above list are not necessarily exclusive, and the relative
importance of a given aspect may vary from one subsystem to the next;
for example, synaptic memory maintenance might be more important for along-term memory system than for primary sensory cortices. Nevertheless,all of the above aspects seem to be important features of synaptic plasticity. However, the development of theoretical learning rules that exhibit
all of the above properties has posed problems in the past. For example,
traditional learning rules that have been proposed as an explanation ofreceptive field development (Bienenstock et al., 1982; Miller et al., 1989)
exhibit a spontaneous separation of synaptic weights into two groups, even
if the input shows no or only weak correlations. This is difficult to reconcile with experimental results in visual cortex of young rats, where aOptimality Model of Unsupervised Spike-Timing-Dependent Plasticity

641

unimodal distribution was found (Sjöström, Turrigiano, & Nelson, 2001).
Moreover model neurons that specialize early in development on one subset of features cannot readily readapt later. Other learning rules, however,that exhibit a unimodal distribution of synaptic weights (Gütig et al., 2003)
do not lead to a long-term stability of synaptic changes.
In this article, we show that all of features 1 to 6 emerge naturally in atheoretical model where we require only a limited number of objectives thatwill be formulated as postulates. In particular, we study how the conflicting
demands on synaptic memory maintenance, plasticity, and distribution ofsynaptic synapses could be satisfied by our model. Although the postulates
are rather general and could be adapted to arbitrary neural systems, we
had in mind excitatory synapses in neocortex or hippocampus and exclude
inhibitory synapses and synapses in specialized systems such as the calix ofHeld in the auditory pathway. Our arguments are based on three postulates:
A. Synapses adapt their weights so as to allow neurons to efficiently
transmit information. More precisely, we impose a theoretical postulate that the mutual information I between presynaptic spike trains
and postsynaptic firing be optimized. Such a postulate stands in the
tradition of earlier theoretical work (Linsker, 1989; Bell & Sejnowski,1995), but is formulated here on the level of spike trains rather thanrates.
B. Homeostatic processes act on synapses to ensure that the long-term
average of the neuronal firing rate becomes close to a target ratethat is characteristic for each neuron. Synaptic rescaling and related
mechanism could be a biophysical implementation of homeostatis(Turrigiano & Nelson, 2004). The theoretical reason for such a postulate is that sustained high firing rates are costly from an energetic
point of view (Laughlin, de Ruyter van Steveninck, and Anderson,
1998; Levy & Baxter, 2002).
C. C1: Maintenance of strong synapses is costly in terms of biophysical machinery, in particular, in view of continued protein synthesis(Fonseca, Nägerl, Morris, & Bonhoeffer, 2004). C2: Synaptic plasticity
is slowed for very weak synapses in order to avoid a (unplausible)
transition from excitatory to inhibitory synapses.
Optimality approaches have a long tradition in the theoretical neurosciences and have been utilized in two different ways. Firstly, optimality
approaches allow deriving strict theoretical bounds against which performance of real neural systems can be compared (Barlow, 1956; Laughlin, 1981; Britten, Shadlen, Newsome, & Movshon, 1992; de Ruyter van
Steveninck & Biale, 1995). Second, they have been used as a conceptual
framework since they allow connecting functional objectives (e.g., “be reliable!”) and constraints (e.g., “don’t use too much energy!”) with electrophysiological properties of single neurons and synapses or neuronal populations (Barlow, 1961; Linsker, 1989; Atick & Redlich, 1990; Levy & Baxter,642

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

2002; Seung, 2003). Our study, a derivation of synaptic update rules froman optimality viewpoint, follows this second, conceptual approach.
2 The Model
2.1 Neuron Model. We simulated a single stochastic point neuronmodel with N = 100 input synapses. Presynaptic spikes at synapse j are def
noted by their arrival time t j and evoke excitatory postsynaptic potentials
ff(EPSPs) with time course exp[−(t − t j )/τm ] for t ≥ t j , where τm = 20 ms isthe membrane time constant. Recent experiments have shown that action
potentials propagating back into the dendrite can partially suppress EPSPs
measured at the soma (Froemke, Poo, & Dan, 2005). Since our model neuron
has no spatial structure, we included EPSP suppression by a phenomenologf
ical amplitude factor a (t j − t̂) that depends on the time difference between
presynaptic spike arrival and the spike trigger time t̂ of the last (somatic)
action potential of the postsynaptic neuron.
In the absence of EPSP suppression, the amplitude of a single EPSP atsynapse j is characterized by its weight w j and its duration by the membrane time constant τm . Summation of the EPSPs caused by presynaptic
spike arrival at all 100 explicitly modeled synapses gives the total postsynaptic potential
u(t) = ur +

N

j=1

wj





f

exp −

f

t j <tt − tjτmfa (t j − t̂),(2.1)

where ur = −70 mV is the resting potential and the sum runs over all spikeffarrival times t j in the recent past, t̂ < t j ≤ t. The EPSP suppression factorfftakes a value of zero if t j < t̂ and is modeled for t j ≥ t̂ as exponentialff

recovery a (t j − t̂) = 1 − exp[−(t j − t̂)/τa ] with time constant τa = 50 ms(see Figure 1a) unless stated otherwise. The parameters w j for 1 ≤ j ≤ Ndenote the synaptic weight of the N = 100 synapses and are updated using
a learning rule discussed below.
In order to account for unspecific background input that was not modeled explicitly, spikes were generated probabilistically with density
ρ(t) = ρr + [u(t) − ur ] · g,

(2.2)

where ρr = 1 Hz is the spontaneous firing rate (in the absence of spike
input at the 100 explicitly modeled synapses) and g = 12.5 Hz/mV is a gain
factor. Thus, the instantaneous spike density increases linearly with the total
postsynaptic potential u(t). Note, however, that due to EPSP suppression,
the total postsynaptic potential increases sublinearly with the number ofOptimality Model of Unsupervised Spike-Timing-Dependent Plasticity

a

b
50

−40
−45

40

−50

ρ [Hz]−55−60−

u [mV]

643

30

20

−65

10
−70
−75
−50

0

50

100

0
0

150

10

20

40

50

ν [Hz]

Time [ms]
d
10Hz
20Hz
30Hz

0.03

0.02

0.01

0

0

20

40

60

ISI [ms]

80

100

Auto correlation [Hz]c
ISI distribution

30

100

10Hz
20Hz
30Hz

80

60

40

20

0
−100

−50

0

50

100

Time [ms]Figure 1: Properties of the stochastically spiking neuron model. (a) Neuronmodel and EPSP suppression. EPSPs arriving just after a postsynaptic spike at
t = 0 are attenuated by a factor (1 − e −t/τa ) and recover exponentially to theirmaximal amplitude w j . (b) The output firing rate ρ̄ of a neuron that receives
stochastic spike input at rate ν at all 100 synapses. Each presynaptic spike
evokes an EPSP with maximal amplitude w j = 0.4 mV. (c) Interspike interval(ISI) distribution with input frequency ν = 10 Hz (solid line), 20 Hz (dashedline), and 30 Hz (dotted line) at all 100 synapses. (d) Autocorrelation functionof postsynaptic action potentials at an input frequency of 10 Hz (solid line), 20
Hz (dashed line), and 30 Hz (dotted line).

input spikes, and so does the mean firing rate of the postsynaptic neuron(see Figure 1b).
The neuron model is simulated in discrete time with time steps of t =1 ms on a standard personal computer using custom-made software written
in Matlab.
2.2 Objective Function. Postulates A and B have been used previously(Toyoizumi, Pfister, Aihara, & Gerstner, 2005a) and lead to an optimality
criterion L = I − γ D, where I is the mutual information between presynaptic input and postsynaptic output and D a measure of the distance of644

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

the mean firing rate of the neuron from its target rate. The parameter γ
scales the importance of the information term I (postulate A) compared tothe homeostatic term D (postulate B). It was shown that optimization ofL by gradient ascent yields a synaptic update rule that shows sensitivity
to correlations (see point 1 above), input selectivity (see point 2 above),and depends on presynaptic firing rates (see point 5 above) (Toyoizumi
et al., 2005a). However, while the learning rule in Toyoizumi et al. (2005a)
showed some dependence on spike timing, it did not (without additional
assumptions) have the typical features of Spike-timing-dependent plasticity (STDP) as measured in vitro (point 6 above), and it exhibited, like earlier
models (Bienenstock et al., 1982), spontaneous synaptic specialization, evenfor very weak input features, which is in contrast to point 3 above.
In this earlier theoretical study, synaptic potentiation was artificiallystopped at some upper bound w max (and synaptic depression was stopped
at weight w = 0), so as to ensure that weights w stayed in a regime 0 ≤w ≤ w max . In this letter, we take the more realistic assumption that strong
weights are more likely to show depression than weaker ones but do not
impose a hard upper bound. Similarly, we require that adaptation speed isslowed for very weak synapses but do not impose a hard bound at zeroweight. We will show that with these assumptions, the resulting synaptic
update rule shows properties of STDP (see point 6 above), is suitable formemory retention (see point 4 above), and leads to synaptic specialization
when driven by strong input (see point 3 above), while keeping properties
1, 2, and 5, which were found in Toyoizumi et al. (2005a).
To avoid hard upper bounds for the synapses, we use postulate C1 andadd a term  to the optimality criterion L that is proportional to w 2 (i.e., thesquare of the synaptic weight) and proportional to the presynaptic firingrate. This term comes with a negative sign, since a cost is associated with
big weights. Hence, from our optimality viewpoint, synapses change so as
to maximize a quantity
L = I − γ D − λ,(2.3)

where I is the information to be maximized, D is a measure of the firing rate mismatch to be minimized, and  is the cost induced by strong
synapses to be minimized. Factors γ and λ control the relative importanceof the three terms. In other words, synapses adjust their weights so as tobe able to transmit information while keeping the mean firing rate andsynaptic weights at low values. Thus, our three postulates A, B, and C
give rise to one unified optimality criterion L. We hypothesize that a significant part of findings regarding synaptic potentiation and depression
can be conceptually understood as the synapse’s attempt to optimize thecriterion L.

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

645

The learning rule used for the update of the synaptic weights w j is
derived from the objective function 2.3, L = I − γ D − λ, which contains
three terms.
The first term is the mutual information between the ensemble of 100
input spike trains (spanning the interval of a single trial from 0 to T; theensemble of all presynaptic trains is formally denoted by X(T) = {x j (t) =fft j δ(t − t j )| j = 1, . . . , 100, 0 ≤ t < T}) and the output spike train of thepostsynaptic neuron over the same interval (denoted by Y(T) = {y(t) =fff
tpost δ(t − t post )|0 ≤ t < T}, where t post represent output spike timing), thatis,P(Y|X)
I = log,
P(Y) Y,X(2.4)

where angular brackets  · Y,X denote averaging over all combinations of
input and output spike trains.1 Here, P(Y|X) is the conditional probability
density of our stochastic neuron model to generate a specific spike train Yfwith (one or several) spike times {tpost } during a trial of duration T given
100 known input spike trains X. This conditional probability density isfgiven as a product of the instantaneous probabilities ρ(tpost ) of firing at thefpostsynaptic spike times {tpost } and the probability of not firing elsewhere,that is,





  T

f
ρ(tpost ) exp −ρ(t)dt .

P(Y|X) = 
f
tpost

(2.5)

0

Similarly, P(Y) is the probability of generating the same output spiketrain Y not knowing the input. Here “not knowing the input” implies that
we have to average over all possible inputs so as to get the expected instantaneous firing density ρ̄(t) at time t. However, because of the EPSP
suppression factor, the expected firing density will also depend on the last
output spike before t. We therefore defineρ̄(t) = ρ(t) X(t)|Y(t) ,(2.6)

that is, we average over the inputs but keep the knowledge of the previousfoutput spikes tpost < t. P(Y) is then given by a formula analogous to equation 2.5, but with ρ replaced by ρ̄. Hence, given our neuron model, both

1 From now on, when X or Y are without argument, we take implicitly X ≡ X(T) andY ≡ Y(T), i.e., the spike trains over the full interval T.

646

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

P(Y|X) and P(Y) in equation 2.4 are well defined. The information term I
of equation 2.4 is the formal instantiation of postulate A.
The second term is the homeostatic termP(Y)D = log,P̃(Y) Y

(2.7)

which compares the actual distribution of output spike trains P(Y) withthat of an ideal distribution P̃(Y) generated by the same neuron firing at atarget rate of ρ̃ = 5 Hz—formula (2.5) with ρ replaced by ρ̃. Mathematicallyspeaking, D is the Kullback-Leibler distance between two distributions(Cover & Thomas, 1991), but in practice, we may think of D simply as ameasure of the difference between actual and target firing rates (Toyoizumi,
et al., 2005a). The term D is our mathematical formulation of postulate B.
The third term is the cost associated with strong synapses. We assume
that the cost increases quadratically with the synaptic weights but that only
synapses that have been activated in the past contribute to the cost. Hence,
the mathematical formulation of postulate C1 yields a cost=1  2 
wj nj X ,
2 j

(2.8)

where n j is the number of presynaptic spikes that have arrived at synapse jduring the duration T of the interval under consideration. Cost terms thatare quadratic in the synaptic weights are common in the theoretical literature (Miller & MacKay, 1994), but the specific dependence on presynaptic
spiking induced by the factor n j in equation 2.8 is not. The dependenceof  on presynaptic spike arrival means that in our model, only activatedsynapses contribute to the cost. The specific formulation of  is mainly due
to theoretical reasons to be discussed below. The intuition is that activation of a synapse in the absence of any postsynaptic activity can weaken
the synapse if the factor λ is sufficiently positive (see also Figure 3d). The
restriction of the cost to previously activated synapses is reminiscent ofsynaptic tagging (Frey & Morris, 1997; Fonseca et al., 2004), although any
relation must be seen as purely hypothetical.
The three terms are given a relative importance by choosing γ = 0.1 (fora discussion of this parameter, see Toyoizumi et al., 2005a) and λ = 0.026 so
as to achieve a baseline of zero in the STDP function (see appendix A).
2.3 Synaptic Update Rule. We optimize the synaptic weights by gradient ascentw j = α(w j )

∂L
,
∂w j

(2.9)

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

647

with a weight-dependent update rate α(w j ). According to postulate C2,plasticity is reduced for very small weights. For simplicity, we chosew4j
α(w j ) = 4 · 10−2 w4 +w
4 , where ws = 0.2 mV, that is, learning slows for weakjssynapses with EPSP amplitudes around or less than 0.2 mV. Note that updates according to Equation 2.9 are always uphill; however, because of thew j dependence of α, the ascent is not necessarily along the steepest gradient.
Using the same mathematical arguments as in Toyoizumi et al. (2005a),
we can transform the optimization by gradient ascent into a synaptic update
rule. First, differentiating each term, we find∂I1
∂ P(Y|X)
P(Y|X)=log,∂w j
P(Y|X) ∂w j
P(Y) Y,X∂ P(Y|X)
1
P(Y)
∂D
=
log
,
∂w j
P(Y|X) ∂w j
P̃(Y) Y,X ∂= wj nj X .
∂w j(2.10)(2.11)(2.12)

We will rewrite the terms appearing in equations 2.10 and 2.11, by introducing the auxiliary variables
c j (t) = ∞dρ/du|u=u(t)ds  (s  )x j (t − s  )[y(t) − ρ(t)]ρ(t)
0(2.13)

and




ρ̄(t)
ρ(t)
− (ρ(t) − ρ̄(t)) − γ y(t) log− (ρ̄(t) − ρ̃) ,B post (t) = y(t) log
ρ̄(t)ρ̃(2.14)
Using the definitions in equations 2.13 and 2.14, we find the derivative of
the conditional probability density that appears in equations 2.10 and 2.11:
 T∂ P(Y|X)= P(Y|X)c j (t  )dt ∂w j0(2.15)andlog

P(Y|X)P(Y)=− γ log
P(Y)
P̃(Y)

 T

B post (t)dt.

(2.16)

0

As a first interpretation, we may say that c j represents the causal correlation
between input and output spikes (corrected for the expected correlation),

648

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstnerand B post is a function of postsynaptic quantities, namely, the output spikes
y, current firing rate ρ via the membrane potential u, average firing rate ρ̄,and the target firing rate ρ̃. More precisely, B post compares the actual output
with the expected output and, modulated by a factor γ , the expected output
with the target.
Hence, with the results from equations 2.10 to 2.16, the derivative of theobjective function is written in terms of averaged quantities ·Y,X as T  T∂L=dt
c j (t  )dt  B post (t) − λw j x j (t)
.
∂w j
0
0
Y,X

(2.17)

An important property of c j is that its average c j Y|X vanishes. On the other
hand, the correlations between c j (t  ) and B post (t) are limited by the timescaleτ AC of the autocorrelation function of the output spike train. Hence, we canlimit the integration to the relevant timescales without loss of generality
and introduce an exponential cut-off factor with time constant τC > τ AC :C j (t) = lim

 t+ε

ε→+0 0c j (t  )e −(t−t )/τC dt  .

(2.18)

With this factor C j , we find a batch learning rule (i.e., with expectations
over the input and output statistics on the right-hand side) of the form
 T


∂L
≈
dt C j (t)B post (t) − λw j x j (t) Y,X .
∂w j
0(2.19)

Finally, for slow learning rate α and stationary input statistics, the system
becomes self-averaging (i.e., expectations can be dropped due to automatic
temporal averaging; Gerstner & Kistler, 2002) so that we arrive at the online
gradient learning rule,dw j= α(w j ) C j (t)B post (t) − λw j x j (t) .
dt(2.20)

See Figure 2 for an illustration of the dynamics of C j and B post . The lastterm has the form of a “weight decay” term common in artificial neuralnetworks (Hertz, Krogh, & Palmer, 1991) and arises from the derivativeof the weight-dependent cost term . The parameter λ is set such thatdw j /dt = 0 for large enough |t pr e − t post | in the STDP in vitro paradigm. A
few steps of calculation (see appendix A) yield λ = 0.026. In our simulations,
we take τC = 100 ms for the cut-off factor in Equation 2.18.
For a better understanding of the learning dynamics defined in equation
2.20, let us look more closely at Figure 2a. The time course w/w of thepotentiation has three components: first, a negative jump at the moment

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

b

−2
−50

0

50

100

B post

−100
0.1

1
0
−1
−2
−50
0.05

Cj

2

0
−0.1
−100

−50

0

50

100

0.4
0.2
0
−0.2
−100

−50

0

50

100

Time [ms]

∆ w/w

∆ w/w

B post

Cj

a

649

0

50

100

150

0

50

100

150

0

50

100

150

0

−0.05
−50
0.05
0
− 0.05− 0.1−50Time [ms]Figure 2: Illustration of the dynamics of C j (top), B post (middle), and w/w =(w − winit )/winit (bottom) during a single pair of pre- and postsynaptic spikes(indicated schematically at the very top). We always start from an initial weight
winit = 4 mV and induce postsynaptic firing at time t post = 0. (a) Pre-before-post
timing with t pr e = −10 ms (solid line) induces a large potentiation, whereas
t pr e = −50 ms (dashed line) induces almost no potentiation. (b) Due to the EPSP
suppression factor, a post-before-pre timing with t pr e = 10 ms (solid line) induces a large depression, whereas t pr e = 50 ms (dashed line) induces a smaller
depression. The marks (circle, cross, square, and diamond) correspond to theweight change due to a single pair of pre- and postsynaptic spike after weights
have converged to their new values. Note that in Figure 3, the corresponding
marks indicate the weight change after 60 pairs of pre- and postsynaptic spikes.

of the presynaptic spike induced by the weight decay term; second, a slow
increase in the interval between pre- and postsynaptic spike times induced
by C j (t)B post (t) > 0; and third, a positive jump immediately after the postsynaptic spike induced by the singularity in B post combined with a positive
Cj.
As it is practically difficult to calculate ρ̄(t) = ρ(t) X(t)|Y(t) , we estimate ρ̄
by the running average of the output spikes, that is,τρ̄d ρ̄ est= −ρ̄ est (t) + y(t),dt(2.21)

with τρ̄ = 1 min. This approximation is valid if the characteristics of the stimulus and output spike trains are stationary and uncorrelated. To simulate invitro STDP experiments, the initial value of ρ̄ est is set equal to the injected
pulse frequency. Other, more accurate estimates of ρ̄ est are possible but lead650

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

to no qualitative change of results (data not shown). In the simulations of
Figures 4 to 6, the initial synaptic strength is set to be winit = 0.4 ± 0.04 mV.
2.4 Stimulation Paradigm
2.4.1 Simulated STDP In Vitro Paradigm. For the simulations of Figure 3,
f
spike timings t pr e = t j at synapse j and postsynaptic spike times t post are
imposed with a given relative timing t pr e − t post . For the calculation of thetotal STDP effect according to a typical in vitro stimulation paradigm, thepairing of pre- and postsynaptic spikes is repeated until 60 spike pairs havebeen accumulated. Spike pairs are triggered at a frequency of 1 Hz except
for Figure 3c, where the stimulation frequency was varied.
2.4.2 Simulated Stochastic Spike Arrival. In most simulations, presynaptic
spike arrival was modeled as Poisson spike input at either a fixed rate (homogeneous Poisson process) or a modulated rate (inhomogeneous Poissonprocess).
For example, for the simulations in Figure 6, with a gaussian profile, spikearrival at synapse j is generated by an inhomogeneous Poisson processwith the following characteristics. During a segment of 200 ms, the rate is
fixed at ν j = (ν max − ν0 ) exp[−0.01 ∗ d( j − k)2 ] + ν0 , where ν0 = 1 Hz is thebaseline firing rate and d( j − k) is the difference between index j and k.
The value of k denotes the location of the maximum. The value of k wasreset every 200 ms to a value chosen stochastically between 1 and 100. (Asindicated in the main text, presynaptic neurons in Figure 6 were considered
to have a ring topology, which has been implemented by evaluating thedifference d( j − k) as d( j − k) = min{| j − k|, 100 − | j − k|}.)However, in the simulations for Figure 5, input spike trains were notindependent Poisson, but we included spike-spike correlations. A correlation index of c = 0.2 implies that between a given pair of synapses, 20% of
spikes have identical timing. More generally, for a given value of c withina group of synaptic inputs, 100 c percent of spike arrival times are identicalat an arbitrary pair of synapse within the group.
3 Results
The mathematical formulation of postulates A, B, and C1 led to an optimality criterion L that was optimized by changing synaptic weights in an uphill
direction. In order to include postulate C2, the adaptation speed was madeto depend on the current value of the synaptic weight so that plasticity wassignificantly slowed for synapses with excitatory postsynaptic potentials(EPSPs) of amplitude less than 0.2 mV (see section 2 for details).
As in a previous study based on postulates A and B (Toyoizumi et al.,2005a), the optimization of synaptic weights can be understood as a synaptic

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

a

b
100

50

4mV
6mV

4000
80

∆w /w [%]

∆w /w [%]

40
30
20
10

3000
2000

60

1000
40

0
0.1

20

1

10

0

0
−10
−100

−50

0

50

−20
2

100

tpre − tpost [ms]

4

6

8

10

Initial EPSP amplitude [mV]

c

d
50

2

0.5Hz
1Hz
2Hz

30
20
10
0
−10
−100

−50

pre −

t

0

post

t

50

[ms]

100

25ms
50ms
No cost1.5∆ w [mV]40

w / w [%]

651

1

0.5

0

−0.5
− 100

− 50

0

50

100

tpre − tpost [ms]Figure 3: The synaptic update rule of the model shares features with STDP.
(a) STDP function (percentage change of EPSP amplitude as a function of t pr e −t post ) determined using 60 pre-and-post spike pairs injected at 1 Hz. The initial
EPSP amplitudes are 4 mV (dashed line) and 6 mV (dotted). Marks (circle, cross,square and diamond) correspond respectively to t pr e = −50 ms, t pr e = −10 ms,
t pr e = 10 ms, and t pr e = 50 ms, also depicted on Figure 2. (b) The percentagechange in EPSP amplitude after 60 pre-and-post spike pairs injected at 1 Hzfor pre-before-post timing ((t pr e − t post ) = −10 ms, solid line) and post-beforepre timing ((t pr e − t post ) = +10 ms, dashed line) as a function of initial EPSP
amplitude. Our model results qualitatively resemble experimental data (seeFigure 5 in Bi & Poo, 1998). (c) Frequency dependence of the STDP function:
spike pairs are presented at frequencies of 0.5 Hz (dashed line), 1 Hz (dottedline), and 2 Hz (dot-dashed line). The STDP function exhibits only a weaksensitivity to the change in stimulation frequency. (d) STDP function for different
choices of model parameters. The extension of the synaptic depression zone
for pre-after-post timing (t pr e − t post > 0) depends on the timescale τa of EPSP
suppression (dot-dashed line, τa = 50 ms; dashed line, τa = 25 ms). The dottedline shows the STDP function in the absence of a weight-dependent cost term
. The STDP function exhibits a positive offset indicating that without the costterm , unpaired presynaptic spikes would lead to potentiation, a non-Hebbianform of plasticity.

652

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

update rule that depends on presynaptic spike arrival, postsynaptic spike
firing, the postsynaptic membrane potential, and the mean firing rate of
the postsynaptic neuron. In addition, the synaptic update rule in this study
included a term that decreases the synaptic weight upon presynaptic spikearrival by a small amount proportional to the EPSP amplitude (see section 2for details). This term can be traced back to the additional weight-dependent
cost term  in equation 2.3 that accounts for postulate C1.
In order to study the consequences of the synaptic update rule derived
from postulates A, B, and C, we used computer simulations of a modelneuron that received presynaptic spike trains at 100 synapses. Each presynaptic spike evoked an EPSP with exponential time course (time constantτm = 20 ms). In order to account for dendritic interaction between somatic
action potentials and postsynaptic potentials, the amplitude of EPSPs was
suppressed immediately after postsynaptic spike firing (Froemke et al.,2005) and recovered with time constant τa = 50 ms (see Figure 1a). As ameasure of the weight of a synapse j, we used the EPSP amplitude w j at
this synapse in the absence of EPSP suppression. With all synaptic parameters w j set to a fixed value, the model neuron fired stochastically with amean firing rate ρ̄ that increases with the presynaptic spike arrival rate (seeFigure 1b), has a broad interspike interval distribution (see Figure 1c), andhas an autocorrelation function with a trough of 10 to 50 ms that is due toreduced excitability immediately after a spike because of EPSP suppression(see Figure 1d).
3.1 The Learning Rule Exhibits STDP. In a first set of plasticity experiments, we explored the behavior of the model system under a simulated
in vitro paradigm as used in typical STDP experiments (Bi & Poo, 1998).
In order to study the influence of the pre- and postsynaptic activity on thechanges of weights as predicted by our online learning rule in equation
2.20, we plotted in Figure 2 the postsynaptic factor B post and the correlation
term C j that both appear on the right-hand side of equation 2.20, together
with the induced weight change w/w as a function of time. Indeed, thelearning rule predicts positive weight changes when the presynaptic spike
occurs 10 ms before the postsynaptic one and negative weight changes
under reversed timing.
For a comparison with experimental results, we used 60 pairs of preand postsynaptic spikes applied at a frequency of 1 Hz and recorded thetotal change w in EPSP amplitude. The experiment is repeated with different spike timings, and the result is plotted as a function of spike timing
difference t pr e − t post . As in experiments (Markram et al., 1997; Zhang, Tao,
Holt, Harris, & Poo, 1998; Bi & Poo, 1998, 2001; Sjöström et al., 2001), we
find that synapses are potentiated if presynaptic spikes occur about 10 ms
before a postsynaptic action potential, but are depressed if the timing isreversed. Compared to synapses with amplitudes in the range of 1 or 2 mV,
synapses that are exceptionally strong show a reduced effect of potentiation

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

653

for pre-before-post timing, or even depression (see Figures 3a and 3b), inagreement with experiments on cultured hippocampal neurons (Bi & Poo,
1998). The shape of the STDP function depends only weakly on the stimulation frequency (see Figure 3c), even though a significant reduction of thepotentiation amplitude with increasing frequency can be observed.
In a recent experimental study (Froemke et al., 2005), a strong correlation
between the timescale of EPSP suppression (which was found to depend on
dendritic location) and the duration of the long-term depression (LTD) partin the STDP function was observed. Since our model neuron had no spatial
structure, we artificially changed the time constant of EPSP suppression
in the model equations. We found that indeed only the LTD part of the
STDP function was affected, whereas the LTP part remained unchanged(see Figure 3d).
In order to study the influence of the weight-dependent cost term  inour optimality criterion L, we systematically changed the parameter λ inequation 2.3. For λ = 0, the weight-dependent cost term has no influence,and because the postsynaptic firing rate is close to the desired rate, synaptic
plasticity in our model is mainly controlled by information maximization.
In this case, synapses with a reasonable EPSP amplitude of one or a few
millivolts are always strengthened, even for post-before-pre timing (seeFigure 3d, dashed line). This can be intuitively understood since an increase of synaptic weight is always beneficial for information transmission
except if spike arrival occurs immediately after a postsynaptic spike. In thiscase, the postsynaptic neuron is insensitive, so that no information can betransmitted. Nevertheless, information transmission is maximal in a situation where the presynaptic spike occurs just before the postsynaptic one.
The weight-dependent cost term derived from postulate C is essential toshift the dashed line in Figure 3d to negative values so as to induce synaptic
depression in our STDP paradigm. The optimal value of λ = 0.026, which
ensures that for large spike timing differences |t pr e − t post | neither potentiation nor depression occurs, has been estimated from a simple analyticalargument (see appendix A).
3.2 Both Unimodal and Bimodal Synapse Distributions Are Stable.
Under random spike arrival with a rate of 10 Hz at all 100 synapses, synaptic
weights show little variability with a typical EPSP amplitude in the range of
0.4 mV. This unspecific pattern of synapses stays stable even if 20 out of the
100 synapses are subject to a common rate modulation between 1 and 30 Hz(see Figure 4a). However, if modulation of presynaptic firing rates becomesstrong, the synapses rapidly develop a specific pattern with large values ofweights at synapses with rate-modulated spike input and weak weights at
those synapses that received input at fixed rates (synaptic specialization; see
Figures 4a and 4b), making the neuron highly selective to input at one groupof synapses (input selectivity). Thus, our synaptic update rule is capable of
selecting strong features in the input, but also allows a stable, unspecific654

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstnerpattern of synapses in case of weak input. This is in contrast to most otherHebbian learning rules, where unspecific patterns of synapses are unstable,
so that synaptic weights move spontaneously toward their upper or lower
bounds (Miller et al., 1989; Miller & MacKay, 1994; Gerstner, Kempter, vanHemmen, & Wagner, 1996; Kempter, Gerstner, & van Hemmen, 1999; Song,
Miller, & Abbott, 2000).
After induction of synaptic specialization by strong modulation of presynaptic input, we reduced the rate modulation back to the value that previously led to an unspecific pattern of synapses. We found that the strong
synapses remained strong and weak synapses remained weak; synaptic

Figure 4: Potentiation and depression depend on the presynaptic firing rate.
Twenty synapses of group 1 receive input with common rate modulation, while
the 80 synapses of group 2 (synapse index 21–100) receive Poisson input at aconstant rate of 10 Hz. The spike arrival rate in group 1 switches stochasticallyevery 200 ms between a low rate νl = 1 Hz and a high rate νh taken as aparameter. (a) Evolution of synaptic weights as a function of time for different
amplitudes of rate modulation, that is, νh changes from 10 Hz during the firsthour to 30 Hz, then to 50 Hz, 30 Hz, and back to 10 Hz. During the first 2 hours of
stimulation, an unspecific distribution of synapses remains stable even thougha slight decrease of weights in group 2 can be observed when the stimulus
switches to νh = 30 Hz. A specialization of the synaptic pattern with large
weights for synapses in group 1 is induced during the third hour of stimulation
and remains stable thereafter. (b) Top: Mean synaptic weights (same data as in a)of group 1 (w̄1 , blue line) and group 2 (w̄2 , green line). Bottom: The stimulation
paradigm νh as a function of time. Note that at νh = 30 Hz (second and fourthhour of stimulation), both an unspecific pattern of synapses with little difference
between w̄1 and w̄2 (second hour, top) and a highly specialized pattern (fourthhour, top, large difference between solid and dashed lines) are possible. (c) The
value of the objective function L (average value per second of time) in color code
as a function of the mean synaptic weight in group 1 (y-axis, w̄1 ) and group 2(x-axis, w̄2 ) during stimulation with νh = 30 Hz. Two maxima can be perceived:a broad maximum for the unspecific synapse pattern (w̄2 ≈ w̄1 ≈ 0.4 mV) anda pronounced elongated maximum for the specialized synapses pattern (w̄2 ≈0; w̄1 ≈ 0.8 mV). The dashed line indicates a one-dimensional section (see d)through the two-dimensional surface. (d) Objective function L as a function of
the difference w̄2 − w̄1 between the mean synaptic weights in groups 2 and 1
along the line indicated in c. For νh = 30 Hz (dashed, same data as in c) two
maxima are visible: a broad maximum at w̄2 − w̄1 ≈ 0 and a narrow but higher
maximum corresponding to the specialized synapse pattern to the very leftof the graph. For νh = 50 Hz (dotted line), the broad maximum at w̄2 − w̄1 ≈0 disappears, and only the specialized synapse pattern remains, whereas forνh = 10 Hz (solid line), the broad maximum of the unspecific synapse pattern
dominates.

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticityb
1.6

Synaptic index

90

1.4

80
70

1.2

60

1

50

0.8

40

0.6

30
0.4

20

0.2

10
0

60

120

180

240

1.5

¯ [mV]
w

100

1
0.5
0
0

60

120

180

240

300

60

120

180

240

300

60

νh [Hz]

a

40
20
0
0

300

Time [min]Time [min]
c

d
1

8

10Hz
30Hz
50Hz

6
0.8

4

0.6

3
0.4
2
0.2

0
0

7

5

1
0
0.2

0.4

¯ 2 [mV]w

0.6

L [s −1]w¯1 [mV]

655

6
5
4

3
−0.8

−0.6

−0.4

−0.2

0

0.2

0.4

¯2 − w¯1 [mV]w

specialization was stable against a change in the input (see Figure 4b). This
result shows that synaptic dynamics exhibits hysteresis, which is an indication of bistability: for the same input, both an unspecific pattern of synapses
and synaptic specialization are stable solutions of synaptic plasticity under
our learning rule. Indeed under rate modulation between 1 and 30 Hz for
20 out of the 100 synapses, the objective function L shows two local maxima(see Figure 4c): a sharp maximum corresponding to synaptic specialization(mean EPSP amplitude about 0.8 mV for synapses receiving rate-modulated
input and less than 0.1 mV for synapses receiving constant input) and abroader but slightly lower maximum where both groups of synapses havea mean EPSP amplitude in the range of 0.3 to 0.5 mV (see appendix B for
details on the method). In additional simulations, we confirmed that boththe unspecific pattern of synapses and the selective pattern representing
synaptic specialization remained stable over several hours of continued
stimulation with rate-modulated input (data not shown). Bistability of selective and unspecific synapse patterns was consistently observed for rates
modulated between 1 and 10 Hz or between 1 and 30 Hz, but the unspecificsynapse pattern was unstable if the rate was modulated between 1 and
50 Hz consistent with the weight dependence of our objective function L(see Figure 4d).

656

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

3.3 Retention of Synaptic Memories. In order to study synaptic memory retention with our learning rule, we induced synaptic specialization
by stimulating 20 out of the 100 synapses by correlated spike input (spikespike correlation index c = 0.2; see section 2.4 for details). The remaining
80 synapses received uncorrelated Poisson spike input. The mean firing rate
(10 Hz) was identical at all synapses. After 60 minutes of correlated input
at the group of 20 synapses, the stimulus was switched to uncorrelated
spike input at the same rate. We studied how well synaptic specialization was maintained as a function of time after induction (see Figures 5aand 5b).
Synaptic specialization was defined by a bimodality index that compared
the distribution of EPSP amplitudes at synapses that received correlated
input with those receiving uncorrelated input. For each of the two groupsof synapses, we calculated the mean w̄ = w and the variance σ 2 = [w j −w̄]2 : w̄ A and σ A2 for the group of synapses receiving correlated input and w̄ B
and σ B2 for those receiving uncorrelated input. We then approximated thetwo distributions by gaussian functions. The bimodality index depends on
the overlap between the two gaussians and is given by b = 0.5 [erf( w̄√A2σ−ŝ ) +A
x√ w̄ B )], where erf(x) = √2
erf( ŝ−exp(−t 2 )dt is the error function and ŝ is
π 0
2σ B
one of the two crossing points of the two gaussians such that w̄ B < ŝ < w̄ A.
The two distributions (strong and weak synapses) started to separate
within the first 5 minutes and remained well separated even after the correlated memory-inducing stimulus was replaced by a random stimulus (seeFigure 5b).
In order to study how synaptic memory retention depended on the induction paradigm, the experiment was repeated with different values of thecorrelation index c that characterizes the spike-spike correlations during
the induction period. For correlations c < 0.1, the two synaptic distributions are not well separated at the end of the induction period (bimodality
index < 0.9), but they are well separated for c ≥ 0.15 (see Figure 5c). Good
separation at the end of the induction period alone is, however, not sufficient to guarantee retention of synaptic memory, since for the synaptic
distribution induced by stimulation with c = 0.15, specialization breaks
down at t = 80 min., that is, after only 20 minutes of memory retention (seeFigure 5d). On the other hand, for c = 0.2 and larger, synaptic memory is
retained over several hours (only the first hour is plotted in Figure 5d).
The long duration of synaptic memory in our model can be explainedby the reduced adaptation speed of synapses with weights close to zero
(postulate C2). If weak synapses change only slowly because of reduced
adaptation speed, strong synapses must stay strong because of homeostatic processes that keep the mean activity of the postsynaptic neuron
close to a target value. Moreover, the terms in the online learning rule derived from information maximization favor the bimodal distribution. Reduced adaptation speed of weak synapses could be caused by a cascade ofOptimality Model of Unsupervised Spike-Timing-Dependent Plasticity

a

657

b
1

100
80

0.8

1.2
1

60

0.8
40

0.6
0.4

20

0.2
0

30

60

90

Distribution

Synaptic index

1.4

0.6
0.4
0.2
0
0

120

0.5

1

1.5

w [mV]

Time [min]
c

d

w¯

1

0.5

0

0

0.05

0.1

0.15

0.2

0.25

Spike-spike correlation

0.3

Specialization index

1.5
1c = 0.05
c = 0.1
c = 0.15
c = 0.2

0.8
0.6
0.4
0.2
0
0

30

60

90

120

Time [min]Figure 5: Synaptic memory induced by input with spike-spike correlations.
(a) Evolution of 100 synaptic weights (vertical axis) as a function of time. During the first 60 minutes, synapses of group A ( j = 1, . . . , 20) receive a Poisson
stimulation of 10 Hz with correlated spike input (c = 0.2), and those of groupB ( j = 21, . . . , 100) receive uncorrelated spike input at 10 Hz. (b) Distributionof the EPSP amplitudes across the 100 synapses after t = 5 min (dotted line),t = 60 min (solid line) and t = 120 min (dot-dashed line). The red lines denotegroup A while the blue ones group B. (c) Mean EPSP amplitude of group A(solid red line) and B (dashed blue line) at t = 60 min for different values ofcorrelation c of the input applied to group A. (d) Bimodality index b of the twogroups of weights as a function of time. The memory induction by correlated
spike input to group B stops at t = 60 min. Memory retention is studied during
the following 60 minutes. A bimodality index close to one implies, as for thecase with c = 0.2, that synaptic memory is well retained.

intracellular biochemical processing stages with different time constants,
as suggested by Fusi, Drew, & Abbott (2005). Thus our synaptic update
rule allows for retention of synaptic memories over timescales that are
significantly longer than the memory induction time, as necessary for658

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

any memory system. Nevertheless, synaptic memory in our model will
eventually decay if random firing of pre- and postsynaptic neurons persists, in agreement with experimental results (Abraham, Logan, Greenwood, & Dragunow, 2002; Zhou, Tao, & Poo, 2003). We note that in theabsence of presynaptic activity, the weights remain unchanged, since thedecay of synaptic weights is conditioned on presynaptic spike arrival (seeequation 2.20).
3.4 Receptive Field Development. Synaptic plasticity is thought to be
involved not only in memory (Hebb, 1949), but also in the development ofcortical circuits (Hubel & Wiesel, 1962; Katz & Shatz, 1996; von der Malsburg, 1973; Bienenstock et al., 1982; Miller et al., 1989) and, possibly, corticalreorganization (Merzenich, Nelson, Stryker, Cynader, Schoppmann, & Zook
1984; Buonomano & Merzenich, 1998). To study how our synaptic update
rule would behave during development, we used a standard paradigm of
input selectivity (Yeung, Shouval, Blais, & Cooper, 2004), which is considered to be a simplified scenario of receptive field development. Our modelneuron was stimulated by a gaussian firing rate profile spanned across the
100 input synapses (see Figure 6a). The center of the gaussian was shifted
every 200 ms to an arbitrarily chosen presynaptic neuron. In order to avoid
border effects, neuron number 100 was considered a neighbor of neuron
number 1, that is, we can visualize the presynaptic neurons as being located on a ring.
Nine postsynaptic neurons with slightly different initial values of synaptic weights received identical input from the same set of 100 presynaptic
neurons. During 1 hour of stimulus presentation, six of the nine neurons
developed synaptic specialization, leading to input selectivity. The optimal stimulus for these six neurons varies (see Figures 6b and 6c), so thatany gaussian stimulus at an arbitrary location excites at least one of thepostsynaptic neurons. In other words, the six postsynaptic neurons havedeveloped input selectivity with different but partially overlapping receptive fields. The distribution of synaptic weights for the selective neurons isbimodal with a first peak for very weak synapses (EPSP amplitudes lessthan 0.1 mV) and a second peak around EPSP amplitudes of 0.6 mV (seeFigure 6d); the amplitude distribution of the unselective neurons is broader,with a single peak at around 0.4 mV.
The number of postsynaptic neurons showing synaptic specialization
depends on the total stimulation time and the strength of the stimulus. If
the stimulation time is extended to 3 hours instead of 1 hour, all postsynaptic
neurons become selective using the same stimulation parameters as before.
However, if the maximal presynaptic firing rate at the center of the gaussian is reduced to 40 Hz instead of 50 Hz, only six of nine are selective after
3 hours of stimulation, and with a further reduction of the maximal rate to
30 Hz, only a single neuron is selective after 3 hours of stimulation (data notshown). We hypothesize that the coexistence of unselective and selective

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

a

659

b
60
50

ν [Hz]

40
30
20
10
0
0

20

40

60

80

100

Synaptic index
c

d
100

0.4

Synaptic index

90

0.2

80

0
0
0.2

70
60
50

0.2

0.4

0.6

0.8

0.2

0.4

0.6

0.8

0.4

0.6

0.8

0.1

40

0
0
0.2

30
20

0.1

10
1

3

5

7

Output neuron

9

0
0

0.2

w [mV]

Figure 6: The synaptic update rule leads to input selectivity of the postsynaptic neuron. (a) gaussian firing rate profile across the 100 presynaptic neurons.
The center of the gaussian is shifted randomly every 200 ms. Presynaptic neurons fire stochastically and send their spikes to nine postsynaptic neurons.
(b) Evolution of synaptic weights of the nine postsynaptic neurons. Some neurons become specialized for a certain input pattern at the early phase of learning,others become specialized later, and the last three neurons have not yet become
specialized. Since the input spike trains are identical for all the nine neurons,
the specialization is due to noise in the spike generator of the postsynaptic neurons. (c) Final synaptic weight values of the nine output neurons after 1 hourof stimulus presentation. (d) The distribution of EPSP amplitudes after 1 hourof stimulation for (top) the specialized output neurons 1, 2, . . . , 6; (middle) for
nonspecialized neurons 7, 8, 9; (bottom) for all nine output neurons.

neurons during development could explain the broad distribution of EPSP
amplitudes seen in some experiments (e.g., Sjöström et al., 2001, in rat visual cortex). For example, if we sample the synaptic distribution across allnine postsynaptic cells, we find the distribution shown at the bottom of
Figure 6d. If the number of unspecific neurons were higher, the relative importance of synapses with EPSP amplitudes of less than 0.1 mV would
diminish. If the number of specialized neurons increased, the distribution would turn into a clear-cut bimodal one, which would be akin to660

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

sampling an ensemble of two-state synapses with all-or-none potentiation
on a synapse-by-synapse basis (Petersen, Malenka, Nicoll, & Hopfield, 1998,
in rat hippocampus).
4 Discussion
4.1 What Can We and What Can We Not Expect from Optimality
Models? Optimality models can be used to clarify concepts, but they areunable to make specific predictions about molecular implementations. In
fact, the synaptic update rule derived in this article shares functional features with STDP and classical long-term potentiation (LTP), but it is blindwith respect to interesting questions such as the role of NMDA, Kainate,
endocannabinoid, or CaMKII in the induction and maintenance of potentiation and depression (Bliss & Collingridge, 1993; Bortolotto, Lauri, Isaac, &Collingridge, 2003; Frey & Morris, 1997; Lisman, 2003; Malenka & Nicoll,
1993; Sjöström, Turrigiano, & Nelson, 2004). If molecular mechanisms are
the focus of interest, detailed mechanistic models of synaptic plasticity(Senn, Tsodyks, & Markram, 2001; Yeung et al., 2004) should be preferred.
On the other hand, the mere fact that similar forms of LTP or LTD seem to be
implemented across various neural systems by different molecular mechanisms leads us to speculate that common functional roles of synapses arepotentially more important for understanding synaptic dynamics than thespecific way that these functions are implemented.
Ideally, optimality approaches such as the one developed in this article
should be helpful to put seemingly diverse experimental or theoreticalresults into a coherent framework. We have listed in section 1 a coupleof points, partially linked to experimental results and partially linked toearlier theoretical investigations. Our aim has been to connect these points
and trace them back to a small number of basic principles. Let us returnto our initial list and discuss the points in the light of the results of thepreceding section.
4.2 Correlations. Hebb (1949) postulated an increase in synaptic coupling in case of repeated coactivation of pre- and postsynaptic neurons asa useful concept for memory storage in recurrent networks. In our optimality framework defined by equation 2.3, correlation-dependent learning
is not imposed explicitly but arises from the maximization of informationtransmission between pre- and postsynaptic neurons. Indeed, informationtransmission is possible only if there are correlations between pre- andpostsynaptic neurons. Information transmission is maximized if these correlations are increased. Gradient ascent of the information term hence leads
to a synaptic update rule that is sensitive to correlations between pre- andpostsynaptic neurons (see equation 2.20). An increase of synaptic weights
enhances these correlations and maximizes information transmission. Weemphasize that in contrast to Hebb (1949), we do not invoke memory

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

661

formation and recall as a reason for correlation dependence, but information transmission. Similar to other learning rules (Linsker, 1986; Oja, 1982),
the sensitivity of our update rule to correlations between pre- and postsynaptic neurons gives the synaptic dynamics a sensitivity to correlations in
the input as demonstrated in Figures 4 to 6.
4.3 Input Selectivity. During cortical development, cortical neurons develop input selectivity typically quantified as the width of receptive fields(Hubel & Wiesel, 1962). As illustrated in the scenario of Figure 6, our synaptic update rule shows input selectivity and hence stands in the researchtradition of many other studies (von der Malsburg, 1973; Bienenstock et al.,
1982; Miller et al., 1989). Input selectivity in our model arises through the
combination of the correlation sensitivity of synapses discussed above withthe homeostatic term D in equation 2.3 (Toyoizumi et al., 2005a). Since thehomeostatic term keeps the mean rate of the postsynaptic neuron close toa target rate, it leads effectively to a normalization of the total synaptic input similar to the sliding threshold mechanism in the Bienenstock-CooperMunro rule (Bienenstock et al., 1982). Normalization of the total synaptic
input through firing rate stabilization has also been seen in previous STDP
models (Song et al., 2000; Kempter et al., 1999; Kempter, Gerstner, & vanHemmen, 2001). Its effect is similar to explicit normalization of synaptic
weights, a well-known mechanism to induce input selectivity (von der
Malsburg, 1973; Miller & MacKay, 1994), but our model does not needan explicit normalization step. So far we have studied only a single or asmall number of independent postsynaptic neurons, but we expect that as
in many other studies (e.g., Erwin, Obermayer, & Schulten, 1995; Song &Abbott, 2001; Cooper et al., 2004), our synaptic update rule would yieldfeature maps if applied to a network of many weakly interacting corticalneurons.
4.4 Unimodal versus Bimodal Synapse Distributions. In several previous models of rate-based or spike-timing-based synaptic dynamics, synaptic weights evolved always toward a bimodal distribution, with somesynapses close to zero and others close to maximal weight (Miller et al.,
1989; Miller & MacKay, 1994; Kempter et al., 1999; Song & Abbott, 2001;
Toyoizumi et al., 2005a). Thus, synapses specialize on certain features of the
input, even if the input has weak or no correlation at all, which seems questionable from a functional point of view and disagrees with experimentallyfound distributions of EPSP amplitudes in rat visual cortex (Sjöström et al.,2001). As shown in Figures 4 to 6, an unspecific pattern of synapses witha broad distribution of EPSP amplitudes is stable with our synaptic update rule if correlations in the input are weak. Synaptic specialization (bimodal distribution) develops only if synaptic inputs show a high degreeof correlations on the level of spikes or firing rates. Thus, neurons specialize only on highly significant input features, not in response to noise. A

662

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

similar behavior was noted in two recent models on spike-timingdependent plasticity (Gütig, Aharonov, Rotter, & Sompolinsky, 2003; Yeunget al., 2004). Going beyond those studies, we also demonstrated stability ofthe specialized synapse distribution over some time, even if the amount ofcorrelation through rate modulation is reduced after induction of synaptic
specialization stimulation (see Figure 4). Thus, for the same input characteristics, our model can show unimodal or bimodal distribution of synaptic
weight, depending on the stimulation history. Moreover, we note that in
the state of bimodal distribution, the EPSP amplitudes at the depressed
synapses are so small that in an experimental setting, they could easily
remain undetected or classified as a “silent” synapse (Kullmann, 2003).
A large proportion of silent synapses has been previously shown to beconsistent with optimal memory storage (Brunel, Hakim, Isope, Nadal, &Barbour, 2004). Furthermore, our results show that in a given populationof cells, neurons with specialized synapses can coexist with others thathave a broad and unspecific synapse distribution. We speculate that thesenonspecialized neurons could then be recruited later for new stimuli, as
hypothesized in earlier models of neural networks (Grossberg, 1987). In
passing, we note that in models of synaptic plasticity with multiplicative
weight dependence, the distribution of weights is always unimodal (Rubin,
Lee, & Sompolinsky, 2001; van Rossum, Bi, & Turrigiano, 2000), but in thiscase, retention of synaptic memories is difficult.
4.5 Synaptic Memory. In the absence of presynaptic input, synapses
in our model do not change significantly. Furthermore, our results showthat even in the presence of random pre- and postsynaptic firing activity,
synaptic memories can be retained over several hours, although a slow
decay occurs. Essential for long-term memory maintenance under random
spike arrival is the reduced adaptation speed for small values of synaptic
weights as formulated in postulate C2. This postulate is similar in spirit to a
recent theory of Fusi et al. (2005), with two important differences. First, wework with a continuum of synaptic states (characterized by the value of w),
whereas Fusi et al. assume a cascade of a finite number of discrete internal
synaptic states where transitions are unidirectional and characterized bydifferent time constants. The scheme of Fusi et al. guarantees slow (i.e.,
not exponential) decay of memories, whereas in our case, decay is always
exponential even if the time constant is weight dependent. Second, Fusi
et al. use a range of different time constants for both the potentiated and theunpotentiated state, whereas in our model, it is sufficient to have a slower
time constant for the unpotentiated state only. If the change of unpotentiatedsynapses is slow compared to homeostatic regulation of the mean firing rate,then the maintenance of the strong synapses is given by homeostasis (andalso supported by information maximization).
The fact that our model has been formulated in terms of a continuum ofsynaptic weights was for convenience only. Alternatively it is conceivable

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

663

to define a number of internal synaptic states that give rise to binary, or asmall number of, synaptic weight values (Fusi, 2002; Fusi et al., 2005; Abarbanel, Talathi, Gibb, & Rabinovich, 2005). The actual number of synaptic
states is unknown with conflicting evidence (Petersen et al., 1998; Lisman,2003).
4.6 Rate Dependence. Our results show that common rate modulationin one group of synapses strengthens these synapses if the modulation
amplitude is strong enough. In contrast, an increase of rates to a fixed valueof 40 Hz (without modulation) in one group of synapses while anothergroup of synapses receives background firing at 10 Hz does not lead to asynaptic specialization, only to a minor readjustment of weights (data notshown). For a comparison with experimental results, it is important to notethat rate dependence is typically measured with extracellular stimulation ofpresynaptic pathways. We assimilate repeated extracelluar stimulation witha strong and common modulation of spike arrival probability at one groupof synapses (as opposed to an increased rate of a homogeneous Poissonprocess). Under this interpretation, our results are qualitatively consistent
with experiments.
4.7 STDP. Our results show that our synaptic update rule shares several
features with STDP as found in experiments (Markram et al., 1997; Bi & Poo,
1998, 2001; Sjöström et al., 2001). The timescale of the potentiation part of
the STDP function depends in our model on the duration of EPSPs. The
timescale of the depression part is determined by the duration of EPSP
suppression in agreement with experiments (Froemke et al., 2005). Ourmodel shows that the relative importance of LTP and LTD depends on
the initial value of the synaptic weight in a way similar to that found inexperiments (Bi & Poo, 1998). However, for EPSP amplitudes between 0.5
and 2 mV, the LTP part clearly dominates over LTD in our model, which
seems to have less experimental support. Also, the frequency dependenceof STDP in our model is less pronounced than in experiments.
Models of STDP have previously been formulated on a phenomenological level (Gerstner et al., 1996; Song et al., 2000; Gerstner & Kistler, 2002)
or on a molecular level (Lisman, 2003; Senn et al., 2001; Yeung et al., 2004).
Only recently have models derived from optimality concepts moved to thecenter of interest (Toyoizumi et al., 2005a; Toyoizumi, Pfister, Aihara, &Gerstner, 2005b; Bohte & Mozer, 2005; Chechik, 2003). There are importantdifferences of the present model to the existing “optimal” models. Chechik(2003) used information maximization, but limited his approach to static
input patterns, while we consider arbitrary inputs. Bell and Parra (2005)
minimize output entropy, and Bohte and Mozer (2005) maximize spikereliability, whereas we maximize the information between full input andoutput spike trains. None of these studies considered optimization underhomeostatic and maintenance cost constraints.

664

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

After we introduced the homeostatic constraint in a previous study,
which gave rise to a learning rule with several interesting properties(Toyoizumi et al., 2005a), we realized that this model did not exhibit properties of STDP in an in vitro situation without some additional assumptions. Indeed, as shown in Figure 3d, the STDP function derived frominformation maximization alone exhibits no depression in the absenceof an additional weight-dependent cost term in the optimality function.
The weight-dependent cost term introduced in this letter hence plays acrucial role in STDP since it shifts the STDP function to more negativevalues.
4.8 How Realistic Is a Weight-Dependent Cost Term? The weightdependent cost term  in the optimality criterion L depends quadraticallyon the value of the weights of all synapses converging onto the same postsynaptic neuron. This turns out to be equivalent to a “decay” term in thesynaptic update (see the last term on the right-hand side of equation 2.20).
Such decay terms are common in the theoretical literature (Bienenstock
et al., 1982; Oja, 1982), but the question is whether such a decay term (leading to a slow depression of synapses) is realistic from a biological point ofview.
We emphasize that the decay term in our synaptic update rule is proportional to presynaptic activity. Thus, in contrast to existing models in
the theoretical literature (Bienenstock et al., 1982; Oja, 1982), a synapsethat receives no input is protected against slow decay. The specific form of
the “decay” term considered in this article was such that synaptic weights
decreased with each spike arrival, but presynaptic activity could also be
represented in the decay term by the mean firing rate rather than spikearrival, with no qualitative changes to the results.
An important aspect of our cost term is that only synapses that haverecently been activated are at risk regarding weight decay. We speculatethat the weight-dependent cost term could, in a loose sense, be related
to the number of plasticity factors that synapses require and compete forduring the first hours of synaptic maintenance (Fonseca et al., 2004). According to the synaptic tagging hypothesis (Frey & Morris, 1997), only
synapses that have been activated in the recent past compete for plasticity
factors, while unpotentiated synapses do not suffer from decay (Fonseca
et al., 2004). We emphasize that such a link of our cost term to the competition for plasticity factors is purely hypothetical. Many relevant detailsof tagging and competition for synaptic maintenance are omitted in our
approach.
4.9 Predictions and Experimental Tests. In order to achieve synaptic
memory that is stable over several hours, the reduced adaptation speed forweak synapses formulated in postulate C2 turns out to be essential. Thus, anessential assumption of our model is testable: for synapses with extremely

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

665

small EPSP amplitudes, in particular silent synapses, the induction of both
LTP and LTD should require stronger stimulation or stimulation sustained
over longer times, compared to synapses that are of average strength. This
aspect is distinct from models (Gütig et al., 2003) that postulate for weak
synapses a reduced adaptation speed for depression only, but maximal
effect of potentiation. Thus, comparison of LTP induction for silent synapses(Kullmann, 2003) with that for average synapses should allow differentiating between the two models. In an alternative formulation, synaptic memory could also be achieved by making strong synapses resistant to further
changes. As an aside, we note that the model of Fusi et al. (2005) assumes
a reduced speed of transition between several internal synaptic states sothat transition would not necessarily be visible as a change in synaptic
weight.
A second test of our model concerns the pattern of synaptic weights
converging on the same postsynaptic neuron. Our results suggest that earlyin development, most neurons would show an unspecific synapse pattern,that is, a distribution of EPSP amplitudes with a single but broad peak,whereas later, a sizable fraction of neurons would show a pattern of synaptic
specialization with some strong synapses and many silent ones, that is, abimodal distribution of EPSP amplitudes. Ideally the effect would be seen
by scanning all the synapses of individual postsynaptic neurons; it remains
to be seen if modern imaging and staining methods will allow doing this.
Alternatively, by electrophysiological methods, distributions of synaptic
strengths could be built up by averaging over many synapses on different
neurons (Sjöström et al., 2001). In this case, our model would predict thatduring development, the histogram of EPSP amplitudes would change in
two ways (see Figure 6d): (1) the number of silent synapses increases sothat the amplitude of the sharp peak at small EPSP amplitude grows, and(2) the location of the second, broader peak shifts to larger values of the
EPSP amplitudes. Furthermore, and in contrast to other models where theunimodal distribution is unstable, the transition to a bimodal distribution
depends in our model on the stimulation paradigm.
4.10 Limitations and Extensions. We emphasize that properties of oursynaptic update rules have so far been tested only for single neurons in
an unsupervised learning paradigm. Extensions are possible in several directions. First, instead of single neurons, a large, recurrent network could
be considered. This could, on one side, further our understanding of themodel properties in the context of cortical map development (Erwin et al.,1995), and on the other side, scrutinize the properties of the synaptic update
rule as a functional memory in recurrent networks (Amit & Fusi, 1994). Second, instead of unsupervised learning where the synaptic update rule treats
all stimuli alike whether they are behaviorally relevant or not, a rewardbased learning scheme could be considered (Dayan & Abbott, 2001; Seung,
2003; Pfister et al., 2006). Behaviorally relevant situations can be taken into

666

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstneraccount by optimizing reward instead of information transmission (Schultz,
Dayan, & Montague, 1997).
Appendix A: Determination of the Parameter λThe parameter λ is set to give w j = 0 for large enough |t pr e − t post | in asimulated STDP in vitro paradigm. In order to find the appropriate valueof λ, we separately consider the effects of a presynaptic spike and a postsynaptic one, which is possible since they are assumed to occur at a large
temporal distance. Since a postsynaptic spike alone does not change synaptic strength (C j (t) = 0, always), we choose a λ that gives no synaptic change
when a presynaptic spike alone is induced. For a given presynaptic spikef
at t j , we have
C j (t) = −g

 t
0



f
t − t j e −(t−t )/τm dt  f f
τC τm− t−t j /τm−(t−t j )/τC
e
.
= −g
−e
τC − τm(A.1)

Since ρ̃ ≈ ρr , and ρ̄ ≈ ρr in this in vitro setting, the factor B post in equation
2.14 is approximated as
 f
B post (t) ≈ −w j ge − t−t j /τm .

(A.2)

Hence, we find the effect of a presynaptic spike as
w =

 T
0

dw j
τC τm
dt = w j g 2dtτC − τm




τC τm
τm
−
− λw j .
τC + τm
2

τC
The condition of no synaptic change gives λ = g 2 ττCm−τmused this λ in the numerical code.



τm τC
− τ2m
τm +τC

(A.3)

. We

Appendix B: Weight-Dependent Evaluation of the Optimality
Criterion
In Figures 4c and 4d the optimality criterion has been evaluated as a function of some artifical weight distribution. Specifically, values of synaptic
weights have been chosen stochastically from two gaussian distributions
with mean w̄1 and standard deviation σ1 for group 1 and w̄2 and σ2 forgroup 2. In order to account for differences in standard deviations due to theOptimality Model of Unsupervised Spike-Timing-Dependent Plasticity

667

weight-dependent update rate α(w), we chose σ (w̄) = 0.1 mV ·w̄ 4 /(ws4 +w̄ 4 ), which gives a variance of synaptic weights in both groups, which isconsistent with the variance seen in Figure 4A.
For a fixed set of synaptic weight values, the network is simulated during
a trial time of 30 minutes, while the synaptic updated rule has been turned
off and the objective function L defined in equation 2.3 is evaluated usingρ̄ est from equation 2.21. The result is divided by the trial time t and plotted in
Figures 4c and 4d in units of s −1 . The mesh size of mean synaptic strength is
0.04 mV. The one-dimensional plot in Figure 4d is taken along the directionw̄2 = 0.8 mV −w̄1 .
Acknowledgments
This research is supported by the Japan Society for the Promotion of Science and Grant-in-Aid for JSPS fellows 03J11691 (T.T.), the Swiss National
Science Foundation 200020-108097/1 (J.P.P), and Grant-in-Aid for Scientific
Research on Priority Areas 17022012 from MEXT of Japan (K.A.).
References
Abarbanel, H., Talathi, S., Gibb, L., & Rabinovich, M. (2005). Synaptic plasticity withdiscrete state synapses. Physical Review E, 72, art. no. 031917 Part 1.
Abraham, W., Logan, B., Greenwood, J., & Dragunow, M. (2002). Induction andexperience-dependent consolidation of stable long-term potentiation lasting
months in the hippocampus. J. Neuroscience, 22, 9626–9634.
Amit, D., & Fusi, S. (1994). Learning in neural networks with material synapses.
Neural Computation, 6, 957–982.
Atick, J., & Redlich, A. (1990). Towards a theory of early visual processing. NeuralComputation, 4, 559–572.
Barlow, H. (1956). Retinal noise and absolute threshold. J. Opt. Soc. Am., 46, 634–639.
Barlow, H. B. (1961). Possible principles underlying the transformation of sensory messages. In W. A. Rosenbluth (Ed.), Sensory communication (pp. 217–234).
Cambridge, MA: MIT Press.
Bell, A. J., & Parra, L. C. (2005). Maximising sensitivity in a spiking network. In L. K.
Saul, Y. Weiss, & L. Bottou (Eds.), Advances in neural information processing systems,
17 (pp. 121–128). Cambridge, MA: MIT Press.
Bell, A., & Sejnowski, T. (1995). An information maximization approach to blindseparation and blind deconvolution. Neural Computation, 7, 1129–1159.
Bi, G., & Poo, M. (1998). Synaptic modifications in cultured hippocampal neurons:
Dependence on spike timing, synaptic strength, and postsynaptic cell type. J.
Neurosci., 18, 10464–10472.
Bi, G., & Poo, M. (2001). Synaptic modification of correlated activity: Hebb’s postulate
revisited. Ann. Rev. Neurosci., 24, 139–166.
Bienenstock, E., Cooper, L., & Munroe, P. (1982). Theory of the development of neuron selectivity: Orientation specificity and binocular interaction in visual cortex.
Journal of Neuroscience, 2, 32–48.

668

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. Gerstner

Bliss, T. V. P., & Collingridge, G. L. (1993). A synaptic model of memory: Long-term
potentiation in the hippocampus. Nature, 361, 31–39.
Bliss, T., & Lomo, T. (1973). Long-lasting potentation of synaptic transmission in thedendate area of anaesthetized rabbit following stimulation of the perforant path.
J. Physiol., 232, 351–356.
Bohte, S. M., & Mozer, M. C. (2005). Reducing spike train variability: A computational theory of spike-timing dependent plasticity. In L. K. Saul, Y. Weiss, & L.
Bottou (Eds.), Advances in neural information processing systems, 17 (pp. 201–208).
Cambridge, MA: MIT Press.
Bortolotto, Z. A., Lauri, S., Isaac, J. T. R., & Collingridge, G. L. (2003). Kainate receptors and the induction of mossy fibre long-term potentiation. Phil. Trans. R. Soc.
Lond B: Biological Sciences, 358, 657–666.
Britten, K., Shadlen, M., Newsome, W., & Movshon, J. (1992). The analysis of visualmotion: A comparison of neuronal and psychophysical performance. J. Neuroscience, 12, 4745–4765.
Brunel, N., Hakim, V., Isope, P., Nadal, J.-P., & Barbour, B. (2004). Optimal information
storage and the distribution of synaptic weights: Perceptron versus Purkinje cell.
Neuron, 43, 745–757.
Buonomano, D., & Merzenich, M. (1998). Cortical plasticity: From synapses to maps.
Annual Review of Neuroscience, 21, 149–186.
Chechik, G. (2003). Spike-timing-dependent plasticity and relevant mututal information maximization. Neural computation, 15, 1481–1510.
Cooper, L., Intrator, N., Blais, B., & Shouval, H. Z. (2004). Theory of cortical plasticity.
Singapore: World Scientific.
Cover, T., & Thomas, J. (1991). Elements of information theory. New York: Wiley.
Dayan, P., & Abbott, L. F. (2001). Theoretical neuroscience. Cambridge, MA: MIT Press.
de Ruyter van Steveninck, R. R., & Bialek, W. (1995). Reliability and statistical efficiency of a blowfly movement-sensitive neuron. Phil. Trans. R. Soc. Lond. Ser. B,
348, 321–340.
Dudek, S. M., & Bear, M. F. (1992). Homosynaptic long-term depression in area CA1of hippocampus and effects of N-methyl-D-aspartate receptor blockade. Proc.
Natl. Acad. Sci. USA, 89, 4363–4367.
Erwin, E., Obermayer, K., & Schulten, K. (1995). Models of orientation and oculardominance columns in the visual cortex: A critical comparison. Neural Comput.,
7, 425–468.
Fonseca, R., Nägerl, U., Morris, R., & Bonhoeffer, T. (2004). Competition for memory:
Hippocampal LTP under the regimes of reduced protein synthesis. Neuron, 44,
1011–1020.
Frey, U., & Morris, R. (1997). Synaptic tagging and long-term potentiation. Nature,
385, 533–536.
Froemke, R. C., Poo, M.-M., & Dan, Y. (2005). Spike-timing-dependent synaptic
plasticity depends on dendritic location. Nature, 434, 221–225.
Fusi, S. (2002). Hebbian spike-driven synaptic plasticity for learning patterns of mean
firing rates. Biol. Cybern., 87, 459–470.
Fusi, S., Annunziato, M., Badoni, D., Salamon, A., & Amit, D. J. (2000). Spike-driven
synaptic plasticity: Theory, simulation, VLSI implementation. Neural Computation,
12, 2227–2258.

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

669

Fusi, S., Drew, P., & Abbott, L. (2005). Cascade models of synaptically stored memories. Neuron, 45, 599–611.
Gerstner, W., Kempter, R., van Hemmen, J. L., & Wagner, H. (1996). A neuronallearning rule for sub-millisecond temporal coding. Nature, 383, 76–78.
Gerstner, W., & Kistler, W. K. (2002). Spiking neuron models. Cambridge: Cambridge
University Press.
Grossberg, S. (1987). The adaptive brain I. Dordrecht: Elsevier.
Gütig, R., Aharonov, R., Rotter, S., & Sompolinsky, H. (2003). Learning input correlations through non-linear temporally asymmetry Hebbian plasticity. J. Neuroscience, 23, 3697–3714.
Hebb, D. O. (1949). The organization of behavior. New York: Wiley.
Hertz, J., Krogh, A., & Palmer, R. G. (1991). Introduction to the theory of neural computation. Redwood City, CA: Addison-Wesley.
Hubel, D. H., & Wiesel, T. N. (1962). Receptive fields, binocular interaction andfunctional architecture in the cat’s visual cortex. J. Physiol. (London), 160, 106–
154.
Katz, L. C., & Shatz, C. J. (1996). Synaptic activity and the construction of corticalcircuits. Science, 274, 1133–1138.
Kelso, S. R., Ganong, A. H., & Brown, T. H. (1986). Hebbian synapses in hippocampus.
Proc. Natl. Acad. Sci. USA, 83, 5326–5330.
Kempter, R., Gerstner, W., & van Hemmen, J. L. (1999). Hebbian learning and spiking
neurons. Phys. Rev. E, 59, 4498–4514.
Kempter, R., Gerstner, W., & van Hemmen, J. L. (2001). Intrinsic stabilization of
output rates by spike-based Hebbian learning. Neural Computation, 13, 2709–
2741.
Kullmann, D. M. (2003). Silent synapses: What are they telling us about long-term
potentiation? Phil. Trans. R. Soc. Lond. B: Biological Sciences, 358, 727–733.
Laughlin, S. (1981). A simple coding procedure enhances a neuron’s information
capacity. Z. Naturforschung, 36, 910–912.
Laughlin, S., de Ruyter van Steveninck, R. R., & Anderson, J. (1998). The metabolic
cost of neural information. Nature Neuroscience, 1, 36–41.
Levy, W., & Baxter, R. (2002). Energy-efficient neuronal computation via quantal
synaptic failures. J. Neuroscience, 22, 4746–4755.
Linsker, R. (1986). From basic network principles to neural architecture: Emergenceof spatial-opponent cells. Proc. Natl. Acad. Sci. USA, 83, 7508–7512.
Linsker, R. (1989). How to generate ordered maps by maximizing the mutual information between input and output signals. Neural Computation, 1(3), 402–
411.
Lisman, J. (1989). A mechanism for Hebb and anti-Hebb processes underlying learning and memory. Proc. Natl. Acad. Sci. USA, 86, 9574–9578.
Lisman, J. (2003). Long-term potentiation: Outstanding questions and attemptedsynthesis. Phil. Trans. R. Soc. Lond B: Biological Sciences, 358, 829–842.
Malenka, R. C., Kauer, J., Zucker, R., & Nicoll, R. A. (1988). Postsynaptic calcium
is sufficient for potentiation of hippocampal synaptic transmission. Science, 242,
81–84.
Malenka, R. C., & Nicoll, R. A. (1993). NMDA-receptor-dependent plasticity: Multiple forms and mechanisms. Trends Neurosci., 16, 480–487.

670

T. Toyoizumi, J.-P. Pfister, K. Aihara, and W. GerstnerMarkram, H., Lübke, J. L, Frotscher, M., & Sakmann, B. (1997). Regulation of synaptic
efficacy by coincidence of postysnaptic AP and EPSP. Science, 275, 213–215.
Merzenich, M., Nelson, R., Stryker, M., Cynader, M., Schoppmann, A., & Zook, J.
(1984). Somatosensory cortical map changes following digit amputation in adultmonkeys. J. Comparative Neurology, 224, 591–605.
Miller, K., Keller, J. B., & Stryker, M. P. (1989). Ocular dominance column development: Analysis and simulation. Science, 245, 605–615.
Miller, K. D., & MacKay, D. J. C. (1994). The role of constraints in Hebbian learning.
Neural Computation, 6, 100–126.
Oja, E. (1982). A simplified neuron model as a principal component analyzer. J.
Mathematical Biology, 15, 267–273.
Petersen, C., Malenka, R., Nicoll, R., & Hopfield, J. (1998). All-or-none potentiationof CA3-CA1 synapses. Proc. Natl. Acad. Sci. USA, 95, 4732–4737.
Pfister, J.-P., Toyoizumi, T., Barber, D., & Gerstner, W. (2006). Optimal spike-timing
dependent plasticity for precise action potential firing in supervised learning.
Neural Computation, 18, 1309–1339.
Rubin, J., Lee, D. D., & Sompolinsky, H. (2001). Equilibrium properties of temporally
asymmetric Hebbian plasticity. Physical Review Letters, 86, 364–367.
Schultz, W., Dayan, P., & Montague, R. (1997). A neural substrate for prediction andreward. Science, 275, 1593–1599.
Senn, W., Tsodyks, M., & Markram, H. (2001). An algorithm for modifying neurotransmitter release probability based on pre- and postsynaptic spike timing.
Neural Computation, 13, 35–67.
Seung, H. (2003). Learning in spiking neural networks by reinforcement of stochastic
synaptic transmission. Neuron, 40, 1063–1073.
Sjöström, P., Turrigiano, G., & Nelson, S. (2001). Rate, timing, and cooperativity
jointly determine cortical synaptic plasticity. Neuron, 32, 1149–1164.
Sjöström, P., Turrigiano, G., & Nelson, S. (2004). Endocannabinoid-dependent neocortical layer-5 LTD in the absence of postsynaptic spiking. J. Neurophysiol., 92,
3338–3343.
Song, S., & Abbott, L. (2001). Column and map development and cortical re-mapping
through spike-timing dependent plasticity. Neuron, 32, 339–350.
Song, S., Miller, K., & Abbott, L. (2000). Competitive Hebbian learning throughspike-time-dependent synaptic plasticity. Nature Neuroscience, 3, 919–926.
Toyoizumi, T., Pfister, J.-P., Aihara, K., & Gerstner, W. (2005a). Generalized
Bienenstock-Cooper-Munro rule for spiking neurons that maximizes information transmission. Proc. National Academy Sciences (USA), 102, 5239–5244.
Toyoizumi, T., Pfister, J.-P., Aihara, K., & Gerstner, W. (2005b). Spike-timing dependent plasticity and mutual information maximization for a spiking neuron model.
In L. K. Saul, Y. Weiss, & L. Bottou (Eds.), Advances in neural information processing
systems, 17 (pp. 1409–1416). Cambridge, MA: MIT Press.
Turrigiano, G., & Nelson, S. (2004). Homeostatic plasticity in the developing nervous
system. Nature Reviews Neuroscience, 5, 97–107.
van Rossum, M. C. W., Bi, G. Q., & Turrigiano, G. G. (2000). Stable Hebbian learningfrom spike timing-dependent plasticity. J. Neuroscience, 20, 8812–8821.
von der Malsburg, C. (1973). Self-organization of orientation selective cells in the
striate cortex. Kybernetik, 14, 85–100.

Optimality Model of Unsupervised Spike-Timing-Dependent Plasticity

671

Yeung, L. C., Shouval, H., Blais, B., & Cooper, L. (2004). Synaptic homeostasis andinput selectivity follow from a model of calcium dependent plasticity. Proc. Nat.
Acad. Sci. USA, 101, 14943–14948.
Zhang, L., Tao, H., Holt, C., Harris, W. A., & Poo, M.-M. (1998). A critical window forcooperation and competition among developing retinotectal synapses. Nature,395, 37–44.
Zhou, Q., Tao, H., & Poo, M. (2003). Reversal and stabilization of synaptic modifications in the developing visual system. Science, 300, 1953–1957.

Received January 13, 2006; accepted June 22, 2006.