// Maurice Petroccine, SUNY Albany (mpetroccione@albany.edu OR ...)
// Last updated on 11/19/2021
// This code applies excitatory and inhibitory inputs onto the reconstructed morphology of a biocytin filled WT MSN.
// It reports the number of action potentials recorded at the soma.
// The simulation is run for different frequencies of excitatory and inhibitory input.
// All inputs are located throughout the dendrites with a random uniform distribution. These inputs are active with a frequency described by a gaussian distribution around the specified mean value.
// The timing of the inputs is generated prior to running this code (see premade.hoc)
// Time step: 0.025 ms
// Window duration: 50 ms
// 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)
// Necessary files needed to run:
load_file ("nrngui.hoc")
load_file("import3d.hoc")
load_file ("ranstream.hoc")
model_type = 1 // choose whether the model should be WT (1) or KO (2)
dt = 0.025 // time step of integration
tstop = 1050 // ms, each simulation stops after it reaches 1050 ms
repetitions = 10 // number of simulations for each frequncy pairing to run
strdef preface, dirstr
preface = "."
sprint(dirstr, "%s/all_tau_vecs.hoc", preface)
xopen(dirstr)
/***** TOPOLOGY *****/
strdef cellToLoad
{
cellToLoad = "WT.swc" // chooses the cell to load
}
begintemplate Cell // This code is based on the implementation of http://www.neuron.yale.edu/phpbb/viewtopic.php?f=13&t=2272;
public soma, axon, dend, apic // creates a dummy cell to facilitate importing the morphology - will be overwritten
create soma[1],axon[1],dend[1],apic[1]
public all,somatic,axonal,basal,apical
objref all,somatic,axonal,basal,apical
proc init() {
all = new SectionList()
somatic = new SectionList()
axonal = new SectionList()
basal = new SectionList()
apical = new SectionList()
}
endtemplate Cell
obfunc mkcell() { localobj import,morph,cell // loads the cell
cell = new Cell()
morph = new Import3d_SWC_read()
morph.input($s1) // $s1 is the morphology name - returns an object "Cell" with the specified morphology
import = new Import3d_GUI(morph,0)
execute("forall delete_section()",cell)
import.instantiate(cell)
return cell
}
objref cell
{
cell = mkcell(cellToLoad)
}
celsius = 23 // experiments were performed at room temperature
forall if (issection(".*soma.*")) { // inserts conductances into all somatic sections
print secname()
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
celsius = 23
ek = -80
ena = 65
ecal = 118
Ra = 100 // Ohm-cm
cm = 1 // uF/cm2
}
forall if (issection(".*dend.*")) { // inserts conductances into all dendritic sections
print secname()
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
celsius = 23
Ra = 100 // Ohm-cm
cm = 2 // uF/cm2, increased cm (membrane capacitance) to account for spines
}
cai0_ca_ion = 0.001 // mM, Churchill & Macvicar (1998)
cao0_ca_ion = 1.2 // mM, Ca2+ concentration in external solution
cali0_cal_ion = 0.001 // mM, Churchill & Macvicar (1998)
calo0_cal_ion = 1.2 // mM, Ca2+ concentration in external solution
CAINF = 1e-5 // mM, steady state intracellular Ca2+ conc.
TAUR = 43 // ms, time const of Ca2+ diffusion - Jackson & Redman (2003)
CA_DRIVE = 10000
CA_PUMP = 0.02
proc set_cainf() { NEW_CAINF = $1 // Ca2+ dynamics, Mattioni & Le Novère (2013)
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}
}
celsius = 23 // recordings took place at room temperature
ena = 65 // sodium reversal potential
ek = -80 // potassium reversal potential
v_init = -71.5 // initial voltage
ISI = 10 // default NetStim parameters - will be overwritten
NUM = 220
START = 50
NOISE = 0
/***** INSTRUMENTATION *****/
objref syn_ampa[NUM], syn_nmda[NUM], ns, syn_gaba[NUM]
objref r2, vec // creates objects for assigning random placement of the input location
r2 = new Random(11021986)
vec = new Vector(1000)
proc setlocation() { local i // procedure for generating a random location for the inputs
r2 = new Random(11021986*(repetitions*l+z)) // change the seed of the RNG for each run
r2.discunif(2,40) // generate values between 0.001 and 0.999 with a uniform distribution
vec.setrand(r2)
for i = 0, NUM-1 { // set the dendritic location of the excitatory inputs based on these generated values
Cell[0].dend[vec.x[i]] syn_nmda[i] = new NMDA(0.5)
Cell[0].dend[vec.x[i]] syn_ampa[i] = new AMPA(0.5)
}
r2.repick // repicks from the set of random numbers
vec.setrand(r2)
for i = 0, NUM-1 { // set the dendritic location of the inhibitory inputs based on these generated values
Cell[0].dend[vec.x[i]] syn_gaba[i] = new GABA(0.5)
}
}
objref nc[NUM], nc2[NUM], nc3[NUM]
objref tvec, tvec2, nil
objref tvec, fvec
fvec = new File() // object for opening pre-generated lists of input times
tvec = new Vector() // vector for transfering excitatory stimulation times
tvec2 = new Vector () // vector for transfering ihibitory stimulation times
objref nsa, nsa2
objref nsalist, nsa2list
nsalist = new List() // creates list for storing excitatory stimulation times
nsa2list = new List() // creates list for storing inhibitory stimulation times
proc makenetstims2() { local i // procedure to produce netstims to trigger the excitatory inputs
nsalist.remove_all() // clears previously stored values
for i = 0, tvec.size - 1 {
nsa = new NetStim()
nsalist.append(nsa)
nsalist.o(i).interval = ISI // ISI doesn't matter because there is only 1 stimulus per netstim
nsalist.o(i).number = 1 // # of stims per netstim
nsalist.o(i).start = tvec.x[i] // start time is read from tvec
nsalist.o(i).noise = 0 // no noise
}
}
proc makenetstims3() { local i // procedure to produce netstims to trigger the inhibitory inputs
nsa2list.remove_all() // clears the inputs from previous runs
for i = 0, tvec2.size - 1 {
nsa2 = new NetStim()
nsa2list.append(nsa2)
nsa2list.o(i).interval = ISI // ISI doesn't matter because there is only 1 stimulus per netstim
nsa2list.o(i).number = 1 // # of stims per netstim
nsa2list.o(i).start = tvec2.x[i] // start time is read from tvec2
nsa2list.o(i).noise = 0 // no noise
}
}
proc instrumentation() { local i // creates netcons to link AMPA and NMDA point processes and sets their weight
for i = 0,nsalist.count()-1 {
if (model_type == 1){ // changes the weights for WT and KO models
nc[i] = new NetCon(nsalist.o(i), syn_ampa[i])
nc[i].weight = 10.84744
nc[i].delay = 0
nc2[i] = new NetCon(nsalist.o(i), syn_nmda[i])
nc2[i].weight = 10.56
nc2[i].delay = 0
}
if (model_type == 2){ // changes the weights for WT and KO models
nc[i] = new NetCon(nsalist.o(i), syn_ampa[i])
nc[i].weight = 10.84744
nc[i].delay = 0
nc2[i] = new NetCon(nsalist.o(i), syn_nmda[i])
nc2[i].weight = 10.56
nc2[i].delay = 0
}
}
}
proc instrumentation3() { local i // creates netcons to link GABA point processes and sets their weight
for i = 0,nsa2list.count()-1 {
nc3[i] = new NetCon(nsa2list.o(i), syn_gaba[i])
nc3[i].weight = 10.171 // should vary for WT = , KO =
nc3[i].delay = 0 // no delay
}
}
objref g // creates a new graph that displays the voltage at the soma
g = new Graph()
addplot(g,0)
g.exec_menu("Keep Lines")
g.size(0,2050,-80,50)
g.addvar("Cell[0].soma[0].v(0.5)", 1, 1, 0.6, 0.9, 2)
objref apc // creates an object for counting the action potentials at soma compartment 0 (threshold = 0 mV)
Cell[0].soma[0] apc = new APCount(0.5)
apc.thresh = 0
objref tfil // creates tfil for writing to output file
tfil = new File()
objref ap_vec // object for storing AP count
strdef tmpstr // creates a temporary string for writing the file name
strdef tmpstr2
strdef tmpstr3
strdef ap_name
proc initialize() { // Sets model to initial state (time, voltage)
t = 0
finitialize(v_init)
fcurrent()
}
u = 0
proc integrate() { // moves the simulation forward
g.begin() // updates the graph
while (t<tstop) {
fadvance()
g.plot(t)
}
g.flush()
}
ap_vec = new Vector()
proc go() { // the run command - note: does not work without the access set to the dendrite
for m = 1,2 {
for l = 1,2 {
for z = 0,repetitions-1 {
sprint(tmpstr2, "WT%dHZ_ex_times_%d.dat", (m*5), z+1) // selects the pregenerated file of random excitatory onset times
fvec.ropen(tmpstr2) // opens file of random onset times
tvec.scanf(fvec) // scans these values into tvec
sprint(tmpstr3, "WT%dHZ_inhib_times_%d.dat", (l*5), z+1) // selects the pregenerated file of random inhibitory onset times
fvec.ropen(tmpstr3) // opens file of random onset times
tvec2.scanf(fvec) // scans these values into tvec2
setlocation() // procedure to place inputs randomly throughout the dendrite
makenetstims2(100) // procedure to set the timing of excitatory inputs
makenetstims3(100) // procedure to set the timing of inhibitory inputs
instrumentation() // procedure to set the weight of excitatory inputs
instrumentation3() // procedure to set the weight of inhibitory inputs
start_sim() // starts the simulation
ap_vec.append(apc.n) // records AP count to ap_vec
print u
u = u+1
}
sprint(tmpstr, "%s%dHZ_ex_%dHZ_inhib.dat", $s1, (m*5),(l*5)) // creates a string to name the file
tfil.wopen(tmpstr) // opens the output file
xytofile(ap_vec, tmpstr) // procedure to write the number of APs to the output file
tfil.close() // closes the ouput file
z = 0
ap_vec.resize(0) // resizes ap_vec to 0 removing previous values
}
}
}
proc start_sim() {
initialize() // initialize the individual simulation
integrate() // run the individual simulation
}
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
ap_label(nsalist.o(1).interval)
tfil.printf("\n%s\n", ap_name)
for b=0,$o1.size()-1 {tfil.printf("%g\n", ap_vec.x[b])}
tfil.printf("MEAN AP COUNT\n%g", ap_vec.mean())
}
proc ap_label (){ // procedure for creating the file label for the current iteration
sprint(ap_name, "APCount\t%g", $1)
}