Non-linear GHK leak for Potassium

[Kadhim]

Hi all, I've just started using neuron for my research project. I'm trying to introduce a non-linear potassium leak using the GHK model but I'm having some problems with NMODL. My code (pasted below) returns the following error:

Code: Select all

Translating ghk.mod into ghk.c
val does not have a proper declaration at line 38 in file ghk.mod
^
make: *** [ghk.c] Error 1

Any help would be really appreciated.

Code [As you can see, I had trouble defining 'coulomb' and so went for a horrible bodge]

Code: Select all

: GHK Leak

NEURON {
	SUFFIX kd
	USEION k READ ek WRITE ik
	GLOBAL cin, cout, val
}

UNITS {
	(S)=(siemens)
	(mV)=(millivolt)
	(mA)=(milliamp)
	(mM)=(millimoles)
	(mol)=(moles)
	FARADAY = (faraday)  (1 coulomb)
	(kel)=(kelvin)
	(j)=(joules)
}

ASSIGNED {
	v (mV)
	:ek (mV)
	:ik (mA/cm2)
	i (mA/cm2)
}

PARAMETER {
	cin = 145 (mM)
	cout = 3 (mM)
	ghk = 1 (S/cm2)
	f = 96400 (faraday/mol)
	r = 8.314 (j/kel/mol)
	temp = 310 (kel)
}

BREAKPOINT {
	i = -v*(ghk/cout)*((cin-cout*exp(-v*f/r/temp))/(1-exp(-v*f/r/temp)))
}

[ted]

First let me thank you for using NEURON in your research.

NMODL gives one a great deal of rope with which it is all too easy to hang oneself. Does your model include a potassium accumulation mechanism (a mechanism with a
USEION k ...
statement that WRITEs ki or ko)?

[ted]

I suppose you already know that this
i = -v*(ghk/cout)*((cin-cout*exp(-v*f/r/temp))/(1-exp(-v*f/r/temp)))
deviates from the GHK equation for current (in particular the ghk/cout bit). The correct formula appears in many places, including
eqn. 5.3 on p 123 of the book "Ion Channels" by Aidley & Stanfield
and
eqn. 5.3.1 on page 123 of the book "Foundations of Cellular Neurophysiology" by Johnston & Wu.
(what a coincidence of equation and page numbers!)

Following the principle of

whenever possible, steal [great] code (Anon.)

I resorted to examining the NMODL examples that accompany NEURON; these are in
nrn/share/examples/nrniv/nmodl
of the expanded nrn*tgz source code, and MSWin users will find them in
c:\nrnxx\examples\nrniv\nmodl
The most relevant files for this case are called fh.mod and standard.inc. I pilfered the essentials from them to create my own kghk.mod:

Code: Select all

: kghk.mod
: potassium leak governed by Goldman-Hodgkin-Katz equation

NEURON {
  SUFFIX kghk
  USEION k READ ki, ko WRITE ik
  RANGE p, i
}

UNITS {
  (mV) = (millivolt)
  (mA) = (milliamp)
  (molar) = (/liter)
  (mM) = (millimolar)
  FARADAY = (faraday) (coulomb)
  R = (k-mole) (joule/degC)
}

PARAMETER {
  v (mV)
  celsius (degC)
  p = 1.2e-3 (cm/s)
}

ASSIGNED {
  i (mA/cm2)
  ik (mA/cm2)
  ko (mM)
  ki (mM)
}

INITIAL {
  i = p*ghk(v, ki, ko)
  ik = i
}

BREAKPOINT {
  i = p*ghk(v, ki, ko)
  ik = i
}

FUNCTION ghk(v(mV), ci(mM), co(mM)) (0.001 coul/cm3) {
  :assume a single charge
  LOCAL z, eci, eco
  z = (1e-3)*FARADAY*v/(R*(celsius+273.15))
UNITSOFF
  eco = co*efun(z)
  eci = ci*efun(-z)
UNITSON
  ghk = (0.001)*FARADAY*(eci - eco)
}

FUNCTION efun(z) {
  if (fabs(z) < 1e-4) {
    efun = 1 - z/2
  }else{
    efun = z/(exp(z) - 1)
  }
}

The advantages of this particular code, aside from the fact that it expects one to specify an actual potassium permeability, are (1) it recognizes celsius as the temperature parameter, and (2) it will "play well with other mechanisms." These are important considerations if you ever want to use it in a model that has any other temperature-dependent mechanisms, or that has some other mechanism that WRITEs ki or ko. Otherwise you'll end up with code that has "magic numbers" (to specify temperature and concentrations) buried in various files, which is a code management headache that can lead to absurdities such as multiple mechanisms in the same sections having different and mutually incompatible values for temperature, equilibrium potentials, and concentrations.

"But my model is very simple and lacks any mechanism that WRITEs ki or ko. Can I specify those from hoc?"

Yes, because if a section lacks any mechanism that WRITEs an ionic concentration, that concentration is simply a range variable parameter, not a state variable, and as such it can be specified by a simple assignment statement. Here's an example of such a model:

Code: Select all

load_file("nrngui.hoc")

create soma
access soma

L = 10
diam = 10/PI // area = 100 um2
// so total current in nA = current density in mA/cm2

insert kghk

celsius = 37 // deg C

// note:  the following assignment statements
// will suffice for setting intra and extracellular
// k concentrations as long as there is no mechanism
// that WRITEs ki or ko
ki = 145 // mM
ko = 3

And here is a bit of "instrumentation and control code" that will exercise this model, sweeping membrane potential from -100 mV to +100 mV in 2 mV increments and plotting i_kghk at each potential. It also sets up a control panel so you can interactively assign new values to temperature, concentrations, and permeability, and generate new plots to see what happens:

Code: Select all

// argument is membrane potential
// returned value is current at that potential
func ikghk() {
  finitialize($1)
  return i_kghk
}

objref g
g = new Graph(0)
g.size(-100,100,0,65)
g.view(-100, 0, 200, 65, 0, 105, 300.48, 200.32)
g.addexpr("ikghk(v)", 1, 1, 0.2, 1, 2)

proc doplot() {
  v = -100
  g.begin()
  for (v=-100; v<=100; v=v+2) g.plot(v)
  g.flush()
  g.exec_menu("View = plot")
}

doplot()

{
xpanel("GHK leak",0)
xvalue("celsius","celsius", 1,"", 0, 1 )
xvalue("ki","soma.ki", 1,"", 0, 0 )
xvalue("ko","soma.ko", 1,"", 0, 0 )
xvalue("p_kghk","soma.p_kghk", 1,"", 0, 0 )
xbutton("New plot","doplot()")
xpanel(327,105)
}

[Kadhim]

Sorry for not replying sooner, Ted.
Thank you so much for your help, I'm sincerely very grateful.