// Author: Etay Hay, 2011
// Models of Neocortical Layer 5b Pyramidal Cells Capturing a Wide Range of
// Dendritic and Perisomatic Active Properties
// (Hay et al., PLoS Computational Biology, 2011)
//
// Modified by Tuomo Maki-Marttunen, 2015, basing on TTC.hoc by Etay Hay, Dendritic excitability and gain control in recurrent cortical microcircuits (Hay and Segev, 2014, Cerebral Cortex)
// Template for models of L5 Pyramidal Cell
begintemplate L5PCtemplate
public init
public synlist
public locateSites, getLongestBranch, distributeSyn, distributeSynGivenDist, distributeGABASynGivenDist, setparameters, initPreSynTrain, queuePreTrains, setpretrains, setpretrainsbranch, initRand
public soma, dend, apic, axon, getAbsSecIndex
public all, somatic, apical, axonal, basal, nSecSoma, nSecApical, nSecBasal, nSecAxonal, nSecAll, nSecAxonalOrig, SecSyn, distribute_channels
objref SecSyn, this
objref all, somatic, apical, axonal, basal
strdef tstr
objref preconlist,synLocList, rList, preTrainList
objref rd1
objref sref, fih
objref synlist
//$s1 - morphology file name
proc init() {localobj nl,import
all = new SectionList()
somatic = new SectionList()
basal = new SectionList()
apical = new SectionList()
axonal = new SectionList()
forall delete_section()
nl = new Import3d_Neurolucida3()
nl.quiet = 0
nl.input($s1)
import = new Import3d_GUI(nl, 0)
import.instantiate(this)
geom_nseg()
biophys()
forsec this.all {
if(diam == 0){
diam = 1
printf("Error : Morphology problem with section [%s] 0 diam \n", secname())
}
}
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)
}
create soma[1], dend[1], apic[1], axon[1]
proc geom() {
}
proc geom_nseg() {local nSec, L1, L2, D1, D2, nSeg1, nSeg2
soma area(.5) // make sure diam reflects 3d points
nSec = 0
forsec all {
nseg = 1 + 2*int(L/40)
nSec = nSec + 1
}
nSecAll = nSec
print nSecAll
nSec = 0
forsec somatic { nSec = nSec + 1}
nSecSoma = nSec
print nSecSoma
nSec = 0
forsec apical { nSec = nSec + 1}
nSecApical= nSec
print nSecApical
nSec = 0
forsec basal { nSec = nSec + 1}
nSecBasal = nSec
print nSecBasal
nSec = 0
forsec axonal { nSec = nSec + 1}
nSecAxonalOrig = nSecAxonal = nSec
print nSecAxonal
}
proc biophys() {localobj bp
delete_axon()
area(0.5)
distance()
access soma
bp = new L5PCbiophys()
bp.biophys(this)
}
// deleting axon, keeping only first 60 micrometers
proc delete_axon(){
forsec axonal{delete_section()}
create axon[2]
access axon[0]{
L= 30
diam = 1
nseg = 1+2*int(L/40)
all.append()
axonal.append()
}
access axon[1]{
L= 30
diam = 1
nseg = 1+2*int(L/40)
all.append()
axonal.append()
}
nSecAxonal = 2
connect axon(0), soma(0.5)
connect axon[1](0), axon[0](1)
access soma
}
proc distribute_channels() {local dist,val,base,maxLength
base = $8
soma distance()
maxLength = getLongestBranch($s1)
forsec $s1 {
if(0==strcmp($s2,"Ra")){
Ra = $8
} else {
for(x) {
if ($3==3) {
dist = distance(x)
} else {
dist = distance(x)/maxLength
}
val = calculate_distribution($3,dist,$4,$5,$6,$7,$8)
sprint(tstr,"%s(%-5.10f) = %-5.10f",$s2,x,val)
execute(tstr)
}
}
}
}
// $1 is the distribution type:
// 0 linear, 1 sigmoid, 2 exponential
// 3 step for absolute distance (in microns)
func calculate_distribution() {local value
if ($1==0) {value = $3 + $2*$4}
if ($1==1) {value = $3 + ($4/(1+exp(($2-$5)/$6)))}
if ($1==2) {value = $3 + $6*exp($4*($2-$5))}
if ($1==3) {
if (($2 > $5) && ($2 < $6)) {
value = $3
} else {
value = $4
}
}
value = value*$7
return value
}
// $s1 section
func getLongestBranch(){local maxL,d localobj distallist,sref
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
}
// $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
}
func getAbsSecIndex(){ local nAbsInd, index localobj str,strObj
strObj = new StringFunctions()
str = new String()
nAbsInd = 0
index = 0
if(strObj.substr($s1, "soma") > 0) {
strObj.tail($s1, "soma", str.s)
if(sscanf(str.s, "%*c%d", &index) < 0) {
index = 0
}
nAbsInd = index
}else if (strObj.substr($s1, "axon") >0) {
strObj.tail($s1, "axon", str.s)
if(sscanf(str.s, "%*c%d", &index) < 0) {
index = 0
}
nAbsInd = nSecSoma + index
}else if (strObj.substr($s1, "dend") >0) {
strObj.tail($s1, "dend", str.s)
if(sscanf(str.s, "%*c%d", &index) < 0) {
index = 0
}
nAbsInd = nSecSoma + nSecAxonalOrig + index
}else if (strObj.substr($s1, "apic") > 0) {
strObj.tail($s1, "apic", str.s)
if(sscanf(str.s, "%*c%d", &index) < 0) {
index = 0
}
nAbsInd = nSecSoma + nSecAxonalOrig + nSecBasal + index
}
return nAbsInd
}
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);i+=1){
if (i==(NsynE)/4-1) { print "25% of synapses done" }
if (i==(NsynE)/2-1) { print "50% of synapses done" }
if (i==3*(NsynE)/4-1) { print "75% of synapses done" }
if (i==(NsynE)-1) { print "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,this)
sprint(cmd2,"%s[siteVec[0]] sref = new SectionRef()",treename)
execute(cmd2,this)
if (i<NsynE){
sref {
synlist.append(new AMPANMDA(siteVec[1]))
syni = synlist.count()-1 //synapse index
synlist.o[syni].E_Glu = 0
synlist.o[syni].gAMPAmax = EsynConductance
synlist.o[syni].gNMDAmax = EsynConductance
synlist.o[syni].MgCon = 0.7
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 distributeSynGivenDist() {local sitenum,syni,preconi,i localobj sl,nilstim
strdef treename,cmd2
treename = $s1
if($3 > getLongestBranch(treename)) {
print "Warning: dist2 = ",$3," is larger than longest branch"
}
if($2 > $3) {
print "Warning: dist1 = ",$2," > dist2 = ",$3
}
for(i=0;i<(NsynE);i+=1){
if (i==(NsynE)/4-1) { print "25% of synapses done" }
if (i==(NsynE)/2-1) { print "50% of synapses done" }
if (i==3*(NsynE)/4-1) { print "75% of synapses done" }
if (i==(NsynE)-1) { print "synapses done" }
sl = locateSites(treename,$2 + ($3 - $2)*rd1.repick())
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,this)
sprint(cmd2,"%s[siteVec[0]] sref = new SectionRef()",treename)
execute(cmd2,this)
if (i<NsynE){
sref {
synlist.append(new AMPANMDA(siteVec[1]))
syni = synlist.count()-1 //synapse index
synlist.o[syni].E_Glu = 0
synlist.o[syni].gAMPAmax = EsynConductance
synlist.o[syni].gNMDAmax = EsynConductance
synlist.o[syni].MgCon = 0.7
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 distributeGABASynGivenDist() {local sitenum,syni,preconi,i localobj sl,nilstim
strdef treename,cmd2
treename = $s1
if($3 > getLongestBranch(treename)) {
print "Warning: dist2 = ",$3," is larger than longest branch"
}
if($2 > $3) {
print "Warning: dist1 = ",$2," > dist2 = ",$3
}
for(i=0;i<(NsynI);i+=1){
if (i==(NsynI)/4-1) { print "25% of synapses done" }
if (i==(NsynI)/2-1) { print "50% of synapses done" }
if (i==3*(NsynI)/4-1) { print "75% of synapses done" }
if (i==(NsynI)-1) { print "synapses done" }
sl = locateSites(treename,$2 + ($3 - $2)*rd1.repick())
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,this)
sprint(cmd2,"%s[siteVec[0]] sref = new SectionRef()",treename)
execute(cmd2,this)
if (i<NsynI){
sref {
synlist.append(new GABA(siteVec[1]))
syni = synlist.count()-1 //synapse index
synlist.o[syni].E_Cl = -80
synlist.o[syni].gGABAmax = 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
}
}
}
}
//$o1 list of vectors (assume synlist has first NsynE exc. synapses and then NsynI inh. synapses)
proc setpretrains(){local j
preTrainList = new List()
for(j=0;j<(NsynE);j+=1){
preTrainList.append($o1.o[j])
}
for(j=0;j<(NsynI);j+=1){
preTrainList.append($o1.o[NsynE+j])
}
}
proc queuePreTrains(){
fih = new FInitializeHandler("initPreSynTrain()",this)
}
// 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])
}
}
}
endtemplate L5PCtemplate