Extracellular Mechanism and Hodgkin-Huxley

[ajpc6d]

What does NEURON do with extracellular potential 'under the hood'?
From the documentation, I took it that NEURON implicitly references both membrane voltage 'v' and extracellular potential 'e' when solving the cable equation, but this gave me unexpected results. When I added the extracellular potential to my NMODL equations explicitly, the results were as expected - but maybe my expectation is flawed, and I'm actually misusing NEURON.
Here's the context:

I have a static, simple cell (a sphere, a cylinder) subjected to a pulsed electric field of very short duration (say, one microsecond). My intuition says the cell membrane should briefly depolarize (in some manner controlled by the membrane time constant), and that if this depolarization is sufficient, the induced intracellular charge gradient will, after a perceptible latency, trigger an action potential.

Solving the Hodgkin-Huxley equations with SciPy yields this response. To incorporate the stimulus, a term is added to the transmembrane current density equation, the difference of the extracellular potential at the membrane and the membrane potential, divided by the product of the resistivity and the thickness of the cell on the axis of the propagation of the electric field. Plotting the membrane voltage shows a latency of several milliseconds between the microsecond depolarization and the triggered action potential.

When I tried to enact the same situation in NEURON, there was no latency (or the latency was imperceptible). As I said, the cell in NEURON is governed by the same ODEs, but without the added stimulus term. Adding the extracellular potential explicitly in the form of g_x*(E_x - V_mem + V_ext) brings the results into much better alignment with the SciPy calculation. But if this is redundant, I must need to reassess my intuition of the physics at work here and question my construction of the SciPy equations.

Appreciate any thoughts--

[ted]

When I tried to enact the same situation in NEURON

Ah, but it turned out not to be "the same situation" because there is a significant mismatch between your conceptual model of what is happening in the physical system that you (partially) describe and the computational model that NEURON is dealing with. Let's try to resolve that.

I have a static, simple cell (a sphere, a cylinder) subjected to a pulsed electric field of very short duration (say, one microsecond).

The discussion will be most useful if it begins with the simplest case, avoids ambiguity, and makes all assumptions explicit. To start with, let's talk about the physical system itself, not a computational model of that system. Also, not sphere OR cylinder--for now, just sphere. And what is the geometry of the electric field? Or better, how does extracellular potential vary as a function of position in space?

[ajpc6d]

I have assumed the cell to experience the propagating electric field as a plane wave (i.e., the potential difference between any points on the plane to which the axis of propagation is orthogonal will be zero).

The extracellular potential, I derived using the the formula V_ecm(t) = -E.l_c (from Boshuo Wang et al 2022 J. Neural Eng. 19 066047, DOI:10.1088/1741-2552/acab30), where V_ecm is the extracellular voltage as a function of time, E is the magnitude of the electric field, and l_c is the displacement vector of the considered segment with respect to the soma center.

So far, so good?

[ted]

I'd like to avoid having to look at any of your hoc or Python code, if possible.

Just a few more questions about your conceptual model.

What is the size of the spherical cell that you are considering?

pulsed electric field

Meaning what, exactly? A spatially and temporally complex stimulus like those used by Wang et al. 2023, or a simple square pulse (zero amplitude except when the field is "on") which produces a linear extracellular potential gradient (and what is the magnitude of that gradient?)?

But in addition to these questions about your conceptual model, maybe it is time to bring up a couple of things that pertain to the computational implementation of that model.

First, if you are using one of NEURON's fixed time step integration methods (the default implicit Euler or one of the Crank-Nicholson methods), do you know that the duration of the stimulus must be >= dt? If it isn't, the stimulus will have no effect. That's just a consequence of fixed time step integration.

