forked from nest/nestml
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathtraub_cond_multisyn.nestml
187 lines (144 loc) · 7.03 KB
/
traub_cond_multisyn.nestml
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
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
"""
Name: traub_cond_multisyn - Traub model according to Borgers 2017.
Reduced Traub-Miles Model of a Pyramidal Neuron in Rat Hippocampus[1].
parameters got from reference [2] chapter 5.
Spike Detection
Spike detection is done by a combined threshold-and-local-maximum search: if
there is a local maximum above a certain threshold of the membrane potential,
it is considered a spike.
- AMPA, NMDA, GABA_A, and GABA_B conductance-based synapses with
beta-function (difference of two exponentials) time course corresponding
to "hill_tononi" model.
References:
[1] R. D. Traub and R. Miles, Neuronal Networks of the Hippocampus,Cam- bridge University Press, Cambridge, UK, 1991.
[2] Borgers, C., 2017. An introduction to modeling neuronal dynamics (Vol. 66). Cham: Springer.
SeeAlso: hh_cond_exp_traub
"""
neuron traub_cond_multisyn:
state:
r integer # number of steps in the current refractory phase
end
initial_values:
V_m mV = -70. mV # Membrane potential
Act_m real = alpha_m_init / ( alpha_m_init + beta_m_init ) # Activation variable m for Na
Inact_h real = alpha_h_init / ( alpha_h_init + beta_h_init ) # Inactivation variable h for Na
Act_n real = alpha_n_init / (alpha_n_init + beta_n_init) # Activation variable n for K
g_AMPA real = 0
g_NMDA real = 0
g_GABAA real = 0
g_GABAB real = 0
g_AMPA$ real = AMPAInitialValue
g_NMDA$ real = NMDAInitialValue
g_GABAA$ real = GABA_AInitialValue
g_GABAB$ real = GABA_BInitialValue
end
equations:
recordable inline I_syn_ampa pA = -convolve(g_AMPA, AMPA) * ( V_m - AMPA_E_rev )
recordable inline I_syn_nmda pA = -convolve(g_NMDA, NMDA) * ( V_m - NMDA_E_rev ) / ( 1 + exp( ( NMDA_Vact - V_m ) / NMDA_Sact ) )
recordable inline I_syn_gaba_a pA = -convolve(g_GABAA, GABA_A) * ( V_m - GABA_A_E_rev )
recordable inline I_syn_gaba_b pA = -convolve(g_GABAB, GABA_B) * ( V_m - GABA_B_E_rev )
recordable inline I_syn pA = I_syn_ampa + I_syn_nmda + I_syn_gaba_a + I_syn_gaba_b
inline I_Na pA = g_Na * Act_m * Act_m * Act_m * Inact_h * ( V_m - E_Na )
inline I_K pA = g_K * Act_n * Act_n * Act_n * Act_n * ( V_m - E_K )
inline I_L pA = g_L * ( V_m - E_L )
V_m' = ( -( I_Na + I_K + I_L ) + I_e + I_stim + I_syn ) / C_m
# Act_n
inline alpha_n real = 0.032 * (V_m / mV + 52.) / (1. - exp(-(V_m / mV + 52.) / 5.))
inline beta_n real = 0.5 * exp(-(V_m / mV + 57.) / 40.)
Act_n' = ( alpha_n * ( 1 - Act_n ) - beta_n * Act_n ) / ms # n-variable
# Act_m
inline alpha_m real = 0.32 * (V_m / mV + 54.) / (1.0 - exp(-(V_m / mV + 54.) / 4.))
inline beta_m real = 0.28 * (V_m / mV + 27.) / (exp((V_m / mV + 27.) / 5.) - 1.)
Act_m' = ( alpha_m * ( 1 - Act_m ) - beta_m * Act_m ) / ms # m-variable
# Inact_h'
inline alpha_h real = 0.128 * exp(-(V_m / mV + 50.0) / 18.0)
inline beta_h real = 4.0 / (1.0 + exp(-(V_m / mV + 27.) / 5.))
Inact_h' = ( alpha_h * ( 1 - Inact_h ) - beta_h * Inact_h ) / ms # h-variable
#############
# Synapses
#############
kernel g_AMPA' = g_AMPA$ - g_AMPA / tau_AMPA_2,
g_AMPA$' = -g_AMPA$ / tau_AMPA_1
kernel g_NMDA' = g_NMDA$ - g_NMDA / tau_NMDA_2,
g_NMDA$' = -g_NMDA$ / tau_NMDA_1
kernel g_GABAA' = g_GABAA$ - g_GABAA / tau_GABAA_2,
g_GABAA$' = -g_GABAA$ / tau_GABAA_1
kernel g_GABAB' = g_GABAB$ - g_GABAB / tau_GABAB_2,
g_GABAB$' = -g_GABAB$ / tau_GABAB_1
end
parameters:
t_ref ms = 2.0 ms # Refractory period 2.0
g_Na nS = 10000.0 nS # Sodium peak conductance
g_K nS = 8000.0 nS # Potassium peak conductance
g_L nS = 10 nS # Leak conductance
C_m pF = 100.0 pF # Membrane Capacitance
E_Na mV = 50. mV # Sodium reversal potential
E_K mV = -100. mV # Potassium reversal potentia
E_L mV = -67. mV # Leak reversal Potential (aka resting potential)
V_Tr mV = -20. mV # Spike Threshold
# Parameters for synapse of type AMPA, GABA_A, GABA_B and NMDA
AMPA_g_peak nS = 0.1 nS # peak conductance
AMPA_E_rev mV = 0.0 mV # reversal potential
tau_AMPA_1 ms = 0.5 ms # rise time
tau_AMPA_2 ms = 2.4 ms # decay time, Tau_1 < Tau_2
NMDA_g_peak nS = 0.075 nS # peak conductance
tau_NMDA_1 ms = 4.0 ms # rise time
tau_NMDA_2 ms = 40.0 ms # decay time, Tau_1 < Tau_2
NMDA_E_rev mV = 0.0 mV # reversal potential
NMDA_Vact mV = -58.0 mV # inactive for V << Vact, inflection of sigmoid
NMDA_Sact mV = 2.5 mV # scale of inactivation
GABA_A_g_peak nS = 0.33 nS # peak conductance
tau_GABAA_1 ms = 1.0 ms # rise time
tau_GABAA_2 ms = 7.0 ms # decay time, Tau_1 < Tau_2
GABA_A_E_rev mV = -70.0 mV # reversal potential
GABA_B_g_peak nS = 0.0132 nS # peak conductance
tau_GABAB_1 ms = 60.0 ms # rise time
tau_GABAB_2 ms = 200.0 ms # decay time, Tau_1 < Tau_2
GABA_B_E_rev mV = -90.0 mV # reversal potential for intrinsic current
# constant external input current
I_e pA = 0 pA
end
internals:
AMPAInitialValue real = compute_synapse_constant( tau_AMPA_1, tau_AMPA_2, AMPA_g_peak )
NMDAInitialValue real = compute_synapse_constant( tau_NMDA_1, tau_NMDA_2, NMDA_g_peak )
GABA_AInitialValue real = compute_synapse_constant( tau_GABAA_1, tau_GABAA_2, GABA_A_g_peak )
GABA_BInitialValue real = compute_synapse_constant( tau_GABAB_1, tau_GABAB_2, GABA_B_g_peak )
RefractoryCounts integer = steps(t_ref) # refractory time in steps
alpha_n_init real = 0.032 * (V_m / mV + 52.) / (1. - exp(-(V_m / mV + 52.) / 5.))
beta_n_init real = 0.5 * exp(-(V_m / mV + 57.) / 40.)
alpha_m_init real = 0.32 * (V_m / mV + 54.) / (1.0 - exp(-(V_m / mV + 54.) / 4.))
beta_m_init real = 0.28 * (V_m / mV + 27.) / (exp((V_m / mV + 27.) / 5.) - 1.)
alpha_h_init real = 0.128 * exp(-(V_m / mV + 50.0) / 18.0)
beta_h_init real = 4.0 / (1.0 + exp(-(V_m / mV + 27.) / 5.))
end
input:
AMPA nS <- spike
NMDA nS <- spike
GABA_A nS <- spike
GABA_B nS <- spike
I_stim pA <- current
end
output: spike
update:
U_old mV = V_m
integrate_odes()
# sending spikes:
if r > 0: # is refractory?
r -= 1
elif V_m > V_Tr and U_old > V_Tr: # threshold && maximum
r = RefractoryCounts
emit_spike()
end
end
function compute_synapse_constant(Tau_1 ms, Tau_2 ms, g_peak real) real:
# Factor used to account for the missing 1/((1/Tau_2)-(1/Tau_1)) term
# in the ht_neuron_dynamics integration of the synapse terms.
# See: Exact digital simulation of time-invariant linear systems
# with applications to neuronal modeling, Rotter and Diesmann,
# section 3.1.2.
exact_integration_adjustment real = ( ( 1 / Tau_2 ) - ( 1 / Tau_1 ) ) * ms
t_peak real = ( Tau_2 * Tau_1 ) * ln( Tau_2 / Tau_1 ) / ( Tau_2 - Tau_1 ) / ms
normalisation_factor real = 1 / ( exp( -t_peak / Tau_1 ) - exp( -t_peak / Tau_2 ) )
return g_peak * normalisation_factor * exact_integration_adjustment
end
end