Rattay F (1986): Analysis of models for external stimulation of axons.
IEEE Trans. Biomed. Eng. BME-33:(10) 974-977
It is the approach taken by McIntyre and Grill in their paper "Extracellular stimulation of central neurons: influence of stimulus waveform and frequency on neuronal output" JNP 88:1592-1604, 2002.
Another method is to assume that the extracellular medium is linear, in which case the coupling between stimulating electrode(s) and locations in space can be represented by transfer impedances. The extracellular potential at a point x,y,z produced by stimulus currents applied through any number of electrodes is then
Code: Select all
Vo(x,y,z) = SUMMA Z_j(x,y,z) I_j Eq. 1 j
Of course, linearity also implies that the extracellular potential recorded by electrode j is composed of the weighted integral of the membrane current over the surface of the cell (weighted by local membrane area and transfer impedance between each point on the cell and the location of the recording electrode). For a spatially discretized model cell, this reduces to the weighted sum
Code: Select all
V_j = SUMMA Z_k I_k Eq. 2 k
If the reactive (capacitive) component of Z can be ignored, then Eq. 1 and 2 are easily implemented via simple multiplications.
The approach outlined by Equations 1 and 2 is illustrated in this file
http://www.neuron.yale.edu/ftp/ted/neur ... nd_rec.zip
which contains complete working implementations of extracellular stimulation and recording of a single neuron model, and could serve as a starting point for a model that involves multiple cells.
To implement either of these approaches in NEURON, one must:
--determine the locations of the nonzero area nodes in a model (i.e. the segment centers). This is easy enough, if somewhat tedious, for a stylized (L, diam) specification of model geometry. For geometries defined by the pt3d method, interpolating node locations from the pt3d data is straightforward (see the zip file).
--force some variable (an activation function current, or e_extracellular) in each segment to follow an appropriate time course. This can be done with an analytical function computed by code in a Point Process, an event handler, or with the Vector class's play() method. See the above links for examples.