Introduction to Neural Networks

Programming Lab 4: Spike-time dependent synaptic plasticity (STDP)


Competitive Hebbian learning through STDP

Today, we will implement the model introduced by Song, Miller, and Abbott in the reading for this week. We will be using the simulation code below, in addition to the introduction to synapses and STDP on the Brian website.

from brian2 import *

N = 1000
taum = 10*ms
taupre = 20*ms
taupost = taupre
Ee = 0*mV
vt = -54*mV
vr = -60*mV
El = -74*mV
taue = 5*ms
F = 15*Hz
wmax = .01
dApre = .01
dApost = -dApre * taupre / taupost * 1.05
dApost *= wmax
dApre *= wmax

eqs_neurons = '''
dv/dt = (ge * (Ee-vr) + El - v) / taum : volt
dge/dt = -ge / taue : 1
'''

input = PoissonGroup(N, rates=F)
neurons = NeuronGroup(1, eqs_neurons, threshold='v>vt', reset='v = vr',
                      method='linear')
S = Synapses(input, neurons,
             '''w : 1
                dApre/dt = -Apre / taupre : 1 (event-driven)
                dApost/dt = -Apost / taupost : 1 (event-driven)''',
             on_pre='''ge += w
                    Apre += dApre
                    w = clip(w + Apost, 0, wmax)''',
             on_post='''Apost += dApost
                     w = clip(w + Apre, 0, wmax)''',
             )
S.connect()
S.w = 'rand() * wmax'

mon = StateMonitor(S, 'w', record=True) 
# the resulting ndarray will now have the shape (N,timesteps),
# which will allow accessing the time evolution of individual 
# synapses (e.g. mon.w[0] will be the time evolution of the first
# input synapse)

s_mon = SpikeMonitor(input)

run(100*second, report='text')

subplot(311)
plot(S.w / wmax, '.k')
ylabel('Weight / wmax')
xlabel('Synapse index')
subplot(312)
hist(S.w / wmax, 20)
xlabel('Weight / wmax')
subplot(313)
plot(mon.t/second, mon.w[0]/wmax) # timecourse of synaptic weight for first input synapse
xlabel('Time (s)')
ylabel('Weight / wmax')

Experiments: (1) Add inhibitory synaptic inputs to the implementation above. (2) As introduced in the paper, start all synaptic weights at wmax and plot the coefficient of variation (\(C_V = \sigma_\tau / \langle\tau\rangle \)) as a function of time (e.g. divided into short windows of 1 second). (3) Replicate Figure 1a,b from the paper. (4) Replicate Figure 1c,d,e from the paper (in separate simulations).