// This is a fully wired network that functions with 500 Granule Cells, 15 Mossy cells, 6 Basket cells and 6 HIPP Cells
// modified version of
// Dentate gyrus network model
// Santhakumar V, Aradi I, Soltesz I (2005) J Neurophysiol 93:437-53
// https://senselab.med.yale.edu/ModelDB/showModel.cshtml?model=51781&file=\dentategyrusnet2005\DG500_M7.hoc
// ModelDB file along with publication:
// Yim MY, Hanuschkin A, Wolfart J (2015) Hippocampus 25:297-308.
// http://onlinelibrary.wiley.com/doi/10.1002/hipo.22373/abstract
// modified by
// Man Yi Yim / 2015
// Alexander Hanuschkin / 2011
// Which figure (1-5) do you want to simulate?
// Please use the original parameters to produce the figures.
// If you want to run the simulation with your own parameters, set fig = 0
fig = 1
// SIMULATION: parameters
dt = 0.1 // simulation time [ms]
tstop = 200 // total simulation time [ms]
trial = 1 // trialnumber = seed of simulation
trial_noise = 1 // seed for noise to PP
// NOTE: // all NetStims share the same random number generator
// used solution given in http://www.neuron.yale.edu/neuron/node/61
err_ = load_file("ranstream.hoc")
objref cvode
cvode = new CVode()
using_cvode_ = 1
debug_ = 2 // print output
// everything: debug==1
// some: debug==0
// FLAGS: Input
p_sprouted_ = 0 // 10% sprouting (i.e. each GC connects to 50 other GCs randomly...) // original p = 0.02
PP_nstimulus_ = 1 // number of stimulated GC (100 in original file), set to 1 to have one input per GC (e.g. uncorrelated Poisson or correlated input (Myers and Scharfman, 2009)
PP_input2MC_ = 0 // stimulate MC with PP
PP_rate_ = 10. // rate of PP Poisson input
PP_rate_noise_ = 0. // rate of PP Poisson noise
PP_box_nr_ = 6 // number of active PPs
PP_box_stop_ = 35. // time of box [ms]
PP_box_start_ = 5. // shift of box [ms]
PP_box_nspk_ = 3 // number of spike per active PP
// FLAGS: Output
print_Vtrace = 1 // print voltage trace to file
print_Raster = 1 // print spike raster to file
print_template = 1 // write out spike sequence into template file
print_GC_sprouting_input_ = 0 // calculates and print out GC sprouting input number distribution
print_stim_ = 1 // print stimulus to GCs for debug purpose..
print_stim_noise = 1 // print noise to GCs for debug purpose.. (same output file as Poisson input to GCs!)
// FLAGS: Scaling
scale_gpas_dg_ = 100 // scaling [%] of gpas of th GC model (100% in the original file)
scale_sk_dg_ = 100 // scaling of Ca-dependent K (SK) of th GC model (100% in the original file)
scale_kir_ = 100 // scaling of KIR conductance in GC model
scale_gabaa_ = 100 // scaling of GABAA in GC model
scale_PP_strength = 10 // scaling of synaptic weight of PP
scale_PP2BC_strength = 100 // scaling of synaptic weight PP->BC (100% for original value which was 50% of PP->GC strength..)
scale_PP2HC_strength = 0 // scaling of synaptic weight PP->HIPP (set to 0% for Santhakumar because there no PP->HIPP connections...)
scale_HC2GC_strength = 100 // scaling of synaptic weight HC->GC (beta_HIPP in Myers and Scharfman, 2009)
scale_MC2GC_strength = 100 // scaling of synaptic weight MC->GC (beta_HIPP in Myers and Scharfman, 2009)
scale_GC2BC_strength = 300 // scaling of synaptic weight GC->BC
scale_BC2GC_strength = 300 // scaling of synaptic weight BC->GC
init_from_file_ = 0 // read in initial potential of all cells/compartments from file
init_to_file_ = 0 // save potential of all cells/compartments to file @ end of run (to be used for init the network later again)
strdef suffix, idname // define file names for simulation results output
suffix = "txt"
idname = "-pp10-gaba1-kir1-st0"
strdef init_state_fname_
init_state_fname_ = "NetworkState_init_10Hz.dat"
// Reset parameters to generate the figures in Yim et al (2015)
if (fig == 1) {
scale_PP_strength = 10
scale_kir_ = 100
scale_gabaa_ = 100
p_sprouted_ = 0
idname = "-pp10-gaba1-kir1-st0"
}
if (fig == 2) {
scale_PP_strength = 10
scale_kir_ = 100
scale_gabaa_ = 100
p_sprouted_ = 30
idname = "-pp10-gaba1-kir1-st30"
}
if (fig == 3) {
scale_PP_strength = 10
scale_kir_ = 400
scale_gabaa_ = 400
p_sprouted_ = 30
idname = "-pp10-gaba4-kir4-st30"
}
if (fig == 4) {
scale_PP_strength = 16
scale_kir_ = 100
scale_gabaa_ = 100
p_sprouted_ = 10
idname = "-pp16-gaba1-kir1-st10"
}
if (fig == 5) {
scale_PP_strength = 16
scale_kir_ = 400
scale_gabaa_ = 400
p_sprouted_ = 10
idname = "-pp16-gaba4-kir4-st10"
}
// End of parameter reset
// NETWORK: parameters
ngcell = 500 // number of GCs
nbcell = 6 // number of BCs (6 for original Samathakumar 2005)
nmcell = 15 // number of MCs
nhcell = 6 // number of HCs (6 for original Samathakumar 2005)
npp = 100 // 100 enorhinal layer II Neurons (Myers and Scharfman, 2009)
ntotal = ngcell + nbcell + nmcell + nhcell + npp + npp // two times npp because of PP input and PP noise...
if (PP_input2MC_ != 0) {
print "PP stimulation not yet implemented if stimulation of also MCs...."
quit()
}
objref PPSt[npp], PPSt_noise[npp]
create aacell // artificial cell if windows (boxed) activation
objref pp_start
objref netcon_d[npp]
// for debug purpose only...
objref netcon, nil
objref vec_debug_
objref PP_init_rnd // random number generator to select an PP stimulus site..
objref nslist, rslist, rs_noiselist
objref ns, rs
//Generate spike matrix and save to file
objref Spike[ntotal-1], Spike_times[ntotal-1]
for i=0, (ngcell+nbcell +nmcell +nhcell -1) {
Spike[i] = new Vector() // vector of spikes for each cell
Spike_times[i] = new Vector() // vector of spikes for each cell
}
strdef Spkstr
objref dfile
dfile = new File()
//********************** GRANULE CELL ****************************************
err_ = load_file("GC.hoc") // moved the granule cell definition to an external file => can be substituted by own reconstructed cell/ different model etc..
// ********* BASKET CELL **************************************
err_ = load_file("BC.hoc")
//**************** MOSSY CELL ***********************************************************
err_ = load_file("MC.hoc")
//*************** HIPP CELL ****************************
err_ = load_file("HIPP.hoc")
//********** Perforant Path Stimulus ***********************************************
err_ = load_file("PP.hoc")
//###############################################################################################################
//############## CONNECTING THE CELLS #############################
// NETWORK SPECIFICATION INTERFACE
for i=0, ngcell-1 {Gcell[i] = new GranuleCell(i)}
for i=0, nbcell-1 {Bcell[i] = new BasketCell(i)}
for i=0, nmcell-1 {Mcell[i] = new MossyCell(i)}
for i=0, nhcell-1 {Hcell[i] = new HIPPCell(i)}
for i=0, npp-1 {PPSt[i] = new PPstim(i)}
for i=0, npp-1 {PPSt_noise[i] = new PPstim(i)}
objref nclist, netcon, cells, net_c, net_d, net_gr, net_bc, net_mc, net_hc, vbc2gc, vmc2gc, vhc2gc
{ cells = new List()
nclist = new List()
}
// Include all cells in cells
func cell_append() {
cells.append($o1)
return cells.count -1
}
func nc_append() { // neuron connect $1 with $2.pre_list.object($3), weight $4, delay $5, threshold $6
// connects:
// cells.object($1) with
// $o1 = cells.object($2).pre_list.object($3) and
// returns:
// netcon = $o2
if ($3 >= 0 ) {
cells.object($1).connect_pre(cells.object($2).pre_list.object($3),netcon) // connect_pre is function in the respective cell definition
netcon.weight = $4 netcon.delay = $5 netcon.threshold = $6
}
nclist.append(netcon)
return nclist.count-1
}
func nc_append_rec() { // neuron connect $1 with $2.pre_list.object($3), weight $4, delay $5, threshold $6
// connects:
// cells.object($1) with
// $o1 = cells.object($2).pre_list.object($3) and
// returns:
// netcon = $o2
// record events to $o7
if ($3 >= 0 ) {
cells.object($1).connect_pre(cells.object($2).pre_list.object($3),netcon) // connect_pre is function in the respective cell definition
netcon.weight = $4 netcon.delay = $5 netcon.threshold = $6
netcon.record($o7) // the difference is here!
}
nclist.append(netcon)
return nclist.count-1
}
// To check for preexisting connections between randomly selected cells
// to avoid multiple contacts between same 2 cells
func is_connected() {local i, c
c=0
for i=0, nclist.count-1 {
net_c= nclist.object(i)
if (($o1 == net_c.postcell()) && ($o2 == net_c.precell())) {c=1}
}
return c
}
objref vbc2gc, vmc2gc, vhc2gc, vgc2bc, vbc2bc, vmc2bc, vhc2bc, vgc2mc, vbc2mc, vmc2mc, vhc2mc, vgc2hc, vmc2hc,vgc2gc
//To delete certain randomly selected cells from net - in this case 8 of 15 MCs
objref killMC
{
vgc2bc = new Vector(nbcell, 0) // defines vectors for each "pre2post" pair, vector length is same as number of post cells fills with 0
vbc2bc = new Vector(nbcell, 0) // increment x[i] for each connection made!
vmc2bc = new Vector(nbcell, 0)
vhc2bc = new Vector(nbcell, 0)
killMC = new Vector(8, 0) //To delete certain randomly selected cells from net - in this case 8 of 15 MCs
vgc2mc = new Vector(nmcell, 0)
vbc2mc = new Vector(nmcell, 0)
vmc2mc = new Vector(nmcell, 0)
vhc2mc = new Vector(nmcell, 0)
vgc2hc = new Vector(nhcell, 0)
vmc2hc = new Vector(nhcell, 0)
vbc2gc = new Vector(ngcell, 0)
vmc2gc = new Vector(ngcell, 0)
vhc2gc = new Vector(ngcell, 0)
vgc2gc = new Vector(ngcell, 0)
}
//initiating random number generator for each pre2post pair
//also randomly select MC to kill "deadMC"
objref rdsynb, rdsyna, rdgc2hc, rdgc2bc, rdgc2mc, rdbc2gc, rdbc2bc, rdbc2mc, deadMC
objref rdmc2gc1, rdmc2gc2, rdmc2bc, rdmc2mc, rdmc2mc1, rdmc2hc, rdhc2gc, rdhc2bc, rdhc2mc, rdgc2gc
objref pp2gc // connectivity vector PP -> GC
objref rnd_pp2gc // new random number generator for drawing connections
objref pp2bc // connectivity vector PP -> BC
objref rnd_pp2bc // new random number generator for drawing connections
objref pp2hc // connectivity vector PP -> HC
objref rnd_pp2hc // new random number generator for drawing connections
err_ = ropen("/proc/uptime") // get a seed that is changing based on the processing time
{
//rseed = fscan() // so simulation will not start with the same seed (original) // Fscan() reads the next number from the file opened with the ropen() command
rseed = trial // set the seed to trialnumber, to get reproduceable connectivity and responses
ropen()
}
//************************************ PP ***********************************************
rnd_pp2gc = new Random(rseed)
proc new_rnd_pp2gc() {rnd_pp2gc.discunif(0,npp-1)}
new_rnd_pp2gc()
rnd_pp2bc = new Random(rseed)
proc new_rnd_pp2bc() {rnd_pp2bc.discunif(0,npp-1)}
new_rnd_pp2bc()
rnd_pp2hc = new Random(rseed)
proc new_rnd_pp2hc() {rnd_pp2hc.discunif(0,npp-1)}
new_rnd_pp2hc()
//************************************ GC ***********************************************
rdgc2bc = new Random(rseed) // use for syn.connections
proc new_rdgc2bc() {rdgc2bc.discunif(-1,1)} // range is based on spread of connections on either side of RING- precell loc =0
new_rdgc2bc()
rdgc2mc = new Random(rseed) // use for syn.connections
proc new_rdgc2mc() {rdgc2mc.discunif(0,2)}
new_rdgc2mc()
rdgc2hc = new Random(rseed) // use for syn.connections
proc new_rdgc2hc() {rdgc2hc.discunif(-2,2)}
new_rdgc2hc()
rdgc2gc = new Random(rseed) // use for syn.connections
//proc new_rdgc2gc() {rdgc2gc.discunif(-50, 50)} // search around a GC ... However this has to be changed for large sprounting values
// n_spout_ = int (p_sprouted_ ) -1 // NOTE: 100% Sprouting = 100 new synapses!!! (compare p. 443 in Santhakumar paper)
proc new_rdgc2gc() {rdgc2gc.discunif(-50, 50)}
new_rdgc2gc()
//************************************ BC ***********************************************
rdbc2gc = new Random(rseed) // use for syn.connections
proc new_rdbc2gc() {rdbc2gc.discunif(-70, 70)} // range is based on spread of connections on either side of RING- precell loc =0
new_rdbc2gc()
rdbc2bc = new Random(rseed) // use for syn.connections
proc new_rdbc2bc() {rdbc2bc.discunif(-1, 1)}
new_rdbc2bc()
rdbc2mc = new Random(rseed) // use for syn.connections
proc new_rdbc2mc() {rdbc2mc.discunif(-3, 3)}
new_rdbc2mc()
//************************************* MC ********************************************
deadMC = new Random(rseed) // randomly select MC to kill "deadMC"
proc new_deadMC() {
deadMC.discunif(ngcell+nbcell, ngcell+nbcell+nmcell-1)
}
new_deadMC()
rdmc2gc1 = new Random(rseed) // use for syn.connections
proc new_rdmc2gc1() {rdmc2gc1.discunif(25, 175)} // range is based on spread of connections on either side of RING- precell loc =0
new_rdmc2gc1()
rdmc2gc2 = new Random(rseed) // use for syn.connections
proc new_rdmc2gc2() {rdmc2gc2.discunif(-175, -25)}
new_rdmc2gc2()
rdmc2bc = new Random(rseed) // use for syn.connections
proc new_rdmc2bc() {rdmc2bc.discunif(-3,3)}
new_rdmc2bc()
rdmc2mc = new Random(rseed) // use for syn.connections
proc new_rdmc2mc() {rdmc2mc.discunif(ngcell+nbcell, ngcell+nbcell+nmcell-1)}
new_rdmc2mc()
rdmc2mc1 = new Random(rseed) // use for syn.connections
proc new_rdmc2mc1() {rdmc2mc1.discunif(-3, 3)}
new_rdmc2mc1()
rdmc2hc = new Random(rseed) // use for syn.connections
proc new_rdmc2hc() {rdmc2hc.discunif(-2, 2)}
new_rdmc2hc()
//************************************* HC ********************************************
rdhc2gc = new Random(rseed) // use for syn.connections
proc new_rdhc2gc() {rdhc2gc.discunif(-130, 130)}
new_rdhc2gc()
rdhc2bc = new Random(rseed) // use for syn.connections
proc new_rdhc2bc() {rdhc2bc.discunif(-2, 2)}
new_rdhc2bc()
rdhc2mc = new Random(rseed) // use for syn.connections
proc new_rdhc2mc() {rdhc2mc.discunif(-2, 2)}
new_rdhc2mc()
//**************** randomizer for the dendritic location of synapse **************************************
rdsyna = new Random(rseed) // initialize random distr.
proc new_rdsyna() {rdsyna.discunif(0, 1)} // randomize among 2 dendrites
new_rdsyna()
rdsynb = new Random(rseed) // initialize random distr.
proc new_rdsynb() {rdsynb.discunif(0, 3)} // randomize among 4 dendrites
new_rdsynb()
// NETWORK INITIATION
for i = 0, ngcell-1 {cell_append(Gcell[i])} // cells 0-499 GCs
for i = 0, nbcell-1 {cell_append(Bcell[i])} // cells 500-505 BC
for i = 0, nmcell-1 {cell_append(Mcell[i])} // cells 506-520 MC
for i = 0, nhcell-1 {cell_append(Hcell[i])} // cells 521-526 HC
for i = 0, npp-1 {cell_append(PPSt[i])} // 527 - xxx PP artificial cell
for i = 0, npp-1 {cell_append(PPSt_noise[i])} // 527 - xxx PP artificial cell for noise
// STIM
// record input stimulus for debug purpose
objref vec_stim[npp]
for i = 0, npp-1 { vec_stim[i] = new Vector() }
// record input noise for debug purpose
objref vec_stim_noise[npp]
for i = 0, npp-1 { vec_stim_noise[i] = new Vector() }
objref pprec // vector if recorded from PP
// record Poisson in for debug purpose
objref distri_gc_input_
distri_gc_input_ = new Vector()
proc initNet_GC() { local i,j
//**************Preforant Path synaptic connections ******************************
if (debug_ == 2 ) {
print "Preforant Path synaptic connections"
print "Correlated input from 100 EC Layer 2 cells"
}
pprec = new Vector(npp, 0)
// ################ PP -> GC
for i=0, ngcell-1 {
pp2gc = new Vector(npp, 0) // init connectivity vector to prevent multiple stimulation connections from the same PP cell
for nr_conn = 0, int(npp/5.)-1 { // select randomly 20% of the PP for input
j = rnd_pp2gc.repick() // id of PP cell
while (pp2gc.x[j] == 1) { // random drawing until unconnected PP was found
j = rnd_pp2gc.repick() // id of PP cell
}
pp2gc.x[j] = 1
if ((print_stim_==1) && (i<npp)) {
if (pprec.x[j] == 0) {
nc_append_rec (j+ngcell+nbcell+nmcell+nhcell, i, 0, 2e-2*scale_PP_strength/100., 3, 10,vec_stim[j]) // record input sequence to first neuron for debug purpose only....
// nc_append_rec (j+ngcell+nbcell+nmcell+nhcell+npp, i, 1, 2e-2, 3, 10, vec_stim_noise[j]) // ADD NOISE TO THE GCs
pprec.x[j] = 1
} else {
nc_append (j+ngcell+nbcell+nmcell+nhcell, i, 0, 2e-2*scale_PP_strength/100., 3, 10) // record input sequence to first neuron for debug purpose only....
// nc_append (j+ngcell+nbcell+nmcell+nhcell+npp, i, 1, 2e-2, 3, 10) // ADD NOISE TO THE GCs
}
} else {
nc_append(j+ngcell+nbcell+nmcell+nhcell, i, 0, 2e-2*scale_PP_strength/100., 3, 10) // connect PP to GC[j],syn[0],wt,del,threshold <AH> NOTE both synapses are equal, ie. weight delay (except position) not important
}
}
if (debug_ == 3) {
print "\nGC: ",i
for ii=0,pp2gc.size()-1{ printf ("%d, ",pp2gc.x[ii])}
}
}
// ################ PP -> BC (input to two BC cells as in Santhakumar 2005)
if (scale_PP2BC_strength != 0) { // if PP2BC connections
for i=ngcell, ngcell+5 {
pp2bc = new Vector(npp, 0) // init connectivity vector to prevent multiple stimulation connections from the same PP cell
for nr_conn = 0, int(npp/5.)-1 { // select randomly 20% of the PP for input
j = rnd_pp2bc.repick() // id of PP cell
// print j
while (pp2bc.x[j] == 1) { // random drawing until unconnected PP was found
j = rnd_pp2bc.repick() // id of PP cell
}
pp2bc.x[j] = 1
nc_append(j+ngcell+nbcell+nmcell+nhcell, i, 0, 1e-2*scale_PP_strength/100.*scale_PP2BC_strength/100., 3, 10) // connect PP to GC[j],syn[0],wt,del,threshold <AH> NOTE both synapses are equal
// nc_append(j+ngcell+nbcell+nmcell+nhcell+npp, i, 0, 1e-2*scale_PP_strength/100.*scale_PP2BC_strength/100., 3, 10) // add noise to the GCs
}
if (debug_ == 3) {
print "\nBC: ",i
for ii=0,pp2bc.size()-1{ printf ("%d, ",pp2gc.x[ii])}
}
}
}
// ################ PP -> HIPP (input to 20% HIPP cells as in Myers and Scharfman 2009)
if (scale_PP2HC_strength != 0) { // if PP2HC connections
if (debug_ == 2) {print "PP -> HIPP (input to two HIPP cells as in Myers and Scharfman 2009)"}
for i=0, nhcell-1 {
pp2hc = new Vector(npp, 0) // init connectivity vector to prevent multiple stimulation connections from the same PP cell
for nr_conn = 0, int(npp/5.)-1 { // select randomly 20% of the PP for input
j = rnd_pp2hc.repick() // id of PP cell
// print j
while (pp2hc.x[j] == 1) { // random drawing until unconnected PP was found
j = rnd_pp2hc.repick() // id of PP cell
}
pp2hc.x[j] = 1
nc_append(j+ngcell+nbcell+nmcell+nhcell, i+ngcell+nbcell+nmcell, 0, 0.5e-3*scale_PP_strength/100.*scale_PP2HC_strength/100., 3, 10) // used here synaptic strength of GC->HIPP / connections PP-> HIPP not in Santhakumar 2005
// nc_append(j+ngcell+nbcell+nmcell+nhcell+npp,i+ngcell+nbcell+nmcell, 0, 0.5e-3*scale_PP_strength/100.*scale_PP2HC_strength/100., 3, 10) // add noise to the PP -> HIPP
}
if (debug_ == 3) {
print "\nHIPP: ",i
for ii=0,pp2bc.size()-1{ printf ("%d, ",pp2hc.x[ii])}
}
}
if (debug_ == 3) { print "PP -> HIPP -DONE-" }
}
//******************************************************************************************
//**************Granule Cell post synaptic connections ******************************
if (debug_ == 2 ) { print "Granule Cell post synaptic connections"}
if (debug_ == 0) {print "n_sprout_ = ",n_sprout_}
if (p_sprouted_>0 && debug_ == 2) {print "NOTE: we have sprouting connections..."}
for i=0, ngcell-1 {
for j=0, 0 {
// Based on the lamellar distribution of the GCs to BCs - 500 GCs were divided into 6 groups proximal to each of the 6(!) BCs
if (i < ngcell/6) { a=0}
if ((i > ngcell/6-1) && (i < ngcell*2/6)) { a=1}
if ((i > ngcell*2/6-1) && (i < ngcell*3/6)) { a=2}
if ((i > ngcell*3/6-1) && (i < ngcell*4/6)) { a=3}
if ((i > ngcell*4/6-1) && (i < ngcell*5/6)) { a=4}
if ((i > ngcell*5/6-1) && (i < ngcell)) { a=5}
Gauz3 = rdgc2bc.repick() // randomly pick location of post synaptic Bcell from distribution [-1:1]
if (a+Gauz3 > nbcell-1) {npost = a+Gauz3-nbcell } // determine appropriate post syn BC
if (a+Gauz3 < 0) {npost = a+Gauz3+nbcell}
if ((a+Gauz3 > -1) && (a+Gauz3 < nbcell)) {npost = a+Gauz3}
dbr = rdsynb.repick() // randomly pick the dendrite to connect to from [0:3] (i.e. // randomize among 4 dendrites)
if (debug_ == 1 ) { print "GC \t",i,"\t to BC\t\t",npost, a}
if (vgc2bc.x[npost] < 90) { // check to make sure that post syn BC does not get more than 90 GC synapses
nc_append(i, ngcell+npost, dbr+2, 4.7e-3*scale_GC2BC_strength/100., .8, 10) // connect GC[i] to BC[j],syn[2]+dendritic_var,wt,del,threshold
if (debug_ == 0 ) { print "GC \t",i,"\t to BC\t\t",npost,"\t",dbr+2}
vgc2bc.x[npost] +=1 // increment the no of synapses to the post cell
} else {
j -= 1
if (debug_ == 1 ) {print "nogc2bc"}
} // for connection that is not made reconnect axon to another cell
}
for j=0, 0 {
// Based on the lamellar distribution of the GCs to MCs - 500 GCs were divided into 5 groups, 3 MCs were distributed in each lamella
// print "Based on the lamellar distribution of the GCs to MCs..."
if (i < ngcell/5) { a=0}
if ((i > ngcell/5-1) && (i < ngcell*2/5)) { a=1}
if ((i > ngcell*2/5-1) && (i < ngcell*3/5)) { a=2}
if ((i > ngcell*3/5-1) && (i < ngcell*4/5)) { a=3}
if ((i > ngcell*4/5-1) && (i < ngcell)) { a=4}
b=a*3 // from [0:12]
npost = rdgc2mc.repick() // from [0:2]
dbr = rdsynb.repick() // from [0:2]
//print npost, b
if (vgc2mc.x[npost+b] < 38){
nc_append(i, ngcell+nbcell+npost+b, dbr+4, 0.2e-3, 1.5, 10) // Gcell[3] to Bcell[1]
if (debug_ == 1 ) {print "GC \t",i,"\t to MC\t\t", npost+b, "\t", dbr+4}
vgc2mc.x[npost+b] +=1
} else {
j -= 1
if (debug_ == 1 ) {print "nogc2mc"}
}
}
for j=0, 2 {
// <AH> comment added:
// Based on the lamellar distribution of the GCs to HIPPs - 500 GCs were divided into 6 groups
// print "Based on the lamellar distribution of the GCs to HIPPs..."
if (i < ngcell/6) { a=0}
if ((i > ngcell/6-1) && (i < ngcell*2/6)) { a=1}
if ((i > ngcell*2/6-1) && (i < ngcell*3/6)) { a=2}
if ((i > ngcell*3/6-1) && (i < ngcell*4/6)) { a=3}
if ((i > ngcell*4/6-1) && (i < ngcell*5/6)) { a=4}
if ((i > ngcell*5/6-1) && (i < ngcell)) { a=5}
Gauz3 = rdgc2hc.repick()
if (a+Gauz3 > nhcell-1) {npost = a+Gauz3-nhcell } // change here to allow more then 6 HIPP cells [-2:2]
if (a+Gauz3 < 0) {npost = a+Gauz3+nhcell}
if ((a+Gauz3 > -1) && (a+Gauz3 < nhcell)) {npost = a+Gauz3}
dbr = rdsynb.repick()
if ((is_connected(HIPPCell[npost], GranuleCell[i]) == 0) && (vgc2hc.x[npost] < 275)) {
nc_append(i, ngcell+nbcell+nmcell+npost, dbr, 0.5e-3, 1.5, 10) // Bcell -> Hcell
if (debug_ == 0 ) {print "GC \t",i,"\t to HC\t\t",npost, "\t", dbr}
vgc2hc.x[npost] +=1
} else {
j -= 1
}
}
// NOTE: THIS IS FOR SPROUTED SYNAPSES
n_sprout_ = p_sprouted_ -1 // NOTE: 100% Sprouting = 100 new synapses! (compare p. 443 in Santhakumar paper)
for j=0,n_sprout_ { // 9 in original file // each GC is recurrent connected to 10 GC (i.e. 10/500 => p = 0.02) but sprouting is diff different -> see above
Gauz3 = rdgc2gc.repick()
//print Gauz3
if (i+Gauz3 > 499) {npost = i+Gauz3-500 }
if (i+Gauz3 < 0) {npost = i+Gauz3+500}
if ((i+Gauz3 > -1) && (i+Gauz3 < 500)) {npost = i+Gauz3}
//print npost
dbr = rdsyna.repick() // [0:1]
if ((is_connected(GranuleCell[npost], GranuleCell[i]) == 0) && (vgc2gc.x[npost] < (n_sprout_*1.5+2) )) { // if is connected AND not more than 14 incoming connections...
// (original file < 15) (assume 1.5 times average connections for upper limit...)
nc_append(i, npost, dbr+7, 2e-3, .8, 10) // Gcell[3] to Bcell[1]
// print npost, dbr+8
if (debug_ == 0 ) {print "GC \t",i,"\t to GC\t\t",npost, "\t", dbr+7}
vgc2gc.x[npost] +=1
} else {
j -= 1
if (debug_ == 0) {print "gc2gc"}
}
}
}
if (print_GC_sprouting_input_ == 1) {
// objref distri_gc_input_
// distri_gc_input_ = new Vector()
max_gc_input_ = 0
if (debug_ ==2) { print "Calculate GC-GC Input Distribution"}
for zz = 0, int(n_sprout_*1.5+2) {distri_gc_input_.append(0)}
for npost=0,ngcell-1 {
distri_gc_input_.x[vgc2gc.x[npost]]+=1
if (vgc2gc.x[npost]>max_gc_input_) { max_gc_input_ = vgc2gc.x[npost]} // find max input number
}
for zz = 0, int(n_sprout_*1.5+2) {print zz,"\t",distri_gc_input_.x[zz]}
print "maximum input number is:\t",max_gc_input_
}
//******************************************************************************************
}
proc initNet_BC() { local i,j
//**************Basket Cell post synaptic connections ******************************
if (debug_ ==2) {print "Basket Cell post synaptic connections ... "}
for i=0, nbcell-1 {
for j=0, 99 {
Gauz3 = rdbc2gc.repick() // [-70:70]
if (debug_ == 1 ) {print Gauz3}
if (i*83+41+Gauz3 > 499) {npost = i*83+41+Gauz3-500 }
if (i*83+41+Gauz3 < 0) {npost = i*83+41+Gauz3+500}
if ((i*83+41+Gauz3 > -1) && (i*83+41+Gauz3 < 500)) {npost = i*83+41+Gauz3}
if (debug_ == 1 ) {print i, npost}
if (nbcell != 6) {max_conn = 4} else {max_conn = 2} // if not original setup use more spread...
if ((is_connected(GranuleCell[npost], BasketCell[i]) == 0) && (vbc2gc.x[npost] < max_conn)) { // change < 2 to < 4
nc_append(i+ngcell, npost, 6, 1.6e-3*scale_BC2GC_strength/100, .85, -10)
vbc2gc.x[npost] +=1
if (debug_ == 1 ) {print i, npost, 6}
} else {
j -= 1
if (debug_ == 1 ) {print "BC2GC"}
}
}
for j=0, 1 {
Gauz3 = rdbc2bc.repick() // [-1,0,1] (postsyn spread around single id...)
//print Gauz3
if (i+Gauz3 > nbcell-1) {npost = i+Gauz3-nbcell } // periodic boundary conditions
if (i+Gauz3 < 0) {npost = i+Gauz3+nbcell}
if ((i+Gauz3 >-1) && (i+Gauz3 < nbcell)) {npost = i+Gauz3}
dbr = rdsyna.repick() // [0:1]
if (nbcell != 6) {max_conn = 4} else {max_conn = 3}
if ((is_connected(BasketCell[npost], BasketCell[i]) == 0) && (vbc2bc.x[npost] < max_conn)) { // change < 3 to < 4
nc_append(i+ngcell, npost+ngcell, dbr+8, 7.6e-3, .8, -10)
if (debug_ == 1 ) {print npost, dbr+8}
vbc2bc.x[npost] +=1
} else {
j -= 1
if (debug_ == 1 ) {print "bc2bc"}
}
}
for j=0, 2 {
Gauz3 = rdbc2mc.repick() // [-3:3]
if (i*2+2+Gauz3 > 14) {npost = i*2+2+Gauz3-15 }
if (i*2+2+Gauz3 < 0) {npost = i*2+2+Gauz3+15}
if ((i*2+2+Gauz3 >-1) && (i*2+2+Gauz3 < 15)) {npost = i*2+2+Gauz3}
// if ((is_connected(MossyCell[npost], BasketCell[i]) == 0) && (vbc2mc.x[npost] < 3) && (killMC.contains(ngcell+nbcell+npost) == 0)) { // use if killing MC
if (nbcell != 6) {max_conn = 4} else {max_conn = 3}
if ((is_connected(MossyCell[npost], BasketCell[i]) == 0) && (vbc2mc.x[npost] < max_conn)) { // change < 3 to < 4
nc_append(i+ngcell, npost+ngcell+nbcell, 12, 1.5e-3, 1.5, -10) // Gcell[3] to Bcell[1]
if (debug_ == 1 ) {print npost, 12}
vbc2mc.x[npost] +=1
} else {
j -= 1
if (debug_ == 1 ) {print "bc2mc"}
}
// if (killMC.contains(ngcell+nbcell+npost) == 1) {j +=1 if (debug_ == 1 ) {print "dead MC"}} // use if killing MC
}
}
//******************************************************************************************
}
proc initNet_rest() { local i,j
//**************Mossy Cell post synaptic connections ******************************
if (debug_ ==2 ) { print "Mossy Cell post synaptic connections"}
for i=0, nmcell-1 {
// MC -> GC
//if (killMC.contains(ngcell+nbcell+i) == 0) // use if killing MC
if (i < 3) { y=0}
if ((i > 2) && (i < 6)) { y=1}
if ((i > 5) && (i < 9)) { y=2}
if ((i > 8) && (i < 12)) { y=3}
if ((i > 11) && (i < 15)) { y=4}
for j=0, 99 {
Gauz1 = rdmc2gc1.repick() // [25:175]
//if (debug_ == 1 ) {print Gauz1}
if (i*33+17+Gauz1 > 499) {
npost1 = i*33+17+Gauz1-500
} else {npost1 =i*33+17+Gauz1}
//if (debug_ == 1 ) {print npost1}
dbr = rdsyna.repick() // [0:1]
if ((is_connected(GranuleCell[npost1], MossyCell[i]) == 0) && (vmc2gc.x[npost1] < 7)) {
nc_append(i+ngcell+nbcell, npost1, dbr+2, 0.3e-3*scale_MC2GC_strength/100., 3, 10) // Gcell[3] to Bcell[1]
vmc2gc.x[npost1] +=1
} else { j -= 1 } // if (debug_ == 1 ) {print "MC2GC1"}
}
for j=0, 99 {
Gauz2 = rdmc2gc2.repick() // [-175:25]
if (i*33+17+Gauz2 < 0) {
npost2 =i*33+17+Gauz2+500
} else {npost2 =i*33+17+Gauz2}
dbr = rdsyna.repick() // [0:1]
if ((is_connected(GranuleCell[npost2], MossyCell[i]) == 0) && (vmc2gc.x[npost2] < 7)) {
nc_append(i+ngcell+nbcell, npost2, dbr+2, 0.3e-3*scale_MC2GC_strength/100., 3, 10) // Gcell[3] to Bcell[1]
vmc2gc.x[npost2] +=1
} else { j -= 1 }
// if (debug_ == 1 ) {print "MC2GC2"}
}
// MC -> BC
for j=0, 0 {
Gauz3 = rdmc2bc.repick() // Gauz3 = [-3:3]
if (y+Gauz3 > nbcell-1) {npost = y+Gauz3-nbcell} // y = [0:4] => y+Gaus3 = [-3:7]
if (y+Gauz3 < 0) {npost = y+Gauz3+nbcell}
if ((y+Gauz3 > -1) && (y+Gauz3 < nbcell)) {npost = y+Gauz3}
dbr = rdsyna.repick()
if ((vmc2bc.x[npost] < 4) && (Gauz3 !=0)) {
nc_append(i+ngcell+nbcell, ngcell+npost, dbr+6, 0.3e-3, 3, 10) // Gcell[3] to Bcell[1]
vmc2bc.x[npost] += 1
} else {
j -= 1
// if (debug_ == 1 ) {print "mc2bc"}
}
}
// MC -> MC
for j=0, 2 {
Gauz3 = rdmc2mc1.repick() //[-3:3]
if (i+Gauz3 > 14) {npost = i+Gauz3-15 }
if (i+Gauz3 < 0) {npost = i+Gauz3+15}
if ((i+Gauz3 >-1) && (i+Gauz3 < 15)) {npost = i+Gauz3}
dbr = rdsynb.repick()
// if ((is_connected(MossyCell[npost], MossyCell[i]) == 0) && (vmc2mc.x[npost] < 4) && (Gauz3 != 0) && (killMC.contains(ngcell+nbcell+npost) == 0)) // use if killing MC
if ((is_connected(MossyCell[npost], MossyCell[i]) == 0) && (vmc2mc.x[npost] < 4) && (Gauz3 != 0)) {
nc_append(i+ngcell+nbcell, npost+ngcell+nbcell, dbr+8, 0.5e-3, 2, 10) // Gcell[3] to Bcell[1]
// print npost, dbr+8
vmc2mc.x[npost] +=1
} else {j -= 1 }// if (debug_ == 1 ) {print "mc2mc"}
// if (killMC.contains(ngcell+nbcell+npost) == 1){ j += 1 print "dead MC"} // use if killing MC
}
// MC -> HC
for j=0, 1 {
Gauz3 = rdmc2hc.repick() // [-2:2]
if (y+Gauz3 > nhcell-1) {npost = y+Gauz3-nhcell} // changed code here to allow > 6 HCells
if (y+Gauz3 < 0) {npost = y+Gauz3+nhcell}
if ((y+Gauz3 > -1) && (y+Gauz3 < nhcell)) {npost = y+Gauz3}
dbr = rdsynb.repick()
if ((is_connected(HIPPCell[npost], MossyCell[i]) == 0) && (vmc2hc.x[npost] < 7) && (Gauz3 != 0)) {
nc_append(i+ngcell+nbcell, ngcell+nbcell+nmcell+npost, dbr+4, 0.2e-3, 3, 10) // Gcell[3] to Bcell[1]
// print npost, dbr+4
vmc2hc.x[npost] +=1
} else {
j -= 1
// if (debug_ == 1 ) {print y, Gauz3, "mc2hc"}
}
}
}
//******************************************************************************************
//**************HIPP Cell post synaptic connections ******************************
if (debug_ ==2 ) { print "HIPP Cell post synaptic connections"}
for i=0, nhcell-1 {
for j=0, 159 { // NOTE number of connections explicitly coded here!
Gauz3 = rdhc2gc.repick() // [-130:130]
//print Gauz3
if (i*83+41+Gauz3 > 499) {npost = i*83+41+Gauz3-500 } // NOTE: 500 is explicitly coded here!
if (i*83+41+Gauz3 < 0) {npost = i*83+41+Gauz3+500}
if ((i*83+41+Gauz3 > -1) && (i*83+41+Gauz3 < 500)) {npost = i*83+41+Gauz3}
//print npost
dbr = rdsyna.repick()
if (nhcell != 6) {max_conn = 5} else {max_conn = 3}
if ((is_connected(GranuleCell[npost], HIPPCell[i]) == 0) && (vhc2gc.x[npost] < max_conn)) { // NOTE < 3 coded explicitly here! -> change to 5
nc_append(i+ngcell+nbcell+nmcell, npost, dbr+4, 0.5e-3*scale_HC2GC_strength/100., 1.6, 10) // Hcell to Gcell
vhc2gc.x[npost] +=1
if (debug_ == 1 ) {print i, npost, dbr+4}
} else {
j -= 1
if (debug_ == 1 ) {print "HC2GC"}
}
}
for j=0, 3 {
Gauz3 = rdhc2bc.repick() // [-2:2]
if (i+Gauz3 > nbcell-1) {npost = i+Gauz3-nbcell}
if (i+Gauz3 < 0) {npost = i+Gauz3+nbcell}
if ((i+Gauz3 > -1) && (i+Gauz3 < nbcell)) {npost = i+Gauz3}
dbr = rdsyna.repick()
if ((is_connected(BasketCell[npost], HIPPCell[i]) == 0) && (vhc2bc.x[npost] < 5)) { // NOTE < 5 coded explicitly here!
nc_append(i+ngcell+nbcell+nmcell, npost+ngcell, dbr+10, 0.5e-3, 1.6, 10) // Hcell to Bcell
if (debug_ == 1 ) {print npost, dbr+10}
vhc2bc.x[npost] += 1
} else {
j -= 1
if (debug_ == 1 ) {print "HC2BC"}
}
}
for j=0, 3 {
Gauz3 = rdhc2mc.repick() // [-2:2]
//print Gauz3
if (i*2+2+Gauz3 > 14) {npost = i*2+2+Gauz3-15 }
if (i*2+2+Gauz3 < 0) {npost = i*2+2+Gauz3+15}
if ((i*2+2+Gauz3 >-1) && (i*2+2+Gauz3 < 15)) {npost = i*2+2+Gauz3}
//print npost
dbr = rdsynb.repick()
// if ((is_connected(MossyCell[npost], HIPPCell[i]) == 0) && (vhc2mc.x[npost] < 2) && (killMC.contains(ngcell+nbcell+npost) == 0)) //use if killing MC
if (nhcell != 6) {max_conn = 4} else {max_conn = 2}
if ((is_connected(MossyCell[npost], HIPPCell[i]) == 0) && (vhc2mc.x[npost] < max_conn)) { // NOTE < 2 coded explicitly here!! => changed to 4
nc_append(i+ngcell+nbcell+nmcell, npost+ngcell+nbcell, dbr+13, 1.5e-3, 1, 10) // Hcell to Mcell
//if (debug_ == 1 ) {print npost, dbr+13}
vhc2mc.x[npost] += 1
} else {
j -= 1
if (debug_ == 1 ) {print "HC2MC"}
}
// if (killMC.contains(ngcell+nbcell+npost) == 1){ j += 1 print "dead MC"} //use if killing MC
}
}
if (debug_ == 2 ) {print "Connections drawn -> Start Simulation"}
}
objref VmT // vector of time points
objref VmMat[cells.count] // for each cell vector of Membrane potential
VmT = new Vector()
for i=0, cells.count-1 {
VmMat[i] = new Vector()
}
proc VecMx() { local i // procedure called in every time step to add voltage of recorded cell at time t
VmT.append(t)
for i=0, (ngcell+nbcell +nmcell +nhcell -1) {
VmMat[i].append( cells.object[i].soma.v(0.5))
}
}
proc ConvVT2Sp() { local i, j, max_id
max_id = npp
if ( print_template ==1) { max_id = npp }
if ( print_Raster ==1) { max_id = ngcell+nbcell +nmcell +nhcell-1 }
for i=0, (max_id) {
Spike[i].spikebin(VmMat[i], 0) // Used to make a binary version of a spike train.
// <data> is a vector of membrane potential.
// <thresh> is the voltage threshold for spike detection.
// <v> is set to all zeros except at the onset of spikes (the first dt which the spike crosses threshold)
}
}
// ********* Print files **************************
err_ = load_file("printfile.hoc")
//################################################################################################
// Define Random Generator for each PPStim object...
print "Define Random Generator for each PPStim object..."
// nslist = new List()
rslist = new List()
rs_noiselist = new List()
// ####### Input Window (box)
random_stream_offset_ = PP_rate_ * 1000
for i=0, npp-1 {
rs = new RandomStream(i) // (stream nr)
rslist.append(rs)
PPSt[i].pp.noiseFromRandom(rs.r)
rs.r.uniform(0,1) // use uniform distribution
rs.start()
}
for i=0, npp-1 {
rs = new RandomStream(i+npp) // (stream nr)
rs_noiselist.append(rs)
PPSt_noise[i].pp.noiseFromRandom(rs.r)
rs.r.uniform(0,1) // use uniform distribution
rs.start()
}
aacell pp_start = new NetStim125(.5) // artificial puls to PPSt[i].pp in order to become active...
pp_start.interval = 1e-5
pp_start.number = 1
pp_start.start = 0
pp_start.forcestop = 1.
pp_start.noise = 1
rs = new RandomStream(npp) // assign Random Generator to init..
return_val_ = rslist.append(rs) // save returned value in return_val_ to suppress screen output
return_val_ = pp_start.noiseFromRandom(rs.r)
return_val_ = rs.r.negexp(1)
return_val_ = rs.start()
for i=0, npp-1 {
netcon_d[i] = new NetCon( pp_start,PPSt[i].pp) // each PPSt[i].pp needs to receive a netevent to become active...
netcon_d[i].weight = 10.
netcon_d[i].delay = 0.001
netcon_d[i].threshold = 10.
}
print "Finished Random Generator for each PPStim object"
objref dataset, data_
proc read_custominit() { local i
data_ = new Vector()
dataset = new File()
err_ = dataset.ropen(init_state_fname_)
if (err_ ==0) {
print "IO code:",err_
print "Initial State File not found!"
quit()
}
err_ = data_.scanf(dataset)
err_ = dataset.close()
cnt = 0
for i = 0, ngcell -1 {
forsec Gcell[i].all { // iterate over all sections of this cell
for (x,0) { // iterate over internal nodes
v = data_.x[cnt] // f() returns the desired initial potential
if (v>-30.) { v = -30. } // cut off initial potential -> no spike artifacts ...
cnt += 1
}
}
}
for i = 0, nbcell -1 {
forsec Bcell[i].all { // iterate over all sections of this cell
for (x,0) { // iterate over internal nodes
v = data_.x[cnt]
cnt += 1
}
}
}
for i = 0, nmcell -1 {
forsec Mcell[i].all { // iterate over all sections of this cell
for (x,0) { // iterate over internal nodes
v = data_.x[cnt]
cnt += 1
}
}
}
for i = 0, nhcell -1 {
forsec Hcell[i].all { // iterate over all sections of this cell
for (x,0) { // iterate over internal nodes
v = data_.x[cnt]
cnt += 1
}
}
}
finitialize()
}
proc init() { local dtsav, temp, secsav //localobj ns, rs // init with a int value =! 100 reset not all random number generators equally...
v_init = -77.71 //-70.45
finitialize(v_init)
t = -1000 // negative t step to initialize to steady state
dtsav = dt // original simulation step size
dt= 10 // larger step size
// if cvode is on, turn it off to do large fixed step
temp= cvode.active()
if (temp!=0) {cvode.active(0)}
while(t<-100) {
fadvance()
}
//restore cvode if reqd
if (temp!=0) {cvode.active(1)}
dt = dtsav
t = 0
if (cvode.active()){
cvode.re_init()
}else{
fcurrent()
}
// restart number generators for PPStim
t = 0
finitialize(v_init)
trial_old = trial
if (debug_ ==2 ) { print "number of Poisson inputs:\t",rslist.count()}
if ($1 == 100) {
for i = 0, rslist.count()-1 rslist.o(i).start()
// for j=0, 5 print j, rslist.o(0).r.repick
if (debug_ ==2 ) {print "reset all Poisson inputs with the correct seed"}
} else {
trial = trial + 1
for i = 0, int((rslist.count()-1)*(1-$1/100.)) rslist.o(i).start()
if (debug_ ==2 ) {print "reset ",int((rslist.count()-1)*(1-$1/100.))," Poisson inputs with a different seed"}
trial = trial_old
for i = int((rslist.count()-1)*(1-$1/100.)), rslist.count()-1 rslist.o(i).start()
if (debug_ ==2 ) {print "reset ", rslist.count()-int((rslist.count()-1)*(1-$1/100.))," Poisson inputs with the correct seed"}
}
// init noise Generators (with seed trial_noise ...)
trial = trial_noise
for i = 0, rs_noiselist.count()-1 rs_noiselist.o(i).start()
trial = trial_old
if (debug_ ==2 ) {print "reset all Poisson noise inputs with the trial_noise seed"}
VmT = new Vector()
for i=0, cells.count-1 {
VmMat[i] = new Vector()
}
for i=0, (ngcell+nbcell +nmcell +nhcell -1) {
Spike[i] = new Vector() // vector of spikes for each cell
Spike_times[i] = new Vector() // vector of spikes for each cell
}
t = 0
finitialize(v_init) // finalize at the end
if (init_from_file_ ==1 ) {
if (debug_ == 2) {print "Read in Init State of Network"}
read_custominit() // read in initial states from file
}
// custominit()
frecord_init()
// ### select input pattern....
for i=PP_box_startid_,PP_box_startid_+PP_box_nr_-1 {
PPSt[i].pp.status = 1
PPSt[i].pp.start = PP_box_start_
PPSt[i].pp.forcestop = PP_box_stop_
PPSt[i].pp.nspk = PP_box_nspk_
}
}
proc run(){
initNet_GC()
initNet_BC()
initNet_rest()
for (PP_box_startid_=0; PP_box_startid_<=12;PP_box_startid_+=1) {
saveNet()
init(100)
if (print_Vtrace == 1) { sMatrix_init() } // file header voltage output file
while (t<tstop) {
fadvance()
if (print_Vtrace == 1) { sMatrix() } // print Voltage trace
if ((print_Raster == 1) || (print_template == 1)) { VecMx() } // record voltage trace for all runs needed
}
if ((print_Raster == 1) || (print_template == 1)) { ConvVT2Sp() }
if (print_Raster == 1) {SpkMx()} // Output Spike Times to file!
if (print_template == 1) {SpkMx_template()} // Output Spike Times to file!!
if (print_stim_==1) {write_stimIn()}
}
if (init_to_file_==1) {write_customstate() } // save state of network for re-init
}
run()