A Computer Model of Electrical Stimulation of Peripheral Nerves in Regional Anesthesia

The primary goal of this investigation was to model the responses of an active nerve fiber to various extraneural current stimuli delivered through a needle. In general, the model consists of two broad descriptions: (1) one of the needle, the extent of its stimulating electrode(s), and the resulting potential field; and (2) one of the active nerve fiber, including a mathematical model of the membrane ionic currents. Tracking the relative orientation of the needle and nerve is essential and requires a well-defined three-dimensional model space and coordinate system.

Figure 1shows the three-dimensional model space. The outer skin surface is defined as a fixed plane at z = 0. The needle insertion point is defined as an x, y  coordinate on the skin plane. In its default position, the needle is oriented perpendicular to the skin. From this vertical position, the needle depth and orientation, through the insertion and rotation angles, can be precisely controlled. A positive insertion angle swings the tip in the positive x  direction, whereas a positive rotation angle produces a clockwise rotation from the perspective of the clinician. The plane containing the nerve (nerve plane) is parallel to the skin plane; this plane can be adjusted to account for the depth of the nerve. Within this plane, the nerve fiber orientation is fixed along the line y = 0. However, the freedom to orient the needle allows for every possible spatial relationship between the needle and nerve. Although the model geometry is obviously a simplification of real tissue, it provides an acceptable approximation for many cases of interest in which, over meaningful distances (i.e. , several centimeters) at the point of needle insertion, nerves run parallel to the skin surface.

The model needle can have varying amounts of insulation, from being fully insulated with just the tip exposed to a complete absence of insulation (uninsulated). An electrode is any contiguous segment of exposed needle, which may act as either the anode or the cathode. An electrode is defined through a distal value and a proximal value relative to the needle tip. For example, identical values of 0.0 cm would represent a needle with just its tip exposed.

The needle electrodes, depending on their length, are modeled by one or more discrete current points equally spaced along the line that defines the needle core. The initial spacing of the points is set equal to the needle radius then slightly adjusted to ensure points at the electrode edges. Electrodes shorter than the needle radius are modeled as single current points centered between the edges. The condition that the conductive needle shaft must be isopotential allows the current amplitude at each point to be determined using available mathematical algorithms. The result is a well-behaved current distribution that closely approximates the continuous properties of real needles. (A related problem, in the context of the inverse problem in electrocardiology, is discussed in a text on bioelectricity by Gulrajani.10)

The discrete current points that define the needle electrode(s) are used to calculate the potential field in the region below the skin plane (i.e. , the tissue space). The presence of the skin imposes a boundary condition in which current introduced into the tissue space does not flow across the skin plane (i.e. , the conductivity of air is assumed to be zero). This condition is satisfied through an application of the image principles from electrostatics, in which a projection (or image) of the needle electrode across the skin plane is included in the computation of the potential field. In other words, each of S current points included along the needle beneath the skin plane will have an associated point above the skin, and the resulting extraneural potential field V^e(x, y, z ), in millivolts, can be found from

where (x  ^s, y  ^s, z  ^s) is the location of the s  th current point, ρ^eis the volume (or extraneural) resistivity (Ω· cm) of the model space, and I^stimis the stimulus current (in mA). Consequently, an electrode, or any portion of an electrode, defined along the needle but lying above the skin plane is not included in the above calculation.

The typical clinical procedure attaches the cathode of a stimulator to the needle and the anode to a patch electrode on the skin surface. This latter electrode, relatively far from the site of stimulation along the nerve, is not explicitly defined in the model. Instead, this remote electrode is assumed to be an infinite distance away from the current points in the model space. This follows from the assumptions implicit in Equation 1, which ensure that the total charge in the model is conserved (i.e. , the sum of all the sources and sinks is zero).

