Extracellular Resistance as a 3D Grid - Bidomain Model

[shailesh]

Hi,

I assume it best to explain my objective first so that you could suggest a better approach to tackle my problem. I am attempting to develop a bidomain model (one grid representing cells, coupled via gap junctions, and the other grid their adjacent extracellular space). This would be a passive system and 3-dimensional with num(cells)>>>100. I have developed the first grid with gap junctions implemented as point processes via NMODL and it seems to be working fine. My task is now to build the second grid.

To start off, I am trying to keep it simple and just develop a 2 cell model (by a method which I can easily and efficiently scale for my larger 3D models). From what I have read on the forum, it means I should rule out using LinearMechanism. Do correct me if I am wrong here.

The most natural approach seemed to be to utilize 'extracellular' mechanism in NEURON. I would by default obtain the connection between the two grids via the cell membrane. I believe, I would only need one layer of extracellular, i.e. nlayer=1. The documentation stated that:

If other than 2 extracellular layers is desired, you may recompile the program by changing the nrn/src/nrnoc/options.h line #define EXTRACELLULAR 2 to the number of layers desired. Be sure to recompile both nrnoc and nrniv as well as any user defined .mod files that use the ELECTRODE_CURRENT statement.

I made the change to options.h file. But I am not clear on what was meant by "recompile both nrnoc and nrniv". I assumed it meant to run mknrndll on 'nrn71\src\nrnoc' folder. I found 'nrniv' in 'nrn71\bin' and ran mknrndll on that as well (I always assumed mknrndll was only for mod files, does it also operate on dll files?). But even then accessing Tools>Distributed Mechanisms>Viewers>Shape Name>(selecting section) I could still see parameters for the 2nd extracellular layer. How can I verify whether nlayer = 1 or 2?

