Axial currents as a function of (x,t)

[biswa]

I wanted to compute the time course of axial currents as a function of space and time in multi-compartment models. To this effect Mike had proposed that I write a mod file whereby the axial current is implemented as a density mechanism. I cooked up the following mod file to achieve this (note that presently I do not have any code to handle parent and child sections connected at 0 or the 1 end respectively).

Code: Select all

COMMENT
This file is used to calculate the axial currents between the defined node(N) and the previous (N-1) node
ENDCOMMENT

NEURON	{
	SUFFIX iax
	RANGE x, i
	POINTER vprev, vnext
}

ASSIGNED	{
	i		(milliamp)
	vprev	(millivolt)
	vnext	(millivolt)
}

PARAMETER {
   x = 0 (1)
}

AFTER INITIAL	{: during a call to finitialize()
	i = (1e-6) * (vprev - vnext)/ri(x)
}

AFTER SOLVE	{ : after each solution step compute
	i = (1e-6) * (vprev - vnext)/ri(x)
}

Using setpointers I tend to attach the pointer reference in the following hoc file:

Code: Select all

load_file("nrngui.hoc")
tstart 		= 0
tstop		= 500
dt 			= 0.01
v_init 		= -70 // mV
celsius 	= 37

create  axon
access axon

// axon
axon {
    cm 			= 0.75
    diam 		= 1
    L 			= 1000
    Ra 			= 150
    nseg 		= int((L/(0.1*lambda_f(100))+.999)/2)*2 + 1
    insert hh
    insert iax
}


// Caution: Such a treatment will not take joined sections into account i.e., if there are two 
// sections (soma & axon) then for nseg > 1, only axial currents as a function of time, in internal nodes'll be computed
forall {
if (ismembrane("iax")) {
xold = 0
for (x, 0)	{
	setpointer vprev_iax(x), secname().v(xold)
	setpointer vnext_iax(x), secname().v(x)
	xold = x
}
}
}

// code for accounting for axial currents through the edges of various sections
run()

I believe that I have grave errors in setting the pointers and the construction of the mod file. Where am I falling short?

[ted]

Try debugging it yourself. Assume that all of your hoc and NMODL code is broken. Then

Code: Select all

for each hoc or NMODL statement
  repeat
    compare what it does with what you wanted it to do
    revise if necessary
  until it does what you wanted

Start with a very small model--no more than 3 segments--and use a very short run time (a couple of msec at most) with dt no smaller than 0.025 so there aren't a lot of garbage results to sort through. printf statements will be particularly helpful; NMODL's printf() uses the same syntax as C does.

[biswa]

With some help I have been able to tackle the mod file. The question that I now have is: say that a soma (nseg=1) is connected to an axon (nseg=3), does ri at axon position 0 mean that the path resistance is calculated between soma(position = 1) and axon(position = 0) or is it between soma(position=0.5) and axon(position=0)? I did corner the source code at

Code: Select all

nrn/src/nrnoc/treeset.c

whereby ri is coded as

Code: Select all

rleft = 1e-2*ra*(dx/2)/(PI*diam*diam/4.); /*left half segment Megohms*/
NODERINV(nd) = 1./(rleft + rright); /*uS*/
rright = rleft;

If the section is defined as root, the 0 end has infinite resistance. Nonetheless, I was not able to follow what happens at the 0/1 end.

[ted]

Sometimes the easiest way to gain understanding is to try a very simple test case.

Code: Select all

oc>create soma,dend
oc>access soma
oc>connect dend(0), soma(1)
oc>soma { L = 10  diam = 10 }
oc>dend { L = 100 diam = 1 }
oc>forall for (x) print secname(), " ", x, ri(x)
first instance of x
soma 0 1e+30 
soma 0.5 0.02253634 
soma 1 0.02253634 
dend 0 0.02253634 
dend 0.5 22.53634 
dend 1 22.53634

which makes sense, given that dend(0) and soma(1) are the same node and therefore the "parent node" of one is the parent node of the other.
Is this phrase in your question

say that a soma (nseg=1) is connected to an axon (nseg=3)

to be taken literally, i.e. was the topology set up with a statement like
connect soma(something), axon(something)
or was it done more "conventionally" i.e. more like
connect axon(0), soma(0)
?

[biswa]

Thanks Ted. What I meant was: I use connect axon(0), soma(1) which is basically a ball and stick model, similar to your toy example.

So in summary (using your toy example), the path resistance at the connection between soma and dend is between soma(1) and dend(0) and not between soma(0.5) and dend(0.5). Therefore, the axial current flowing out of the soma should be (soma.v(1) - dend.v(0))/(axial resistance at dend(0)) (note that this is more of a pseudo code than the actual syntax) and not (soma.v(0.5) - dend.v(0.5))/(axial resistance at dend(0.5)) i.e., the axial current flowing out at the connection between the soma and dend is not between the centre of the compartments? Pardon me for the confusion.

[ted]

connect dend(0), soma(1)
attaches the 0 end of dend to the 1 end of soma. These two nodes are now identical.
The resistance between them is 0.
dend.v(0) will be identical to soma.v(1)
Assuming that no other sections and no point processes are attached to that node, the current that enters it from dend will equal the current drawn out of it by soma. The resistance from dend(0.5) to dend(0) (aka soma(1)) is found by
dend print ri(0.5)
and the resistance from dend(0) (aka soma(1)) to soma(0.5) is found by
dend print ri(0)
or
soma print ri(1)

[biswa]

Thanks a bunch, this is indeed very helpful.

[ted]

It's not immediately obvious nor is it explicitly spelled out in the documentation. All the hints are there, but one has to be an ace at drawing inferences in order to grasp all this just by reading the docs a few times. For most (including me), constructing and examining toy models is often essential to understanding some of this stuff.