Coherence detection in a spiking neuron via Hebbian learning ...
Coherence detection in a spiking neuron via Hebbian learning Laurent Perrinet and Manuel Samuelides ONERA-DTIM, 2 Av. E. Belin, BP 4025, F-31055 Toulouse CEDEX, France email: {perrinet, samuelid}@cert.fr Abstract It is generally assumed that neurons in the central nervous system communicate through temporal £ring patterns. As a £rst step, we will study the learning of a layer of realistic neurons in the particular case where the relevant messages are formed by temporally correlated patterns, or syn£re patterns. The model is a layer of Integrate-and-Fire (IF) neurons with synaptic current dynamics that adapts by minimizing a cost ac-
1.1 Coding scheme We will represent (as in [Gerstner99]) the signal Si at synapse i by the sum of Dirac pulses located at the spiking times tki drawn from the lists of spikes Γi (see Figure 1-left). X Si = δ(t − tki ) (1)
1 t − ti ) ci(t) = exp(− τg τg 1 t − ti t − ti pi(t) = (exp(− ) − exp(− )) τV − τ g τg τv t − tj η(t) = (Vreset − Vthreshold)exp(− ) τV
k∈Γi
Syn£re patterns are generated in analogy with the response of a retina to ¤ashed binary images. The input of the synapses is characterized as the output of single-synapse IF neurons responding to a speci£c binary input.This response may be described as the sum of two random point processes with different time scales. At a narrow time scale, the input is the spontaneous activity, i.e. a background noise independent of time and synapses that may be described by a Poisson point process of rate 1/τnoise. At a larger time scale, the syn£re pattern activates a given subset M of synapses once per ¤ash with a correlation de£ned by its jitter τjitter . (see £gure 1-left)
45
−10
40
potential (mV)
synapse #
35
25 20 15
−40
−70 5 80
100
−80 time (ms)
20
40 60 time (ms)
80
Figure 1:
Neuron Model : (left) Input spikes (with a syn£re pattern at t = 25ms), are (middle) modulated in time and amplitude forming postsynaptic current pulses and are £nally (right) integrated at the soma. When the potential (blue line) reaches the threshold (dotted green line), a spike is emitted and the potential is decreased. A sample PSP is shown in red.
2
The learning mechanism
A possible cost function may therefore be the squared distance to the potentials of neurons at the £ring time t n added to the total sum of the squared weights: 2.E = λn(
We will consider N1 synapses (indexed by i) connected to a layer of N2 neurons j. Those are generalized version of IF neurons with synaptic current dynamics, a one compartment model with no delay and the synapses have contacts characterized by their weight wij . The state variables are the N1.N2 synaptic driving conductances gij and the N2 membrane potentials Vj . Incoming spikes trigger those conductances by opening the driving gates with time constant τg : dgij τg = −gij + wij .Si dt
1.3 Reduced equations We introduce reduced equations for this IF Layer to study its dynamical behavior and simplify its implementation. In fact, this reduction follows the concept of the Spike-Response Model (SRM) which was extensively studied in [Gerstner99, Kempter99]. It is similar to [Perrinet01] which aimed at reducing the STDHP equations to a set of £rst order equations. ½ ci τg dt = −ci + Si (4) pi τV dt = −pi + ci P then, Vrest + ( 1≤i≤N1 wij .pi) veri£es the equation system (2, 3), with gij = wij .ci. To account for the threshold mechanism at a time tkj , we may then add a resetting value to Vj by setting ηj (tkj ) = Vreset − Vthreshold and then: (5)
So that £nally, an equivalent version of the IF Layer consists of (4, 5) and : X Vj (t) = Vrest + ( wij .pi(t)) + ηj (t) (6) 1≤i≤N 1
3.1 Response to syn£re patterns To achieve this experiment we presented syn£re patterns to the layer. The weights were set at random so that the network could £re to all the inputs. The patterns were presented at random times that were suf£ciently distant. This unsupervised learning converges quickly, and as may be observed in neuromuscular connectivity, the synapses tend to sparsify and the neurons tend to respond to only one input (see Figure 2). 100 90 0
target 2 ) )+(
(Vj (tn) − Vj
wij 2) (10)
1≤i≤N1
(11) (12)
−20
60 50 40 30
−40 −60 −80 6
20
2000
4
10
1000
2 200
400 600 time (ms)
800
1000
neuron #
0 0
time (ms)
Figure 2:
Coherence detection: (left) different input patterns (t = 100ms, 300ms, 500ms, 700ms, 900ms) are (right) learnt by the system : only one neuron per input £res (100 learning steps)
3.2 Response to oriented bars The next experiment consisted in applying those results to a basic retina which is fed with centered rotated lines. A £xed analogical contrast layer (ON and OFF) sends then spikes to the learning layer that learns with the rule we presented. We observe unsupervised emergence of V1like receptors £elds sensitive to the orientation (see Figure 3). Further experiments with lateral interactions and accounting for dendritic delay show even more realistic £lters.
Where λn is a scaling parameter and ∆V ≈ 5 mV .
Figure 3:
Oriented bars detection: after learning, the weights show sensitivity to orientation (black:OFF; white:ON; gray:neutral)
2.2 Gradient descent It follows from equations 10and 6 :
Conclusion ∂E target ∂Vj = λn.(Vj (tn) − Vj ) + wij ∂wij ∂wij target = λn.(Vj (tn) − Vj ).pi + wij
We may therefore formulate the gradient descent algorithm in our model as :
(3)
When Vj reaches the threshold’s potential Vthreshold ≈ −54 mV from below, the target neuron £res (see Figure 1-right) and is shunted (V j is set e.g. to Vreset ≈ −75 mV ).
ηj (t) τV = −ηj (t) dt
1≤j≤N X 2
target Vj = Vrest for j 6= wn target = Vthreshold + ∆V Vwn
(2)
and the potential Vj at the soma integrates with time constant τV the driving currents and the leaking current gleak (with a potential Vrest ≈ −70 mV ): (see Figure 1-middle): X dVj τV = gleak (Vrest − Vj ) + ( gij ) 1≤i≤N dt 2
X
1≤j≤N2
1.2 Integrate-and-Fire Layer
We implemented this model using discrete versions of the differential equations (forward Euler method) on a MATLAB system.
70
3. economy of the total synaptic ef£cacy should be respected.
100
Numerical results
80
2. to discrimante between the different input patterns, the output voltage should be close to a target value : the potential of the winning neuron (which we index j = wn) should be above threshold whereas other neurons should be hyperpolarized,
−50
10
40 60 time (ms)
Those equations (6, 7, 8 and 9) are the equivalent SRM version of our IF model. More precisely eq. 7 represents the PostSynaptic Current (PSC), see £g. 1-middle, and eq. 8 the PostSynaptic Potential (PSP), see £g. 1-right. .
1. the learning is associated with a response : the nth learning step occcurs at the nth output £ring time t n,
−30
−60
20
(9)
Based on neurophysiological studies, we set the following principles :
−20
30
(8)
2.1 De£nition of the cost function
0
50
(7)
3
potential (mV)
Description of the model
This formulation depends only on the present state and not on the past values. It is therefore biologically more plausible and computationally cheaper. Integrating these equations after emission of a presynaptic spike at t i or a postsynaptic spike at tj leads to :
synapse #
1
cording to a gradient descent scheme. This leads to a rule similar to Spike-Time Dependent Hebbian Plasticity (STDHP). The different units compete in an unsupervised manner Our results show that the rule that we derive is biologically plausible and leads to the detection of the coherence in the input. An application to shape recognition is shown as an illustration.
∂E n+1 n wij = wij − γn ∂wij n + γ λ .(V target − V (t )).p = (1 − γn)wij n n j j n i
Acknowledgements
with factor γn which satis£es n=1...∞ γn → ∞ and P a learning 2 < ∞. We may then set the parameter λ in accordance γ n n=1...∞ with the input. For instance, it may be set to zero except for the winning neuron which maximizes the chance to £re again P with the similar pattern pi(tn). This means that for j = wn, ( 1≤i≤N1 wij .pi) ≈ Vthreshold − Vrest + ∆V : X X X target n+1 n ( wij .pi) = (1 − γn)( wij .pi) + γn λ(Vj − Vj ).pi2 i X i i λ∆V. pi2 = Vthreshold − Vrest + ∆V P
i
£nally,
λ = Vthreshold−V P rest2+∆V ∆V
i pi
n + λ(V target − V ).p wn+1 = (1 − γn)wij j i ij j
We have presented an original gradient descent method to £nd a learning rule for a layer of spiking neurons. The simplicity of the rule gives a new insight into the comprehension of the mechanism behind the observed STDHP. Further work is done for the detection of asynchronous patterns. However, this study should be extended to more realistic spike trains (e.g. bursts), account for more complex behavior (e.g. facilitation and depression) and may be extended to population of neurons and recurrent systems.
(13)
2.3 Spike-time Dependent Plasticity
A closer look at equation 13 shows that the change of wij is proportiontarget nal to (Vj − Vj ).pi. This is a hebbian type of rule : when a neuron j £res after the £ring of synapse j, there is a mechanism that strengthen the connection. The strengthening depends therefore on the relative time of the pre- and post-synaptic spikes as is observed in biological systems [Bi98].
This work has been initiated during the EU Advanced Course in Computational Neuroscience. LP wish to thank its organizers, the teachers, the course-mates and my tutor, S. Panzeri.
Online simulations http://laurent.perrinet.free.fr/app/cns.html
References [Bi98] Bi G.-q., Poo M.-m. (1998). Synaptic modi£cations in cultured hippocampal neurons: dependence on spike timing, synaptic strength, and postsynaptic cell type, J. Neurosci., 18:10464–10472. [Gerstner99] Gerstner, W. (1999). Spiking neurons. In Pulsed Coupled Neural Networks, W. Maass, and C. M. Bishop, chapter 1, pages 3–54. MIT Press, Cambridge. [Kempter99] Kempter, R., Gerstner, W. and van Hemmen, J. (1999). Hebbian learning and spiking neurons. Physical Review E, 59(4):4498-4514. [Perrinet01] Perrinet L., Delorme A. and Thorpe S. (2001). Network of integrate-and-£re neurons using Rank Order Coding A: How to implement spike timing dependent plasticity. Neurocomputing, in press.
1.1 Coding scheme We will represent (as in [Gerstner99]) the signal Si at synapse i by the sum of Dirac pulses located at the spiking times tki drawn from the lists of spikes Γi (see Figure 1-left). X Si = δ(t − tki ) (1)
1 t − ti ) ci(t) = exp(− τg τg 1 t − ti t − ti pi(t) = (exp(− ) − exp(− )) τV − τ g τg τv t − tj η(t) = (Vreset − Vthreshold)exp(− ) τV
k∈Γi
Syn£re patterns are generated in analogy with the response of a retina to ¤ashed binary images. The input of the synapses is characterized as the output of single-synapse IF neurons responding to a speci£c binary input.This response may be described as the sum of two random point processes with different time scales. At a narrow time scale, the input is the spontaneous activity, i.e. a background noise independent of time and synapses that may be described by a Poisson point process of rate 1/τnoise. At a larger time scale, the syn£re pattern activates a given subset M of synapses once per ¤ash with a correlation de£ned by its jitter τjitter . (see £gure 1-left)
45
−10
40
potential (mV)
synapse #
35
25 20 15
−40
−70 5 80
100
−80 time (ms)
20
40 60 time (ms)
80
Figure 1:
Neuron Model : (left) Input spikes (with a syn£re pattern at t = 25ms), are (middle) modulated in time and amplitude forming postsynaptic current pulses and are £nally (right) integrated at the soma. When the potential (blue line) reaches the threshold (dotted green line), a spike is emitted and the potential is decreased. A sample PSP is shown in red.
2
The learning mechanism
A possible cost function may therefore be the squared distance to the potentials of neurons at the £ring time t n added to the total sum of the squared weights: 2.E = λn(
We will consider N1 synapses (indexed by i) connected to a layer of N2 neurons j. Those are generalized version of IF neurons with synaptic current dynamics, a one compartment model with no delay and the synapses have contacts characterized by their weight wij . The state variables are the N1.N2 synaptic driving conductances gij and the N2 membrane potentials Vj . Incoming spikes trigger those conductances by opening the driving gates with time constant τg : dgij τg = −gij + wij .Si dt
1.3 Reduced equations We introduce reduced equations for this IF Layer to study its dynamical behavior and simplify its implementation. In fact, this reduction follows the concept of the Spike-Response Model (SRM) which was extensively studied in [Gerstner99, Kempter99]. It is similar to [Perrinet01] which aimed at reducing the STDHP equations to a set of £rst order equations. ½ ci τg dt = −ci + Si (4) pi τV dt = −pi + ci P then, Vrest + ( 1≤i≤N1 wij .pi) veri£es the equation system (2, 3), with gij = wij .ci. To account for the threshold mechanism at a time tkj , we may then add a resetting value to Vj by setting ηj (tkj ) = Vreset − Vthreshold and then: (5)
So that £nally, an equivalent version of the IF Layer consists of (4, 5) and : X Vj (t) = Vrest + ( wij .pi(t)) + ηj (t) (6) 1≤i≤N 1
3.1 Response to syn£re patterns To achieve this experiment we presented syn£re patterns to the layer. The weights were set at random so that the network could £re to all the inputs. The patterns were presented at random times that were suf£ciently distant. This unsupervised learning converges quickly, and as may be observed in neuromuscular connectivity, the synapses tend to sparsify and the neurons tend to respond to only one input (see Figure 2). 100 90 0
target 2 ) )+(
(Vj (tn) − Vj
wij 2) (10)
1≤i≤N1
(11) (12)
−20
60 50 40 30
−40 −60 −80 6
20
2000
4
10
1000
2 200
400 600 time (ms)
800
1000
neuron #
0 0
time (ms)
Figure 2:
Coherence detection: (left) different input patterns (t = 100ms, 300ms, 500ms, 700ms, 900ms) are (right) learnt by the system : only one neuron per input £res (100 learning steps)
3.2 Response to oriented bars The next experiment consisted in applying those results to a basic retina which is fed with centered rotated lines. A £xed analogical contrast layer (ON and OFF) sends then spikes to the learning layer that learns with the rule we presented. We observe unsupervised emergence of V1like receptors £elds sensitive to the orientation (see Figure 3). Further experiments with lateral interactions and accounting for dendritic delay show even more realistic £lters.
Where λn is a scaling parameter and ∆V ≈ 5 mV .
Figure 3:
Oriented bars detection: after learning, the weights show sensitivity to orientation (black:OFF; white:ON; gray:neutral)
2.2 Gradient descent It follows from equations 10and 6 :
Conclusion ∂E target ∂Vj = λn.(Vj (tn) − Vj ) + wij ∂wij ∂wij target = λn.(Vj (tn) − Vj ).pi + wij
We may therefore formulate the gradient descent algorithm in our model as :
(3)
When Vj reaches the threshold’s potential Vthreshold ≈ −54 mV from below, the target neuron £res (see Figure 1-right) and is shunted (V j is set e.g. to Vreset ≈ −75 mV ).
ηj (t) τV = −ηj (t) dt
1≤j≤N X 2
target Vj = Vrest for j 6= wn target = Vthreshold + ∆V Vwn
(2)
and the potential Vj at the soma integrates with time constant τV the driving currents and the leaking current gleak (with a potential Vrest ≈ −70 mV ): (see Figure 1-middle): X dVj τV = gleak (Vrest − Vj ) + ( gij ) 1≤i≤N dt 2
X
1≤j≤N2
1.2 Integrate-and-Fire Layer
We implemented this model using discrete versions of the differential equations (forward Euler method) on a MATLAB system.
70
3. economy of the total synaptic ef£cacy should be respected.
100
Numerical results
80
2. to discrimante between the different input patterns, the output voltage should be close to a target value : the potential of the winning neuron (which we index j = wn) should be above threshold whereas other neurons should be hyperpolarized,
−50
10
40 60 time (ms)
Those equations (6, 7, 8 and 9) are the equivalent SRM version of our IF model. More precisely eq. 7 represents the PostSynaptic Current (PSC), see £g. 1-middle, and eq. 8 the PostSynaptic Potential (PSP), see £g. 1-right. .
1. the learning is associated with a response : the nth learning step occcurs at the nth output £ring time t n,
−30
−60
20
(9)
Based on neurophysiological studies, we set the following principles :
−20
30
(8)
2.1 De£nition of the cost function
0
50
(7)
3
potential (mV)
Description of the model
This formulation depends only on the present state and not on the past values. It is therefore biologically more plausible and computationally cheaper. Integrating these equations after emission of a presynaptic spike at t i or a postsynaptic spike at tj leads to :
synapse #
1
cording to a gradient descent scheme. This leads to a rule similar to Spike-Time Dependent Hebbian Plasticity (STDHP). The different units compete in an unsupervised manner Our results show that the rule that we derive is biologically plausible and leads to the detection of the coherence in the input. An application to shape recognition is shown as an illustration.
∂E n+1 n wij = wij − γn ∂wij n + γ λ .(V target − V (t )).p = (1 − γn)wij n n j j n i
Acknowledgements
with factor γn which satis£es n=1...∞ γn → ∞ and P a learning 2 < ∞. We may then set the parameter λ in accordance γ n n=1...∞ with the input. For instance, it may be set to zero except for the winning neuron which maximizes the chance to £re again P with the similar pattern pi(tn). This means that for j = wn, ( 1≤i≤N1 wij .pi) ≈ Vthreshold − Vrest + ∆V : X X X target n+1 n ( wij .pi) = (1 − γn)( wij .pi) + γn λ(Vj − Vj ).pi2 i X i i λ∆V. pi2 = Vthreshold − Vrest + ∆V P
i
£nally,
λ = Vthreshold−V P rest2+∆V ∆V
i pi
n + λ(V target − V ).p wn+1 = (1 − γn)wij j i ij j
We have presented an original gradient descent method to £nd a learning rule for a layer of spiking neurons. The simplicity of the rule gives a new insight into the comprehension of the mechanism behind the observed STDHP. Further work is done for the detection of asynchronous patterns. However, this study should be extended to more realistic spike trains (e.g. bursts), account for more complex behavior (e.g. facilitation and depression) and may be extended to population of neurons and recurrent systems.
(13)
2.3 Spike-time Dependent Plasticity
A closer look at equation 13 shows that the change of wij is proportiontarget nal to (Vj − Vj ).pi. This is a hebbian type of rule : when a neuron j £res after the £ring of synapse j, there is a mechanism that strengthen the connection. The strengthening depends therefore on the relative time of the pre- and post-synaptic spikes as is observed in biological systems [Bi98].
This work has been initiated during the EU Advanced Course in Computational Neuroscience. LP wish to thank its organizers, the teachers, the course-mates and my tutor, S. Panzeri.
Online simulations http://laurent.perrinet.free.fr/app/cns.html
References [Bi98] Bi G.-q., Poo M.-m. (1998). Synaptic modi£cations in cultured hippocampal neurons: dependence on spike timing, synaptic strength, and postsynaptic cell type, J. Neurosci., 18:10464–10472. [Gerstner99] Gerstner, W. (1999). Spiking neurons. In Pulsed Coupled Neural Networks, W. Maass, and C. M. Bishop, chapter 1, pages 3–54. MIT Press, Cambridge. [Kempter99] Kempter, R., Gerstner, W. and van Hemmen, J. (1999). Hebbian learning and spiking neurons. Physical Review E, 59(4):4498-4514. [Perrinet01] Perrinet L., Delorme A. and Thorpe S. (2001). Network of integrate-and-£re neurons using Rank Order Coding A: How to implement spike timing dependent plasticity. Neurocomputing, in press.