The geometry of the one-dimensional model nerve fiber is diagrammed in figure 1(bottom). The fiber is of finite length with an outer diameter δ^fiber. To model excitation, the fiber is discretized into N segments of length L^node, the spacing between the nodes of Ranvier. Each segment contains a single node of length x  ^nodewith diameter δ^axon(set to 0.7δ^fiber).11Nodes are flanked by myelinated sheath, which is assumed to be a perfect insulator. Therefore, current crossing the fiber membrane is confined to the nodes. The nodal membranes are modeled as a capacitance c^node(μF/cm²) in parallel with a series of ionic channels; current flows between nodes via  the axoplasmic space, defined with a resistivity ρⁱ.

Making use of the electrical description above, the response of the fiber to the needle current stimulus is tracked through changes in the transmembrane potential V^m(mV) at each node, where V^mis defined as the difference between the intracellular potential Vⁱand the extracellular potential V^e. Of particular interest are stimuli that generate propagating action potentials along the model nerve; in the context of PNS, these are stimuli that produce a muscle twitch. Solving for the time course of V^mat each node requires equating two expressions for the total transmembrane current. The system of equations to solve can be expressed in matrix form as

where variables in bold are vectors, i.e. , V^m= {V^1m,V^2m,…,V^Nm}.

The left-hand side of Equation 2is the sum of the capacitive and ionic currents at each node. The total ionic current i^ion  (mA/cm²) can be expressed as

the superposition of the currents from each individual ion (e.g. , sodium, Na⁺, and potassium, K⁺), including any passive, leakage currents. The driving force for each local ionic current is the difference between V^mand the constant Nernst potential E^xassociated with that ion. The conductances (g^x , S/cm²) require a mathematical description of the gating mechanisms for the ionic channels—a set of nonlinear, ordinary differential equations that depend on V^m.

A good review of the mathematical details of nerve membrane models can be found elsewhere.12In general, the data required to describe these models involve extensive in vitro  work.13In this article, results were generated by using a modified version of the Frankenhaeuser-Huxley model based on voltage clamp data from myelinated frog axons.14The Frankenhaeuser-Huxley model was modified by adjusting the Nernst potentials, the resting membrane potential, and the rate constants governing the gating mechanisms of the ionic channels for body temperature. (The data underlying the original Frankenhaeuser-Huxley equations were recorded at 20°C.) These adjustments and their effects on current thresholds for nerve activation are discussed in detail by Rattay et al.  15 

The right-hand side of Equation 2expresses the transmembrane current at each fiber segment arising from the spatial differences in the transmembrane potential (GV  ^m) and the stimulus current delivered through the needle (GV  ^e), where V  ^eis the distribution of extracellular potentials along the fiber (measured at the center of each segment) calculated from Equation 1. G is an N × N symmetric, tridiagonal conductance matrix that defines the passive electrical connections between each segment and its neighboring segments, or

Special attention must be paid to the fiber ends. The fiber length, which can be varied, should be long enough to ensure that the current thresholds for activation are not affected by the terminal segments of the fiber. (The boundary condition we used set the axial current at the terminal segments to zero.) End effects are also minimized by positioning the fiber such that its central segment is always the closest point to the needle tip. The results we present were generated with a 14-cm fiber, found to be a conservative value.

The mathematical framework outlined above was implemented in the Matlab® computing environment (The MathWorks, Natick, MA). Users interact with the “PNS Simulation System” through a graphical user interface organized around a central main menu. From this menu, additional windows can be opened. Some of these windows group related model parameters that can be manipulated by the user (i.e. , model inputs), including: the position and orientation of the needle; the location and extent of the needle electrode(s); the electrical properties of the nerve fiber; and the pulse width, amplitude, and frequency of the stimulus. Out of range or otherwise invalid inputs are prevented, and appropriate explanations are provided as feedback to the main menu. Other windows provide model output or display information about the current state of the system.

The software provides two primary results: 1) two-dimensional maps of the potential field around the needle; and 2) computation of the nerve fiber response to the stimulus. The latter result requires numerical integration of Equation 2. Given an initial condition (the entire fiber at a resting potential), the transmembrane potential V^malong the fiber is tracked over a series of discrete time steps, and the results are animated as the computation progresses. Of course, a given stimulus may or may not produce a suprathreshold response. If it does, impulses (action potentials) may propagate in both directions along the fiber, or they may be blocked in one or both directions. Model feedback (to the main menu upon completion of a simulation) includes the distance each action potential propagated along the fiber and its conduction velocity.

