Mathematical model for windup (Aguiar et al. 2010)

 Download zip file   Auto-launch 
Help downloading and running models
Accession:128559
"Windup is characterized as a frequency-dependent increase in the number of evoked action potentials in dorsal horn neurons in response to electrical stimulation of afferent C-fibers. ... The approach presented here relies on mathematical and computational analysis to study the mechanism(s) underlying windup. From experimentally obtained windup profiles, we extract the time scale of the facilitation mechanisms that may support the characteristics of windup. Guided by these values and using simulations of a biologically realistic compartmental model of a wide dynamic range (WDR) neuron, we are able to assess the contribution of each mechanism for the generation of action potentials windup. ..."
Reference:
1 . Aguiar P, Sousa M, Lima D (2010) NMDA channels together with L-type calcium currents and calcium-activated nonspecific cationic currents are sufficient to generate windup in WDR neurons. J Neurophysiol 104:1155-66 [PubMed]
Model Information (Click on a link to find other models with that property)
Model Type: Neuron or other electrically excitable cell;
Brain Region(s)/Organism:
Cell Type(s): Wide dynamic range neuron;
Channel(s): I Na,p; I Na,t; I L high threshold; I N; I K; I K,Ca;
Gap Junctions:
Receptor(s): GabaA; AMPA; NMDA;
Gene(s):
Transmitter(s):
Simulation Environment: NEURON; MATLAB;
Model Concept(s): Activity Patterns; Action Potentials;
Implementer(s):
Search NeuronDB for information about:  GabaA; AMPA; NMDA; I Na,p; I Na,t; I L high threshold; I N; I K; I K,Ca;
/
WDR-Model
readme.html
AMPA_DynSyn.mod
CaIntraCellDyn.mod
GABAa_DynSyn.mod *
GABAb_DynSyn.mod *
HH2.mod *
iCaAN.mod *
iCaL.mod *
iKCa.mod *
iNaP.mod *
mGluR_DynSyn.mod
NK1_DynSyn.mod *
NMDA_DynSyn.mod *
herreroscatter.m
interneuron.hoc *
loadsynapticcurrents.m
mosinit.hoc
screenshot.jpg
WDR.hoc
wdr_spike_times.dat *
wdr-complete-model.hoc
wdr-complete-model.ses
wdr-complete-model-exportsyns.hoc
                            
//Created by Paulo Aguiar [pauloaguiar@fc.up.pt]

// Steps:
// 1 - CREATE WDR neuron and interneuron from the respective template cells
// 2 - SET INPUT SIGNALS/STIMULI
// 3 - SET SYNAPSES/RECEPTORS
// 4 - CONNECT STIMULI TO SYNAPSES
// 5 - PROVIDE NOISE TO WDR
// 6 - FINAL COMMANDS





// 1 *****************************************************************************************************

// CREATE WDR neuron and interneuron from the respective template cells

load_file("WDR.hoc")
load_file("interneuron.hoc")

objref wdr, interneuron
wdr = new WDR()
interneuron = new Interneuron()





// 2 *****************************************************************************************************

// SET INPUT SIGNALS/STIMULI
// C-fibres and Ad-fibres signals are driven by spike-times stored in vectors
// noise signals are driven by NetStim objects

// signal periods; all times are in ms
n_burst_spikes = 1		 	// define the number of spikes created by each stim in the fibres
T0 = 10.0					// if each stim produces a short burst of spikes, this is their ISI [ms]
n_close_stims = 1			// define the number of stims in each set (ex: 1=singlet; 3=triplets)
T1 = 300.0					// in heterogeneous stimulation n_close_stims are separated by T1 [ms]
n_stim_sets = 15			// define the number of stim sets
T2 = 1000.0					// period between stimuli sets [ms]
start_time = 1000.0			// time for first stim [ms]

//tstop is defined after loading the session (close to the end of this file)

// NOTES:	average stimuli period = T2 / n_close_stims
//			total stims = n_close_stims * n_stim_sets

// contruct the stimulus times vector
objref stim_times
stim_times  = new Vector(n_stim_sets*n_close_stims*n_burst_spikes, 0.0)
index = 0
for i=0,n_stim_sets-1 {
    for j=0,n_close_stims-1 {
		for k=0,n_burst_spikes-1 {
		    stim_times.x[index] = start_time + i * T2 + j * T1 + k * T0
		    index = index + 1
		}
    }
}

