Where to find squid.hoc associated with hh.mod

NMODL and the Channel Builder.
Post Reply
pascal
Posts: 71
Joined: Thu Apr 26, 2012 11:51 am

Where to find squid.hoc associated with hh.mod

Post by pascal » Tue Oct 09, 2018 12:13 pm

I am learning to write my own MOD files, and I'm starting with hh.mod as a learning tool. The comment section in hh.mod says to see squid.hoc as an example of the hh mechanism's implementation...but I cannot find squid.hoc anywhere. It's not in my installation of NEURON, and I can't find it anywhere online. A pointer in the right direction would be appreciated. Thanks.

pascal
Posts: 71
Joined: Thu Apr 26, 2012 11:51 am

Re: Where to find squid.hoc associated with hh.mod

Post by pascal » Tue Oct 09, 2018 5:43 pm

Okay, I believe I found it here: https://github.com/neurodroid/neurodroi ... /squid.hoc

Here is the code:

Code: Select all

/* Squid giant axon ap */

load_file("stdrun.hoc")

/* Print v at every time step */
proc advance() {
    fadvance()
    print "ND ", v
}

create axon
access axon
insert hh

gnabar_hh = gna
gkbar_hh = gk

objref ic
ic = new IClamp(0)
ic.dur = 0.5 /* ms */
ic.amp = 50  /* nA */

print "ND dt ", dt 
run()
But this raises a new question: hh.mod includes the lines USEION na READ ena WRITE ina and USEION k READ ek WRITE ik. Where are 'ena' and 'ek' initialized? It doesn't happen in hh.mod, nor does it happen in squid.hoc. In reading the documentation on ion_style(), it seems that the initialization of these reversal potentials is done somewhere under the hood.

Are there standard values for intracellular/extracellular concentrations of these ions that are used to initialize the reversal potentials? And what would be the best way to use my own custom values for these reversal potentials?

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

Re: Where to find squid.hoc associated with hh.mod

Post by ted » Tue Oct 09, 2018 10:56 pm

Okay, I believe I found it here:
Nope. The COMMENT block at the start of hh.mod was apparently written by Stuart Jaslove in 1992, a time at which NEURON was known only to a small circle of people. Stuart was at Duke then, where he interacted with Michael Hines and John Moore. That was well before the time of Christoph Schmidt-Hieber, who set up the github archive you discovered.
Are there standard values for intracellular/extracellular concentrations of these ions that are used to initialize the reversal potentials? And what would be the best way to use my own custom values for these reversal potentials?
Good questions. The answers are in chapters 8 and 9 of The NEURON Book. Chapter 9 is a deep introduction to NMODL. Chapter 8 discusses initialization in detail. If you don't have the book, get these preprints of those chapters:
https://www.neuron.yale.edu/ftp/ted/boo ... xedref.pdf
https://www.neuron.yale.edu/ftp/ted/boo ... xedref.pdf

It'll take some time to read that material. You may find it faster to start by reading
https://www.neuron.yale.edu/neuron/stat ... modl2.html
especially the stuff about USEION
and
https://www.neuron.yale.edu/neuron/stat ... /ions.html

Basically, NEURON knows nothing about any ions or their equilibrium potentials unless your model has at least one mechanism with a USEION statement that names the ion. The ion names na and k are special, in that NEURON "knows" default values of initial values for their intra- and extracellular concentrations and their equilibrium potentials--but those default values are appropriate only for squid axon at 6.3 degC.

If your model knows about some ion foo but doesn't have an ion accumulation mechanism for foo (that is, it doesn't have a mechanism that WRITEs the intra- and/or extracellular concentration of foo), then foo's equilibrium potential efoo is simply a constant parameter, and you can give efoo whatever value you like with a simple assignment statement as in these examples:

// hoc example
efoo = -80

# Python example
h.efoo = -80

If your model DOES have a mechanism that WRITEs fooi and/or fooo (intracellular and extracellular concentrations of foo, respectively), then by default efoo will be calculated from the Nernst equation at every fadvance().