To evaluate the model, a series of simulations were performed for an insulated needle and then for an uninsulated needle approaching a nerve fiber positioned 2.0 cm below the skin. For each needle design, data were collected for the needle attached to the cathode and then repeated with the needle attached to the anode. The needle was oriented perpendicular to the skin plane and inserted at (0, 0.05). The 0.05-cm offset (y  ^insert) allowed the needle, with an outer diameter of 0.1 cm (19- to 20-gauge), to be advanced while avoiding direct contact with the nerve. The other model parameter values where chosen to match in vivo  conditions as closely as possible (table 1).11The protocol for each simulation was as follows:

1. The needle depth that just produced a propagated response along the fiber at a maximal current amplitude of 5 mA (7 mA for the uninsulated anode case) was found.

2. The needle was advanced 0.1 cm from this initial depth toward the nerve.

3. At the new depth, the lower boundary (i.e. , threshold) for generating a propagated impulse was found at a current resolution of 0.001 mA, as well as an upper boundary if it existed.

4. Steps 2 and 3 were repeated until the needle tip was beyond the nerve plane to where the maximal current amplitude no longer activated the nerve (or, for the uninsulated cathode, to approximately 8 mm beyond the nerve).

The data were used to construct current-distance maps that highlight combinations of current amplitude and tip-to-nerve distance that evoked a propagated response.

Two additional features of the PNS Simulation System were used to speed up data collection. First, simulations were run as batch processes. This allowed the above protocol to be performed automatically (i.e. , with minimal user interaction) by defining a range of needle depths and current amplitudes and having the resulting propagation boundaries at each depth output to a file. Second, an algorithm was incorporated to locate the boundaries with a minimal number of simulations. In short, the results from previous steps in the protocol were used to determine the next current amplitude. For the results presented herein, an average of approximately 17 simulations was required at each needle depth to find the propagation boundaries.

The current-distance maps for the four needle permutations (insulated cathode/anode, uninsulated cathode/anode) are presented in figure 2. The lower boundary of the outlined and shaded region represents the current thresholds typically reported in the literature—the line that clinicians follow when using the PNS technique. The insulated needle attached to the cathode (fig. 2A) activated the nerve fiber at the smallest current thresholds for a given tip-to-nerve distance. For example, when the needle tip was 2 mm proximal to the nerve plane, the insulated cathode required a current greater than 0.457 mA to produce a propagated response, compared with 2.395 mA for the uninsulated cathode (fig. 2C) and 2.354 mA for the insulated anode (fig. 2B). When the uninsulated needle was attached to the anode (fig. 2D), currents up to 7 mA were inadequate to produce a response at this tip-to-nerve distance.

Interestingly, for cathodal stimulation (fig. 2, A and C), there was also an upper boundary that existed for close tip-to-nerve distances. (±1.1 mm relative to the nerve plane for the insulated needle; 0.3 to -1.0 mm for the uninsulated needle). Therefore, despite delivering what would seem to be an adequate current with a needle attached to the cathode (the typical clinical scenario) and in close proximity to the nerve, interactions between the stimulus current and the nerve fiber blocked the propagation of the generated action potential. For example, when the insulated needle was 1 mm proximal to the nerve plane, the lower threshold was 0.177 mA but propagation was blocked for currents greater than 2.37 mA. Viewed another way, if a simulation were performed with the current held constant at 2 mA, the same needle would first produce propagation when the needle tip was approximately 4.6 mm from the nerve plane; however, the response would cease when the needle was advanced to within 1.0 mm of the nerve plane.