The next concern was to have the extracellular space of adjacent cells connected/coupled (the ECF at last segment of one cell connected to the ECF of first segment at the other cell). I believe, by default, the ECFs for different sections in NEURON are disconnected. I searched through the forum for any related articles, and one of the likely ones seemed to be regarding ephaptic coupling (https://www.neuron.yale.edu/phpBB/viewtopic.php?f=12&t=1812&p=6503&hilit=ephaptic#p6503). But it seemed to suggest using LinearMechanism which again seems to be infeasible for the larger model here.

A seemingly weird thought that came up, considering that I had already coupled the cells via gap junctions, was to have a similar mechanism for coupling/connecting the extracellular resistances. Though I had my own doubts, I gave it a try. I made another NMODL file (pretty much identical to gap.mod):

Code: Select all

NEURON {
   POINT_PROCESS GapE
   POINTER vext,vgapext
   RANGE r, i
   NONSPECIFIC_CURRENT i
}

PARAMETER {
   r = 1e10 (megohm)
}

ASSIGNED {
   vext (millivolt)
   vgapext (millivolt)
   i (nanoamp)
}

BREAKPOINT {
   i = (vext - vgapext)/r
}


The relevant hoc code being:

Code: Select all

//Coupling Extracellular
objref ge[2]
for i = 0 ,1 {
   ge[i] = new GapE(0.5)
   ge[i].r = 10
}
cella ge[0].loc(0.9999)
cellb ge[1].loc(0.0001)
setpointer ge[0].vext, cella.vext[0](0.9999)
setpointer ge[1].vext, cellb.vext[0](0.0001)
setpointer ge[0].vgapext, cellb.vext[0](0.0001)
setpointer ge[1].vgapext, cella.vext[0](0.9999)


But I realized (should have thought of it in the first place) that the NONSPECIFIC_CURRENT i would only be injected in the corresponding intracellular space and not at the extracellular nodes. Is there any way to overcome the same?

Another tame effort was to change the BREAKPOINT block to:

Code: Select all

BREAKPOINT {
   vext = vext - ((vext - vgapext)/r)
}


But obviously this was also doomed to fail as vext would just get overwritten with values calculated here vs those for extracellular mechanism rather than each contributing to the final value. I am overdoing it by this stage, but can setting e_extracellular of first segment of cell2 = vext[0] of last segment of cell1 (and vice versa) offer any possibility? Apologies in advance if it doesn't make sense, I think I need to take my mind off it for a while.

I realize I am thinking along a single track, so do advice me on what would be a better approach to develop the extracellular grid of the bidomain model.

Thanks,
Shailesh Appukuttan

[ted]

To start off, I am trying to keep it simple and just develop a 2 cell model (by a method which I can easily and efficiently scale for my larger 3D models).

Always a good way to start.

From what I have read on the forum, it means I should rule out using LinearMechanism. Do correct me if I am wrong here.

Done in a different post.

The most natural approach seemed to be to utilize 'extracellular' mechanism in NEURON. I would by default obtain the connection between the two grids via the cell membrane.

Yes. Use extracellular's i_membrane as current sources to drive the nodes in your extracellular grid, and couple the grid's local potentials back to the extracellular membrane via e_extracellular. The default values of extracellular's parameters are such that current through the xraxial and xc elements is vanishingly small, and voltage drop across the xg elements is negligible. Consequently, membrane current flows "radially" to ground through the e_extracellular voltage sources.

I believe, I would only need one layer of extracellular, i.e. nlayer=1.

True in theory, but why divert your effort to compiling from source code unless it's necessary--and you won't know if it's necessary until you run into space or speed issues. I'd advise sticking with plain vanilla NEURON for now.

The next concern was to have the extracellular space of adjacent cells connected/coupled

That is the purpose of your extracellular grid.

A seemingly weird thought that came up, considering that I had already coupled the cells via gap junctions, was to have a similar mechanism for coupling/connecting the extracellular resistances.

If you're thinking "how can I join the extracellular layers of one cell to those of another cell at the site of a gap junction," forget about it. Coupling between cells is to be done by your extracellular grid.

Some investigators have literally split the task into two parts, with NEURON simulating a cell and some other program simulating the surrounding conductive medium. Simulation execution then proceeds like so

Code: Select all

Initialize the model cell and model bulk conductor to t=0
REPEAT
  given the extracellular field and membrane currents at time t
    use NEURON to calculate the membrane currents at t+dt
    use a Laplace solver to calculate the field in the bulk conductor to t+dt
    exchange data between NEURON and the Laplace solver
      i.e. NEURON sends membrane currents at t+dt to the Laplace solver,
      and the Laplace solver sends the extracellular potentials at t+dt to NEURON
  update t to t+dt
UNTIL done

Could all of this be done in NEURON? In principle, yes. The chief problem is how to deal with extracellular space. The task is easiest if the conductive medium is resistive and isotropic, in which case a network model whose cells contain a total of N intracellular nodes (segments) requires an NxN matrix to represent the extracellular interactions between segments. The simplest implementation of a simulation run would then proceed like this:

Code: Select all

Initialize the cells and extracellular potentials at t=0
REPEAT
  given the e_extracellular values and membrane currents at time t
    advance all cells by dt
    from the new i_membrane values, calculate new e_extracellular values
      i.e. by e_ex = Rx * i_m where i_m is the vector of membrane currents,
        Rx is the matrix of transfer resistances between the model's segments,
        and e_ex is the vector of extracellular potentials
          that result from the membrane currents
  update t to t+dt
  update the e_extracellular values from the elements of e_ex
UNTIL done

[shailesh]

A BIG thanks for the point by point response. I did see several other posts where extracellular was being modeled in another software (Matlab/Comsol/Maxwell etc). But I did prefer to stick to NEURON for the same and hence will try the second approach that you have suggested.

And yes I will go ahead with vanilla NEURON (nlayer=2) for the time being to keep things simple.

Thanks again!
Shailesh Appukuttan