You might wonder how to initialize ionic concentrations to desired values. The appropriate answer depends on exactly what you are trying to achieve; you'll find chapter 8 to be very informative about this.

pascal
Posts: 71
Joined: Thu Apr 26, 2012 11:51 am

Re: Where to find squid.hoc associated with hh.mod

Post by pascal » Wed Oct 10, 2018 11:04 am

Thanks for the resources, Ted. I'll read through them and get back to you if I have any more questions.

pascal
Posts: 71
Joined: Thu Apr 26, 2012 11:51 am

Re: Where to find squid.hoc associated with hh.mod

Post by pascal » Thu May 09, 2019 1:49 pm

Okay, I think I am getting a handle on ionic mechanisms in mod files. Those two links were very helpful. I still have a few questions:

1) When writing mod files for custom ionic currents, is there ever a situation where I would not want to include "THREADSAFE" within the NEURON block?

2) I noticed that in the rates procedure within hh.mod, UNITSOFF is used. Is this generally acceptable practice, or should I avoid this when writing my own mod files? It doesn't seem like a good idea in general...

3) In hh.mod, what is the point of the question mark lines (? currents, ? states, ? rates). They do not appear to be comments, but I'm not sure what their purpose is.

Thanks for the help!

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

Re: Where to find squid.hoc associated with hh.mod

Post by ted » Fri May 10, 2019 6:58 pm

Okay, I believe I found it here:
Nope. FYI anybody can put anything on the web and call it anything they like, but that doesn't mean that it is the same as something else with the same name. This is a corollary of "On the Internet, nobody knows you're a dog."

pascal
Posts: 71
Joined: Thu Apr 26, 2012 11:51 am

Re: Where to find squid.hoc associated with hh.mod

Post by pascal » Mon May 20, 2019 1:08 pm

Right, I get that the squid.hoc file I found was wrong, but I'm talking about the following hh.mod file (which I got from the NEURON 7.4 source code):

Code: Select all

TITLE hh.mod   squid sodium, potassium, and leak channels
 
COMMENT
 This is the original Hodgkin-Huxley treatment for the set of sodium, 
  potassium, and leakage channels found in the squid giant axon membrane.
  ("A quantitative description of membrane current and its application 
  conduction and excitation in nerve" J.Physiol. (Lond.) 117:500-544 (1952).)
 Membrane voltage is in absolute mV and has been reversed in polarity
  from the original HH convention and shifted to reflect a resting potential
  of -65 mV.
 Remember to set celsius=6.3 (or whatever) in your HOC file.
 See squid.hoc for an example of a simulation using this model.
 SW Jaslove  6 March, 1992
ENDCOMMENT
 
UNITS {
        (mA) = (milliamp)
        (mV) = (millivolt)
	(S) = (siemens)
}
 
? interface
NEURON {
        SUFFIX hh
        USEION na READ ena WRITE ina
        USEION k READ ek WRITE ik
        NONSPECIFIC_CURRENT il
        RANGE gnabar, gkbar, gl, el, gna, gk
        :GLOBAL minf, hinf, ninf, mtau, htau, ntau
        RANGE minf, hinf, ninf, mtau, htau, ntau
	THREADSAFE : assigned GLOBALs will be per thread
}
 
PARAMETER {
        gnabar = .12 (S/cm2)	<0,1e9>
        gkbar = .036 (S/cm2)	<0,1e9>
        gl = .0003 (S/cm2)	<0,1e9>
        el = -54.3 (mV)
}
 
STATE {
        m h n
}
 
ASSIGNED {
        v (mV)
        celsius (degC)
        ena (mV)
        ek (mV)

	gna (S/cm2)
	gk (S/cm2)
        ina (mA/cm2)
        ik (mA/cm2)
        il (mA/cm2)
        minf hinf ninf
	mtau (ms) htau (ms) ntau (ms)
}
 
? currents
BREAKPOINT {
        SOLVE states METHOD cnexp
        gna = gnabar*m*m*m*h
	ina = gna*(v - ena)
        gk = gkbar*n*n*n*n
	ik = gk*(v - ek)      
        il = gl*(v - el)
}
 
 
INITIAL {
	rates(v)
	m = minf
	h = hinf
	n = ninf
}

? states
DERIVATIVE states {  
        rates(v)
        m' =  (minf-m)/mtau
        h' = (hinf-h)/htau
        n' = (ninf-n)/ntau
}
 
:LOCAL q10
? rates
PROCEDURE rates(v(mV)) {  :Computes rate and other constants at current v.
                      :Call once from HOC to initialize inf at resting v.
        LOCAL  alpha, beta, sum, q10
:        TABLE minf, mtau, hinf, htau, ninf, ntau DEPEND celsius FROM -100 TO 100 WITH 200
UNITSOFF
        q10 = 3^((celsius - 6.3)/10)
                :"m" sodium activation system
        alpha = .1 * vtrap(-(v+40),10)
        beta =  4 * exp(-(v+65)/18)
        sum = alpha + beta
	mtau = 1/(q10*sum)
        minf = alpha/sum
                :"h" sodium inactivation system
        alpha = .07 * exp(-(v+65)/20)
        beta = 1 / (exp(-(v+35)/10) + 1)
        sum = alpha + beta
	htau = 1/(q10*sum)
        hinf = alpha/sum
                :"n" potassium activation system
        alpha = .01*vtrap(-(v+55),10) 
        beta = .125*exp(-(v+65)/80)
	sum = alpha + beta
        ntau = 1/(q10*sum)
        ninf = alpha/sum
}
 
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/(exp(x/y) - 1)
        }
}
 
