na_ion mechanism not inserted in section axon

Anything that doesn't fit elsewhere.
Wosen

na_ion mechanism not inserted in section axon

Post by Wosen »

Dear Neuron support,

I have the following puzzle.

I have defined the following point mechanism to get an ion current (note that it is not a distributed mechanism as usual but rather a point process).

Code: Select all

UNITS {
        (mA) = (milliamp)
        (mV) = (millivolt)
        (S) = (siemens)
}

NEURON {
        POINT_PROCESS sodium_brette
        USEION na READ ena WRITE ina
        RANGE gnabar, gna, ena, v_resting, trigger
        GLOBAL minf, mtau, v1_2
}

PARAMETER {
        gnabar = 0 (S)  <0,1e9>  : no density as this is a point process
        trigger = 0
        v1_2 = -40 (mV)
        mtau = 0.1 (ms)   : fixed to 100 microsec in Brette 2013
        v_resting = -75 (mV)
}

STATE {
        m
}

ASSIGNED {
        v (mV)
        ena (mV)

        gna (S)
        ina (nA)
        minf 
}

LOCAL mexp, nexp

BREAKPOINT {
        SOLVE states METHOD cnexp
        gna = gnabar*m
        ina = gna*(v - ena)
}

INITIAL {
        rates(v)
        m = minf
}

DERIVATIVE states {  
        rates(v)
        m' =  (minf-m)/mtau
}

PROCEDURE rates(v(mV)) {  :Computes rate and other constants at current v.
                      :Call once from HOC to initialize inf at resting v.
        TABLE minf FROM -100 TO 100 WITH 200

UNITSOFF
        
        minf = 1 / (1 + exp((v1_2 - v)/6))
}

FUNCTION vtrap(x,y) {  :Traps for 0 in denominator of rate eqns.
        if (fabs(x/y) < 1e-6) {
                vtrap = -y*(1 - x/y/2)
        }else{
                vtrap = x/(1 - exp(x/y))
        }
}

NET_RECEIVE (w) {
        trigger = 1
        m = 1 / (1 + exp((v1_2 - v_resting)/6))
}

UNITSON
I would like to set ena as is done with distributed mechanisms. However, if I load my neuron model with this mechanism in the axon after writing ena in the command line I get:

Code: Select all

na_ion mechanism not inserted in section axon
/home/david/Programme/Neuron/nrn/x86_64/bin/nrniv: 
 near line 20
 ena
If I try to access it by sodiumchan.ena instead (which is an instance of the point process sodium_brette), I get:

Code: Select all

sodiumchan.ena
ena not a public member of sodium_brette
/home/david/Programme/Neuron/nrn/x86_64/bin/nrniv: sodium_brette ena
 near line 21
 sodiumchan.ena
My question now is: how can I tell neuron that my point process is an na_ion mechanism? Should the statement

Code: Select all

USEION na READ ena WRITE ina
not be enough?

Thank you very much for your help!

Best regards,
David
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: na_ion mechanism not inserted in section axon

Post by ted »

hoc doesn't know anything about ena unless a section contains a mechanism that uses the na ion. The properties of such a mechanism may be specified by a Channel Builder or NMODL code.

If a mechanism's properties are specified by NMODL code, and that code's NEURON block contains a
USEION foo . . .
statement, then
1. inserting that mechanism into a section will automatically create a mechanism called foo_ion.
and
2. foo's equilibrium potential will be known to hoc as efoo, and the membrane current carried by foo will be known as ifoo

Note that foo_ion is a "distributed" or "density" mechanism. This means that
1. you must specify which section you mean when you refer to efoo and ifoo
2. both efoo and ifoo are range variables, so they can have different values in different segments.

Attach your point process to a section, and then you will be able to specify ena in that section.


Comments on the NMODL code you provided:
It should be checked with modlunit, which will reveal problems that should be fixed.. First you will find that nA must be declared in the UNITS block.
(nA) = (nanoamp)
will do the job.

