// Maurice Petroccine, SUNY Albany (mpetroccione@albany.edu)
// Last updated on 11/15/2021
// This is a ball and stick model in which a 15 um diameter soma is connected to a 180 um (WT) or 266 um (KO) long dendrite.
// The dendrite receives two inputs.
// The first input occurs at increasing distances from the soma (0.001<dist<0.999).
// The second input is either (1) located halfway between the first stimulus and the most distal end of the dendrite or (2) halfway between the first stimulus and the soma.
// Time step: 0.025 ms
// This code should record the maximum voltage reached during each iteration of the simulation and ouput these values in a dat file
// to run type go("FileNameGoesHere") where you substitute what you want to name the file for FileNameGoesHere
// File will appear in the same folder as the hoc file after it is done (make sure simulation is finished before opening)
load_file ("nrngui.hoc")
model_type = 1 // chose whether the model should be (1) WT 180 um (2) EAAC1 KO 266 um
stim_order = 1 // chose whether the stimulation order is (1) proximal to the soma first or (2)distal first
stimtime = 20 // time between stimulation
num_t_step = 50 // the number times to increment the inter-input interval
num_d_step = 50 // the number of points along the axon
timeoffset = 5 // the amount to increment the interstumulus interval by between each simulation
dt = 0.025 // time step of integration
tstop = 1000 // ms, each simulation stops after it reaches 1000 ms
strdef preface, dirstr
preface = "."
sprint(dirstr, "%s/all_tau_vecs.hoc", preface)
xopen(dirstr)
// Topology
create soma, dend
connect dend(0), soma(1)
access dend
soma { // creates soma and a dendrite
L = 15 // um, length of the soma
diam = 15 // um, diameter of the some
nseg = 1
}
if (model_type == WT){
dend {
L = 180 // um, length of the dendrite, WT = 180 um
diam = 2.75 // um, starting diameter of the dendrite
nseg = 51 // # of segments
diam(0:1) = 2.75:1 // diameter of dendrite tapers from 2.75 um to 1 um
}
}
if (model_type == WT){
dend {
L = 266 // um, length of the dendrite, KO = 266 um
diam = 2.75 // um, starting diameter of the dendrite
nseg = 51 // # of segments
diam(0:1) = 2.75:1 // diameter of dendrite tapers from 2.75 um to 1 um
}
}
// Biophysics
soma {
insert caL
pbar_caL = 0.0001
insert caL13
pcaLbar_caL13 = 0.0001
insert kas
gkbar_kas = 0.00025
insert kir
gkbar_kir = 0.00025
qfact_kir = 1
insert krp
gkbar_krp = 0.002
insert nap
gnabar_nap = 0.0001325
insert caldyn
insert na3
gbar_na3 = 0.035
insert kdr
gkdrbar_kdr = 0.03
insert cadyn
insert pas
g_pas = 0.95e-4
e_pas = -65
ek = -80
ena = 65
ecal =118
Ra = 100 // Ohm-cm
cm = 1 // uF/cm2
}
dend {
insert kdr
gkdrbar_kdr = 0.003
insert caL
pbar_caL = 1e-5
insert caL13
pcaLbar_caL13 = 1e-5
insert kas
gkbar_kas = 0.000025
insert kir
gkbar_kir = 0.00025
qfact_kir = 1
insert krp
gkbar_krp = 0.0002
insert pas
g_pas = 0.95e-4
e_pas = -65
ek = -80
ecal =118
Ra = 100 // Ohm-cm
cm = 2 // uF/cm2
}
celsius = 23 // recordings took place at room temperature
ena = 65 // sodium reversal potential
ek = -80 // potassium reversal potential
v_init = -71.5 // initial voltage
// Ca2+ concentrations and dynamics, from Mattioni
cai0_ca_ion = 0.001 // mM, Churchill 1998
cao0_ca_ion = 1.2 // mM, Ca2+ concentration in external solution
cali0_cal_ion = 0.001 // mM, Churchill 1998
calo0_cal_ion = 1.2 // mM, Ca2+ concentration in external solution
CAINF = 1e-5 // mM, steady state intracell ca conc.
TAUR = 43 // ms, time const of ca diffusion - Jackson 2003
CA_DRIVE = 10000
CA_PUMP = 0.02
proc set_cainf() { NEW_CAINF = $1 //Ca2+ dynamics, from Mattioni
nCA_INF = NEW_CAINF
forall if (ismembrane("cadyn")) {cainf_cadyn = NEW_CAINF}
forall if (ismembrane("caldyn")) {cainf_caldyn = NEW_CAINF}
}
proc set_taur() { NEW_TAUR = $1
nCA_TAUR = NEW_TAUR
forall if (ismembrane("cadyn")) {taur_cadyn = NEW_TAUR}
forall if (ismembrane("caldyn")) {taur_caldyn = NEW_TAUR}
}
proc set_cadrive() { NEW_DRIVE = $1
nCA_DRIVE = NEW_DRIVE
forall if (ismembrane("cadyn")) {drive_cadyn = NEW_DRIVE}
forall if (ismembrane("caldyn")) {drive_caldyn = NEW_DRIVE}
}
proc set_pump() { NEW_PUMP = $1
nCA_PUMP = NEW_PUMP
forall if (ismembrane("cadyn")) {pump_cadyn = NEW_PUMP}
forall if (ismembrane("caldyn")) {pump_caldyn = NEW_PUMP}
}
// Instrumentation
objref syn_ampa, syn_nmda, syn2_ampa, syn2_nmda, ns, ns2, nc, nc2, nc3, nc4
ns = new NetStim(0.5) // creates a netstim to trigger the first input
ns.interval = 1
ns.start = stimtime
ns.noise = 0
ns.number = 1
ns2 = new NetStim(0.5) // creates a netstim to trigger the second input
ns2.interval = 1
ns2.start = stimtime
ns2.noise = 0
ns2.number = 1
dend syn_ampa = new AMPA(0.1) //Creates four point-processes on the dendrite (2 AMPA and 2 NMDA inputs)
dend syn_nmda = new NMDA(0.1)
dend syn2_ampa = new AMPA(0.1)
dend syn2_nmda = new NMDA(0.1)
nc = new NetCon(ns, syn_ampa) // creates a netcon linking the AMPA input to the to the first netstim
nc.weight = 0.66
nc.delay = 0
nc2 = new NetCon(ns, syn_nmda) // creates a netcon linking the NMDA input to the to the first netstim
nc2.weight = 0.56
nc2.delay = 0
nc3 = new NetCon(ns, syn2_ampa) // creates a netcon linking the AMPA input to the to the second netstim
nc3.weight = 0.66
nc3.delay = 0
nc4 = new NetCon(ns, syn2_nmda) // creates a netcon linking the NMDA input to the to the second netstim
nc4.weight = 0.56
nc4.delay = 0
// Graphical display
objref g // create a new graph that displays the voltage at the soma
g = new Graph()
addplot(g,0)
g.exec_menu("Keep Lines")
g.size(0,50,-80,40)
g.addvar("soma.v(0.5)", 1, 1, 0.6, 0.9, 2)
//g.addvar("dend.v(0.7)", 2, 1, 0.6, 0.9, 2)
//g.addvar("dend.v(0.1)", 3, 1, 0.6, 0.9, 2)
// Sim control
v_init = -71.5 // intital voltage
objref recVm // creates vector object, records V at soma
recVm = new Vector()
recVm.record(&soma.v(0.5))
objref vChange // creates a vector one for keeping track of max V at soma per sweep
vChange = new Vector(101,0) // vector length is set to 100
objref diststep
diststep = new Vector(101) // used for plotting vMax vs the distance
diststep.indgen(0,1,0.01)
objref distdrop // creates a new graph - later used to plot V vs Dist
distdrop = new Graph()
distdrop.size(0,1,-70,50) // SCALE OF GRAPH IS WRONG ****FIX******
distdrop.exec_menu("Keep Lines")
objref vMax_list [num_t_step] // creates an object for storing the list of Vmax from of each sweep
objref vMax [num_t_step] // creates an object for calculating Vmax
objref tfil // creates tfil for writing to file
tfil = new File()
strdef tmpstr // creates a temporary string for writing the file name
strdef tstepname // creates a temporary string for writing the tstep name
proc initialize() { // Sets model to initial state (time, voltage)
t = 0
finitialize(v_init)
fcurrent()
}
proc integrate() { // moves the simulation forward
g.begin() // and updates the graph
while (t<tstop) {
fadvance()
g.plot(t)
}
g.flush()
}
proc go() { // the run command - note: does not work without have the access set to the dendrite
sprint(tmpstr, "%s.dat", $s1) // creates a string for the output file name based on what is written when the simulation is run
tfil.wopen(tmpstr) // opens that file
for z = 0,(num_t_step-1) { // runs simulations for the num_t_steps
vMax_list[z] = new Vector() // sets vMax_list and vMax to be vectors
vMax[z] = new Vector (101,0)
j=0
nc3.delay= (timeoffset*z) // sets the offset of the second set of synapses (AMPA component) to increment by "timeoffset"
nc4.delay= (timeoffset*z) // sets the offset of the second set of synapses (NMDA component) to increment by "timeoffset"
change_dist() // runs "change_dist" procedure to run the simulation at each point along the dendrite
vMax[z].line(distdrop, diststep) // prints the vMax to the graph
print z // this prints an update every 100 loops - comment out to run faster simulations
}
xytofile(vMax_list[w], tmpstr) // procedure to write vMax to a file
tfil.close() // closes the file
}
proc change_dist() { // loop to run the simulation at each point along the dendrite
for i = 0,(num_d_step-1) { // loops for the specified number of distance steps
if (j==0){ // avoids placing the inputs at postion 0 on the dendrite
j=0.001
}
if (j==1){ // avoids placing the inputs at postion 0 on the dendrite
j=0.999
}
syn_ampa.loc(j) // changes the first syn location (AMPA component)
syn_nmda.loc(j) // changes the first syn location (NMDA component)
if (stim_order == 1) { // places the second input halfway between the first input and the distal tip of the dendrite
syn2_ampa.loc(((1-j)/2)+j)
syn2_nmda.loc(((1-j)/2)+j)
} else if (stim_order == 2) { // places the second input halfway between the first input and the soma
syn2_nmda.loc(j/2)
syn2_ampa.loc(j/2)
} else { // error message if incorrect settings are selected
print "no stimulation order selected - please enter and run again"
}
initialize() // initialize the individual simulation
integrate() // run the individual simulation
vMax[z].x[i] = recVm.max() // records the most depolarized voltage experienced at the soma to vMax
vMax_list[z].append(vMax[z].x[i]) // appends the peak potential to vector "vMax_list"
j=j+(1/num_d_step) // moves the input down the axon by an amount = 1/dstep
}
}
proc xytofile() { local b // procedure to write the number of action potentials to a file and name the file
w=0
print "writing to ", $s2 // notifies when writing to file
for w=0,num_t_step-1 {
tsteplabel(w+1)
tfil.printf("\n%s\n", tstepname) // prints the specific time step and then the list of vMax values ordered by distance from the soma
for b=0,$o1.size()-1 tfil.printf("%g\n", vMax_list[w].x[b])}
}
proc tsteplabel (){ // procedure for creating the file lable for the tstep iteration
sprint(tstepname, "tstep\t%g", $1)
}