/*
for the paper:
Yu Y.G., Shu Y., Duque A., Haider B., and McCormick D.A.
Cortical Action Potential Back-propagation Explains Spike Threshold Variability and Rapid-Onset Kinetics.
Journal of Neuroscience, 28: 7260-7272, (2008).
McCormick D.A., Shu Y., and Yu Y.G.
Neurophysiology: Hodgkin and Huxley model - still standing?
Nature. 445: E1-E2, (2007).
Create a single cell containing a dendrite, soma, an axon with 60 micrometer and a axon bleb
to investigate generation of spike kink at soma
As observed in experiment, spike usually initiated at the axon site 40~70 micrometer away from soma.
There is a rapid rising phase in somatic spike, due to the initial segment traveling spike.
*/
load_file("nrngui.hoc")
// --------------------------------------------------------------
// redefine some things in stdrun.hoc
// --------------------------------------------------------------
objref st, st1
tstart = 0
tstop = 400
dt = 0.01/1
steps_per_ms = 1/dt
// Create the neuron
rm = 30000
v_init = -70
celsius = 37
Ek = -80 //-90 might affect the axon spike phase slope value low, and less noisy, increase to be -85?
// maybe not, still need check
Ena = 60
ra = 150 //the larger, the soma kink slope large is
c_m = 0.75 // the optimal value, both increase and decrease will decrease slope value
//however, small C value will have high dv/dt value
xn = 1
gna = 8000 //7500
gk = 1600 //1800
ndend=1 // increase this will increase kink slope
//create dend[ndend], soma, axon
create dend, soma, axon, bleb
gkm = .3// .3 //.3~ 0.5
gca = .3//.3 //.3~ 0.5
gkca = .3// 3 //.3~ 0.5
gca_soma = gca
gkm_soma = gkm
gkca_soma = gkca
dend {
L=3000 //there some optimal length make stronger kink?
nseg=L/50 //don't why increase or decrease nseg will decrease slope
diam = 5 //sure this is also optimal value; increase this strenth the kink and biphase
}
soma {
L=30 //increase size make stronger kink
nseg=L/5
diam = 20
}
axon {
L=50
nseg=L/5
//diam=1
// diam(0:0.1) = 2:1 // weeken the biphase, but looks closer to the real data
// diam(0.1:1) = 1:1 // weeken the biphase, but looks closer to the real data
diam=1
}
bleb {
nseg=1
L=5 //remove the bleb doesn't changing sth much except increase axon spike
// phase slope from 2.5 to 4---which is good!
diam = 5
}
soma connect dend(0), 0
soma connect axon(0), 1
axon connect bleb(0), 1
// connect dend(1), soma(0)
// connect axon(0), soma(1)
// connect bleb(0), axon(1)
proc init_cell() {
forall {
insert pas
Ra = ra
cm = c_m
g_pas = 1/rm
e_pas = v_init
}
bleb g_pas = 1/rm // 5/rm
axon cm=c_m*0.75
soma cm=c_m
bleb cm = c_m*0.75
forall insert na
soma.gbar_na = 0.1*gna //750 //gna/8 // increase this will increase dv/dt amplitude, but weak kink value
axon.gbar_na = gna //7500
bleb.gbar_na = gna/3 //3000 //gna/3
dend.gbar_na = 20 //100 //10 //increase this will make cell has spontaneous spikes.
forall insert kv
soma.gbar_kv = gk/5 //250 //gk/2 // increase this will increase dv/dt amplitude and weak biphase
axon.gbar_kv = gk
bleb.gbar_kv = gk/3 //gk/2
dend.gbar_kv=10
forall {
insert km gbar_km = 0 //gkm
insert kca gbar_kca = 0 //gkca
insert ca gbar_ca = 0 //gca
insert cad
}
axon gbar_km=0
bleb gbar_km=0
axon gbar_kca=0
bleb gbar_kca=0
axon gbar_ca=0
bleb gbar_ca=0
forall if(ismembrane("ca_ion")) {
eca = 140
ion_style("ca_ion",0,1,0,0,0)
vshift_ca = 0
}
forall if(ismembrane("k_ion")) ek = Ek
forall if(ismembrane("na_ion")) {
ena = Ena
// seems to be necessary for 3d cells to shift Na kinetics -5 mV
vshift_na = -6 //-5
}
finitialize(v_init)
if (cvode.active()) {
cvode.re_init()
} else {
fcurrent()
}
frecord_init()
}
init_cell()
nrnmainmenu()
nrncontrolmenu()
access soma
proc newplot() { local wd,ht
graphItem = new Graph(0)
wd=300 ht=200
graphItem.save_name("graphList[0].")
graphList[0].append(graphItem)
// graphItem.view(100,-1,106,90,300,200,wd,ht)
graphItem.view(0,-80,tstop,140,300,200,wd,ht)
}
proc newplot1() { local wd,ht
graphItem = new Graph(0)
wd=300 ht=200
graphItem.save_name("graphList[0].")
graphList[0].append(graphItem)
graphItem.view(0,-3.5,tstop,4,300,200,wd,ht)
}
newplot1()
graphItem.addvar("soma.ina(0.5)",1,1)
graphItem.addvar("axon.ina(0.9)",2,1)
//graphItem.addvar("bleb.ina(0.5)",3,1)
newplot1()
graphItem.addvar("soma.i_cap(.5)",1,1)
graphItem.addvar("axon.i_cap(.9)",2,1)
//graphItem.addvar("bleb.i_cap(.5)",3,1)
graphItem.addvar("soma.ina(.5)",5,1)
graphItem.addvar("axon.ina(.9)",6,1)
objectvar scene_vector_[14]
objectvar save_window_, rvp_
save_window_ = new Graph(0)
save_window_.size(-100,200,-80,40)
scene_vector_[5] = save_window_
{save_window_.view(-100, -80, 200, 120, 33, 410, 300, 200)}
flush_list.append(save_window_)
save_window_.save_name("flush_list.")
objectvar rvp_
//rvp_ = new RangeVarPlot("gbar_na")
//rvp_ = new RangeVarPlot("ina")
rvp_ = new RangeVarPlot("v")
dend rvp_.begin(0)
//bleb rvp_.end(1)
//soma rvp_.begin(0)
//bleb rvp_.end(1)
save_window_.addobject(rvp_, 2, 1, 0.8, 0.9) //addobject: Adds the RangeVarPlot to the list of items to be plotted on flush
newplot()
graphItem.addvar("soma.v(.5)",1,1)
graphItem.addvar("axon.v(.9)",2,1)
//graphItem.addvar("bleb.v(.5)",3,1)
/*
//------controlled by Gluct.mod---------------------------------------------------------------
access soma
objref fl
fl = new Gfluct2(0.5)
fl.g_e0 = 0.0012 // 5 times larger
fl.g_i0 = 0.0012*5
//fl.std_e = 0.012*0.8 // 4 times larger
//fl.std_i = 0.0264*0.8
fl.std_e = 0 //0.012/3 // 4 times larger
fl.std_i = 0 //0.0264/3
*/
// --------------------------------------------------------------
// stimulus
// --------------------------------------------------------------
st=new IClamp(.5)
st.dur = 5000
st.del = 5
st.amp = .0
st1=new IClamp(.5)
st1.dur = 5000
st1.del = 5
st1.amp = .0
proc set_stim() {
st.loc($1) st.amp = $2 st.del = $3 st.dur = $4
}
proc set_stim1() {
st1.loc($1) st1.amp = $2 st1.del = $3 st1.dur = $4
}
objref rect, recv, rec_somaic, rec_axonv, rec_blebv
objref rec_axon1, rec_axon2, rec_axon3, rec_axon4, rec_axon5 , rec_axon6, rec_axon7, rec_axon8, rec_axon9, rec_axon10
objref rec_somaic, rec_axonic, rec_blebic, rec_somaina, rec_axonina//used to record values of variables
objref rec_somaik, recaxonik, recblebina, rec_blebik
objref csna1, cana1, cana2, cana3, cana4, cana5, cana6, cana7, cana8, cana9, cana10
rect = new Vector()
recv = new Vector()
rec_axonv=new Vector()
rec_axon1=new Vector()
rec_axon2=new Vector()
rec_axon3=new Vector()
rec_axon4=new Vector()
rec_axon5=new Vector()
rec_axon6=new Vector()
rec_axon7=new Vector()
rec_axon8=new Vector()
rec_axon9=new Vector()
rec_axon10=new Vector()
rec_somaic=new Vector()
rec_axonic=new Vector()
rec_somaina=new Vector()
rec_axonina=new Vector()
rec_somaik=new Vector()
recaxonik=new Vector()
objref savdata //used for saving file
savdata = new File() //savdata is a file pointer
savdata.wopen("neuron_soma.dat") //open a file with a defined name to save data
//objref savdata1 //used for saving file
//savdata1 = new File() //savdata is a file pointer
//savdata1.wopen("ina.dat") //open a file with a defined name to save data
soma set_stim1(.5,0,0,50000)
proc soma_inj() {
//soma set_stim(0.5,0.,50,200)
soma set_stim(0.5,0.1,50,250)
//recv.record(&axon.v(0)) //put value of soma.v to recv
recv.record(&soma.v(0.5)) //put value of soma.v to recv
rect.record(&t) //put value of t to rect
rec_axon1.record(&axon.v(0.1))
rec_axon2.record(&axon.v(0.2))
rec_axon3.record(&axon.v(0.3))
rec_axon4.record(&axon.v(0.4))
rec_axon5.record(&axon.v(0.5))
rec_axon6.record(&axon.v(0.6))
rec_axon7.record(&axon.v(0.7))
rec_axon8.record(&axon.v(0.8))
rec_axon9.record(&axon.v(0.9))
rec_axon10.record(&axon.v(1))
rec_somaic.record(&soma.i_cap(0.5))
rec_axonic.record(&axon.i_cap(0.9))
rec_somaina.record(&soma.ina(0.5))
rec_axonina.record(&axon.ina(0.9))
rec_somaik.record(&soma.ik(0.5))
recaxonik.record(&axon.ik(0.9))
run() // this is the right position
for i=0,rect.size()-1 {
savdata.printf("%g %g %g %g %g %g %g %g %g %g %g %g\n", rect.x(i), recv.x(i), rec_axon1.x(i), rec_axon2.x(i), rec_axon3.x(i), rec_axon4.x(i), rec_axon5.x(i), rec_axon6.x(i), rec_axon7.x(i), rec_axon8.x(i), rec_axon9.x(i), rec_axon10.x(i))
//savdata1.printf("%g %g %g %g %g %g %g %g %g\n", rect.x(i), recv.x(i), rec_axonv.x(i),rec_somaic.x(i),rec_axonic.x(i),rec_somaina.x(i),rec_axonina.x(i),rec_somaik.x(i),recaxonik.x(i))
}
savdata.close()
//savdata1.close()
}
xpanel("Stuart & Sakmann")
xvalue("stim amp","st.amp")
xbutton("inject soma","soma_inj()")
xpanel()