Permalink
Newer
Older
100644 253 lines (218 sloc) 9.87 KB
Apr 24, 2014 @miscco Code cleanups.
1 /*
Feb 1, 2016 @miscco Updated the naming convention and switched to the Random_stream header
2 * Copyright (c) 2015 University of Lübeck
3 *
4 * Permission is hereby granted, free of charge, to any person obtaining a copy
5 * of this software and associated documentation files (the "Software"), to deal
6 * in the Software without restriction, including without limitation the rights
7 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
8 * copies of the Software, and to permit persons to whom the Software is
9 * furnished to do so, subject to the following conditions:
10 *
11 * The above copyright notice and this permission notice shall be included in
12 * all copies or substantial portions of the Software.
13 *
14 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
15 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
16 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
17 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
18 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
19 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
20 * THE SOFTWARE.
21 *
22 * AUTHORS: Michael Schellenberger Costa: mschellenbergercosta@gmail.com
23 *
24 * Based on: A thalamocortical neural mass model of the EEG during NREM sleep and its response
25 * to auditory stimulation.
26 * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz,
27 * JC Claussen.
Sep 13, 2016 @miscco Major cleanup and code modernization
28 * PLoS Computational Biology http://dx.doi.org/10.1371/journal.pcbi.1005022
Feb 1, 2016 @miscco Updated the naming convention and switched to the Random_stream header
29 */
Apr 24, 2014 @miscco Code cleanups.
30
Sep 13, 2016 @miscco Major cleanup and code modernization
31 /******************************************************************************/
32 /* Functions of the thalamic module */
33 /******************************************************************************/
Apr 18, 2013 @miscco renamed some files, fixed some typos and some descriptions.
34 #include "Thalamic_Column.h"
Sep 13, 2016 @miscco Major cleanup and code modernization
35
36 /******************************************************************************/
37 /* Initialization of RNG */
38 /******************************************************************************/
Mar 10, 2014 @miscco Code cleanups.
39 void Thalamic_Column::set_RNG(void) {
Sep 13, 2016 @miscco Major cleanup and code modernization
40 extern const double dt;
41 unsigned numRandomVariables = 1;
42
43 MTRands.reserve(2*numRandomVariables);
44 Rand_vars.reserve(2*numRandomVariables);
45 for (unsigned i=0; i < numRandomVariables; ++i){
46 /* Add the RNG for I_{l}*/
Sep 26, 2016 @miscco Fixed minor casting bug in renadom stream and renamed them to camel case
47 MTRands.push_back(randomStreamNormal(0.0, dphi*dt));
Sep 13, 2016 @miscco Major cleanup and code modernization
48
49 /* Add the RNG for I_{l,0} */
Sep 26, 2016 @miscco Fixed minor casting bug in renadom stream and renamed them to camel case
50 MTRands.push_back(randomStreamNormal(0.0, dt));
Sep 13, 2016 @miscco Major cleanup and code modernization
51
52 /* Get the random number for the first iteration */
53 Rand_vars.push_back(MTRands[2*i]());
54 Rand_vars.push_back(MTRands[2*i+1]());
55 }
56 }
57 /******************************************************************************/
58 /* RK noise scaling */
59 /******************************************************************************/
60 double Thalamic_Column::noise_xRK(int N, int M) const{
61 return gamma_e * gamma_e * (Rand_vars[2*M] + Rand_vars[2*M+1]/std::sqrt(3))*B[N];
Mar 10, 2014 @miscco Code cleanups.
62 }
63
Sep 13, 2016 @miscco Major cleanup and code modernization
64 double Thalamic_Column::noise_aRK(int M) const{
65 return gamma_e * gamma_e * (Rand_vars[2*M] - Rand_vars[2*M+1]*std::sqrt(3))/4;
66 }
Mar 10, 2014 @miscco Code cleanups.
67
Sep 13, 2016 @miscco Major cleanup and code modernization
68 /******************************************************************************/
69 /* Firing Rate functions */
70 /******************************************************************************/
Apr 18, 2013 @miscco renamed some files, fixed some typos and some descriptions.
71 double Thalamic_Column::get_Qt (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
72 return Qt_max/ (1 + exp(-C1 * (Vt[N] - theta_t) / sigma_t));
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
73 }
74
Apr 18, 2013 @miscco renamed some files, fixed some typos and some descriptions.
75 double Thalamic_Column::get_Qr (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
76 return Qr_max / (1 + exp(-C1 * (Vr[N]- theta_r) / sigma_r));
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
77 }
Apr 29, 2013 @miscco Cleaned up the code and structurized it
78
Sep 13, 2016 @miscco Major cleanup and code modernization
79 /******************************************************************************/
80 /* Synaptic currents */
81 /******************************************************************************/
Apr 24, 2014 @miscco Code cleanups.
82 /* Excitatory input to TC population */
Jul 3, 2013 @miscco Changed some parameters and formulas to match Arnes Formulation. Shows
83 double Thalamic_Column::I_et (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
84 return g_AMPA * s_et[N]* (Vt[N]- E_AMPA);
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
85 }
86
Apr 24, 2014 @miscco Code cleanups.
87 /* Inhibitory input to TC population */
Feb 1, 2016 @miscco Updated the naming convention and switched to the Random_stream header
88 double Thalamic_Column::I_gt (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
89 return g_GABA * s_gt[N]* (Vt[N]- E_GABA);
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
90 }
Apr 24, 2014 @miscco Code cleanups.
91 /* Excitatory input to RE population */
Jul 3, 2013 @miscco Changed some parameters and formulas to match Arnes Formulation. Shows
92 double Thalamic_Column::I_er (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
93 return g_AMPA * s_er[N]* (Vr[N]- E_AMPA);
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
94 }
95
Apr 24, 2014 @miscco Code cleanups.
96 /* Inhibitory input to RE population */
Feb 1, 2016 @miscco Updated the naming convention and switched to the Random_stream header
97 double Thalamic_Column::I_gr (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
98 return g_GABA * s_gr[N]* (Vr[N]- E_GABA);
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
99 }
100
Sep 13, 2016 @miscco Major cleanup and code modernization
101 /******************************************************************************/
102 /* I_T gating functions */
103 /******************************************************************************/
Apr 24, 2014 @miscco Code cleanups.
104 /* Activation in TC population after Destexhe 1996 */
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
105 double Thalamic_Column::m_inf_T_t (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
106 return 1/(1+exp(-(Vt[N]+59)/6.2));
Jan 31, 2013 @miscco Added the possibility to insert connectivities as parameters of the
107 }
108
Apr 24, 2014 @miscco Code cleanups.
109 /* Activation in RE population after Destexhe 1996 */
Nov 6, 2013 @miscco Added references.
110 double Thalamic_Column::m_inf_T_r (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
111 return 1/(1+exp(-(Vr[N]+52)/7.4));
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
112 }
113
Apr 24, 2014 @miscco Code cleanups.
114 /* Deactivation in TC population after Destexhe 1996 */
Jul 3, 2013 @miscco Changed some parameters and formulas to match Arnes Formulation. Shows
115 double Thalamic_Column::h_inf_T_t (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
116 return 1/(1+exp( (Vt[N]+81)/4));
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
117 }
118
Apr 24, 2014 @miscco Code cleanups.
119 /* Deactivation in RE population after Destexhe 1996 */
Nov 6, 2013 @miscco Added references.
120 double Thalamic_Column::h_inf_T_r (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
121 return 1/(1+exp( (Vr[N]+80)/5));
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
122 }
123
Oct 20, 2014 @miscco Changed some parameters to better fit with previous detailed single
124 /* Deactivation time in RE population after Destexhe 1996 */
Nov 6, 2013 @miscco Added references.
125 double Thalamic_Column::tau_h_T_t (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
126 return (30.8 + (211.4 + exp((Vt[N]+115.2)/5))/(1 + exp((Vt[N]+86)/3.2)))/3.7371928;
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
127 }
128
Apr 24, 2014 @miscco Code cleanups.
129 /* Deactivation time in RE population after Destexhe 1996 */
Apr 29, 2013 @miscco Cleaned up the code and structurized it
130 double Thalamic_Column::tau_h_T_r (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
131 return (85 + 1/(exp((Vr[N]+48)/4) + exp(-(Vr[N]+407)/50)))/3.7371928;
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
132 }
133
Sep 13, 2016 @miscco Major cleanup and code modernization
134 /******************************************************************************/
135 /* I_h gating functions */
136 /******************************************************************************/
Apr 24, 2014 @miscco Code cleanups.
137 /* Activation in TC population after Destexhe 1993 */
Apr 29, 2013 @miscco Cleaned up the code and structurized it
138 double Thalamic_Column::m_inf_h (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
139 return 1/(1+exp( (Vt[N]+75)/5.5));
Apr 29, 2013 @miscco Cleaned up the code and structurized it
140 }
141
Nov 5, 2014 @miscco Cleaned up the code (removed P_h as variable)
142 /* Activation time for slow components in TC population after Chen2012 */
Jul 3, 2013 @miscco Changed some parameters and formulas to match Arnes Formulation. Shows
143 double Thalamic_Column::tau_m_h (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
144 return (20 + 1000/(exp((Vt[N]+ 71.5)/14.2) + exp(-(Vt[N]+ 89)/11.6)));
Apr 29, 2013 @miscco Cleaned up the code and structurized it
145 }
Nov 5, 2014 @miscco Cleaned up the code (removed P_h as variable)
146
147 /* Instantaneous calcium binding onto messenger protein after Chen2012 */
148 double Thalamic_Column::P_h (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
149 //return k1 * pow(Ca[N], n_P)/(k1*pow(Ca[N], n_P)+k2);
150 return k1 * Ca[N] * Ca[N] * Ca[N] * Ca[N]/(k1 * Ca[N] * Ca[N] * Ca[N] * Ca[N]+k2);
Nov 5, 2014 @miscco Cleaned up the code (removed P_h as variable)
151 }
May 5, 2015 @miscco Added a project file for QT Creator.
152
153 /* Return I_h activation */
154 double Thalamic_Column::act_h (void) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
155 return m_h[0] + g_inc * m_h2[0];
May 5, 2015 @miscco Added a project file for QT Creator.
156 }
Apr 29, 2013 @miscco Cleaned up the code and structurized it
157
Sep 13, 2016 @miscco Major cleanup and code modernization
158 /******************************************************************************/
159 /* Intrinsic currents */
160 /******************************************************************************/
Apr 24, 2014 @miscco Code cleanups.
161 /* Leak current of TC population */
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
162 double Thalamic_Column::I_L_t (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
163 return g_L * (Vt[N]- E_L_t);
164
Sep 18, 2013 @miscco Slightly changed the parameters.
165 }
166
Apr 24, 2014 @miscco Code cleanups.
167 /* Potassium leak current of TC population */
Sep 18, 2013 @miscco Slightly changed the parameters.
168 double Thalamic_Column::I_LK_t (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
169 return g_LK * (Vt[N]- E_K);
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
170 }
171
Apr 24, 2014 @miscco Code cleanups.
172 /* Leak current of RE population */
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
173 double Thalamic_Column::I_L_r (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
174 return g_L * (Vr[N]- E_L_r);
Sep 18, 2013 @miscco Slightly changed the parameters.
175 }
176
Apr 24, 2014 @miscco Code cleanups.
177 /* Potassium leak current of RE population */
Sep 18, 2013 @miscco Slightly changed the parameters.
178 double Thalamic_Column::I_LK_r (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
179 return g_LK * (Vr[N]- E_K);
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
180 }
181
Apr 24, 2014 @miscco Code cleanups.
182 /* T-type current of TC population */
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
183 double Thalamic_Column::I_T_t (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
184 return g_T_t * m_inf_T_t(N) * m_inf_T_t(N) * h_T_t[N] * (Vt[N]- E_Ca);
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
185 }
186
Apr 24, 2014 @miscco Code cleanups.
187 /* T-type current of RE population */
Apr 29, 2013 @miscco Adapted the different currents as proposed by Zygierewicz2001. The I_T
188 double Thalamic_Column::I_T_r (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
189 return g_T_r * m_inf_T_r(N) * m_inf_T_r(N) * h_T_r[N] * (Vr[N]- E_Ca);
Feb 13, 2013 @miscco Reworked the cell. New HH-currents and physiological parameters.
190 }
191
Apr 24, 2014 @miscco Code cleanups.
192 /* h-type current of TC population */
Apr 18, 2013 @miscco renamed some files, fixed some typos and some descriptions.
193 double Thalamic_Column::I_h (int N) const{
Sep 13, 2016 @miscco Major cleanup and code modernization
194 return g_h * (m_h[N] + g_inc * m_h2[N]) * (Vt[N]- E_h);
Mar 25, 2013 @miscco Added I_h, I_CAN, I_KCa currents to the populations, I_h only on TC
195 }
Apr 29, 2013 @miscco Cleaned up the code and structurized it
196
Sep 13, 2016 @miscco Major cleanup and code modernization
197 /******************************************************************************/
198 /* SRK iteration */
199 /******************************************************************************/
Mar 10, 2014 @miscco Code cleanups.
200 void Thalamic_Column::set_RK (int N) {
Sep 13, 2016 @miscco Major cleanup and code modernization
201 extern const double dt;
202 Vt [N+1] = Vt [0] + A[N]*dt*(-(I_L_t(N) + I_et(N) + I_gt(N))/tau_t - C_m * (I_LK_t(N) + I_T_t(N) + I_h(N)));
203 Vr [N+1] = Vr [0] + A[N]*dt*(-(I_L_r(N) + I_er(N) + I_gr(N))/tau_r - C_m * (I_LK_r(N) + I_T_r(N)));
204 Ca [N+1] = Ca [0] + A[N]*dt*(alpha_Ca * I_T_t(N) - (Ca[N] - Ca_0)/tau_Ca);
205 h_T_t [N+1] = h_T_t[0] + A[N]*dt*(h_inf_T_t(N) - h_T_t[N])/tau_h_T_t(N);
206 h_T_r [N+1] = h_T_r[0] + A[N]*dt*(h_inf_T_r(N) - h_T_r[N])/tau_h_T_r(N);
207 m_h [N+1] = m_h [0] + A[N]*dt*((m_inf_h(N) * (1 - m_h2[N]) - m_h[N])/tau_m_h(N) - k3 * P_h(N) * m_h[N] + k4 * m_h2[N]);
208 m_h2 [N+1] = m_h2 [0] + A[N]*dt*(k3 * P_h(N) * m_h[N] - k4 * m_h2[N]);
209 s_et [N+1] = s_et [0] + A[N]*dt*(x_et[N]);
210 s_er [N+1] = s_er [0] + A[N]*dt*(x_er[N]);
211 s_gt [N+1] = s_gt [0] + A[N]*dt*(x_gt[N]);
212 s_gr [N+1] = s_gr [0] + A[N]*dt*(x_gr[N]);
213 x_et [N+1] = x_et [0] + A[N]*dt*(gamma_e*gamma_e * ( - s_et[N]) - 2 * gamma_e * x_et[N]) + noise_xRK(N,0);
214 x_er [N+1] = x_er [0] + A[N]*dt*(gamma_e*gamma_e * (N_rt * get_Qt(N) - s_er[N]) - 2 * gamma_e * x_er[N]);
215 x_gt [N+1] = x_gt [0] + A[N]*dt*(gamma_g*gamma_g * (N_tr * get_Qr(N) - s_gt[N]) - 2 * gamma_g * x_gt[N]);
216 x_gr [N+1] = x_gr [0] + A[N]*dt*(gamma_g*gamma_g * (N_rr * get_Qr(N) - s_gr[N]) - 2 * gamma_g * x_gr[N]);
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
217 }
Mar 10, 2014 @miscco Code cleanups.
218 /****************************************************************************************************/
219 /* end */
220 /****************************************************************************************************/
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
221
Mar 10, 2014 @miscco Code cleanups.
222
223 void Thalamic_Column::add_RK(void) {
Sep 13, 2016 @miscco Major cleanup and code modernization
224 add_RK(Vt);
225 add_RK(Vr);
226 add_RK(Ca);
227 add_RK(s_et);
228 add_RK(s_er);
229 add_RK(s_gt);
230 add_RK(s_gr);
231 add_RK_noise(x_et, 0);
232 add_RK(x_er);
233 add_RK(x_gt);
234 add_RK(x_gr);
235 add_RK(h_T_t);
236 add_RK(h_T_r);
237 add_RK(m_h);
238 add_RK(m_h2);
239
240 /* Generate noise for the next iteration */
241 for (unsigned i=0; i<Rand_vars.size(); ++i) {
242 Rand_vars[i] = MTRands[i]() + input;
243 }
Jan 16, 2013 @miscco Basic implementation of the model proposed in Steyn-Ross2009.
244 }
Aug 27, 2014 @miscco Removed unnecessary ODE.h file and made iteration a Class function. Not
245
246 void Thalamic_Column::iterate_ODE(void) {
Sep 13, 2016 @miscco Major cleanup and code modernization
247 /* First calculating every ith RK moment. Has to be in order, 1th moment first */
248 for (unsigned i=0; i<4; ++i) {
249 set_RK(i);
250 }
251 add_RK();
Aug 27, 2014 @miscco Removed unnecessary ODE.h file and made iteration a Class function. Not
252 }