1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
|
#!/usr/bin/env python3
"""
Fig. 3 of
Simple Model of Spiking Neurons
IEEE Transactions on Neural Networks ( Volume: 14, Issue: 6, Nov. 2003)
Eugene M. Izhikevich
based on net.m by Eugene M. Izhikevich (http://izhikevich.org/publications/spikes.htm)
Akif Erdem Sağtekin and Sebastian Schmitt, 2022
"""
import matplotlib.pyplot as plt
import numpy as np
from brian2 import NeuronGroup, Synapses, SpikeMonitor, StateMonitor
from brian2 import ms, mV
from brian2 import defaultclock, run
tfinal = 1000 * ms
Ne = 800
Ni = 200
re = np.random.uniform(size=Ne)
ri = np.random.uniform(size=Ni)
weights = np.hstack(
[
0.5 * np.random.uniform(size=(Ne + Ni, Ne)),
-np.random.uniform(size=(Ne + Ni, Ni)),
]
).T
defaultclock.dt = 1 * ms
eqs = """dv/dt = (0.04*v**2 + 5*v + 140 - u + I + I_noise )/ms : 1
du/dt = (a*(b*v - u))/ms : 1
I : 1
I_noise : 1
a : 1
b : 1
c : 1
d : 1
"""
N = NeuronGroup(Ne + Ni, eqs, threshold="v>=30", reset="v=c; u+=d", method="euler")
N.v = -65
N_exc = N[:Ne]
N_inh = N[Ne:]
spikemon = SpikeMonitor(N)
statemon = StateMonitor(N, 'v', record=0, when='after_thresholds')
N_exc.a = 0.02
N_exc.b = 0.2
N_exc.c = -65 + 15 * re**2
N_exc.d = 8 - 6 * re**2
N_inh.a = 0.02 + 0.08 * ri
N_inh.b = 0.25 - 0.05 * ri
N_inh.c = -65
N_inh.d = 2
N_exc.u = "b*v"
N_inh.u = "b*v"
S = Synapses(
N,
N,
"w : 1",
on_pre={"up": "I += w", "down": "I -= w"},
delay={"up": 0 * ms, "down": 1 * ms},
)
S.connect()
S.w[:] = weights.flatten()
N_exc.run_regularly("I_noise = 5*randn()", dt=1 * ms)
N_inh.run_regularly("I_noise = 2*randn()", dt=1 * ms)
run(tfinal)
fig, (ax, ax_voltage) = plt.subplots(2, 1, sharex=True,
gridspec_kw={'height_ratios': (3, 1)})
ax.scatter(spikemon.t / ms, spikemon.i[:], marker="_", color="k", s=10)
ax.set_xlim(0, tfinal / ms)
ax.set_ylim(0, len(N))
ax.set_ylabel("neuron number")
ax.set_yticks(np.arange(0, len(N), 100))
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
ax.axhline(Ne, color="k")
ax.text(500, 900, 'inhibitory', backgroundcolor='w', color='k', ha='center')
ax.text(500, 400, 'excitatory', backgroundcolor='w', color='k', ha='center')
ax_voltage.plot(statemon.t / ms, np.clip(statemon.v[0], -np.inf, 30),
color='k')
ax_voltage.text(25, 0, 'v₁(t)')
ax_voltage.set_xticks(np.arange(0, tfinal / ms, 100))
ax_voltage.spines['right'].set_visible(False)
ax_voltage.spines['top'].set_visible(False)
ax_voltage.set_xlabel("time, ms")
plt.show()
|
