Cell Howto

How to create the specification of the model cell in a form that is easy to use with the MRF.

Anatomy

The model used for this exercise borrows the anatomy of the pyramidal cell that is built into neurondemo.

Run neurondemo and select the "Pyramidal" model.
Build / CellBuilder
Management / Import
Click on the big "Import" button, then select "Top Level".
A "warning" message pops up.
Click on "Go ahead and import".
Click on "Turn off indexed name display" and "Don't draw short sections as circles".
Save the CellBuilder to a ses file called cell0.ses, then exit NEURON.

Specify the discretization strategy

Start nrngui and load cell0.ses .
On the CellBuilder's Geometry strategy page, specify the d_lambda discretization strategy for the all subset.
On the Geometry parameters page, verify that the d_lambda criterion is 0.1 length constant at 100 Hz.
Save the CellBuilder to a ses file.

Initial specification of biophysical properties

What mechanisms are present?

Every model cell has Ra and cm.

Since the cell has a linear I-V relationship, and a resting potential of -70 mV, this model cell should also have the pas mechanism, with e_pas = -70 mV.

The exact values of Ra, cm, and g_pas are unknown, but plausible starting values are

Ra = 100 ohm cm
cm = 1 uf/cm2
g_pas = 0.0001 S/cm2
It is easy to use the Biophysics strategy and parameters pages to specify these parameters for the all subset.

How to specify a channel density gradient

The question is how to set up the g_pas gradient in the apical tree.

Save the CellBuilder to a ses file for future use. Next export the model specification to a hoc file:

Management / Export / Export to file
Specify "cell.hoc" as the file name, then click on Save.
Now exit NEURON.

Revising cell.hoc so the MRF can take advantage of its parameters and procedures

At this point, cell.hoc contains hoc statements that, when executed, will create a model cell that has the desired anatomy plus our initial guess at biophysical properties. To use this model cell with the MRF, we will have to tell the MRF the names of the parameters it should adjust, and how to change the gradient of g_pas in the apical tree.

We could do this by writing our own procedures in hoc, but why reinvent the wheel? With model specifications that are complex, the more we can reuse working code and avoid writing our own code, the less likely we are to introduce mistakes.

So let's see if we can reuse the procedures that the CellBuilder created for us in cell.hoc

Before you start, make sure to register cell.hoc with your version control system. If you're not using a version control system (why not?), at least make a copy of cell.hoc and put it in a safe place.
Open cell.hoc in a text editor. Note that it begins with
load_file("subiter.hoc")

proc celldef() {
  topol()
  subsets()
  geom()
  biophys()
  geom_nseg()
  biophys_inhomo()
}
celldef() is the "main" procedure in cell.hoc. It is "in charge" of everything that cell.hoc does to create a model cell. Notice that celldef() calls two procedures whose names indicate they may be involved in specifying the biophysical properties of the model: biophys() and biophys_inhomo().

Search for biophys() and find

proc biophys() {
  forsec all {
    Ra = 100
    cm = 1
    insert pas
      g_pas = 0.0001
      e_pas = -70
  }
}
Uh-oh. "Magic numbers" buried deep inside code.
We need to replace the right hand sides of the assignment statements for Ra, cm, and g_pas with symbolic names whose values are declared at the start of cell.hoc. So change biophys() to
proc biophys() {
  forsec all {
//    Ra = 100
//    cm = 1
    Ra = Ra_
    cm = cm_
    insert pas
//      g_pas = 0.0001
      g_pas = g_pas_
      e_pas = -70
  }
}
and insert the following at the start of cell.hoc
Ra_ = 100
cm_ = 1
g_pas_ = 0.0001
Now in order to change Ra to a new value (say, 200), just execute Ra_=200 and then call biophys().

We must also make sure that A0, A, and d can be specified by the MRF. Search for biophys_inhomo() and find

proc biophys_inhomo() {
  // Path Length from root with no translation
  //   and no normalization ranges from 41.2433 to 959.532
  apicals_x = new SubsetDomainIterator(apicals, 0, 0, 0)
  g_pas_apicals_x()
}
Note that the next procedure is
proc g_pas_apicals_x() {local x, p, p0, p1, A0, A, k, d
  apicals_x.update()
  p0 = apicals_x.p0  p1 = apicals_x.p1
  A0 = 0.0001
  A = 0.0009
  k = 0.02
  d = 1000
  for apicals_x.loop() {
    x = apicals_x.x  p = apicals_x.p
    g_pas(x) = A0 + A/(1 + exp(k*(d - p)))
  }
}
More magic numbers, buried deep inside code. And these parameters (A0, A, and d) are all "local" so they aren't visible outside g_pas_apicals_x().

The workaround is identical to what we did to proc biophys(): invent new variables A0_, A_, and d_ that are defined at the beginning of cell.hoc

A0_ = 0.0001
A_ = 0.0009
d_ = 1000
and change the corresponding assignment statements in g_pas_apicals_x()
proc g_pas_apicals_x() {local x, p, p0, p1, A0, A, k, d
  . . .
//  A0 = 0.0001
  A0 = A0_
//  A = 0.0009
  A = A_
  k = 0.02
//  d = 1000
  d = d_
  . . .
}
Finished! cell.hoc now starts with these statements
Ra_ = 100
cm_ = 1
g_pas_ = 0.0001
A0_ = 0.0001
A_ = 0.0009
d_ = 1000

load_file("subiter.hoc")

proc celldef() {
 . . .
so the MultipleRunFitter can control A0_, A_, d_, Ra_, cm_, g_pas_, then call biophys() and biophys_inhomo() to make these parameter changes affect the model.
A final comment: Remember that changing Ra or cm may affect spatial accuracy. After obtaining what looks like a good fit, it is a good idea to verify that the optimized model's spatial grid is sufficiently fine, e.g. execute forall nseg*=3 then run a simulation to see if this affects the results. If it turns out that the grid is not fine enough, it may be necessary to increase nseg in one or more sections, then run another optimization. It may save some time if you do some exploratory simulations in order to discover which sections need larger nseg values.


NEURON hands-on course
Copyright © 1998-2012 by N.T. Carnevale and M.L. Hines, all rights reserved.