Is this right

[davidfourie]

I am developing a channel from a journal and am not sure weather I have done it correctly.

The journal uses the following equations:

I = g *m^3*h*(V-E)

tm(V)*(dm/tm) = m0(V)-m
th(V)*(dh/th) = h0(V)-h

m0(V) = (1/1+exp(-(V+26.8)/8.2)
h0(V) = (1/1+exp(-(V+48.5)/4.8)

tm(V) = 19.8 - (10.7/1+exp[-(V+26.5)/8.6]
th(V) = 666 - (379/1+exp[-(V+33.6)/11.7]

g = 2.7
E = 50

Does the below code represent the above

Code: Select all

COMMENT
  The Sodium P channel
  Based on the INaP channel in (Soto-Trevino 2005)
  
ENDCOMMENT

TITLE Sodium P Channel

NEURON {
  SUFFIX nap                      	:the mechanism is refered to by nap
  USEION na READ ena WRITE ina  	
  RANGE gbar, g, i, e              	:gbar is max conductance, g macroscopic conductance, i g's current, e is the reversal potential
}

UNITS {
  (S) 	= (siemens)
  (mS)  = (millisiemens)
  (uS)  = (microsiemens)
  (mV) 	= (millivolt)
  (mA) 	= (milliamp)
  (uA) 	= (microamp)
}

PARAMETER { 
  gbar 	 = 2.7 	    (uS/cm2)	    :maximum conductance
  e 	   = 50   		(mV)	:reversal potential
  vm	   = 26.8 	  (mV)		  :half maximum potential 
  vh	   = 48.5 	  (mV)		  :half maximum potential 
  sm	   = 8.2 		  (mV) 		  :Step Width
  sh             = 4.8      (mV) 		  :Step Width
  am            = 19.8     (ms)
  am1          = 10.7     (ms)
  km            = 26.5     (mV)
  km1          = 8.6      (mV)
  ch             = 666      (ms)
  ch1           = 379      (ms)
  kh             = 33.6     (mV)
  kh1           = 11.7     (mV)
}

ASSIGNED {
  v    (mV)		        :current membrane potential
  ena (mV)          	:previous equilibrium potential
  ina  (mA/cm2)      	:future sodium current
  g    (mS/cm2)        	:macroscopic conductance
  i     (mA/cm2)        :current passing through g
}

STATE { m h }		      :m = activation h = inactivation

BREAKPOINT {
  SOLVE states METHOD cnexp
      g = ((0.001)*gbar) * m^3 * h   :macroscopic conductance calculation
      i = (0.001)*(g * (v - e))	    :current calculation
      ina = i		                    :assigning the write variable
}
INITIAL {
  m = mStead(v)
  h = hStead(v)
}
DERIVATIVE states {
  m' = (mStead(v) - m) * (1/Tm(v))
  h' = (hStead(v) - h) * (1/Tm(v))
}


FUNCTION mStead(Vm (mV)) {
  mStead = 1/(1 + exp (-(Vm + vm)/sm))
}

FUNCTION hStead(Vm (mV)) {
  hStead = 1/(1 + exp ((Vm + vh)/sh))
}

FUNCTION Tm(Vm (mV)) (ms) {
  Tm = am - (am1/(1+exp(-(Vm + km)/km1)))
}

FUNCTION Th(Vm (mV)) (ms) {
  Th = ch - (ch1/(1+exp(-(Vm + kh)/kh1)))
}

Does the code give the same results as the equations.

[ted]

Check it out yourself. First, fix any syntax errors or units inconsistencies that modlunit may
detect. Then compile it with mknrndll (nrnivmodl under UNIX/Linux) and look for any error
messages. Finally, use it to generate plots of the voltage dependence of the gating variables'
time constants and steady state values. Use the Grapher to do this--see
http://www.neuron.yale.edu/neuron/stati ... ml#Grapher
For more detailed instructions on how to use the Grapher, see this post on the Forum:
plotting taus and steady-state gs
https://www.neuron.yale.edu/phpBB2/viewtopic.php?p=588
Hint: you can have several Graphers open at one time, each one showing its own plot.

Compare the Grapher's plots with plots that you generate straight from the published
equations (e.g. with a spreadsheet program, SigmaPlot, gnuplot, MathCad, MATLAB, or
even by writing a few funcs in hoc and plotting them . . . in yet another Grapher, of course).

Another test is to reproduce published families of responses to current or voltage clamp,
but few articles bother to show such figures.

[davidfourie]

I am getting no syntax errors and the code is comiling correctly.


I would like to check the formulas using a spreadsheet.I put the equations in from the paper. After each iteration I will have a current value for each current. To generate the V value for the subsequent iteration is it :

V=IR

where I is the total current and R is the cytoplasmic resistivity (Ra).

[ted]

I would like to check the formulas using a spreadsheet.I put the equations in from the paper. After each iteration I will have a current value for each current.

Whoa, not so fast. If there are no syntax errors, units inconsistencies, or compilation errors,
the first practical test is to verify that the gating states' time constants and steady state values
are correct. This doesn't involve any numerical integration--it's just a matter of
evaluating the FUNCTIONs mStead, Tm etc. at a series of voltages, and comparing these
values with results you can grind out with a spreadsheet or even some funcs written in hoc.
For this you can use a for loop in hoc, or the Grapher as per the example I mentioned in my
previous message.

Example: make a hoc file called testhh.hoc that contains these lines:

Code: Select all

create soma
access soma
insert hh
for (v=-100; v<=50; v+=10) print v, " ", minf_hh, " ", mtau_hh, " ", ninf_hh, " ", ntau_hh

Use NEURON to execute this file.
Make another file that contains your own hoc funcs called minf, mtau, ninf, and ntau, that
evaluate the corresponding formulas and print the results.
Finally make a third file that, for a series of voltages, calculates and reports the differences
between minf_hh etc. and your minf(v) etc.

The next test is to calculate the steady state conductance as a function of voltage. Easily
done--for a series of v values, do an finitialize() and then print out the value of g calculated
by this mechanism. Again, use a for loop in hoc, or a Grapher. Compare these results
with what you get from a spreadsheet or a hoc func that you write. Here's an example
that prints out steady state gna_hh:

Code: Select all

create soma
access soma
insert hh
for (v=-100; v<=50; v+=10) { finitialize(v)  print gna_hh }

A further test that one might consider is to do a voltage clamp experiment and compare
the time courses of the simulated conductances and current with values computed by
some other method. This is rarely necessary, because the tests described above are
sufficient to ensure proper model specification. Its chief utility is as a means for comparing
or debugging numerical integration engines.