Then you will discover that there are serious inconsistencies of units. Since gna is in Siemens, the calculation of ina in the BREAKPOINT block is incorrect by 6 orders of magnitude. The fix is to change it to
ina = (1e6)*gna*(v - ena)

If this NMODL code adhered to NEURON-specific conventions for naming and units,
1. its name would not include an underscore character.
2. its name would use upper case because this is a class. I'd choose NaB because it's reasonably mnemonic, and easy to type and read.
3. gnabar and gna would both be in microsiemens. This would eliminate the need for a conversion factor in the calculation of ina. It could be done by declaring
(uS) = (microsiemens)
in the NEURON block, then applying this new unit in the declarations of gna and gnabar.

If you ever want to use multithreaded parallel execution to speed up simulation of a model that uses this mechanism, you'll discover that it is not thread safe. The problem is the assignment of values to the global variable minf. This could be fixed by including the keyword
THREADSAFE
in the NEURON block. That would change the global variables to "per-thread globals" i.e. variables that are global only for the mechanism instances that exist on a single thread.
Wosen

Re: na_ion mechanism not inserted in section axon

Post by Wosen »

Dear Ted,

Thank you very much for pointing me out to the mistake in the mechanism!

Still, the error persists but thanks to your comment I could further encompass the problem. Indeed, whether or not na_ion is created in the specific section seems to be dependent on the position I insert my point process. The following command line code shows this:

Code: Select all

oc>access axon
oc>sodiumchan
	NULLobject 
oc>sodiumchan = new sodium_brette(0.4)
oc>ena
na_ion mechanism not inserted in section axon
/home/david/Programme/Neuron/nrn/x86_64/bin/nrniv: 
 near line 8
 ena
    ^
oc>sodiumchan = new sodium_brette(0.5)
oc>ena
	50 
Maybe I can do a workaround by inserting the same point process two times and setting the one at position 0.5 to zero conductivity to have it not influence the current?

I would be very greatful for further recommendations!

Best regards,
David
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: na_ion mechanism not inserted in section axon

Post by ted »

Actually the situation turns out to be more complex than you or I thought, and your suggested workaround is not guaranteed to work.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: na_ion mechanism not inserted in section axon

Post by hines »

For now, work around the problem by explicitly inserting na_ion in any section before you instantiate a sodium_brette instance in that section. ie with the statement
section { insert na_ion }
This will ensure that ena and ina will exist in every non-zero area segment of the section. na_ion is a density mechanism so you should avoid locating your sodium_brette POINT_PROCESS
at 0 or 1.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: na_ion mechanism not inserted in section axon

Post by hines »

This bug is now fixed in the hg repository.
http://www.neuron.yale.edu/hg/neuron/nr ... cd5f896742
aNNe
Posts: 4
Joined: Fri Sep 13, 2013 10:59 am

Re: na_ion mechanism not inserted in section axon

Post by aNNe »

First let me thank you for fixing the bug in a blizz. Since we can change ena, the model performs exactly as we hoped. Going through the results, we realized, that two things are not completely clear to us.
1. the "location" at which the point process injects the current
as you, Ted, pointed out earlier, na_ion is a distributed mechanism. Therefore, the current induced by this particular point process should be mathematically equivalent to a current created by any distributed mechanism that also uses the na_ion in the same node. Is this correct?

2. how to "read out" ina of this point process in the gui
Following your advice to sort that six orders of magnitude error, the units of ina in the definition of the point process are nanoamp. However the units of ina, read out for instance in a phase plane plot, is normally mA/cm2, because na_ion is a distributed process.
So is there a direct way to read out the current injected by this point process, instead of going through the distributed ina that might also be influenced by other processes?

Best regards,
Andreas
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: na_ion mechanism not inserted in section axon

Post by ted »

