From 03541c040ff885335d01deb2e23cd1af5d18c545 Mon Sep 17 00:00:00 2001 From: Michael Schellenberger Costa Date: Wed, 23 Mar 2016 11:42:44 +0100 Subject: [PATCH] Fixed error in license and added option for noise free setting --- Data_Storage.h | 4 ++-- NM_Thalamus.pro | 4 ++-- Thalamic_Column.cpp | 56 ++++++++++++++++++++++++-------------------- Thalamic_Column.h | 67 ++++++++++++++++++++++++++++------------------------- Thalamus.cpp | 2 +- Thalamus_mex.cpp | 5 ++-- 6 files changed, 73 insertions(+), 65 deletions(-) diff --git a/Data_Storage.h b/Data_Storage.h index 6ebed5a..9e5d34b 100644 --- a/Data_Storage.h +++ b/Data_Storage.h @@ -25,7 +25,7 @@ * to auditory stimulation. * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, * JC Claussen. - * PLoS Computational Biology In Review (in review). + * PLoS Computational Biology (in review). */ /****************************************************************************************************/ @@ -44,4 +44,4 @@ void get_data(int counter, Thalamic_Column& Col, double* Vt, double* Vr, double* } /****************************************************************************************************/ /* end */ -/****************************************************************************************************/ \ No newline at end of file +/****************************************************************************************************/ diff --git a/NM_Thalamus.pro b/NM_Thalamus.pro index 40a9042..f6c53f8 100644 --- a/NM_Thalamus.pro +++ b/NM_Thalamus.pro @@ -9,8 +9,8 @@ SOURCES += Thalamic_Column.cpp \ Thalamus.cpp \ Thalamus_mex.cpp -HEADERS += ODE.h \ - Data_Storage.h \ +HEADERS += Data_Storage.h \ + ODE.h \ Random_Stream.h \ Thalamic_Column.h diff --git a/Thalamic_Column.cpp b/Thalamic_Column.cpp index 82c9525..29de11b 100644 --- a/Thalamic_Column.cpp +++ b/Thalamic_Column.cpp @@ -25,7 +25,7 @@ * to auditory stimulation. * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, * JC Claussen. - * PLoS Computational Biology In Review (in review). + * PLoS Computational Biology (in review). */ /****************************************************************************************************/ @@ -45,8 +45,14 @@ void Thalamic_Column::set_RNG(void) { /* Add the RNG for I_{l}*/ MTRands.push_back(random_stream_normal(0.0, dphi*dt)); - /* Add the RNG for I_{l,0} */ - MTRands.push_back(random_stream_normal(0.0, dt)); + /* Add the RNG for I_{l,0} + * In the case of a noise free system (dphi=0) set I_{l,0} + * to zero. + */ + if(dphi==0) + MTRands.push_back(random_stream_normal(0.0, 0)); + else + MTRands.push_back(random_stream_normal(0.0, dt)); /* Get the random number for the first iteration */ Rand_vars.push_back(MTRands[2*i]()); @@ -82,24 +88,24 @@ double Thalamic_Column::get_Qr (int N) const{ /****************************************************************************************************/ /* Excitatory input to TC population */ double Thalamic_Column::I_et (int N) const{ - double I = g_AMPA * s_tt[N]* (Vt[N]- E_AMPA); + double I = g_AMPA * s_et[N]* (Vt[N]- E_AMPA); return I; } /* Inhibitory input to TC population */ double Thalamic_Column::I_gt (int N) const{ - double I = g_GABA * s_rt[N]* (Vt[N]- E_GABA); + double I = g_GABA * s_gt[N]* (Vt[N]- E_GABA); return I; } /* Excitatory input to RE population */ double Thalamic_Column::I_er (int N) const{ - double I = g_AMPA * s_tr[N]* (Vr[N]- E_AMPA); + double I = g_AMPA * s_er[N]* (Vr[N]- E_AMPA); return I; } /* Inhibitory input to RE population */ double Thalamic_Column::I_gr (int N) const{ - double I = g_GABA * s_rr[N]* (Vr[N]- E_GABA); + double I = g_GABA * s_gr[N]* (Vr[N]- E_GABA); return I; } /****************************************************************************************************/ @@ -251,21 +257,21 @@ double Thalamic_Column::noise_aRK(int M) const{ /****************************************************************************************************/ 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 - (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 - (I_LK_r(N) + I_T_r(N))); + 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_tt [N+1] = s_tt [0] + A[N]*dt*(x_tt[N]); - s_tr [N+1] = s_tr [0] + A[N]*dt*(x_tr[N]); - s_rt [N+1] = s_rt [0] + A[N]*dt*(x_rt[N]); - s_rr [N+1] = s_rr [0] + A[N]*dt*(x_rr[N]); - x_tt [N+1] = x_tt [0] + A[N]*dt*(pow(gamma_e, 2) * ( - s_tt[N]) - 2 * gamma_e * x_tt[N]); - x_tr [N+1] = x_tr [0] + A[N]*dt*(pow(gamma_e, 2) * (N_tr * get_Qt(N) - s_tr[N]) - 2 * gamma_e * x_tr[N]); - x_rt [N+1] = x_rt [0] + A[N]*dt*(pow(gamma_g, 2) * (N_rt * get_Qr(N) - s_rt[N]) - 2 * gamma_g * x_rt[N]); - x_rr [N+1] = x_rr [0] + A[N]*dt*(pow(gamma_g, 2) * (N_rr * get_Qr(N) - s_rr[N]) - 2 * gamma_g * x_rr[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]); + x_et [N+1] = x_et [0] + A[N]*dt*(pow(gamma_e, 2) * ( - s_et[N]) - 2 * gamma_e * x_et[N]) + noise_xRK(N,0); + x_er [N+1] = x_er [0] + A[N]*dt*(pow(gamma_e, 2) * (N_tr * get_Qt(N) - s_er[N]) - 2 * gamma_e * x_er[N]); + x_gt [N+1] = x_gt [0] + A[N]*dt*(pow(gamma_g, 2) * (N_rt * get_Qr(N) - s_gt[N]) - 2 * gamma_g * x_gt[N]); + x_gr [N+1] = x_gr [0] + A[N]*dt*(pow(gamma_g, 2) * (N_rr * get_Qr(N) - s_gr[N]) - 2 * gamma_g * x_gr[N]); } /****************************************************************************************************/ /* end */ @@ -279,14 +285,14 @@ void Thalamic_Column::add_RK(void) { Vt [0] =(-3*Vt [0] + 2*Vt [1] + 4*Vt [2] + 2* Vt [3] + Vt [4])/6; Vr [0] =(-3*Vr [0] + 2*Vr [1] + 4*Vr [2] + 2* Vr [3] + Vr [4])/6; Ca [0] =(-3*Ca [0] + 2*Ca [1] + 4*Ca [2] + 2* Ca [3] + Ca [4])/6; - s_tt [0] =(-3*s_tt [0] + 2*s_tt [1] + 4*s_tt [2] + 2* s_tt [3] + s_tt [4])/6; - s_tr [0] =(-3*s_tr [0] + 2*s_tr [1] + 4*s_tr [2] + 2* s_tr [3] + s_tr [4])/6; - s_rt [0] =(-3*s_rt [0] + 2*s_rt [1] + 4*s_rt [2] + 2* s_rt [3] + s_rt [4])/6; - s_rr [0] =(-3*s_rr [0] + 2*s_rr [1] + 4*s_rr [2] + 2* s_rr [3] + s_rr [4])/6; - x_tt [0] =(-3*x_tt [0] + 2*x_tt [1] + 4*x_tt [2] + 2* x_tt [3] + x_tt [4])/6; - x_tr [0] =(-3*x_tr [0] + 2*x_tr [1] + 4*x_tr [2] + 2* x_tr [3] + x_tr [4])/6; - x_rt [0] =(-3*x_rt [0] + 2*x_rt [1] + 4*x_rt [2] + 2* x_rt [3] + x_rt [4])/6; - x_rr [0] =(-3*x_rr [0] + 2*x_rr [1] + 4*x_rr [2] + 2* x_rr [3] + x_rr [4])/6; + s_et [0] =(-3*s_et [0] + 2*s_et [1] + 4*s_et [2] + 2* s_et [3] + s_et [4])/6; + s_er [0] =(-3*s_er [0] + 2*s_er [1] + 4*s_er [2] + 2* s_er [3] + s_er [4])/6; + s_gt [0] =(-3*s_gt [0] + 2*s_gt [1] + 4*s_gt [2] + 2* s_gt [3] + s_gt [4])/6; + s_gr [0] =(-3*s_gr [0] + 2*s_gr [1] + 4*s_gr [2] + 2* s_gr [3] + s_gr [4])/6; + x_et [0] =(-3*x_et [0] + 2*x_et [1] + 4*x_et [2] + 2* x_et [3] + x_et [4])/6 + noise_aRK(0); + x_er [0] =(-3*x_er [0] + 2*x_er [1] + 4*x_er [2] + 2* x_er [3] + x_er [4])/6; + x_gt [0] =(-3*x_gt [0] + 2*x_gt [1] + 4*x_gt [2] + 2* x_gt [3] + x_gt [4])/6; + x_gr [0] =(-3*x_gr [0] + 2*x_gr [1] + 4*x_gr [2] + 2* x_gr [3] + x_gr [4])/6; h_T_t [0] =(-3*h_T_t[0] + 2*h_T_t[1] + 4*h_T_t[2] + 2* h_T_t [3] + h_T_t [4])/6; h_T_r [0] =(-3*h_T_r[0] + 2*h_T_r[1] + 4*h_T_r[2] + 2* h_T_r [3] + h_T_r [4])/6; m_h [0] =(-3*m_h [0] + 2*m_h [1] + 4*m_h [2] + 2* m_h [3] + m_h [4])/6; diff --git a/Thalamic_Column.h b/Thalamic_Column.h index 52b3830..498b2d9 100644 --- a/Thalamic_Column.h +++ b/Thalamic_Column.h @@ -25,7 +25,7 @@ * to auditory stimulation. * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, * JC Claussen. - * PLoS Computational Biology In Review (in review). + * PLoS Computational Biology (in review). */ /************************************************************************************************/ @@ -59,7 +59,7 @@ class Thalamic_Column { /* Constructor for simulation */ Thalamic_Column(double* Par) - : g_LK (Par[1]), g_h (Par[0]) + : g_LK (Par[0]), g_h (Par[1]) {set_RNG();} /* Iterate one time step through SRK4 */ @@ -116,20 +116,20 @@ class Thalamic_Column { /* Declaration and Initialization of parameters */ /* Membrane time in ms */ - const double tau_t = 20; - const double tau_r = 20; + const double tau_t = 20.; + const double tau_r = 20.; /* Maximum firing rate in ms^-1 */ const double Qt_max = 400.E-3; const double Qr_max = 400.E-3; /* Sigmoid threshold in mV */ - const double theta_t = -58.6; - const double theta_r = -58.6; + const double theta_t = -58.5; + const double theta_r = -58.5; /* Sigmoid gain in mV */ - const double sigma_t = 6; - const double sigma_r = 6; + const double sigma_t = 6.; + const double sigma_r = 6.; /* Scaling parameter for sigmoidal mapping (dimensionless) */ const double C1 = (M_PI/sqrt(3)); @@ -138,7 +138,10 @@ class Thalamic_Column { const double gamma_e = 70E-3; const double gamma_g = 100E-3; - /* Conductivities */ + /* Membrane capacitance in muF/cm^2 */ + const double C_m = 1.; + + /* Weights/ conductivities */ /* Leak in aU */ const double g_L = 1.; @@ -147,36 +150,36 @@ class Thalamic_Column { const double g_GABA = 1.; /* Potassium leak current in mS/m^2 */ - const double g_LK = 0.02; + double g_LK = 0.033; /* T current in mS/m^2 */ const double g_T_t = 3; const double g_T_r = 2.3; /* h current in mS/m^2 */ - const double g_h = 0.051; + double g_h = 0.02; /* Reversal potentials in mV */ /* Synaptic */ - const double E_AMPA = 0; - const double E_GABA = -70; + const double E_AMPA = 0.; + const double E_GABA = -70.; /* Leak */ - const double E_L_t = -70; - const double E_L_r = -70; + const double E_L_t = -70.; + const double E_L_r = -70.; /* Potassium */ - const double E_K = -100; + const double E_K = -100.; /* I_T current */ - const double E_Ca = 120; + const double E_Ca = 120.; /* I_h current */ - const double E_h = -40; + const double E_h = -40.; /* Calcium parameters */ const double alpha_Ca = -51.8E-6; /* influx per spike in nmol */ - const double tau_Ca = 10; /* calcium time constant in ms */ + const double tau_Ca = 10.; /* calcium time constant in ms */ const double Ca_0 = 2.4E-4; /* resting concentration */ /* I_h activation parameters */ @@ -184,8 +187,8 @@ class Thalamic_Column { const double k2 = 4E-4; const double k3 = 1E-1; const double k4 = 1E-3; - const double n_P = 4; - const double g_inc = 2; + const double n_P = 4.; + const double g_inc = 2.; /* Noise parameters in ms^-1 */ const double mphi = 0E-3; @@ -193,9 +196,9 @@ class Thalamic_Column { double input = 0.0; /* Connectivities (dimensionless) */ - const double N_tr = 3; - const double N_rt = 5; - const double N_rr = 30; + const double N_tr = 3.; + const double N_rt = 5.; + const double N_rr = 25.; /* Parameters for SRK4 iteration */ const vector A = {0.5, 0.5, 1.0, 1.0}; @@ -211,14 +214,14 @@ class Thalamic_Column { vector Vt = _INIT(E_L_t), /* TC membrane voltage */ Vr = _INIT(E_L_r), /* RE membrane voltage */ Ca = _INIT(Ca_0), /* Calcium concentration of TC population */ - s_tt = _INIT(0.0), /* PostSP from TC population to TC population */ - s_tr = _INIT(0.0), /* PostSP from TC population to RE population */ - s_rt = _INIT(0.0), /* PostSP from RE population to TC population */ - s_rr = _INIT(0.0), /* PostSP from RE population to RE population */ - x_tt = _INIT(0.0), /* derivative of s_tt */ - x_tr = _INIT(0.0), /* derivative of s_tr */ - x_rt = _INIT(0.0), /* derivative of s_rt */ - x_rr = _INIT(0.0), /* derivative of s_rr */ + s_et = _INIT(0.0), /* PostSP from TC population to TC population */ + s_er = _INIT(0.0), /* PostSP from TC population to RE population */ + s_gt = _INIT(0.0), /* PostSP from RE population to TC population */ + s_gr = _INIT(0.0), /* PostSP from RE population to RE population */ + x_et = _INIT(0.0), /* derivative of s_et */ + x_er = _INIT(0.0), /* derivative of s_er */ + x_gt = _INIT(0.0), /* derivative of s_gt */ + x_gr = _INIT(0.0), /* derivative of s_gr */ h_T_t = _INIT(0.0), /* inactivation of T channel */ h_T_r = _INIT(0.0), /* inactivation of T channel */ m_h = _INIT(0.0), /* activation of h channel */ diff --git a/Thalamus.cpp b/Thalamus.cpp index e2156d8..edec717 100644 --- a/Thalamus.cpp +++ b/Thalamus.cpp @@ -25,7 +25,7 @@ * to auditory stimulation. * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, * JC Claussen. - * PLoS Computational Biology In Review (in review). + * PLoS Computational Biology (in review). */ /****************************************************************************************************/ diff --git a/Thalamus_mex.cpp b/Thalamus_mex.cpp index ef858fb..ad1ef58 100644 --- a/Thalamus_mex.cpp +++ b/Thalamus_mex.cpp @@ -25,14 +25,13 @@ * to auditory stimulation. * M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, * JC Claussen. - * PLoS Computational Biology In Review (in review). + * PLoS Computational Biology (in review). */ /****************************************************************************************************/ /* Implementation of the simulation as MATLAB routine (mex compiler) */ /* mex command is given by: */ -/* mex CXXFLAGS="\$CXXFLAGS -std=c++11" Thalamus_mex.cpp Thalamic_Column.cpp */ -/* The Simulation requires the following boost libraries: Random */ +/* mex CXXFLAGS="\$CXXFLAGS -std=c++11 -O3" Thalamus_mex.cpp Thalamic_Column.cpp */ /****************************************************************************************************/ #include "mex.h" #include "matrix.h"