Permalink
Cannot retrieve contributors at this time
Fetching contributors…
| /* | |
| * Copyright (c) 2015 University of Lübeck | |
| * | |
| * Permission is hereby granted, free of charge, to any person obtaining a copy | |
| * of this software and associated documentation files (the "Software"), to deal | |
| * in the Software without restriction, including without limitation the rights | |
| * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell | |
| * copies of the Software, and to permit persons to whom the Software is | |
| * furnished to do so, subject to the following conditions: | |
| * | |
| * The above copyright notice and this permission notice shall be included in | |
| * all copies or substantial portions of the Software. | |
| * | |
| * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR | |
| * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, | |
| * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE | |
| * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER | |
| * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, | |
| * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN | |
| * THE SOFTWARE. | |
| * | |
| * AUTHORS: Michael Schellenberger Costa: mschellenbergercosta@gmail.com | |
| * | |
| * Based on: A thalamocortical neural mass model of the EEG during NREM sleep and its response | |
| * to auditory stimulation. | |
| * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, | |
| * JC Claussen. | |
| * PLoS Computational Biology http://dx.doi.org/10.1371/journal.pcbi.1005022 | |
| */ | |
| /******************************************************************************/ | |
| /* Functions of the thalamic module */ | |
| /******************************************************************************/ | |
| #include "Thalamic_Column.h" | |
| /******************************************************************************/ | |
| /* Initialization of RNG */ | |
| /******************************************************************************/ | |
| void Thalamic_Column::set_RNG(void) { | |
| extern const double dt; | |
| unsigned numRandomVariables = 1; | |
| MTRands.reserve(2*numRandomVariables); | |
| Rand_vars.reserve(2*numRandomVariables); | |
| for (unsigned i=0; i < numRandomVariables; ++i){ | |
| /* Add the RNG for I_{l}*/ | |
| MTRands.push_back(randomStreamNormal(0.0, dphi*dt)); | |
| /* Add the RNG for I_{l,0} */ | |
| MTRands.push_back(randomStreamNormal(0.0, dt)); | |
| /* Get the random number for the first iteration */ | |
| Rand_vars.push_back(MTRands[2*i]()); | |
| Rand_vars.push_back(MTRands[2*i+1]()); | |
| } | |
| } | |
| /******************************************************************************/ | |
| /* RK noise scaling */ | |
| /******************************************************************************/ | |
| double Thalamic_Column::noise_xRK(int N, int M) const{ | |
| return gamma_e * gamma_e * (Rand_vars[2*M] + Rand_vars[2*M+1]/std::sqrt(3))*B[N]; | |
| } | |
| double Thalamic_Column::noise_aRK(int M) const{ | |
| return gamma_e * gamma_e * (Rand_vars[2*M] - Rand_vars[2*M+1]*std::sqrt(3))/4; | |
| } | |
| /******************************************************************************/ | |
| /* Firing Rate functions */ | |
| /******************************************************************************/ | |
| double Thalamic_Column::get_Qt (int N) const{ | |
| return Qt_max/ (1 + exp(-C1 * (Vt[N] - theta_t) / sigma_t)); | |
| } | |
| double Thalamic_Column::get_Qr (int N) const{ | |
| return Qr_max / (1 + exp(-C1 * (Vr[N]- theta_r) / sigma_r)); | |
| } | |
| /******************************************************************************/ | |
| /* Synaptic currents */ | |
| /******************************************************************************/ | |
| /* Excitatory input to TC population */ | |
| double Thalamic_Column::I_et (int N) const{ | |
| return g_AMPA * s_et[N]* (Vt[N]- E_AMPA); | |
| } | |
| /* Inhibitory input to TC population */ | |
| double Thalamic_Column::I_gt (int N) const{ | |
| return g_GABA * s_gt[N]* (Vt[N]- E_GABA); | |
| } | |
| /* Excitatory input to RE population */ | |
| double Thalamic_Column::I_er (int N) const{ | |
| return g_AMPA * s_er[N]* (Vr[N]- E_AMPA); | |
| } | |
| /* Inhibitory input to RE population */ | |
| double Thalamic_Column::I_gr (int N) const{ | |
| return g_GABA * s_gr[N]* (Vr[N]- E_GABA); | |
| } | |
| /******************************************************************************/ | |
| /* I_T gating functions */ | |
| /******************************************************************************/ | |
| /* Activation in TC population after Destexhe 1996 */ | |
| double Thalamic_Column::m_inf_T_t (int N) const{ | |
| return 1/(1+exp(-(Vt[N]+59)/6.2)); | |
| } | |
| /* Activation in RE population after Destexhe 1996 */ | |
| double Thalamic_Column::m_inf_T_r (int N) const{ | |
| return 1/(1+exp(-(Vr[N]+52)/7.4)); | |
| } | |
| /* Deactivation in TC population after Destexhe 1996 */ | |
| double Thalamic_Column::h_inf_T_t (int N) const{ | |
| return 1/(1+exp( (Vt[N]+81)/4)); | |
| } | |
| /* Deactivation in RE population after Destexhe 1996 */ | |
| double Thalamic_Column::h_inf_T_r (int N) const{ | |
| return 1/(1+exp( (Vr[N]+80)/5)); | |
| } | |
| /* Deactivation time in RE population after Destexhe 1996 */ | |
| double Thalamic_Column::tau_h_T_t (int N) const{ | |
| return (30.8 + (211.4 + exp((Vt[N]+115.2)/5))/(1 + exp((Vt[N]+86)/3.2)))/3.7371928; | |
| } | |
| /* Deactivation time in RE population after Destexhe 1996 */ | |
| double Thalamic_Column::tau_h_T_r (int N) const{ | |
| return (85 + 1/(exp((Vr[N]+48)/4) + exp(-(Vr[N]+407)/50)))/3.7371928; | |
| } | |
| /******************************************************************************/ | |
| /* I_h gating functions */ | |
| /******************************************************************************/ | |
| /* Activation in TC population after Destexhe 1993 */ | |
| double Thalamic_Column::m_inf_h (int N) const{ | |
| return 1/(1+exp( (Vt[N]+75)/5.5)); | |
| } | |
| /* Activation time for slow components in TC population after Chen2012 */ | |
| double Thalamic_Column::tau_m_h (int N) const{ | |
| return (20 + 1000/(exp((Vt[N]+ 71.5)/14.2) + exp(-(Vt[N]+ 89)/11.6))); | |
| } | |
| /* Instantaneous calcium binding onto messenger protein after Chen2012 */ | |
| double Thalamic_Column::P_h (int N) const{ | |
| //return k1 * pow(Ca[N], n_P)/(k1*pow(Ca[N], n_P)+k2); | |
| return k1 * Ca[N] * Ca[N] * Ca[N] * Ca[N]/(k1 * Ca[N] * Ca[N] * Ca[N] * Ca[N]+k2); | |
| } | |
| /* Return I_h activation */ | |
| double Thalamic_Column::act_h (void) const{ | |
| return m_h[0] + g_inc * m_h2[0]; | |
| } | |
| /******************************************************************************/ | |
| /* Intrinsic currents */ | |
| /******************************************************************************/ | |
| /* Leak current of TC population */ | |
| double Thalamic_Column::I_L_t (int N) const{ | |
| return g_L * (Vt[N]- E_L_t); | |
| } | |
| /* Potassium leak current of TC population */ | |
| double Thalamic_Column::I_LK_t (int N) const{ | |
| return g_LK * (Vt[N]- E_K); | |
| } | |
| /* Leak current of RE population */ | |
| double Thalamic_Column::I_L_r (int N) const{ | |
| return g_L * (Vr[N]- E_L_r); | |
| } | |
| /* Potassium leak current of RE population */ | |
| double Thalamic_Column::I_LK_r (int N) const{ | |
| return g_LK * (Vr[N]- E_K); | |
| } | |
| /* T-type current of TC population */ | |
| double Thalamic_Column::I_T_t (int N) const{ | |
| return g_T_t * m_inf_T_t(N) * m_inf_T_t(N) * h_T_t[N] * (Vt[N]- E_Ca); | |
| } | |
| /* T-type current of RE population */ | |
| double Thalamic_Column::I_T_r (int N) const{ | |
| return g_T_r * m_inf_T_r(N) * m_inf_T_r(N) * h_T_r[N] * (Vr[N]- E_Ca); | |
| } | |
| /* h-type current of TC population */ | |
| double Thalamic_Column::I_h (int N) const{ | |
| return g_h * (m_h[N] + g_inc * m_h2[N]) * (Vt[N]- E_h); | |
| } | |
| /******************************************************************************/ | |
| /* SRK iteration */ | |
| /******************************************************************************/ | |
| void Thalamic_Column::set_RK (int N) { | |
| extern const double dt; | |
| 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))); | |
| 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))); | |
| Ca [N+1] = Ca [0] + A[N]*dt*(alpha_Ca * I_T_t(N) - (Ca[N] - Ca_0)/tau_Ca); | |
| 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); | |
| 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); | |
| 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]); | |
| m_h2 [N+1] = m_h2 [0] + A[N]*dt*(k3 * P_h(N) * m_h[N] - k4 * m_h2[N]); | |
| s_et [N+1] = s_et [0] + A[N]*dt*(x_et[N]); | |
| s_er [N+1] = s_er [0] + A[N]*dt*(x_er[N]); | |
| s_gt [N+1] = s_gt [0] + A[N]*dt*(x_gt[N]); | |
| s_gr [N+1] = s_gr [0] + A[N]*dt*(x_gr[N]); | |
| y [N+1] = y [0] + A[N]*dt*(x [N]); | |
| x_et [N+1] = x_et [0] + A[N]*dt*(gamma_e*gamma_e * ( + N_tp * Cortex->y[N] - s_et[N]) - 2 * gamma_e * x_et[N]) + noise_xRK(N,0); | |
| x_er [N+1] = x_er [0] + A[N]*dt*(gamma_e*gamma_e * (N_rt * get_Qt(N) + N_rp * Cortex->y[N] - s_er[N]) - 2 * gamma_e * x_er[N]); | |
| 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]); | |
| 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]); | |
| x [N+1] = x [0] + A[N]*dt*(nu * nu * ( get_Qt(N) - y [N]) - 2 * nu * x [N]); | |
| } | |
| void Thalamic_Column::add_RK(void) { | |
| add_RK(Vt); | |
| add_RK(Vr); | |
| add_RK(Ca); | |
| add_RK(s_et); | |
| add_RK(s_er); | |
| add_RK(s_gt); | |
| add_RK(s_gr); | |
| add_RK(y); | |
| add_RK_noise(x_et, 0); | |
| add_RK(x_er); | |
| add_RK(x_gt); | |
| add_RK(x_gr); | |
| add_RK(x); | |
| add_RK(h_T_t); | |
| add_RK(h_T_r); | |
| add_RK(m_h); | |
| add_RK(m_h2); | |
| /* Generate noise for the next iteration */ | |
| for (unsigned i=0; i<Rand_vars.size(); ++i) { | |
| Rand_vars[i] = MTRands[i]() + input; | |
| } | |
| } |