Abstract
Classical experiments on spike timingdependent plasticity (STDP) use a protocol based on pairs of presynaptic and postsynaptic spikes repeated at a given frequency to induce synaptic potentiation or depression. Therefore, standard STDP models have expressed the weight change as a function of pairs of presynaptic and postsynaptic spike. Unfortunately, those pairedbased STDP models cannot account for the dependence on the repetition frequency of the pairs of spike. Moreover, those STDP models cannot reproduce recent triplet and quadruplet experiments. Here, we examine a triplet rule (i.e., a rule which considers sets of three spikes, i.e., two pre and one post or one pre and two post) and compare it to classical pairbased STDP learning rules. With such a triplet rule, it is possible to fit experimental data from visual cortical slices as well as from hippocampal cultures. Moreover, when assuming stochastic spike trains, the triplet learning rule can be mapped to a Bienenstock–Cooper–Munro learning rule.
Introduction
During the last decade, an increasing number of experiments have shown that synaptic strength changes as a function of the precise spike timing of the presynaptic and postsynaptic neurons. In the early experiments (Markram et al., 1997; Bi and Poo, 1998, 2001; Zhang et al., 1998), potentiation has been elicited by a sequence of n pairs of “pre then post” spikes, whereas depression occurred when the timing was reversed (i.e., when each postsynaptic spike precedes a presynaptic one). At this point, it was natural to characterize synaptic plasticity as a function of the time difference Δt = t^{post} − t^{pre} between pairs of spikes. However, performing experiments with pairs of spikes does not mean that pairs of spikes are the elementary building block. There is no a priori reason to think that pairs of spikes are more relevant than three spikes (triplets), four spikes (quadruplets), or even more. It is clear that a lot of other neuronal variables such as calcium concentration (Malenka et al., 1988; Lisman, 1989; Lisman and Zhabotinsky, 2001; Shouval et al., 2002) or postsynaptic membrane potential (Rao and Sejnowski, 2001; Sjöström et al., 2001; Lisman and Spruston, 2005) play an important role in triggering potentiation or depression. The point of this study was to see how far we can explain experiments that only use spike timing as a parameter with models that only use spike timing.
Recent experiments (Bi and Wang, 2002; Froemke and Dan, 2002; Wang et al., 2005; Froemke et al., 2006) have studied the detailed role of spike timing by triggering synaptic plasticity with spike triplets (one presynaptic spike combined with two postsynaptic spikes or one postsynaptic spike with two presynaptic spikes). The results of those experiments indicate that classical spike timingdependent plasticity (STDP) models based on pairs of spikes are not sufficient to explain synaptic changes triggered by triplets or quadruplets of spikes.
In the first part of this study, we review some experimental protocols performed in visual cortex (Sjöström et al., 2001) and hippocampal culture (Wang et al., 2005) and show why the classical pairbased STDP models fail to reproduce those experimental data. In the second part of this study, we show that if we assume that synaptic plasticity is governed by a suitable combination of pairs and triplets of spikes, the results from the above mentioned protocols can be surprisingly well reproduced. Moreover, we show that our triplet learning rule elicits input selectivity analogous to that of the Bienenstock–Cooper–Monro (BCM) theory (Bienenstock et al., 1982).
Claiming that triplet of spikes are more relevant than pairs of spikes is not enough to construct a model of synaptic plasticity. It is also necessary to determine how those pairs or triplets of spikes integrate. For both the pairbased models and the tripletbased models, we consider the case in which a presynaptic spike interacts with all previous postsynaptic ones or vice versa (we call this the AlltoAll interaction) (Gerstner et al., 1996; Kempter et al., 1999; Kistler and van Hemmen, 2000; Song et al., 2000) and the case where only neighboring spikes are taken into account (NearestSpike interaction) (van Rossum et al., 2000; Bi, 2002; Izhikevich and Desai, 2003; Burkitt et al., 2004; Pfister and Gerstner, 2006). We found a slight preference for AlltoAll interactions.
Materials and Methods
We compared a new tripletbased model with experimental data from the hippocampus and visual cortex. The visual cortex data set (Sjöström et al., 2001) used in this study consists of a standard pairing protocol in which the frequency of the pairing has been changed. We also considered a hippocampal culture data set (Wang et al., 2005) which consists of pair, triplet, and quadruplet protocols. Because both data sets disagree on some specific protocols (at low frequency of the pairing protocol, no potentiation is elicited in Sjöström's data, whereas a large amount of potentiation is present in Wang's data) and because the preparations are different, we fitted our models with different parameters for each data set.
Synaptic learning rule.
Our new tripletbased model of STDP is an extension of classical pairbased STDP models. Traditional mechanistic models of STDP involve a small number of variables that are updated by presynaptic and postsynaptic firing events (Kistler and van Hemmen, 2000; Abarbanel et al., 2002; Gerstner and Kistler, 2002; Karmarkar and Buonomano, 2002). The new triplet rule is formulated in such a framework.
To introduce the variables used in our model, we considered the process of synaptic transmission. Whenever a presynaptic spike arrives at an excitatory synapse, glutamate is released into the synaptic cleft and binds to glutamate receptors. Let r_{1} denote the amount of glutamate bound to a postsynaptic receptor. The variable r_{1} increases whenever there is a presynaptic spike and decreases back to zero otherwise with a time constant of τ_{+}. This can be written as follows:
Here, t^{pre} denotes the moment of spike arrival at the presynaptic terminal. The units of r_{1} are chosen such that glutamate binding increases by one unit after spike arrival. We emphasize that r_{1} is an abstract variable. Instead of glutamate binding, it could describe equally well some other quantity that increases after presynaptic spike arrival. We call r_{1} a detector of presynaptic events.
Instead of having only one process triggered by a presynaptic spike, it is possible to consider several different quantities, which increase in the presence of a presynaptic spike. In our model, we considered two different detectors of presynaptic events, namely r_{1} and r_{2}. The dynamics of r_{2} is analogous to that of r_{1} except that its time constant τ_{x} is larger than τ_{+}. Similarly, we assume that each postsynaptic spike t^{post} induces an increase of two different quantities that we denote o_{1} and o_{2}. Potential interpretations of o_{1} and o_{2} are given below. In the absence of postsynaptic spiking, these postsynaptic detectors decrease their value with a time constant τ_{−} and τ_{y}, respectively. Formally, this gives the following:
We do not want to identify the variables r_{1}, r_{2}, o_{1}, and o_{2} with specific biophysical quantities. Candidates of detectors of presynaptic events are, for example, the amount of glutamate bound (Karmarkar and Buonomano, 2002) or the number of NMDA receptors in an activated state (Senn et al., 2001). Postsynaptic detectors o_{1} and o_{2} could represent the influx of calcium concentration through voltagegated Ca^{2+} channels and NMDA channels (Karmarkar and Buonomano, 2002) or the number of secondary messengers in a deactivated state of the NMDA receptor (Senn et al., 2001) or the voltage trace of a backpropagating action potential (Shouval et al., 2002).
Because our present model is formulated as a mechanistic model, it is possible to define changes of synaptic efficacies for our triplet learning rule with AlltoAll interactions as a function of those four detectors without making any assumption on the biophysical quantities they represent. We assume that the weight decreases after presynaptic spike arrival by an amount that is proportional to the value of the postsynaptic variable o_{1} but depends also on the value of the second presynaptic detector r_{2}. Hence, presynaptic spike arrival at time t^{pre} triggers a change given by the following:
Similarly, a postsynaptic spike at time t^{post} triggers a change that depends on the presynaptic variable r_{1} and the second postsynaptic variable o_{2} as follows:
Here, A_{2}^{+} and A_{2}^{−} denote the amplitude of the weight change whenever there is a prepost pair or a postpre pair. Similarly, A_{3}^{+} and A_{3}^{−} denote the amplitude of the triplet term for potentiation and depression, respectively (Fig. 1A). All of the four amplitude parameters are assumed to be greater or equal to zero. ε is a small positive constant to ensure that the weight is updated before the detectors r_{2} or o_{2}. In other words, r_{2} is zero unless a previous presynaptic spike has led to an increase of r_{2}. This ensures the detection of spike triplets.
Figure 1B illustrates how a 1pre2post triplet is detected by the learning rule. At the time of a postsynaptic spike, the learning rule “reads” the value of the second postsynaptic variable o_{2} just before the spike (see the red dot at time t^{post} − ε in Fig. 1B) as well as the value of the presynaptic detector r_{1} (see the blue dot at time t^{post} in Fig. 1B) and increases the weight by an amount A_{2}^{+}r_{1}(t^{post}) + A_{3}^{+}r_{1}(t^{post})o_{2}(t^{post}  ε) (see Eq. 4).
Note that if we set A_{3}^{+} = 0 and A_{3}^{−} = 0, the model becomes a classical pairbased STDP model (Gerstner et al., 1996; Kempter et al., 1999; Kistler and van Hemmen, 2000; Song et al., 2000). This pairbased STDP model was used for the results of Figure 2. It should be further noted that the two extra triplet terms vanish if a single spike pair is presented or if spike pairings are repeated at low frequency. This means that in the limit of low frequency, the classical pairbased learning rule is identical to our triplet learning rule.
The triplet learning rule of Equations 3 and 4 can also describe a NearestSpike interaction scheme if we redefine the update rule of presynaptic and postsynaptic detectors. Instead of simply lowpass filtering the spike trains (i.e., adding the effects of all spikes), the detector variables saturate at 1 (i.e., 0 ≤ r_{1}, r_{2}, o_{1}, o_{2} ≤ 1). This is achieved by updating the variables to the value of 1 instead of updating by at step of 1. In this way, the synapse forgets all other previous spikes and keeps only the memory of the last one (Fig. 1C).
In this paper, we consider first a full triplet model, which takes into account all four terms of Equations 3 and 4. Then, we will see that only some of the terms are really necessary. This is why we define two different minimal models. The first one is intended to fit the visual cortex data and disregards two terms (i.e., A_{2}^{+} = 0 and A_{3}^{−} = 0). For the hippocampal culture data set, we consider a slightly different minimal model, which disregards only one term (i.e., A_{3}^{−} = 0).
In principle, the amplitude parameters A_{2}^{+}, A_{2}^{−}, A_{3}^{+}, and A_{3}^{−} could change on a slow time scale. For example, similar to the threshold in the BCM rule, those parameters could, because of homeostatic processes, depend on the mean postsynaptic firing rate ρ_{y} averaged over a time scale of 10 min or more.
Protocols.
To compare our model to experimental data, we followed three different experimental protocols (see Fig. 2) in which the synaptic weight changes as a function of the presynaptic and postsynaptic spike statistics. The forth protocol is of a more theoretical value in the sense that it can be compared with the BCM learning rule, which has interesting computational properties.
Pairing protocol.
This is the classical STDP protocol (see Fig. 2A) (Markram et al., 1997; Bi and Poo, 1998, 2001; Zhang et al., 1998; Sjöström et al., 2001; Froemke and Dan, 2002). n = 60 pairs of presynaptic and postsynaptic spikes shifted by Δt are elicited at regular intervals of 1/ρ. The interest of the study of Sjöström et al. (2001) is that the authors analyzed, in this pairing protocol, the weight change as a function of the frequency ρ for a fixed Δt. Changing the frequency ρ is a good way to check the validity of a model, especially at high frequency, where many spikes are potentially in the temporal range of interaction.
It should be noted that the amount of potentiation for a prepost (Δt = 10 ms) pair reported by Wang et al. (2005) is significantly lower than the one originally measured (Bi and Poo, 1998) under the same conditions. As mentioned by Wang et al. (2005), this can be accounted for by the difference in initial synaptic strength, which was higher in the study by Bi and Poo (1998). To test our model on a consistent set of data, we took the measurements of Wang et al. (2005) (compare their supplemental Fig. 1) (i.e., Δw ≃ 0.25 ± 0.05 for Δt = +10 ms and Δw ≃ − 0.17 ± 0.05 for Δt = −10 ms. Data from Wang et al. (2005), including error bars, are redrawn in Figure 3. Because the potentiation and depression time constant are not present in the study by Wang et al. (2005), we took τ_{+} = 16.8 ms and τ_{−} = 33.7 from (Bi and Poo, 2001).
Triplet protocol.
The first triplet protocol (see Fig. 2C) consists of n = 60 sets of three spikes repeated at a given frequency ρ = 1 Hz. Each triplet consists of two presynaptic spikes and one postsynaptic spike characterized by Δt_{1} = t^{post} − t_{1}^{pre} and Δt_{2} = t^{post} − t_{2}^{pre} where t_{1}^{pre} and t_{2}^{pre} are the first and second presynaptic spikes of the triplet.
The second triplet protocol (see Fig. 2D) also consists of n = 60 triplets. The only difference is that each triplet consists of one presynaptic and two postsynaptic spikes. In this case, Δt_{1} = t_{1}^{post} − t^{pre} and Δt_{2} = t_{2}^{post} − t^{pre}, where t_{1}^{post} and t_{2}^{post} are, respectively, the first and second postsynaptic spikes of the triplet.
Experiments with such a triplet protocol have been performed by Froemke and Dan (2002) in L2/3 pyramidal neurons of the rat visual cortex and by Wang et al. (2005) in hippocampal cultures. To have a consistent and broad data set (i.e., pair, triplet, and quadruplet experiments), we decided, in the present study, to focus only on the data of Wang et al. (2005), because we did not find enough quantitative information about quadruplets in the study by Froemke and Dan (2002).
Quadruplet protocol.
This protocol consists of n = 60 quadruplets at frequency ρ = 1 Hz (see Fig. 2B). It was used by Wang et al. (2005) and is characterized as follows: a postpre pair with a delay of Δt_{1} = t_{1}^{post} − t_{1}^{pre} < 0 is followed after a time T by a prepost pair with a delay of Δt_{2} = t_{2}^{post} − t_{2}^{pre} > 0. When T is negative, the opposite happens. A prepost pair (Δt_{2} = t_{2}^{post} − t_{2}^{pre} > 0) is followed by a postpre pair (Δt_{1} = t_{1}^{post} − t_{1}^{pre} < 0). Formally, T is defined by T = (t_{2}^{pre} + t_{2}^{post})/2 − (t_{1}^{pre} + t_{1}^{post})/2. Throughout this paper, we took Δt = −Δt_{1} = Δt_{2} = 5 ms.
Poisson protocol.
The presynaptic and postsynaptic spike trains are Poisson spike trains with firing rate ρ_{x} and ρ_{y}, respectively. The interest of such a protocol is that it is possible to establish a link with the BCM learning rule (Bienenstock et al., 1982), which has attractive theoretical properties. Indeed, this learning rule was originally used to explain the emergence of orientation selectivity in the visual cortex. Even if this protocol has less experimental support than the other protocols, some aspects of it have been indirectly measured in the visual cortex (Kirkwood et al., 1996) and in hippocampal slice (Artola et al., 1990; Dudek and Bear, 1992).
Data fitting.
To fit the amplitude parameters A_{2}^{+}, A_{2}^{−}, A_{3}^{+}, and A_{3}^{−} and the time constants τ_{x} and τ_{y} (τ_{+} = 16.8 ms and τ_{−} = 33.7 ms are kept fixed), we calculated the total weight change Δw_{i}^{mod} for a given pairing or triplet protocol and compared it to the experimental value Δw_{i}^{exp}. For the optimization of the parameters, we performed a minimization of the normalized meansquare error E defined by the following;
where Δw_{i}^{exp} and σ_{i} are the experimental mean weight change and SEM weight change for a given data point i. P is the number of data points within a data set; p = 10 for the visual cortex data set (Table 1), and p = 13 for the hippocampal culture data set (Table 2). Δw_{i}^{mod} is the weight change for a given model (pair or triplet model).
Numerical procedures.
The weight change Δw^{mod} for a given model and a given protocol can be either simulated numerically with Equations 1–4 or calculated analytically. See supplemental material (available at www.jneurosci.org) for an example of analytical calculation of the weight change of the triplet model applied to the pairing protocol with NearestSpike interactions.
In the present study, the weight changes predicted by all different models (pairbased models, minimal and full tripletbased models with both NearestSpike and AlltoAll interactions) have been calculated analytically and then evaluated numerically with MatLab (MathWorks, Natick, MA) on a Sun machine. The normalized meansquare error E of Equation 5 has been minimized with the MatLab builtin function Isqnonlin, which uses a reflective Newton method.
Results
Standard pairbased STDP models fail to reproduce frequency effects
In a first series of experiments, we applied a classical pairbased STDP learning rule (compare Eqs. 3 and 4 with A_{3}^{−} = 0 and A_{3}^{+} = 0) to the pairing protocol with 60 pairs of presynaptic and postsynaptic spikes (see Materials and Methods). Obviously, the weight change predicted by the model depends on the precise choice of the parameters A_{2}^{+}, A_{2}^{−}, τ_{+}, and τ_{−}. We therefore set those parameters in such a way that the normalized mean square error E across all experimental protocols is minimal (see Eq. 5). We found that even with the best set of parameters, the classical STDP model fails, for both the AlltoAll interaction and the NearestSpike interaction, to reproduce the experimental data (Fig. 2A). This is attributable to the following reasons.
First, as pointed out by Sjöström et al. (2001), a surprising aspect of their finding is that at low repetition frequency, ρ, there is no potentiation. This cannot be captured by standard pairbased STDP models, because for any choice of the parameter A_{2}^{+} > 0, the pairbased model induces LTP if a presynaptic spike precedes a postsynaptic one by a few milliseconds.
Second, as we can see in Figure 2, for Δt > 0, potentiation increases when frequency increases. This behavior can also not be reproduced by classical STDP models. Indeed, in pairbased STDP models, as soon as the frequency increases, the prepost pairs approach each other and generate an interaction between the postsynaptic spike of one pair and the presynaptic spike of the next pair. The effect of these postpre pairs should increase with frequency and therefore depress the synapse, which is not what is seen in the experiments. Therefore, classical pairbased models fail to reproduce the pairing experiment of Sjöström et al. (2001).
It should be noted that the absence of potentiation at low frequency is in direct conflict with the results of Bi and Poo (1998), Zhang et al. (1998), and Froemke and Dan (2002), where there is a reasonable potentiation at low frequency. Because the preparation of Sjöström et al. (2001) is different from the one of Bi and Poo (1998) and Wang et al. (2005) and the results in conflict, it seems natural to use different parameters in our model for each data set.
Standard pairbased STDP models fail to reproduce triplet and quadruplet experiments
The following is a second set of evidence of the limits of pairbased STDP learning rules. In triplet experiments (Fig. 2C,D), there is a clear asymmetry between a prepostpre and a postprepost experiment. For example, 60 repetitions of a prepostpre triplet with relative timing (Δt_{1}, Δt_{2}) = (5 ms, −5 ms) yields no weight change, whereas the same number of repetitions of a postprepost triplet with (Δt_{1}, Δt_{2}) = (−5 ms, 5 ms) yields a weight change of ∼30%. However, any pairbased model would predict the same result for prepostpre and postprepost experiments, because the same pairs occur. Therefore, triplet results cannot be explained by a sum of a prepost potentiation term and a postpre depression term (Fig. 2C,D).
Finally, the asymmetry present in the quadruplets experiments (Fig. 2B) also causes some problems for pairbased STDP models. A quadruplet consists of a prepostpostpre sequence or a postpreprepost sequence, and T denotes the interval between the first and last pair of spikes within the quadruplet (see Materials and Methods for more details). In a pairbased model with AlltoAll interactions and for a given interval T between the pairs, the weight changes for postpreprepost and prepostpostpre are strictly identical because of the symmetry of the protocol and the symmetry of the AlltoAll interaction. The weight change predicted by a pair model can therefore not explain the asymmetry seen in the data. With NearestSpike interactions, the situation gets even worse: prepostpostpre quadruplets consist of two prepost pairs and one postpre term, whereas for the postpreprepost case, the opposite occurs: two postpre pairs and only one prepost pair. Therefore, the NearestSpike interaction scheme leads to an asymmetry that is opposite to the one found in experiments (Fig. 2B).
Triplet rule
So far, we have shown that standard pairbased STDP models fail to reproduce frequency effects of the pairing protocol as well as triplet and quadruplet experiments. This is mainly because of the fact that pairbased models are intrinsically symmetric, in the sense that they predict the same weight change for a prepost pair followed by a postpre pair with the same delay Δt as for the inverted order [i.e., a postpre pair followed by a prepost pair (with the same delay Δt)]. However, there is no a priori reason to think that a prepostpre and a postprepost triplet should give the same result because they will activate different presynaptic and postsynaptic pathways.
We therefore included extra terms in the learning rule to break the symmetry induced by pairbased models. Specifically, we added a triplet depression term (i.e., a 2pre1post term) as well as a triplet potentiation term (i.e., a 1pre2post term) (see Materials and Methods for more details). We call this model a full triplet model, because it includes both pair terms and triplet terms. The full triplet model is described by eight parameters: four amplitude parameters A_{2}^{−}, A_{2}^{+}, A_{3}^{−}, and A_{3}^{+} and four time constants τ_{+}, τ_{−}, τ_{x}, and τ_{y}. Note that pairbased models are described by four parameters (A_{2}^{−}, A_{2}^{+}, τ_{+}, and τ_{−}).
In analogy to our approach in the previous subsection, we applied our triplet model to the protocols described in Materials and Methods. More precisely, we calculated analytically for each protocol the weight change predicted by our triplet learning rule (see supplemental material, available at www.jneurosci.org, for an example of explicit expression of the weight change). As before, we want our triplet learning rule to fit as best as possible to the experimental data of Sjöström et al. (2001) or Wang et al. (2005). We therefore minimized the normalized mean square error across all data points of a given data set (Table 1 or 2) by adjusting the eight parameters mentioned above. The resulting parameters are summarized in Tables 3 and 4.
As a first test for the triplet learning rule, we checked whether it can reproduce the biphasic learning window observed by Bi and Poo (1998). Our triplet learning rule succeeds to reproduce the classical STDP learning window (Fig. 3), because the triplet terms specific to our model play a minor role at a fixed low frequency.
Triplet learning rules can reproduce frequency effects
In this section, we study the pairing protocol used by Sjöström et al. (2001) in visual cortex (i.e., we apply 60 pairs of presynaptic and postsynaptic spikes at a given frequency ρ). As shown in Figure 4A, our full triplet learning rule succeeds to reproduce frequency effects of the pairing protocol. Indeed, the two main problems of the pairbased STDP models explained in section three for the pairing protocol are solved by the triplet model for the following reasons. First, the absence of potentiation at low frequency is achieved by setting A_{2}^{+} to a low value; second, the increase of potentiation with frequency is implemented via the triplet potentiation term controlled by A_{3}^{+}, which has a stronger effect than the triplet depression term A_{3}^{−} (Table 3). Thus, our model can explain results at different frequencies without an explicit “potentiation wins” mechanism suggested previously (Sjöström et al., 2001; Wang et al., 2005).
Because some of the optimized parameters of the triplet learning rule have values close to zero, we concluded that the terms controlled by these parameters can be neglected. This allowed us to define a minimal triplet model with less parameters. The first parameter we can easily drop is the amplitude A_{2}^{+} of the pair potentiation term, because it is extremely small in both the AlltoAll and NearestSpike interaction scheme (Table 3). The second parameter we neglect is A_{3}^{−}. This is possible for the following reason. In the AlltoAll interaction scheme, we have A_{3}^{−} ≪ A_{2}^{−}. Therefore, the effect of the triplet depression term is negligible compared with the depression induced by spike pairs.
Results with the minimal triplet model show good agreement with experimental data (Fig. 5A). Hence, the minimal model with five parameters can explain the visual cortex data that the classical pairbased STDP model with four parameters fails to explain.
Triplet learning rules can reproduce triplet and quadruplet experiments
By following the same procedure as the one described in the previous paragraph, we applied our full triplet learning rule to the second set of data (i.e., the hippocampal culture data set) (Bi and Poo, 1998, 2001; Wang et al., 2005). The parameters resulting from the minimization of the normalized mean square error across the pair, triplet, and quadruplet data are summarized in Table 3.
Our triplet learning rule does not only reproduce the classical STDP learning window (Fig. 3), but it also captures the results of most of the triplet and the quadruplet experiments. See Figure 4B–D. For example, the asymmetry between the prepostpre [(Δt_{1}, Δt_{2}) = (5 ms, −5 ms)] and the postprepost [(Δt_{1}, Δt_{2}) = (−5 ms, 5 ms)] triplets can be well captured by our model. For those two specific triplet protocols, the predicted weight change of the full triplet learning rule with AlltoAll interactions is within 1.1 σ (SEM) off the experimental mean weight change, whereas the pairbased learning predictions are off by >4 σ. We should, however, note that even if our triplet learning rule captures most of the triplet experiments, the fit is not perfect. For example, the prepostpre [with (Δt_{1}, Δt_{2}) = (5 ms, −15 ms)] triplet experiment is not well reproduced by our triplet learning rule (Figs. 4C, 5C).
With arguments similar to those applied above to model the visual cortex data set, it is possible to reduce the complexity of the model of the hippocampal data set. Specifically, we have set A_{3}^{−} = 0 as done previously. However, in contrast to the above minimal model for visual cortex data, the pair term controlled by A_{2}^{+} is kept as part of the model, because it is necessary to explain the potentiation at 1 Hz repetition frequency. The resulting weight change of the minimal model applied to the triplet and quadruplet experiments is depicted in Figure 5B–D. We emphasize that the minimal model for the hippocampal data is different from the one used for the visual cortex data.
To compare the pair models and the minimal and full triplet models, we plotted the fitting error given by Equation 5 as a function of the number of parameters in the model (Fig. 6). The best types of model are those that can predict the experimental data as well as possible while being as simple as possible (i.e., having as few parameters as possible). In this sense, the minimal models are the best, because they perform almost as well as the full triplet models while having only one extra parameter compared with standard pairbased models (two extra parameters for the hippocampal culture data set).
Finally, for future tests of the triplet models, we propose two new protocols that have not yet been used experimentally. The first protocol consists of prepostpre triplets with relative timing (Δt_{1}, Δt_{2}) = (5 ms, −5 ms), and the second protocol consists of postprepost triplets with relative timing (Δt_{1}, Δt_{2}) = (−5 ms, 5 ms). Triplets are repeated 60 times at different frequencies ρ. Figure 6, C and D, depicts the weight change predicted by the minimal triplet models (with AlltoAll and NearestSpike interactions) for the two triplet protocols. The models predict a frequency dependence with a positive slope. However, the overall level of potentiation predicted by the AlltoAll model is clearly different from that of the NearestSpike interaction model. Thus, the above experimental protocol would allow to test the triplet models and distinguish between its two variants.
Triplet learning rule can be mapped to the BCM learning rule
Functional consequences of our new triplet model can be studied in two different ways (i.e., analytically or by numerical simulations). We used a combination of the two and proceeded as follows. First, we show analytically a close analogy (“mapping”) between our triplet model and the traditional BCM theory. As a result of this mapping, we may conclude that, under random spike arrival with rate ρ_{x}, our triplet model behaves as a BCM model and inherits all of its functional properties. In particular, we expect that our triplet model exhibits synaptic competition leading to input selectivity as required for receptive field development. In a second step, we tested this prediction of input selectivity by numerical simulation.
First, we show that unlike standard pairbased STDP learning rules, our triplet learning rule can be mapped to the BCM learning rule. If we assume that the presynaptic and postsynaptic spike trains have Poisson statistics with ρ_{x} and ρ_{y}, respectively, as firing rate, the expected weight change can be calculated analytically. Intuitively, we may expect that a triplet term with one presynaptic and two postsynaptic spikes leads to a weight change that is proportional to the postsynaptic rate and the square of the presynaptic rate. An analogous argument holds for the other terms. Indeed, a detailed calculation for the AlltoAll triplet learning rule based on Equations 1–4 yields an expected weight change as follows:
Figure 7A depicts the expected weight change of Equation 6 as a function of the postsynaptic frequency ρ_{y}. The above weight dynamics can be written as a BCM learning rule. Indeed, the BCM theory requires first that the weight change can be written as dw/dt = φ (ρ_{y}, θ)ρ_{x}, where φ is such that φ (ρ_{y} < θ, θ) < 0, φ (ρ_{y} > θ, θ) > 0, and φ (0, θ) = 0. Our Equation 6 can satisfy this condition if A_{3}^{−} = 0, as is the case for our minimal triplet models. The second requirement is that the threshold θ between potentiation and depression is proportional to the expectation of the p^{th} power of the postsynaptic firing rate, i.e., θ = α〈ρ_{y}^{p}〉, where p > 1 (Bienenstock et al., 1982; Intrator and Cooper, 1992). This second requirement can be fulfilled if the parameters A_{2}^{−} and A_{2}^{+} depend on the mean firing rate 〈ρ_{y}〉 (or powers thereof) of the postsynaptic neuron. Specifically, we set A_{2}^{−} → A_{2}^{−}〈ρ_{y}^{p}〉/ρ_{0}^{p} as well as A_{2}^{+} → A_{2}^{+}〈ρ_{y}^{p}〉/ρ _{0}^{p}. By doing so, the threshold becomes θ = 〈ρ_{y}^{p}〉(A_{2}^{−}τ_{−} A_{2}^{+}τ_{+})/(ρ_{0}^{p}A_{3}^{+}τ _{+}τ _{y}).
Strictly speaking, 〈ρ_{y}^{p}〉 corresponds to the expectation over the input statistics of the p^{th} power of the postsynaptic firing rate. Practically, this quantity can be evaluated online by lowpass filtering ρ_{y}^{p} with a time constant of the order of 10 min or more. With this range of time scale, 〈ρ_{y}^{p}〉 can be considered as constant (i.e., 〈ρ_{y}^{p}〉 ≃ ρ_{0}^{p}) over the duration of the pairing, triplet, and quadruplet protocols we used in this study. As an aside, we note that with NearestSpike interactions, our triplet learning rule can almost (but not strictly) be mapped to a BCM learning rule.
Because the triplet rule shares properties with BCM theory, we expect that it generates input selectivity if a neuron receives a large number of inputs. Development of input selectivity is thought be an important property to account for receptive fields development.
For a numerical illustration of the input selectivity property of the triplet learning rule, we simulated the following scenario. We assume that our model neuron receives 100 afferents (1 ≤ i ≤ 100), which are stimulated with Gaussian profiles ν_{i} = 1 Hz + 50 Hz exp[− (i − μ)^{2}/(2 × 10^{2})], i = 1, …, 100, the center μ of which is shifted randomly every 200 ms (Fig. 7C) over 10 possible positions. Presynaptic spikes are generated at time t_{i}^{f} with a rate ν_{i}. Each presynaptic spike generates an exponential postsynaptic potential with decay time constant τ = 10 ms, so that the total potential is u = Σ_{i} Σ_{tif<t} w_{i}ε_{0} exp[−(t − t_{i}^{f})/τ], where ε_{0} = 1 mV. The postsynaptic firing rate increases with the membrane potential according to ν^{post} = 1 Hz + g × u, where g = 10Hz/mV. The neuron is stimulated over 60 s, whereas synapses change according to our triplet learning rule. As we can see in Figure 7D, the neuron becomes automatically specialized to one of the 10 input patterns (i.e., the one with μ = 70). In other words, learning leads to input selectivity, a necessary property for receptive field development.
It is interesting to note that the dynamics of Equation 6 can be seen as a gradient ascent of an objective function L (i.e., Δw α ∂L/∂w). Let p = 2 and β = A_{3}^{+}τ _{+}τ _{y}. L can be written as L = (β/3)ρ_{y}^{3} − (β/α)θ^{2}. If the model neuron has only two input afferents and, hence, only two synapses, this objective function L (or energy landscape) (Fig. 7B) elicits two selective points, which correspond to the two maxima of the function L. The first maximum is at w_{1} = 1 and w_{2} = 0, and the second is at w_{1} = 0 and w_{2} = 1. Therefore, the pattern of synaptic weight corresponds to input selectivity (i.e., the neuron is sensitive to only one of two inputs). Thus, the objective function can be used for a mathematical demonstration of input selectivity. Note that the objective function exists only if we assume that φ is a function of 〈ρ_{y}^{p}〉 and not 〈ρ_{y}〉^{p} (See also Cooper et al., 2004).
Discussion
In this paper, we first showed the limitations of the standard pairbased STDP models in terms of predicting the outcome of several spike timingbased protocols. We then showed that a triplet learning rule is more suitable to reproduce those experimental protocols, namely, the frequency dependence of the pairing protocol as well as the triplet and quadruplet protocols. Finally, we showed the link between our triplet learning rule and the BCM learning rule. We found noteworthy and somewhat unexpected that our detailed modeling of frequency dependence of pairbased protocols and asymmetries in triplet protocols should lead under the assumption of Poisson spike trains to a known theoretical rule with well characterized features.
Throughout this paper, we compared the AlltoAll interactions versus the NearestSpike interactions for pairbased models as well as for tripletbased models. Although NearestSpike interactions induce some potentially interesting nonlinearities in pairbased models (van Rossum et al., 2000; Izhikevich and Desai, 2003; Burkitt et al., 2004) (especially in the Poisson protocol), it is not possible to make a strict mapping of NearestSpike interactions models to the BCM rule, and more importantly pairbased models with NearestSpike interactions fail to reproduce the correct frequency dependence in the pairing protocol as well as triplet and quadruplet experiments.
Limitations
Even if our triplet model can capture most of the triplet and quadruplet experiments, it is necessary to keep in mind the kind of experiments this model cannot reproduce. Because our model predicts weight changes as a function of the spike timing only, it fails to make any kind of inference for experiments that trigger explicitly other biophysical parameters such as Ca^{2+} concentration or postsynaptic membrane potential. We nevertheless think that this approach is interesting, because in one way or another, those biophysical parameters depend on the timing of the presynaptic and postsynaptic spikes. For example, the calcium concentration depends on the timing of the postsynaptic spike (via the backpropagating action potential) and the presynaptic spike (via voltagegated calcium channels and NMDA channels).
Recent experiments (Froemke et al., 2005) show that the shape of the depression part of the learning window depends on the position of the synapse on the dendritic tree. Although our model does not include such geometrical properties, it is possible to account for the position of the synapse by changing explicitly the time constant τ_{−} (characterizing the LTD part of the learning window) as a function of the distance between the synapse and the soma.
Even in the context of typical STDP experiments, some aspects are not covered by our model. In most STDP experiments, plasticity is induced after a repetition of a fixed number of pairs of presynaptic and postsynaptic spikes. Clearly, the amount of plasticity depends on the number of pairs. In fact, the amount of potentiation increases with the number of pairs of presynaptic and postsynaptic spikes and saturates at a given value (Senn et al., 2001; Froemke et al., 2006). This saturation is not taken into account in our present model, because the weight dependence is not explicitly mentioned. The dependence on the weights can easily be added in the triplet models (the parameters A_{2}^{+}, A_{2}^{−}, A_{3}^{+}, and A_{3}^{−} could also depend on w). Even if we have some indications (Bi and Poo, 2001; Wang et al., 2005) of how synapses change as a function of w, more experimental data are clearly needed to determine the correct weight dependence. It should be noted that if we add the dependence on the weight, there would not be an unambiguous mapping to the BCM theory.
Alternative interpretations of the experimental data
The goal of this study was to go as far as possible in the prediction of the weight change with only spike timing and no other neuronal variables or mechanism. It is interesting to note that our triplet learning rule can reproduce both the that have been explained by a postsynaptic potential effect (Sjöström et al., 2001) or by a suppression effect (Wang et al., 2005). In Sjöström's experiment, the increased potentiation at high frequency is explained by the increased membrane potential because of the accumulation of presynaptic inputs, whereas in our model the increased potentiation is attributable to the increase of the postsynaptic variable o_{2}. Combined with a suitable neuron model, an increased frequency would of course yield a higher potential on average. Wang et al. (2005) interpreted their triplet and quadruplet experiment as a result of a suppression mechanism (i.e., if a prepost pair is followed by a postpre pair, the latter depression term suppresses the first potentiation term, and not the other way round). This phenomenon is captured in our framework by the extra potentiation attributable to the 1pre2post triplet term.
Expansion perspective
It is possible to see pair terms and triplet terms of Equations 3 and 4 in a more general framework (Gerstner and Kistler, 2002). From the point of view of pure spike timingdependent plasticity, the instantaneous weight change is an unknown functional of the presynaptic spike train X and the postsynaptic spike train Y (i.e., dw/dt = H[X, Y]). In this framework, a pairbased STDP learning rule corresponds to the Volterra expansion of H to the second order. In this paper, we pushed the expansion to the third order and showed that the prediction power increased a lot with only one or two extra parameters.
It is interesting to note that quadruplets data can be fitted with a triplet rule. This suggests that thirdorder terms (triplet terms) are good enough, and therefore there is no need to take into account higherorder terms. It is of course possible that new experiments will show the limitations of a triplet model and force us to consider higherorder terms. Clearly, the relevance of such an approach depends on how far we have to go in the expansion.
The learning rule as it is now does not depend directly on the membrane potential and therefore cannot reproduce the experiments of Sjöström et al. (2001) in which the membrane potential is controlled by the experimentalist. However, we should note that in the framework of the Volterra expansion of the unknown functional H, it is possible to assume that H depends explicitly on the membrane potential (Gerstner and Kistler, 2002) and therefore captures the voltage dependence experiments of Sjöström et al. (2001).
Comparison to other models
It is interesting to note that our triplet learning rule has some similarities with the Senn–Markram–Tsodyks (SMT) model (Senn et al., 2001). In their model, they need a first presynaptic spike to activate a fraction of NMDA channels and then a postsynaptic spike to set some secondary messengers in an upstate and finally a last postsynaptic spike to trigger synaptic potentiation. Thus, their rule essentially consists of a triplet prepostpost term for potentiation and a triplet postprepre term for depression. Even if their model makes implicit use of triplet terms, it should be noted that the order of spikes is different compared with our present triplet model. In our present triplet model, we need one presynaptic and two postsynaptic spikes regardless of the order (i.e., it encompasses prepostpost as well as postprepost triplets of spikes). This difference is of particular importance if we want to fit the postprepost triplet experiments or quadruplet experiments of Wang et al. (2005). For those protocols, the SMT model cannot reproduce the data.
Another model that takes into account a multispike interaction is the Froemke–Dan learning rule (Froemke and Dan, 2002). Their model (which is in fact a quadruplet model) predicts a synaptic behavior that is in direct contrast with the synaptic dynamics given by Equations 3 and 4. In their model, if a postsynaptic spike precedes a prepost pair of spikes, the effective potentiation will decrease as soon as the two postsynaptic spikes get closer to each other, whereas in our model, the opposite occurs. This is the main reason why the Froemke–Dan model, under a Poisson assumption for the presynaptic and postsynaptic spike trains, predicts an increasing depression as the postsynaptic firing rate increases, as reported by Izhikevich and Desai (2003), which seems unplausible in view of the results in Figure 8C of Sjöström et al. (2001). It should be noted that the revised suppression model of Froemke et al. (2006) gets around this problem by setting two different saturation values: one for depression and one for potentiation.
Until now, we have not made any assumption about the cellular processes that are described by our triplet learning rule. It is known that the amount of potentiation or depression expressed by a synapse depends critically on the concentration of calcium. Moreover, we know that there is a supralinear summation of calcium on the postsynaptic site when an EPSP precedes an action potential (Waters et al., 2003). From this point of view, two closely spaced postsynaptic spikes can increase the level of calcium and therefore increase potentiation. This would correspond to the 1pre2post triplet term of our formalism.
In this study, we showed that a minimal triplet model can capture most, but not all, aspects of the pairing, triplet, and quadruplet experiments. A natural extension of this study could be to include explicitly the dependence on biophysical quantities such as the Ca^{2+} concentration, the postsynaptic membrane potential, or other neuronal quantities. A more appealing approach would be to consider existing detailed biophysical models of synaptic plasticity and try to reduce them to a triplet model and therefore identify the underlying biological quantities of our triplet model.
Footnotes

This work was supported by the Swiss National Science Foundation 200020108097. We gratefully acknowledge discussions with Taro Toyoizumi and Magnus Richardson.
 Correspondence should be addressed to Prof. Wulfram Gerstner, Ecole Polytechnique Fédérale de Lausanne, Laboratory of Computational Neuroscience, Station 15, 1015 Lausanne, Switzerland. wulfram.gerstner{at}epfl.ch