An important difference between insulated (fig. 2, A and B) and uninsulated needles (fig. 2, C and D) was how the current-distance relationship changed as the needle tip was advanced beyond the nerve. For the insulated needle, the current-distance data had a high degree of symmetry across the nerve plane; some asymmetry was introduced by the presence of the skin boundary, but at a nerve depth of 2 cm, this effect was negligible at the resolution of the current-distance maps. However, in the case of the uninsulated needle, a pronounced asymmetry was introduced by the distribution of the delivered stimulus current along the length of the submerged shaft (fig. 3). As the tip of the uninsulated needle moved past the nerve plane, the current leaking from the shaft in close proximity to the fiber kept the lower boundary in the current-distance maps relatively flat. For the uninsulated needle, the minimal current at the lower boundary occurred when the needle tip was just past the nerve plane. The minimal current was 0.461 mA for the uninsulated cathode with the tip 0.1 mm past the nerve plane and 4.100 mA for the uninsulated anode with the tip 0.2 mm past the nerve plane. For the insulated needle, the minimal current at the lower boundary occurred with the needle tip positioned in the nerve plane—0.061 mA with the needle attached to the cathode and 0.458 mA with the needle attached to the anode.

We present a computer (mathematical) model of peripheral electrical nerve stimulation that draws extensively from past work in electrophysiology on field stimulation of nerve axons lying in a large volume conductor.11,15–17It extends this work by including needle electrodes of varying lengths and incorporating the ability to track the relative position of the needle and nerve below the skin boundary. The model was implemented in Matlab® to provide a graphical interface for manipulating model parameters and viewing model data.

Mathematical models of nerve fiber stimulation are limited in their ability to account for all the complexities of real tissue. Therefore, it is important to discuss the assumptions inherent in the mathematical framework we present. However, despite model limitations, the current-distance data it produced were largely consistent with the clinical performance of insulated and uninsulated needles.18–20This supports our view that many of the fundamental characteristics of nerve fibers, needle electrodes, and the surrounding tissue have been faithfully represented. Initial results indicate that the model can indeed produce data that mimic relevant experimental and clinical situations.

The two-dimensional potential maps generated by the model (fig. 3A) are similar to some previously published maps in which the potential field around a needle attached to a current stimulator was sampled either through manual measurements in a shallow saline bath or by immersing the needle in an electrophoresic gel.18,21One benefit of the model compared with these experimental approaches is the ability to generate two-dimensional potential maps for a range of needle designs quickly and at finer resolutions. Furthermore, in vitro  measurements do not directly link the electrical footprint of the needle to the actual process of nerve activation, whereas the model allows the potentials generated by a given needle design to be associated with current-distance data. However, the model does not take into account the three-dimensional geometries (including various tip shapes or curvatures) of real needles, but instead represents the needle electrodes as a finite number of discrete current points equally spaced along a single line.

Current-distance data are difficult to obtain in vivo .22One study conducted by Ford et al.  compared insulated and uninsulated needles for locating the sciatic and saphenous nerves in anesthetized cats.23Needles were attached to the cathode of a Grass stimulator, and their movement was controlled with a micromanipulator. Their current-distance plots qualitatively match those of our model, but the thresholds found by Ford et al.  were significantly higher. For example, the insulated needle attached to the cathode activated the sciatic nerve with a minimal current of 0.9 ± 0.4 mA when the needle tip was touching the nerve, or 14.8 times the model value at 0.05 cm. For the uninsulated needle, the average minimal current was 2.1 ± 1.0 mA when the tip was 0.6 ± 0.3 cm from the nerve, or 4.5 times the model value at 0.05 cm. When an insulated needle was attached to the anode, they reported that 4.3 ± 0.4 times the cathodal current was required to reach threshold (i.e. , a minimal current of 3.87 mA, or 8.4 times the model value).

The quantitative differences between the data of Ford et al.  and the model data may reflect differences in the two protocols (exact needle approaches and tip-to-nerve distances), as well as differences between the multi-fascicle nerves present in vivo  and the single myelinated fiber in the model. The model fiber has a uniform geometry with constant inner and outer diameters, equidistant nodes of Ranvier of equal length, homogeneous intracellular resistivity, and homogeneous ionic channel properties. The model fiber lies in a tissue space defined by a planar skin boundary and having a uniform resistivity. In contrast, most peripheral nerves targeted in regional anesthesiology contain a distribution of both myelinated and unmyelinated axons closely packed and surrounded by a fascia layer (or other tissue inhomogeneities) that may alter current flow, and the skin through which the needle is inserted may not be the only significant—or even the closest—boundary to current flow. Furthermore, the return electrode is not explicitly included in the model but is assumed to be an infinite distance away from the needle current points; this may be a reasonable simplification considering how the stimulus-induced potential field is inversely proportional to source-field distance (Equation 1) but may not accurately reproduce some clinical cases in which the target nerve is superficial and the skin electrode is placed close to the needle insertion point.

