// The axon is the venerable Mainen 1996 axon. Not perfect, yet very stable.
// Modified by Tuomo Maki-Marttunen, 2016, basing on TTC.hoc by Etay Hay, Dendritic excitability and gain control in recurrent cortical microcircuits (Hay and Segev, 2014, Cerebral Cortex)
// Modified by Tuomo Maki-Marttunen, 2020, added functions that distribute synapses to particular sections
create a_soma
n_axon_seg = 2
create iseg,hill,myelin[n_axon_seg],node[n_axon_seg]
objref preconlist,synLocList, rList, preTrainList
objref rd1
objref sref, fih
objref synlist
proc create_axon() {
n_axon_seg = 2
create iseg,hill,myelin[n_axon_seg],node[n_axon_seg]
a_soma {
equiv_diam = sqrt(area(.5)/(4*PI))
}
iseg { // initial segment between hillock + myelin
L = 15
nseg = 5
diam = equiv_diam/10 // see Sloper and Powell 1982, Fig.71
}
hill {
L = 20
nseg = 5
diam(0:1) = 2*iseg.diam:iseg.diam
}
// construct myelinated axon with nodes of ranvier
for i=0,n_axon_seg-1 {
myelin[i] { // myelin element
nseg = 2
L = 100
diam = iseg.diam
}
node[i] { // nodes of Ranvier
nseg = 1
L = 1.0
diam = iseg.diam*.75 // nodes are thinner than axon
}
}
a_soma connect hill(0), 0.5
hill connect iseg(0), 1
iseg connect node[0](0), 1
node[0] connect myelin[0](0), 1
for i=0,n_axon_seg-2 {
myelin[i] connect node[i+1](0), 1
node[i+1] connect myelin[i+1](0), 1
}
}
// --------
// Spines
// --------
proc add_spines() { local a
forsec "dend" {
cm*=$1
g_pas*=$1
}
forsec "apic" {
cm*=$1
g_pas*=$1
}
}
proc load_3dcell() {
// $s1 filename
forall delete_section()
xopen($s1)
access a_soma
// make sure no compartments exceed 50 uM length
forall {
diam_save = diam
n = int(L/20)
if (n==0) n=1
nseg = n
if (n3d() == 0) diam = diam_save
}
create_axon()
init_cell()
}
proc init_cell() {
// passive
forall {
insert pas
Ra = ra
cm = c_m
g_pas = 1/rm
e_pas = epas_sim
}
// exceptions along the axon
forsec "myelin" cm = cm_myelin
forsec "node" g_pas = g_pas_node
forall {
insert iA
insert kslow
insert na
insert iH
insert cah
insert car
insert cad
insert bk
insert sk
}
forall if(ismembrane("k_ion")) ek = Ek
forall if(ismembrane("na_ion")) ena = Ena
forall eh=-33
forall if(ismembrane("ca_ion")) {
vshift_ca = 0
}
synlist = new List()
preconlist = new List()
preTrainList = new List()
lengthA = 0
lengthB = 0
forsec "apic" {
lengthA = lengthA + L
}
forsec "dend" {
lengthB = lengthB + L
}
pA = lengthA/(lengthA + lengthB)
}
proc density() {
forall {
Ra = ra
cm = c_m
g_pas = 1/rm
e_pas = epas_sim
vshiftm_na= na_shift1
vshifth_na= na_shift2
taum_scale_na=na_taum_scale
tauh_scale_na=na_tauh_scale
q10_iH=ih_q10
shift_cah=cah_shift
shifth_cah=cah_shifth
shift_car=car_shift
shifth_car=car_shifth
qm_car=car_qm
}
add_spines(2)
// exceptions along the axon
forsec "myelin" cm = cm_myelin
forsec "node" g_pas = g_pas_node
forall vshift_na = -5
forsec "myelin" {
gbar_na = gna_soma
gbar_kslow =gkslow_beta
gbar_iA = gka_beta
}
hill {
gbar_na = gna_node
vshiftm_na=7
vshifth_na=3
gbar_iA = gka_node
gbar_kslow = gkslow_node
}
iseg {
gbar_na = gna_node
vshiftm_na=7
vshifth_na=3
gbar_iA = gka_node
gbar_kslow = gkslow_node
}
forsec "node" {
gbar_na = gna_node
vshiftm_na=7
vshifth_na=3
gbar_iA = gka_node
gbar_kslow = gkslow_node
}
// we put in all of teh dendrites the same values as in the a_soma and then overwrite the appical values
forsec "dend" {
gbar_na = gna_soma
gbar_kslow = gkslow_start+gkslow_beta
gbar_iA = gka_start+gka_beta
gbar_iH=gih_start
pbar_car=pcar_soma
pbar_cah=pcah_soma
gbar_sk=gsk_dend
gbar_bk=gbk_dend
}
a_soma {
gbar_na = gna_soma
gbar_kslow = gkslow_start+gkslow_beta
gbar_iA = gka_start+gka_beta
gbar_iH=gih_start
pbar_car=pcar_soma
pbar_cah=pcah_soma
gbar_sk=gsk_soma
gbar_bk=gbk_soma
}
a_soma distance(0,0.5)
forsec ApicalDendSectionName {
dist1=distance(0)
dist2=distance(1)
gbar_na = gna_api
if (dist1<dist_na){
gbar_na(0:1) = (gna_soma+dist1*(gna_api-gna_soma)/dist_na):(gna_soma+dist2*(gna_api-gna_soma)/dist_na)
}
pbar_cah = pcah_api
if (dist1<dist_cah){
pbar_cah(0:1) = (pcah_soma+dist1*(pcah_api-pcah_soma)/dist_cah):(pcah_soma+dist2*(pcah_api-pcah_soma)/dist_cah)
}
pbar_car = pcar_api
if (dist1<dist_car){
pbar_car(0:1) = (pcar_soma+dist1*(pcar_api-pcar_soma)/dist_car):(pcar_soma+dist2*(pcar_api-pcar_soma)/dist_car)
//(pcar_soma+dist1*(pcar_api-pcar_soma)/dist_car):(pcar_soma+dist2*(pcar_api-pcar_soma)/dist_car)
}
gbar_iA(0:1) = (gka_start+gka_beta*exp(gka_alpha*dist1)):(gka_start+gka_beta*exp(gka_alpha*dist2))
gbar_kslow(0:1) = (gkslow_start+gkslow_beta*exp(gkslow_alpha*dist1)):(gkslow_start+gkslow_beta*exp(gkslow_alpha*dist2))
gbar_iH(0:1) = (gih_start+gih_end/(1+exp(gih_alpha*(dist1-gih_x2)))):(gih_start+gih_end/(1+exp(gih_alpha*(dist2-gih_x2))))
gbar_bk = gbk_dend
if (dist1<dist_bk){
gbar_bk(0:1) = (gbk_soma+dist1*(gbk_dend-gbk_soma)/dist_bk):(gbk_soma+dist2*(gbk_dend-gbk_soma)/dist_bk)
}
gbar_sk = gsk_dend
if (dist1<dist_sk){
gbar_sk(0:1) = (gsk_soma+dist1*(gsk_dend-gsk_soma)/dist_sk):(gsk_soma+dist2*(gsk_dend-gsk_soma)/dist_sk)
}
}
forsec ApicalDendSectionName {
if (pcah_soma<pcah_api){
if (pbar_cah>pcah_api){
pbar_cah=pcah_api
}
}
if (pcah_soma>pcah_api){
if (pbar_cah<pcah_api){
pbar_cah=pcah_api
}
}
if (pcar_soma<pcar_api){
if (pbar_car>pcar_api){
pbar_car=pcar_api
}
}
if (pcar_soma>pcar_api){
if (pbar_car<pcar_api){
pbar_car=pcar_api
}
}
if (gbar_na<gna_api) gbar_na=gna_api
if (gbar_bk<gbk_dend) gbar_bk=gbk_dend
if (gbar_sk<gsk_dend) gbar_sk=gsk_dend
}
}
// Copied from Hay model
// $s1 section
func getLongestBranch(){local maxL,d localobj distallist,sref
//Commented out: Always calculate the distance with respect to soma 0.5
//sprint(tstr,"%s distance()",$s1)
//execute(tstr,this)
//if(0==strcmp($s1,"axon")){
// sprint(tstr,"%s[0] distance(1)",$s1)
// execute(tstr,this)
//}
maxL = 0
d = 0
distallist = new SectionList()
forsec $s1 {
sref = new SectionRef()
if (sref.nchild==0) distallist.append()
}
forsec distallist{
d = distance(1)
if(maxL<d) maxL = d
}
// for the soma case
if (maxL == 0) {
$s1 {
maxL = L
}
}
return maxL
}
//Copied from Hay model
// $s1 section
// $2 distance x in micrometers
// return list of [1,2] vectors - of the appropriate section and the location in each vector
obfunc locateSites() {local maxL,site,d0,d1,siteX,i localobj vv,ll
ll = new List()
//sprint(tstr,"%s distance()",$s1)
//execute(tstr,this)
//if(0==strcmp($s1,"axon")){
// sprint(tstr,"%s[0] distance(1)",$s1)
// execute(tstr,this)
//}
maxL = getLongestBranch($s1)
site = $2
i = 0
forsec $s1 {
if (distance(0) < distance(1)) {
d0 = distance(0)
d1 = distance(1)
} else {
d1 = distance(0)
d0 = distance(1)
}
if (site <= d1 && site >= d0) {
siteX = (site-d0)/(d1-d0)
secNum = i
vv = new Vector()
ll.append(vv.append(secNum,siteX))
}
i = i+1
}
return ll
}
proc setparameters(){
EsynConductance = $1
IsynConductance = $2
NsynE = $3
NsynI = $4
}
proc initRand() {
rList = new List() //for stochastic synapses
rd1 = new Random($1) // unique to this TTC
rd1.uniform(0,1)
}
double siteVec[2]
proc distributeSyn() {local sitenum,syni,preconi,i localobj sl,nilstim
strdef treename,cmd2
for(i=0;i<(NsynE+NsynI);i+=1){
if (i==(NsynE+NsynI)/4-1) { print "25% of synapses done" }
if (i==(NsynE+NsynI)/2-1) { print "50% of synapses done" }
if (i==3*(NsynE+NsynI)/4-1) { print "75% of synapses done" }
if (rd1.repick()<pA){
treename = "apic"
} else {
treename = "dend"
}
sl = locateSites(treename,rd1.repick()*getLongestBranch(treename))
sitenum = int((sl.count()-1)*rd1.repick())
siteVec[0] = sl.o[sitenum].x[0]
siteVec[1] = sl.o[sitenum].x[1]
sprint(cmd2,"access %s[siteVec[0]]",treename)
execute(cmd2)
//print "distributesSyn(): ",cmd2
//print "siteVec = ",siteVec[0],siteVec[1]
sprint(cmd2,"%s[siteVec[0]] sref = new SectionRef()",treename)
//print "distributesSyn(): ",cmd2
execute(cmd2)
if (i<NsynE){
sref {
synlist.append(new ProbAMPANMDA_EMST(siteVec[1]))
syni = synlist.count()-1 //synapse index
rList.append(new Random(int(1000000*rd1.repick())))
rList.o[syni].negexp(1)
synlist.o[syni].setRNG(rList.o[syni])
synlist.o[syni].tau_r_AMPA = 0.3
synlist.o[syni].tau_d_AMPA = 3
synlist.o[syni].tau_r_NMDA = 2
synlist.o[syni].tau_d_NMDA = 65
synlist.o[syni].e = 0
synlist.o[syni].Dep = 800
synlist.o[syni].Fac = 0
synlist.o[syni].Use = 0.6
synlist.o[syni].u0 = 0
synlist.o[syni].gmax = EsynConductance
preconlist.append(new NetCon(nilstim, synlist.o[syni]))
preconi = preconlist.count()-1 //connection index
preconlist.o[preconi].weight = 1
preconlist.o[preconi].delay = 0
}
} else {
sref {
synlist.append(new ProbUDFsyn2(siteVec[1]))
syni = synlist.count()-1 //synapse index
rList.append(new Random(int(1000000*rd1.repick())))
rList.o[syni].negexp(1)
synlist.o[syni].setRNG(rList.o[syni])
synlist.o[syni].tau_r = 1
synlist.o[syni].tau_d = 20
synlist.o[syni].e = -80
synlist.o[syni].Dep = 800
synlist.o[syni].Fac = 0
synlist.o[syni].Use = 0.25
synlist.o[syni].u0 = 0
synlist.o[syni].gmax = IsynConductance
preconlist.append(new NetCon(nilstim, synlist.o[syni]))
preconi = preconlist.count()-1 //connection index
preconlist.o[preconi].weight = 1
preconlist.o[preconi].delay = 0
}
}
if (i==(NsynE+NsynI)-1) { print "synapses done" }
}
}
proc distributeSynBranch() {local syni,preconi,i,itip,NsynsBranch,x localobj nilstim
strdef treename,cmd2
NsynsBranch = $1
treename = $s2
itip = $3
x = rd1.repick()
for(i=0;i<(NsynsBranch);i+=1){
sprint(cmd2,"access %s[%g]",treename,itip)
execute(cmd2)
sprint(cmd2,"%s[%g] sref = new SectionRef()",treename,itip)
execute(cmd2)
sref {
synlist.append(new ProbAMPANMDA_EMST(x))
syni = synlist.count()-1 //synapse index
rList.append(new Random(int(1000000*rd1.repick())))
rList.o[syni].negexp(1)
synlist.o[syni].setRNG(rList.o[syni])
synlist.o[syni].tau_r_AMPA = 0.3
synlist.o[syni].tau_d_AMPA = 3
synlist.o[syni].tau_r_NMDA = 2
synlist.o[syni].tau_d_NMDA = 65
synlist.o[syni].e = 0
synlist.o[syni].Dep = 800
synlist.o[syni].Fac = 0
synlist.o[syni].Use = 0.6
synlist.o[syni].u0 = 0
synlist.o[syni].gmax = EsynConductance
preconlist.append(new NetCon(nilstim, synlist.o[syni]))
preconi = preconlist.count()-1 //connection index
preconlist.o[preconi].weight = 1
preconlist.o[preconi].delay = 0
}
}
}
proc setsynapseconductances(){local syni, i
for(i=0;i<(NsynE+NsynI);i+=1){
if (i<NsynE){
synlist.o[i].gmax = $1
} else {
synlist.o[i].gmax = $2
}
}
}
//$o1 list of vectors
proc setpretrains(){local j
for(j=0;j<(NsynE+NsynI);j+=1){
preTrainList.append($o1.o[j])
}
}
proc setpretrainsbranch(){local j,NsynsBranch
NsynsBranch = $o1.count()
for(j=0;j<(NsynsBranch);j+=1){
preTrainList.append($o1.o[j])
}
}
proc queuePreTrains(){
fih = new FInitializeHandler("initPreSynTrain()")
}
// sets presynaptic spike events
proc initPreSynTrain(){local ti,si
for(ti=0;ti<preTrainList.count();ti+=1){
for(si=0;si<preTrainList.o[ti].size();si+=1){
preconlist.o[ti].event(preTrainList.o[ti].x[si])
}
}
}