load_file("nrngui.hoc")
/* simulation of neuronal model with varying complexity (set nr 'branch_levels').
A constant apical trunk diameter is used and the rall diameter rule is applied
license:
author: Luuk van der Velden, University of Amsterdam
paper: Altered dendritic complexity affects firing properties of cortical layer 2/3 pyramidal neurons in mice lacking the 5-HT3A receptor
journal of neurophysiology, 2012
no warranties :)
run as: 'nrngui filename' (easier on the mechanism loads)
*/
// main experimental parameters and parameter initialization
max_branch_levels = 7 // maximum nr. of branch levels (from 0 to 6)
branch_levels = 1 // nr of branch levels in dendrite (if 0, only apical dendrite is present)
begin_branch_diam = 2.5 // diameter micrometer of apical dendrite
exponent = 3/2 // diameter rule exponent
max_seg_length = 30 // max length of segments (micrometer)
//TotalPump (total density of pump sites on the cell membrane (mol/cm2))
TotalPump_cadp = 0.3e-12
// TotalBuffer (total amount of buffer present)
TotalBuffer_cadp = 0.000
v_init = -70 // initial membrane potential mV (same as rest potential
NDEND = 0 // number of dendrites
NDEND2 = 0 // will be 1 if branch_levels is set to zero
// in that case an unconnected piece of dendrite is added to the model (treat exceptional case)
// compute the number of dendritic segments (apart from the apical segment)
for i=1, branch_levels {
NDEND=NDEND+(2^i)
}
NDEND2 = NDEND
// if branch_levels == 0, (NDEND==0) then treat exceptional case (add an unconnected piece of dendrite)
if (NDEND<1){
NDEND2=1 // one disconnected segment will be created
}
objref leveldiam
leveldiam = new Vector(branch_levels+1,0)
objref levellength
levellength = new Matrix(max_branch_levels+1,max_branch_levels+1)
// Topology
create soma, apical, dend[NDEND2], hillock
access soma
connect soma(0), hillock(1)
connect apical(0), soma(1)
// connect dendritic segments if branch_levels
if (branch_levels>0) {
for i = 0, branch_levels-1{
start_dend = 0
for k = 0, i{
start_dend=start_dend+2^k
}
start_dend=start_dend-1
for p=0, (2^(i+1))-1{
if (i<1){
connect dend[start_dend+p](0),apical(1)
} else{
connect dend[start_dend+p](0),dend[(start_dend-(2^i))+int(p/2)](1)
}
}
}
}
// Geometry
path_length_dendr = 250 // path length along dendrite
sc = 0.5 // scaling factor of length between branch levels
// data array of dendritic compartment lengths (branch level dependent)
levellength.x[0][0] = path_length_dendr
levellength.x[1][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[1][1] = path_length_dendr*sc
levellength.x[2][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[2][1] = path_length_dendr*sc - (path_length_dendr*(sc^2))
levellength.x[2][2] = path_length_dendr*(sc^2)
levellength.x[3][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[3][1] = path_length_dendr*sc - (path_length_dendr*(sc^2))
levellength.x[3][2] = (path_length_dendr*(sc^2))-(path_length_dendr*(sc^3))
levellength.x[3][3] = path_length_dendr*(sc^3)
levellength.x[4][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[4][1] = path_length_dendr*sc - (path_length_dendr*(sc^2))
levellength.x[4][2] = (path_length_dendr*(sc^2)) - (path_length_dendr*(sc^3))
levellength.x[4][3] = (path_length_dendr*(sc^3))- (path_length_dendr*(sc^4))
levellength.x[4][4] = (path_length_dendr*(sc^4))
levellength.x[5][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[5][1] = path_length_dendr*sc - (path_length_dendr*(sc^2))
levellength.x[5][2] = (path_length_dendr*(sc^2)) - (path_length_dendr*(sc^3))
levellength.x[5][3] = (path_length_dendr*(sc^3))- (path_length_dendr*(sc^4))
levellength.x[5][4] = (path_length_dendr*(sc^4))- (path_length_dendr*(sc^5))
levellength.x[5][5] = (path_length_dendr*(sc^5))
levellength.x[6][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[6][1] = path_length_dendr*sc - (path_length_dendr*(sc^2))
levellength.x[6][2] = (path_length_dendr*(sc^2)) - (path_length_dendr*(sc^3))
levellength.x[6][3] = (path_length_dendr*(sc^3))- (path_length_dendr*(sc^4))
levellength.x[6][4] = (path_length_dendr*(sc^4))- (path_length_dendr*(sc^5))
levellength.x[6][5] = (path_length_dendr*(sc^5)) - (path_length_dendr*(sc^6))
levellength.x[6][6] = (path_length_dendr*(sc^6))
levellength.x[7][0] = path_length_dendr - (path_length_dendr*sc)
levellength.x[7][1] = path_length_dendr*sc - (path_length_dendr*(sc^2))
levellength.x[7][2] = (path_length_dendr*(sc^2)) - (path_length_dendr*(sc^3))
levellength.x[7][3] = (path_length_dendr*(sc^3))- (path_length_dendr*(sc^4))
levellength.x[7][4] = (path_length_dendr*(sc^4))- (path_length_dendr*(sc^5))
levellength.x[7][5] = (path_length_dendr*(sc^5)) - (path_length_dendr*(sc^6))
levellength.x[7][6] = (path_length_dendr*(sc^6)) - (path_length_dendr*(sc^7))
levellength.x[7][7] = (path_length_dendr*(sc^7))
nr_end_branches = 2^(branch_levels) // number of end points / branches
// compute total length of dendritic tree
total_length = 0 // initial value for sum of total_length
for i = 0, branch_levels{
total_length = total_length+(levellength.x[branch_levels][i]*(2^i))
}
//compute DCI of dendritic trees
DCI = (nr_end_branches+(nr_end_branches*branch_levels))*total_length
print "dci"
print DCI
print "total length"
print total_length
// the Rall value to keep constant
diameter_constant = begin_branch_diam^exponent
// compute diameters for various branch levels (apply Rall rule)
for i=0, branch_levels {
leveldiam.x[i] = (diameter_constant/(2^i))^(1/exponent)
print leveldiam.x[i]
}
dend_L = path_length_dendr/(branch_levels+1)
//apical geometry
apical{
L = levellength.x[branch_levels][0]
diam = leveldiam.x[0]
nseg = int(levellength.x[branch_levels][0]/max_seg_length)+1
}
//set dendritic geometry (diameter, length, nr of segments
if (branch_levels>0){
for i=1, branch_levels{
start_dend = 0
for k = 0, i-1 {
start_dend=start_dend+2^k
}
start_dend=start_dend-1
for p=0, (2^(i))-1{
dend[start_dend+p]{
diam = leveldiam.x[i]
L = levellength.x[branch_levels][i]
nseg = int(levellength.x[branch_levels][i]/max_seg_length)+1
}
}
}
}
soma_L = 10 // somatic length (micrometer)
soma_diam = 15 // soma diameter (micrometer)
hillock_L = 2 // axon hillock length (micrometer)
hillock_diam = 2 // axon hillock diameter (micrometer)
// set soma geometric properties
soma{
L = soma_L
diam = soma_diam
nseg = 1
}
// set axon hillock geometric properties
hillock{
L = hillock_L
diam = hillock_diam
}
// Biophysics
ra = 150 // axial resistance
global_ra = ra
rm = 30000 // 30 kilo ohm
c_m = 0.75 // micro farad
// passive apical properties
apical {
insert pas
Ra = ra
cm = c_m
g_pas = 1/rm
e_pas = v_init
}
// passive soma properties
soma {
insert pas
Ra = ra
cm = c_m
g_pas = 1/rm
e_pas = v_init
}
celsius = 32 // degrees celsius of experiment
Ek = -90 // reversal potential potassium
Ena = 50 // reversal potential sodium
// dendritic conductances
gna = 1
gca = 2
gkca = 10
gkm = 0.1
// soma conductances
gkv_soma= 800
gna_soma= 300
gca_soma = 0.1
gkca_soma = 0
gkm_soma = 0
//hillock conductances
gkv_hillock = 2000
gna_hillock = 30000
// somatic active channels
soma{
insert na gbar_na = gna_soma
insert kv gbar_kv = gkv_soma
insert ca gbar_ca = gca_soma
insert cadp
}
// apical active channels
apical{
insert na gbar_na = gna
insert km gbar_km = gkm
insert ca gbar_ca = gca
insert cadp
insert kca gbar_kca = gkca
}
// axon hilock active channels
hillock{
insert na gbar_na = gna_hillock
insert kv gbar_kv = gkv_hillock
}
// insert passive and active properties in dendritic compartments
for i=0, NDEND-1{
dend[i]{
insert pas
Ra = ra
cm = c_m
g_pas = 1/rm
e_pas = v_init
insert na gbar_na = gna
insert km gbar_km = gkm
insert ca gbar_ca = gca
insert cadp
insert kca gbar_kca = gkca
}
}
// set general extracellular properties for calcium
forall if(ismembrane("ca_ion")) {
eca = 140
ion_style("ca_ion",0,1,0,0,0)
vshift_ca = 0
}
// Instrumentation
// current input
objref stim
soma stim = new IClamp(0.5)
stim.amp = .02
stim.del = 50
stim.dur = 1000
//Simulation control
tstop = 2000
steps_per_ms = 40
dt = 0.025
nr_cells=int(tstop/dt)+2
// adjust output data size to number of branch levels
if(branch_levels<5){
nvars = 5 // no thin dendrite calcium dynamics
}else{
nvars = 6
}
double data[nvars][nr_cells] // create empty array for data storage
// Plotting
// graph of somatic membrane potential
objref g
g = new Graph()
g.size(0,1050,-70,60)
g.addvar("soma.v(0.5)",1,1,0.6,0.9,2)
graphList[1].append(g)
g.save_name("somatic voltage against time")
// Experimental control
proc initialize() {
t = 0
cnt = 0
finitialize(v_init)
fcurrent()
}
proc integrate() {
g.begin()
while (t<tstop) {
// next section write to the 'data' array
if(branch_levels<5){
data[0][cnt]= t
data[1][cnt] = soma.v(.5)
data[2][cnt] = soma.cai(.5)
data[3][cnt] = apical.cai(0.5)
data[4][cnt] = dend[0].cai(.5) // 1st branch level calcium readout
}else{
data[0][cnt]= t
data[1][cnt] = soma.v(.5)
data[2][cnt] = soma.cai(.5)
data[3][cnt] = apical.cai(0.5)
data[4][cnt] = dend[0].cai(.5)
data[5][cnt] = dend[14].cai(.5) // 4th branch level calcium readout
}
cnt=cnt+1
g.plot(t)
fadvance()
}
g.flush()
}
// function runs experiment
proc go() {
initialize()
integrate()
finalize()
}
// finalize writes data to file
proc finalize(){
strdef filename, prefix, postfix
prefix = "altered_complexity_branch_levels_ratio"
postfix="txt"
num = branch_levels
num2 = exponent
sprint(filename, "%s_%d_%f.%s", prefix, num, num2,postfix)
wopen(filename)
// write data to text file
for i=0,nr_cells-1{
if(branch_levels<5){
fprint("%f\t %f\t %f\t %f\t %f\n",data[0][i],data[1][i],data[2][i],data[3][i],data[4][i])
}else{
fprint("%f\t %f\t %f\t %f\t %f\t %f\n",data[0][i],data[1][i],data[2][i],data[3][i],data[4][i],data[5][i])
}
}
wopen()
}
go() // run experiment