// for Poisson stimulation: contruct the stimulus times vector
/*
objref r, stim_times
r = new Random(666)
r.negexp(T2)
stim_times  = new Vector(n_stim_sets, 0.0)
last_time = start_time
for i=0,n_stim_sets-1 {
    stim_times.x[i] = last_time
    last_time = last_time + r.repick()
}
*/


printf("\nSTIMULATION TIMES")
stim_times.printf()




// 3 *****************************************************************************************************

// SET CONNECTIONS/SYNAPSES
//	Ad-fibres	--> NMDAR, AMPAR	--> wdr
//	C-fibres	--> NMDAR, AMPAR, NK1R	--> interneuron
//	interneuron	--> NMDAR, AMPAR, GABAR	--> wdr
//	interneuron	--> GABAAR		--> wdr
//	noise_inh	--> GABAAR		--> wdr
//	noise_exc	--> AMPAR		--> wdr

// Create post-synaptic receptors for all synapses



// SYNAPSES FROM Ad-FIBRES
// -----------------------

N_A  = 20		// number of A synapses impinging on WDR
Ad_delay = 30.0		// delay in the activation of Ad synapses [ms]
Ad_dispersion = 5.0	//dispersion in the arrival time of Ad signals [ms]

objref nil
objref ampar_Ad[N_A], nmdar_Ad[N_A]
objref ampa_Ad_conn[N_A], nmda_Ad_conn[N_A]

for i=0,N_A-1 {
    
    wdr.dendrite ampar_Ad[i]  = new AMPA_DynSyn(0.5)
    ampar_Ad[i].tau_fac   = 0.1 // NOT SUBJECT TO SHORT-TERM PLASTICITY
    ampar_Ad[i].tau_rec   = 0.1
    ampar_Ad[i].U1        = 1.0
    ampar_Ad[i].tau_rise  = 0.1
    ampar_Ad[i].tau_decay = 5.0
    
    //netcon = new NetCon(source, target, threshold [mV], delay [ms], weight [uS]) -> delay is overwritten by FIinitializeHandler
    ampa_Ad_conn[i] = new NetCon(nil, ampar_Ad[i], -30.0, 0.0 , 0.0008)//0.0008 (0.00085 in NMDA block exp)
    
    
    wdr.dendrite nmdar_Ad[i]  = new NMDA_DynSyn(0.5)
    nmdar_Ad[i].tau_fac   = 0.1 // NOT SUBJECT TO SHORT-TERM PLASTICITY
    nmdar_Ad[i].tau_rec   = 0.1
    nmdar_Ad[i].U1        = 1.0
    nmdar_Ad[i].tau_rise  = 2.0
    nmdar_Ad[i].tau_decay = 100.0
    
    //netcon = new NetCon(source, target, threshold [mV], delay [ms], weight [uS]) -> delay is overwritten by FIinitializeHandler
    nmda_Ad_conn[i] = new NetCon(nil, nmdar_Ad[i], -30.0, 0.0 , 0.0001)//0.0001
    
}

// the activation times of Ad-fibres are set together with C-fibres - below



// SYNAPSES FROM C-FIBRES
// ----------------------

objref nil
objref ampar_C, nmdar_C, nk1r_C, nk1r_C_WDR
objref ampa_C_conn, nmda_C_conn, nk1_C_conn, nk1_C_WDR_conn

interneuron.dendrite ampar_C = new AMPA_DynSyn(0.5)
ampar_C.tau_fac   = 0.1//3000.0
ampar_C.tau_rec   = 0.1//200.0
ampar_C.U1        = 1.0//0.2
ampar_C.tau_rise  = 0.1
ampar_C.tau_decay = 5.0

interneuron.dendrite nmdar_C = new NMDA_DynSyn(0.5)
nmdar_C.tau_fac   = 0.1//3000.0
nmdar_C.tau_rec   = 0.1//200.0
nmdar_C.U1        = 1.0//0.2
nmdar_C.tau_rise  = 2.0
nmdar_C.tau_decay = 100.0

interneuron.dendrite nk1r_C = new NK1_DynSyn(0.5)
nk1r_C.tau_fac    = 0.1//3000.0
nk1r_C.tau_rec    = 0.1//200.0
nk1r_C.U1         = 1.0//0.2 
nk1r_C.tau_rise   = 100.0
nk1r_C.tau_decay  = 3000.0