The limitations inherent in the model must be considered when contemplating one of the more interesting and pregnant results: that relatively low-amplitude currents delivered through an insulated needle attached to the cathode can produce propagation block at close tip-to-nerve distances. Clinically, this means that the loss of muscle contractions may indicate not only too little current at too great a distance from the nerve, but also too much current delivered too close to the nerve. Some current amplitudes that blocked propagation in the model would more likely produce a graded response when stimulating multi-fascicle nerves in vivo  (i.e. , a decrease in the intensity of muscle contractions through propagation block in a subset of the available motor fibers). However, the paradoxical cessation of the evoked muscle response may need to be considered when using the PNS technique. Conceptually, this could explain clinical situations in which the needle is advanced quickly at relatively high current but fails to elicit a response until gradually withdrawn.

A recent study by Urmey and Stanton sought to relate the elicitation of paresthesia with the ability to elicit a motor response with PNS.24In 30 patients, reported paresthesia was immediately followed by an attempt to elicit a motor response. Despite 30 paresthesias leading to 30 clinically successful blocks, there were only nine cases of recorded motor responses after paresthesia, with eight of those patients in the group in which an uninsulated needle was used. They concluded, perhaps correctly, that although the needle was in direct contact with the nerve root, the needle tip was too far from the motor nerve components; therefore, the uninsulated needles produced a higher incidence of motor response through stimulation from their shaft. It is interesting to interpret these data in light of the propagation block seen in the model. Given that the needle is brought into contact with the nerve before stimulation is initiated, the current amplitude may be too high to produce a motor response with the insulated needle tip. The minimal current that resulted in propagation block in the model with an insulated needle was approximately 0.32 mA (fig. 2A) with the needle tip 0.05 cm from the nerve; this value would be even lower at closer tip-to-nerve spacings. In contrast, to block propagation, the uninsulated needle required currents greater than 3.28 mA (fig. 2C).

The mechanisms underlying propagation block in field stimulation of nerve fibers have been well defined in the electrophysiology literature and can be understood by studying how the stimulus affects the transmembrane potential V^malong the fiber for several needle tip to nerve distances (fig. 4).15For cathodal stimulation, the region of the fiber closest to the stimulus current point(s) is depolarized; depolarization is associated with activation. Flanking this central region are hyperpolarized regions that work to inhibit activation. The interaction of these regions and the temporal evolution of the fiber response while the stimulus is being delivered can lead to subtle differences in the current-distance data. The widths of these regions (and their amplitudes) depend directly on the distance between the current point(s) and nerve; as the needle approaches the nerve, the depolarized region narrows and the peak depolarization and hyperpolarization increase. Less current is required for activation at closer spacings, but then the greater inhibitory effects of the hyperpolarized regions may block propagation from the central region. When the needle is attached to the anode, the polarity of the response is reversed. The central region is hyperpolarized and propagation is initiated from the smaller depolarized regions on either side; thus, the higher threshold currents and the absence of block (fig. 2B).

Despite the widespread use of PNS, many technical aspects of the technique have not been adequately investigated. This may be partly due to the difficulty in obtaining detailed current-distance data experimentally or clinically. The model we present provides a tool for generating these data quickly for a range of needle designs or other parameters that may affect block performance (e.g. , tissue and nerve properties). Various blocks and approaches can be mimicked easily, and the typical clinical variables (needle depth, needle orientation, and stimulus current) can be virtually manipulated through the model interface. The model implementation has the potential to enhance clinical and experimental data collection or generate new hypotheses about the PNS technique.