UNITSON

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

Re: Where to find squid.hoc associated with hh.mod

Post by ted » Mon May 20, 2019 2:53 pm

Code: Select all

When writing mod files for custom ionic currents, is there ever a situation where I would not want to include "THREADSAFE" within the NEURON block?
There are many. The only reason to use the THREADSAFE directive is in NMODL code that contains statements that assign values to one or more GLOBAL variables, e.g. to hold the results of intermediate calculations. THREADSAFE tells the NMODL compiler to treat GLOBAL variables as "per-thread GLOBALs", so that each thread has its own instance of each GLOBAL and the compiled mechanism can be safely used in multithreaded simulations. However, THREADSAFE also suppresses warning messages about other things that may make an NMODL-specified mechanism inappropriate for multithreaded simulation, such as VERBATIM blocks. Bottom line: don't insert THREADSAFE unless the only non-threadsafe aspect of a mechanism is the fact that it assigns values to GLOBALs.
I noticed that in the rates procedure within hh.mod, UNITSOFF is used. Is this generally acceptable practice, or should I avoid this when writing my own mod files? It doesn't seem like a good idea in general...
UNITSOFF . . . UNITSON is often used to surround blocks of code that would otherwise gag modlunit, the model units checker. Examples of such code include statements that calculate rate constants steady-state and time constant values, or kinetic schemes. Obviously this largely defeats the purpose of checking NMODL files with modlunit, but it has the advantage of enabling use of modlunit to check the other code in a file (the "easy stuff").

modlunit stops checking a file as soon as it runs into a single statement that has inconsistent units. If you remove UNITSOFF . . . UNITSON and fix every inconsistency that it finds, you can be sure that the resulting code contains no unit errors (e.g. calculating current as microamps/um2 instead of mA/cm2). However, the corrected code may be difficult to read because of all of the unit-specifying statements that must be inserted.
In hh.mod, what is the point of the question mark lines (? currents, ? states, ? rates). They do not appear to be comments
They are comments. The only thing in favor of ? is that, on modern computer displays, it is easier to see than : -- but if you use ? in your mod files, expect questions.

pascal
Posts: 71
Joined: Thu Apr 26, 2012 11:51 am

Re: Where to find squid.hoc associated with hh.mod

Post by pascal » Fri May 24, 2019 12:47 am

Thanks, Ted, that's very helpful.

Post Reply