wdr.dendrite nk1r_C_WDR = new NK1_DynSyn(0.5) // substance P diffusion to deeper laminae reaches WDR which have some NK1R
nk1r_C_WDR.tau_fac    = 0.1//3000.0
nk1r_C_WDR.tau_rec    = 0.1//200.0
nk1r_C_WDR.U1         = 1.0//0.2 
nk1r_C_WDR.tau_rise   = 200.0
nk1r_C_WDR.tau_decay  = 3000.0


//netcon = new NetCon(source, target, threshold [mV], delay [ms], weight [uS]) -> delay is overwritten by FIinitializeHandler
ampa_C_conn     = new NetCon(nil, ampar_C   , -30.0, 0.0, 0.008000) //0.008000
nmda_C_conn     = new NetCon(nil, nmdar_C   , -30.0, 0.0, 0.004000) //0.004000
nk1_C_conn      = new NetCon(nil, nk1r_C    , -30.0, 0.0, 0.000020) //0.000020; if bigger interneuron will show strong windup (but constrained by iKCa)
nk1_C_WDR_conn  = new NetCon(nil, nk1r_C_WDR, -30.0, 0.0, 0.000015) //0.000015; smaller than above to reflect diffusion and less expression


//set activation times
objref fih
fih = new FInitializeHandler("loadqueue()")

proc loadqueue() { local ii, jj, event_time  localobj r
        
    r = new Random(123456789) // use a different seed if you need different random streams; "123456789" is the seed used in the paper figures
    
    //load Ad-fiber spike times
    for ii=0,stim_times.size()-1 {
	
	//distribute through all synapses
	for jj=0,N_A-1 {
	    
	    event_time = stim_times.x[ii] + Ad_delay + Ad_dispersion * r.repick()

	    ampa_Ad_conn[jj].event( event_time )
	    nmda_Ad_conn[jj].event( event_time )
	    
	}
    }        
    
    //load C-fiber spike times
    for ii=0,stim_times.size()-1 {
	
	event_time = stim_times.x[ii] // THE C-FIBRES DELAY IS PLACED AFTER THE INTERNEURON TO ALLOW CONTROLED DISPERSION

	ampa_C_conn.event( event_time )
	nmda_C_conn.event( event_time )
	nk1_C_conn.event( event_time )
	nk1_C_WDR_conn.event( event_time )	    
	
    }
    
}




// SYNAPSES FROM INTERNEURON
// -------------------------

N_IC = 20		// number of interneuron (from C) synapses impinging on WDR
C_delay = 200.0		// conduction in C-fibres is slower; this is the delay
C_dispersion = 20.0	// dispersion in the arrival time of C signals [ms]

objref ampar_interneuron[N_IC], nmdar_interneuron[N_IC], gabaar_interneuron[N_IC]
objref ampa_interneuron_conn[N_IC], nmda_interneuron_conn[N_IC], gabaa_interneuron_conn[N_IC]

objref r
r = new Random(123456789) // use a different seed if you need different random streams; the seed used for the figures in the paper was 123456789

for i=0,N_IC-1 {
    
    wdr.dendrite ampar_interneuron[i]  = new AMPA_DynSyn(0.5)
    ampar_interneuron[i].tau_fac   = 0.1
    ampar_interneuron[i].tau_rec   = 0.1
    ampar_interneuron[i].U1        = 1.0
    ampar_interneuron[i].tau_rise  = 0.1
    ampar_interneuron[i].tau_decay = 5.0
    
    wdr.dendrite nmdar_interneuron[i]  = new NMDA_DynSyn(0.5)
    nmdar_interneuron[i].tau_fac   = 0.1
    nmdar_interneuron[i].tau_rec   = 0.1
    nmdar_interneuron[i].U1        = 1.0
    nmdar_interneuron[i].tau_rise  = 2.0
    nmdar_interneuron[i].tau_decay = 100.0
    
    wdr.dendrite gabaar_interneuron[i]  = new GABAa_DynSyn(0.5)
    gabaar_interneuron[i].tau_fac   = 0.1 // NOT SUBJECT TO SHORT-TERM PLASTICITY
    gabaar_interneuron[i].tau_rec   = 0.1
    gabaar_interneuron[i].U1        = 1.0
    gabaar_interneuron[i].tau_rise  = 0.1
    gabaar_interneuron[i].tau_decay = 10.0
    
    
    event_time = C_delay + C_dispersion * r.repick()  
    printf("\n C-fibres activations: fibre [%d] - delay time = %f ms", i, event_time)
    
    //netcon = new NetCon(source, target, threshold [mV], delay [ms], weight [uS])
    interneuron.hillock ampa_interneuron_conn[i]  = new NetCon(&v(0.5), ampar_interneuron[i],  -30.0, event_time,     0.00020) //0.00020 (0.00047 in NMDA block exp)
    interneuron.hillock nmda_interneuron_conn[i]  = new NetCon(&v(0.5), nmdar_interneuron[i],  -30.0, event_time,     0.00020) //0.00020
    interneuron.hillock gabaa_interneuron_conn[i] = new NetCon(&v(0.5), gabaar_interneuron[i], -30.0, event_time+1.0, 0.00020) //+1.0, 0.00020 (this is not required for windup!; it serves only the purpose of demostrating that GABA_A blockers enhance windup profile responses)
    
}




	
// 5 ****************************************************************************************************