Second, for extracellular stimulation implemented by driving e_extracellular to be effective, the computational model must have two or more compartments, and they can't all have the same value of e_extracellular. This is a consequence of Gauss's law (also of Maxwell's equations). You and I had this discussion a few months ago; interested readers should see viewtopic.php?p=10782#p10782

If I were modeling the effect of an extracellular field on a spherical cell, I'd probably represent the cell with two sections whose axes are parallel to the extracellular potential gradient, so that each section represented a hemisphere. Also I'd assume the gradient is parallel to the x axis. Then I'd use the pt3d method to specify each section's geometry, placing the connection between the two sections at 0,0,0, and basing the x3d and diam3d values on the locations of the segment boundaries.

[ajpc6d]

I'd like to avoid having to look at any of your hoc or Python code, if possible.

Totally understand - my original question was just about the mechanics of how NEURON handles e_extracellular in a custom Hodgkin-Huxley-style .mod file, specifically, about the use of the NMODL 'v'. Is e_extracellular automatically calculated as part of this quantity? Is it double-dipping to include it separately in an ionic current term [ g_x( E_x - v - e_extracellular) ]? If not, what happens to e_extracellular inside the black box?

Everything was else was context (in case relevant). That said, I can keep providing specifics if that question can't be answered (or if it changes the answers..?)

Size of the cell: 5 um.
E-field: currently simple square pulse, will advance to spatiotemporally complex stimulus next.
Extracellular potential gradient: dependent on pulse duration. At rheobase, SciPy finds excitation at about 10mV depolarization, therefore 2kV/m E-field magnitude.
Time-step: variable, determined by a function to never be greater than 2% the stimulus duration.
Compartmentalization: confirmed - the cell has multiple segments (tested with 2-20).

[ted]

v is not the sum of e_extracellular and anything else. v is transmembrane potential, the integral of transmembrane capacitive current, the electrical potential difference across the cell membrane at a particular location in a particular section. Its complete name is
secname.v(range) in hoc
secname(range).v in Python
where
secname is a user-assigned name of a section
range is a user-assigned name of variable whose value lies in the closed interval 0..1 and specifies the normalized distance of a point of interest from the 0 end of the section.
In an NMODL file, it is just called v and its value is the value of transmembrane potential at the midpoint of a segment in which the density mechanism or point process defined by the NMODL file exists.

e_extracellular is used to represent user-specified perturbations of extracellular potential. It is an ideal voltage source, and its default value is 0. NMODL knows nothing about e_extracellular. Suggest you (re)read the Python or hoc documentation of extracellular at nrn.readthedocs.io

[ajpc6d]

v is not the sum of e_extracellular and anything else

This was my deduction both from the documentation at https://nrn.readthedocs.io/en/latest/ho ... /mech.html --

Vext is solved simultaneously with the v. When the extracellular mechanism is present, v refers to the membrane potential and vext (i.e. vext[0]) refers to the extracellular potential just next to the membrane. Thus the internal potential is v+vext

-- and from the different results between methods.

My question then becomes, why does the transmembrane voltage even react if I play a non-zero amplitude pulse through e_extracellular, but never do anything with it? Shouldn't this be a case of a tree falling with no one to hear it?

[ted]

why does the transmembrane voltage even react if I play a non-zero amplitude pulse through e_extracellular, but never do anything with it?

Good question. I bet you can figure this out for yourself with a toy example. Consider a model that consists of two adjacent segments. Draw the equivalent circuit. Simplify the diagram by

omitting all of extracellular's resistances whose default numeric values are 1e9
omitting all of extracellular's capacitances whose default numeric values are 0
replacing all of extracellular's conductances whose default numeric values are 1e9 with "wires"

This means omitting all xraxial and xc elements, and replacing all xg elements with "wires"

Now look at the circuit and predict what happens if one of the e_extracellular elements increases to 1 while the other one stays equal to 0.

[ted]

2kV/m E-field magnitude

For a model of a 5 um diam spherical cell with hh membrane properties, I'm getting 3.14522 kV/m; for a 10 um diam cell it's half as much (as expected). With extracellular potential gradient 0, rheobase to 4.5 places is 0.0017373 nA, and comparing the resulting trace with the response to a stimulus of 0.0017372 nA (which fails to trigger a spike) could lead one to estimate "threshold potential" to be about -55 mV (10 mV above rest).

That said, do two more runs that add one more significant figure to the estimate of rheobase--one with stimulus 0.00173724 nA (subthreshold) and the other with 0.00173725 nA (suprathreshold)--and the apparent "threshold potential" becomes more like -53 mV (12 mV above rest). An observed value that is so sensitive to experimental manipulation (changing stimulus amplitude by 1 part in ~170000 changes the "threshold potential" by 20%!) is a questionable peg from which to hang a hat.

[ted]

https://www.neuron.yale.edu/ftp/ted/neu ... sphere.zip contains source code for demonstrations of extracellular stimulation of a spherical cell with HH membrane properties in an electric field that induces a linear extracellular potential gradient parallel to the x axis. One implementation is entirely in hoc, the other in Python. Expand the zip file, cd into the resulting directory, execute nrnivmodl, then execute either
nrngui stimsphere.hoc
or
python3 -i stimsphere.py

In either case, launch a simulation by clicking on the RunControl panel's Init & Run button. See the comments in the hoc and py files for more details. The soma diameter is 5 um by default, but can be changed by editing sphere.hoc. Stimulus intensity (potential gradient) can be changed by executing makestim() with appropriate arguments at the oc> or >>> prompt. Be sure to read the hoc or Python references at nrn.readthedocs.io as needed when encountering unfamiliar keywords.

[ajpc6d]

Wow, that's incredible - thank you! I've been working to meet a deadline this week, but hoping to return to this issue first thing next week--