aNNe wrote:1. the "location" at which the point process injects the current
as you, Ted, pointed out earlier, na_ion is a distributed mechanism. Therefore, the current induced by this particular point process should be mathematically equivalent to a current created by any distributed mechanism that also uses the na_ion in the same node. Is this correct?
What is your question? I don't understand what you mean by "mathematically equivalent." Are you suggesting that if a segment has more than one mechanism that WRITEs ina, then all of these mechanisms will produce the same current? The answer is no. Each one will produce a current whose magnitude is governed by its own equations.
how to "read out" ina of this point process in the gui
If you want to discover the current produced by an instance of this particular mechanism, the current that it produes must be declared to be a RANGE variable. Otherwise it will be pooled with all sodium currents produced by whatever mechanisms WRITE ina, i.e. it will affect charge balance correctly, and it will affect sodium mass balance correctly, but you won't be able to access its numerical value directly. So just add the statement
RANGE ina
to the NEURON block and you will be able to do this:

Code: Select all

objref foo
soma foo = new sodum_brette(0.5)
run()
print foo.ina
menica
Posts: 71
Joined: Wed Sep 02, 2015 11:02 am

Re: na_ion mechanism not inserted in section axon

Post by menica »

Hi Ted,
I read this thread.
In order to know the current generated by a point process, I declared the current as RANGE. How can I access to the conversion in mA/cm2 which neuron does? I tried to convert the current my myself ...but when I look at the total current generated by the point process specific of a certain ion species I have different values. In other words: I have a synaptic point process which writes icl current (syn.icl) in nA. In the compartment where this point process is there are density mechanisms which generate others icl currents in mA/cm2 (icl_leak, icl_pump). I expect that the icl total will be the sum of all the other currents involved : icl= icl_leak+icl_pump+syn.icl*0.000001/pi*diam*diam, but I cannot see this.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: na_ion mechanism not inserted in section axon

Post by hines »

The conversion of a density current sec.icl(x) which has the units mA/cm2 to absolute current (nA) through the segment containing position x is
in the hoc interpreter
sec {
absolute_icl = icl(x)*area(x)*0.01
}

and in the python interpreter
absolute_icl = sec(x).icl*sec(x).area()*0.01
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: na_ion mechanism not inserted in section axon

Post by ted »

menica wrote: . . . In the compartment where this point process is there are density mechanisms which generate others icl currents in mA/cm2 (icl_leak, icl_pump). I expect that the icl total will be the sum of all the other currents involved : icl= icl_leak+icl_pump+syn.icl*0.000001/pi*diam*diam, but I cannot see this.
NEURON expects point processes to generate currents that are in nA. A properly written point process that WRITEs icl will declare the units of icl to be in nA. Then if foo is an instance of that point process, the hoc or Python interpreter will return the value of foo.icl in nA. The icl that the interpreter returns for any segment is indeed the total chloride current for that segment in mA/cm2, and that total includes all chloride currents generated by whatever density mechanisms and point processes that are in that segment; in calculating that total, NEURON automatically takes into account the segment area and the fact that point process currents are in nA.
menica
Posts: 71
Joined: Wed Sep 02, 2015 11:02 am

Re: na_ion mechanism not inserted in section axon

Post by menica »

Dear Hines and Ted,
thanks for your explanation.
There is a way in order to know the current generated by the point process that neuron will sum up to the other currents and convert in mA/cm2, directly in mA/cm2? Is this synaptic mechanism wrong? :

Code: Select all

NEURON {	
	POINT_PROCESS is   
	USEION cl WRITE icl VALENCE -1	   
        RANGE tau1, tau2,  g,icl,i		   
	POINTER e
}


PARAMETER {
	tau1 = 2   (ms)
	tau2 = 6   (ms)
}

ASSIGNED {
	v    (mV)
	e    (mV)
	ecl  (mV)
	icl  (nanoamp)
        g (microsiemens)
        factor
	area (um2)
	i (mA/cm2)
}

STATE {
        A (microsiemens)
        B (microsiemens)
}