// PROVIDE NOISE TO THE WDR MEMBRANE POTENTIAL

// create source of stimulation for inhibitory and excitatory noise synapses
objref stim_exc, stim_inh

wdr.soma stim_exc = new NetStim(0.5)
stim_exc.interval	= 10.0
stim_exc.start		= 0.0
stim_exc.number		= 1000000
stim_exc.noise		= 1.0

wdr.soma stim_inh = new NetStim(0.5)
stim_inh.interval	= 10.0
stim_inh.start		= 0.0
stim_inh.number		= 1000000
stim_inh.noise		= 1.0


objref ampar_noise, gabaar_noise

wdr.soma ampar_noise = new AMPA_DynSyn(0.5)   // excitatory
ampar_noise.tau_fac    = 0.1
ampar_noise.tau_rec    = 0.1
ampar_noise.U1         = 1.0 
ampar_noise.tau_rise   = 0.1
ampar_noise.tau_decay  = 5.0

wdr.soma gabaar_noise = new GABAa_DynSyn(0.5) // inhibitory
gabaar_noise.tau_fac    = 0.1
gabaar_noise.tau_rec    = 0.1
gabaar_noise.U1         = 1.0 
gabaar_noise.tau_rise   = 0.1
gabaar_noise.tau_decay  = 5.0


//netcon = new NetCon(source, target, threshold [mV], delay [ms], weight [uS])
objref ampa_noise_conn, gabaa_noise_conn
ampa_noise_conn  = new NetCon(stim_exc, ampar_noise,  -30.0, 0.0, 0.0001) //0.0001
gabaa_noise_conn = new NetCon(stim_exc, gabaar_noise, -30.0, 0.0, 0.0001) //0.0001
//the mean of the noise current is roughly 0.001 mA/cm²

//ampa_noise_conn  = new NetCon(stim_exc, ampar_noise,  -30.0, 0.0, 0.000) //0.000100
//gabaa_noise_conn = new NetCon(stim_inh, gabaar_noise, -30.0, 0.0, 0.000) //0.000429
//the mean of this synaptic noise current (0.000429 GABAA + 0.000100 AMPA) is roughly 0.0 mA/cm²;
//the stdev of fluctuations produced in the membrane potential is 0.5226 mV
//this assumes a resting potential of -64.8806 mV (E_GABAA=-80mV; E_AMPA=0mV)




// 6 ****************************************************************************************************

// FINAL COMMANDS

// store spike-times
objref nc, nil, vec
wdr.soma nc = new NetCon(&v(.5), nil, -30.0, 0.0, 1.0)
vec = new Vector()
nc.record(vec)



access wdr.soma

load_file("wdr-complete-model.ses")

tstop = (n_stim_sets + 1) * T2 + start_time

celsius = 36

dt = 0.0125 //has finite representation in binary


init()

//forall cvode.print_event_queue()
xpanel("Aguiar et al. 2010")
  xbutton("Click here to START simulation to reproduce Fig 2A","run()")
  xlabel("(The simulation takes a few minutes)")
xpanel()


// write spike times to file "wdr_spike_times.dat"
objref fileobj
fileobj = new File()
fileobj.wopen("wdr_spike_times.dat")
vec.printf(fileobj)
fileobj.close()





//************************************************************************************
//UNITS	    
//Category						Variable		Units
//Time							t				[ms]
//Voltage						v				[mV]
//Current						i				[mA/cm2] (distributed)	[nA] (point process)
//Concentration					ko, ki, etc.	[mM]
//Specific capacitance			cm				[uf/cm2]
//Length						diam, L			[um]
//Conductance					g				[S/cm2] (distributed)	[uS] (point process)
//Cytoplasmic resistivity		Ra				[ohm cm]
//Resistance					Ri				[10E6 ohm]

Loading data, please wait...