/NM_TC

Permalink
Browse files

Fixed wrong lebeling of connectivities

  • Loading branch information...
1 parent d31e02d commit 4b7a6f688d65396d3121c24be378273a686494e1 @miscco committed Mar 24, 2016
Unified Split
Showing with 105 additions and 105 deletions.
  1. +3 −3 Cortical_Column.cpp
  2. +19 −19 Cortical_Column.h
  3. +1 −1 Data_Storage.h
  4. +21 −21 Plots.m
  5. +9 −12 Stimulation.h
  6. +23 −23 Thalamic_Column.cpp
  7. +29 −26 Thalamic_Column.h
View
6 Cortical_Column.cpp
@@ -178,9 +178,9 @@ void Cortical_Column::set_RK (int N) {
s_gp[N+1] = s_gp[0] + A[N] * dt*(x_gp[N]);
s_gi[N+1] = s_gi[0] + A[N] * dt*(x_gi[N]);
y [N+1] = y [0] + A[N] * dt*(x [N]);
- x_ep[N+1] = x_ep[0] + A[N] * dt*(pow(gamma_e, 2) * (N_ep * get_Qp(N) + N_tp * Thalamus->y[N] - s_ep[N]) - 2 * gamma_e * x_ep[N]) + noise_xRK(N, 0);
- x_ei[N+1] = x_ei[0] + A[N] * dt*(pow(gamma_e, 2) * (N_ei * get_Qp(N) + N_ti * Thalamus->y[N] - s_ei[N]) - 2 * gamma_e * x_ei[N]) + noise_xRK(N, 1) ;
- x_gp[N+1] = x_gp[0] + A[N] * dt*(pow(gamma_g, 2) * (N_ip * get_Qi(N) - s_gp[N]) - 2 * gamma_g * x_gp[N]);
+ x_ep[N+1] = x_ep[0] + A[N] * dt*(pow(gamma_e, 2) * (N_pp * get_Qp(N) + N_pt * Thalamus->y[N] - s_ep[N]) - 2 * gamma_e * x_ep[N]) + noise_xRK(N, 0);
+ x_ei[N+1] = x_ei[0] + A[N] * dt*(pow(gamma_e, 2) * (N_ip * get_Qp(N) + N_it * Thalamus->y[N] - s_ei[N]) - 2 * gamma_e * x_ei[N]) + noise_xRK(N, 1) ;
+ x_gp[N+1] = x_gp[0] + A[N] * dt*(pow(gamma_g, 2) * (N_pi * get_Qi(N) - s_gp[N]) - 2 * gamma_g * x_gp[N]);
x_gi[N+1] = x_gi[0] + A[N] * dt*(pow(gamma_g, 2) * (N_ii * get_Qi(N) - s_gi[N]) - 2 * gamma_g * x_gi[N]);
x [N+1] = x [0] + A[N] * dt*(pow(nu, 2) * ( get_Qp(N) - y [N]) - 2 * nu * x [N]);
}
View
38 Cortical_Column.h
@@ -64,8 +64,8 @@ class Cortical_Column {
{set_RNG();}
Cortical_Column(double* Param, double* Con)
- :sigma_p (Param[0]), g_KNa (Param[1]), dphi (Param[2]),
- N_tp (Con[2]), N_ti (Con[3])
+ :sigma_p (Param[0]), g_KNa (Param[1]), dphi (Param[2]),
+ N_pt (Con[2]), N_it (Con[3])
{set_RNG();}
/* Connect to the thalamic module */
@@ -170,12 +170,12 @@ class Cortical_Column {
double input = 0.0;
/* Connectivities (dimensionless) */
- const double N_ep = 115;
- const double N_ei = 72;
- const double N_ip = 90;
+ const double N_pp = 115;
+ const double N_ip = 72;
+ const double N_pi = 90;
const double N_ii = 90;
- const double N_tp = 2.5;
- const double N_ti = 2.5;
+ const double N_pt = 2.5;
+ const double N_it = 2.5;
/* Pointer to thalamic column */
Thalamic_Column* Thalamus;
@@ -192,18 +192,18 @@ class Cortical_Column {
/* Population variables */
vector<double> Vp = _INIT(E_L_p), /* excitatory membrane voltage */
- Vi = _INIT(E_L_g), /* inhibitory membrane voltage */
- Na = _INIT(Na_eq), /* Na concentration */
- s_ep = _INIT(0.0), /* PostSP from excitatory to excitatory population */
- s_ei = _INIT(0.0), /* PostSP from excitatory to inhibitory population */
- s_gp = _INIT(0.0), /* PostSP from inhibitory to excitatory population */
- s_gi = _INIT(0.0), /* PostSP from inhibitory to inhibitory population */
- y = _INIT(0.0), /* axonal flux */
- x_ep = _INIT(0.0), /* derivative of s_ep */
- x_ei = _INIT(0.0), /* derivative of s_ei */
- x_gp = _INIT(0.0), /* derivative of s_gp */
- x_gi = _INIT(0.0), /* derivative of s_gi */
- x = _INIT(0.0); /* derivative of y */
+ Vi = _INIT(E_L_g), /* inhibitory membrane voltage */
+ Na = _INIT(Na_eq), /* Na concentration */
+ s_ep = _INIT(0.0), /* PostSP from excitatory to excitatory population */
+ s_ei = _INIT(0.0), /* PostSP from excitatory to inhibitory population */
+ s_gp = _INIT(0.0), /* PostSP from inhibitory to excitatory population */
+ s_gi = _INIT(0.0), /* PostSP from inhibitory to inhibitory population */
+ y = _INIT(0.0), /* axonal flux */
+ x_ep = _INIT(0.0), /* derivative of s_ep */
+ x_ei = _INIT(0.0), /* derivative of s_ei */
+ x_gp = _INIT(0.0), /* derivative of s_gp */
+ x_gi = _INIT(0.0), /* derivative of s_gi */
+ x = _INIT(0.0); /* derivative of y */
};
/****************************************************************************************************/
/* end */
View
2 Data_Storage.h
@@ -39,7 +39,7 @@
/* Save data */
/****************************************************************************************************/
void get_data(int counter, Cortical_Column& Cortex, Thalamic_Column& Thalamus, vector<double*> Data) {
- Data[0][counter] = Cortex.Ve [0];
+ Data[0][counter] = Cortex.Vp [0];
Data[1][counter] = Thalamus.Vt [0];
Data[2][counter] = Thalamus.Ca [0];
Data[3][counter] = Thalamus.act_h ();
View
42 Plots.m
@@ -9,53 +9,53 @@ function Plots(type)
%mex CXXFLAGS="\$CXXFLAGS -std=c++11 -O3" TC_mex.cpp Cortical_Column.cpp Thalamic_Column.cpp
if type == 1
- Param_Cortex = [4.7; % sigma_e
+ Param_Cortex = [4.7; % sigma_p
1.33; % g_KNa
- 20E-1]; % dphi
+ 2.]; % dphi
- Param_Thalamus = [0.048; % g_h
- 0.024]; % g_LK
+ Param_Thalamus = [0.034; % g_LK
+ 0.052]; % g_h
else
- Param_Cortex = [6; % sigma_e
- 2.05; % g_KNa
- 20E-1]; % dphi
+ Param_Cortex = [5.8; % sigma_p
+ 1.8; % g_KNa
+ 2.]; % dphi
- Param_Thalamus = [0.048; % g_h
- 0.02]; % g_LK
+ Param_Thalamus = [0.038; % g_LK
+ 0.052]; % g_h
end
-Connectivity = [ 2.6; % N_et
- 2.6; % N_er
- 5; % N_te
- 10]; % N_ti
+Connectivity = [ 2.5; % N_tp
+ 2.5; % N_rp
+ 5; % N_pt
+ 10]; % N_it
% stimulation parameters
% first number is the mode of stimulation
% 0 == none
% 1 == semi-periodic
% 2 == phase dependend
-var_stim = [ 0; % mode of stimulation
- 60; % strength of the stimulus in Hz (spikes per second)
- 120; % duration of the stimulus in ms
+var_stim = [ 2; % mode of stimulation
+ 200; % strength of the stimulus in Hz (spikes per second)
+ 100; % duration of the stimulus in ms
5; % time between stimulation events in s (ISI)
0; % range of ISI in s [ISI-range,ISI+range]
3; % Number of stimuli per event
1050; % time between stimuli within a event in ms
- 450]; % time until stimuli after minimum in ms
+ 400]; % time until stimuli after minimum in ms
T = 30; % duration of the simulation
-[Ve, Vt, Ca, ah, Marker_Stim] = TC_mex(T, Param_Cortex, Param_Thalamus, Connectivity, var_stim);
+[Vp, Vt, Ca, ah, Marker_Stim] = TC_mex(T, Param_Cortex, Param_Thalamus, Connectivity, var_stim);
L = length(Vt);
timeaxis = linspace(0,T,L);
figure(2)
-subplot(411), plot(timeaxis,Ve)
+subplot(411), plot(timeaxis,Vp)
title('Pyramidal membrane voltage'),
xlabel('Time in s'),
-ylabel('V_{e} in mV')
+ylabel('V_{p} in mV')
xlim([0,30]);
ylim([-80, -40]);
% vertical line for markers
@@ -69,7 +69,7 @@ function Plots(type)
xlabel('Time in s'),
ylabel('V_{t} in mV')
xlim([0,30]);
-ylim([-80,-35]);
+%ylim([-80,-35]);
% vertical line for markers
View
21 Stimulation.h
@@ -89,7 +89,7 @@ class Stim {
int time_between_stimuli = 1050E1;
/* Threshold for phase dependent stimulation */
- double threshold = -68.5;
+ double threshold = -68;
/* Internal variables */
/* Simulation on for TRUE and off for FALSE */
@@ -135,7 +135,7 @@ class Stim {
int count_pause = 0;
/* Old voltage value for minimum detection */
- double Ve_old = 0.0;
+ double Vp_old = 0.0;
/* Pointer to columns */
Cortical_Column* Cortex;
@@ -194,11 +194,8 @@ void Stim::setup (double* var_stim) {
/* If ISI is random create RNG */
if (ISI_range != 0){
- /* Create the generator */
- Generator = ENG(rand());
-
- /* Combine RNG with distribution */
- Uniform_Distribution = DIST_int(ISI-ISI_range, ISI+ISI_range);
+ /* Generate uniform distribution */
+ Uniform_Distribution = random_stream_uniform_int(ISI-ISI_range, ISI+ISI_range);
}
} else {
/* In case of phase dependent stimulation, time_to_stim is the time from minimum detection to start of stimulation */
@@ -236,7 +233,7 @@ void Stim::check_stim (int time) {
}
/* After last stimulus in event update the timer with respect to (random) ISI*/
else {
- time_to_stimuli += (ISI_range==0)? ISI : Uniform_Distribution(Generator);
+ time_to_stimuli += (ISI_range==0)? ISI : Uniform_Distribution();
/* Reset the stimulus counter for next stimulation event */
count_stimuli = 1;
@@ -248,19 +245,19 @@ void Stim::check_stim (int time) {
case 2:
/* Search for threshold */
if(!stimulation_started && !minimum_found && !threshold_crossed && time>onset_correction && !stimulation_paused) {
- if(Cortex->Ve[0]<=threshold) {
+ if(Cortex->Vp[0]<=threshold) {
threshold_crossed = true;
}
}
/* Search for minimum */
if(threshold_crossed) {
- if(Cortex->Ve[0]>Ve_old) {
+ if(Cortex->Vp[0]>Vp_old) {
threshold_crossed = false;
minimum_found = true;
- Ve_old = 0;
+ Vp_old = 0;
} else {
- Ve_old = Cortex->Ve[0];
+ Vp_old = Cortex->Vp[0];
}
}
View
46 Thalamic_Column.cpp
@@ -89,18 +89,18 @@ double Thalamic_Column::I_et (int N) const{
/* Inhibitory input to TC population */
double Thalamic_Column::I_gt (int N) const{
- double I = g_AMPA * 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_GABA * s_er[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;
}
/****************************************************************************************************/
@@ -135,15 +135,15 @@ double Thalamic_Column::h_inf_T_r (int N) const{
return h;
}
-/* deactivation time in RE population after Destexhe 1996 */
+/* Deactivation time in RE population after Destexhe 1996 */
double Thalamic_Column::tau_h_T_t (int N) const{
- double tau = (30.8 + (211.4 + exp((Vt[N]+115.2)/5))/(1 + exp((Vt[N]+86)/3.2)))/pow(3, 1.2);
+ double tau = (30.8 + (211.4 + exp((Vt[N]+115.2)/5))/(1 + exp((Vt[N]+86)/3.2)))/3.7371928;
return tau;
}
/* Deactivation time in RE population after Destexhe 1996 */
double Thalamic_Column::tau_h_T_r (int N) const{
- double tau = (85 + 1/(exp((Vr[N]+48)/4) + exp(-(Vr[N]+407)/50)))/pow(3, 1.2);
+ double tau = (85 + 1/(exp((Vr[N]+48)/4) + exp(-(Vr[N]+407)/50)))/3.7371928;
return tau;
}
/****************************************************************************************************/
@@ -169,7 +169,7 @@ double Thalamic_Column::tau_m_h (int N) const{
/* Instantaneous calcium binding onto messenger protein after Chen2012 */
double Thalamic_Column::P_h (int N) const{
//double P_h = k1 * pow(Ca[N], n_P)/(k1*pow(Ca[N], n_P)+k2);
- double P_h = k1 * Ca[N]* Ca[N]* Ca[N]/(k1* Ca[N]* Ca[N]* Ca[N]+k2);
+ double P_h = k1 * Ca[N] * Ca[N] * Ca[N] * Ca[N]/(k1 * Ca[N] * Ca[N] * Ca[N] * Ca[N]+k2);
return P_h;
}
@@ -252,23 +252,23 @@ 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);
- 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_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]);
- y [N+1] = y [0] + A[N]*dt*(x [N]);
- x_et [N+1] = x_et [0] + A[N]*dt*(pow(gamma_e, 2) * ( + N_pt * 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*(pow(gamma_e, 2) * (N_tr * get_Qt(N) + N_pr * Cortex->y[N] - s_er[N]) - 2 * gamma_e * x_er[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]);
- x [N+1] = x [0] + A[N]*dt*(pow(nu, 2) * ( get_Qt(N) - y [N]) - 2 * nu * x [N]);
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*(pow(gamma_e, 2) * ( + 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*(pow(gamma_e, 2) * (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*(pow(gamma_g, 2) * (N_tr * 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]);
+ x [N+1] = x [0] + A[N]*dt*(pow(nu, 2) * ( get_Qt(N) - y [N]) - 2 * nu * x [N]);
}
/****************************************************************************************************/
/* end */
@@ -284,13 +284,13 @@ void Thalamic_Column::add_RK(void) {
Ca [0] =(-3*Ca [0] + 2*Ca [1] + 4*Ca [2] + 2* Ca [3] + Ca [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_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;
+ 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;
y [0] =(-3*y [0] + 2*y [1] + 4*y [2] + 2* y [3] + y [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_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;
+ 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;
x [0] =(-3*x [0] + 2*x [1] + 4*x [2] + 2* x [3] + x [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;
Oops, something went wrong.

0 comments on commit 4b7a6f6

Please sign in to comment.