INITIAL {
 	LOCAL tp
 	if (tau1/tau2 > .9999) {
    	 tau1 = .9999*tau2
 	}
 	A = 0
 	B = 0
 	tp = (tau1*tau2)/(tau2 - tau1) * log(tau2/tau1)
 	factor = -exp(-tp/tau1) + exp(-tp/tau2)
 	factor = 1/factor
}

BREAKPOINT {
 	SOLVE state METHOD cnexp
  	g = B - A
	icl= g*(v - e)
	i=icl/(area*0.01)
}

DERIVATIVE state {
 	A' = -A/tau1
 	B' = -B/tau2
}

 NET_RECEIVE(weight (microsiemens)) {
 	A = A + weight*factor
 	B = B + weight*factor
} 
in this way, i can leave neuron summing up in the correct way the icl currents because icl is in nanoampere, but at the same time, I will have access to the I currents in mA/cm2.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: na_ion mechanism not inserted in section axon

Post by hines »

when I used modlunit on the fragment you sent (copied to temp.mod) I get

Code: Select all

$ modlunit temp.mod
model   1.1.1.1   1994/10/12 17:22:51
Checking units of temp.mod
Need declaration in UNITS block of the form:
	(mV)	(units)
Cannot recognize the units: mV at line 1 in file temp.mod
NEURON<<ERROR>> {	
adding a units block

Code: Select all

UNITS {
  (mV) = (millivolt)
  (um) = (micron)
  (mA) = (milliamp)
}
then gives

Code: Select all

$ modlunit temp.mod
model   1.1.1.1   1994/10/12 17:22:51
Checking units of temp.mod
The previous primary expression with units: 1000 coul/m2-sec
is missing a conversion factor and should read:
  (100)*()
 at line 52 in file temp.mod
	i=icl/(area*0.01)<<ERROR>>
because POINT_PROCESS currents must be nanoamps
Changeing the declaration of i to (nanoamp) and removing the area factor
in the BREAKPOINT block to i = icl causes the file to be accepted by
modlunit.

Note, though, that USEION cl refers to a density mechanism so that icl
has external units of mA/cm2. Whereas the icl units in this mod file
are (nanoamp). So the contribution of an instance of this model's
current at the
interpreter is cl[...].icl whereas the total icl at this location (say
sec(x) in hoc is sec.icl(x) which must be in mA/cm2.

So the question is whether nocmodl produces dimensionally consistent
values for internal icl and external icl. Looking at the translated c
file, I see the external value is translated as

Code: Select all

  _ion_icl += icl * 1.e2/ (_nd_area);
so everything seems to be fine.

By the way, in your mod file the statement ' i = icl ' is useless as it is icl that contributes to the total current and not i. In this modlel the
RANGE icl is not the total icl but only the part contibuted by this POINT_PROCESS instance.
menica
Posts: 71
Joined: Wed Sep 02, 2015 11:02 am

Re: na_ion mechanism not inserted in section axon

Post by menica »

Dear Hines,
thanks. This means that it is not possible to have access to the icl produced by the point process in mA/cm2 in the nmodl code.
Maybe I should convert it in the hoc code in the post-processing. If I save the syn.icl current in a vector, then I can multiply each value by conversion factor: conv_fac = 1/(0.01*PI*diam*diam)? Is it the right conversion? Should I obtain in this way the same value that neuron obtains when it automatically convert in mA/cm2?
Last edited by menica on Fri Dec 15, 2017 10:46 am, edited 1 time in total.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: na_ion mechanism not inserted in section axon

Post by hines »

This means that it is not possible to have access to the icl produced by the point process in mA/cm2 in the nmodl code
No. Only icl mus be nA. My change to i was just one of several ways to establish units consistency. It wold be ok to declare
icl (mA/cm2) but then your original
i=ica/(area*0,01)
would have to be changed to
i=ica/area * (100)
or
i = ica/(area*(0.01))
A number surrounded by () is treated as a units conversion factor as opposed to a numeric factor.
Post Reply