Accurate simulation of cuff electrode stimulation predicting in-vivo strength-duration thresholds

Abstract

Background: In-silico experiments used to optimize and inform how peripheral nerve based electrode designs perform hold the promise of greatly reducing the guesswork with new designs as well as the number of animals used to identify and prove promising designs. Given adequate realism, in-silico experiments offer the promise of identifying putative mechanisms that further inform exploration of novel stimulation and recording techniques and their interactions with bioelectric phenomena. However, despite using validated nerve fiber models, when applied to the more complex case of an implanted extracellular electrode, the in-silico experiments often do not compare quantitatively with the results of experiments conducted in in-vivo experiments. This suggests that the accuracy/realism of the environment and the lamination of the nerve bundle plays an important role in this discrepancy. This paper describes the sensitivity of in-silico models to the electrical parameter estimates and volume conductor type used.

Methods: In-vivo work was performed on rat vagus nerves (N = 2) to characterize the strength-duration curve for various peaks identified in a compound nerve action potential (CAP) measured via a needle electrode. The vagus nerve has several distinct populations of nerve fiber calibers and types. Recruitment of a fiber caliber/type generates distinct peaks that can be identified, and whose conduction delay correlates to a conduction velocity. Peaks were identified by their recruitment thresholds and associated to their conduction velocities by the conduction delays of their peaks. An in-silico analog of the in-vivo experiment was constructed and experiments were run at the two extreme volume conductor cases: (1) The nerve in-saline, and (2) the nerve in-air. The specifically targeted electrical parameters were extraneural environment (in-air versus saline submersion), the resistivity (ρ) of the epineurium and perineurium, and the relative permittivity (ε_r ) of those same tissues. A time varying finite element method (FEM) model of the potential distribution vs time was quantified and projected onto a modified McIntyre, Richardson, and Grill (MRG), myelinated spinal nerve, active fiber model in NEURON to identify the threshold of activation as a function of stimulus pulse amplitude versus pulse width versus fiber diameter. The in-silico results were then compared to the in-vivo results.

Results: The finite element method simulations spanned two macro environments: in-saline and in-air. For these environments, the resistivities for low and high frequencies as well as two different permittivity cases were used. Between these 8 cases unique cases it was found that the most accurate combination of those variables was the in-air environment for low-frequency resistivity (ρ₀ ) and ex-vivo a measured permittivity (ε_r,measured ) from unpublished ex-vivo experiments in canine vagal nerve, achieving a high degree of convergence (r² = 0.96). As the in-vivo work was conducted in in-air, the in-air boundary condition test case was convergent with the in-silico results.

Conclusions: The results of this investigation suggest that increasing realism in simulations begets more accurate predictions. Of particular importance are (ρ) and extraneural environment, with reactive electrical parameters becoming important for input waveforms with energy in higher frequencies.

Keywords: NEURON; finite element analysis; neural simulations; threshold prediction.

FIGURE 1

In‐vivo experimental configuration showing electrode placement on the rat vagus nerve and instrumentation used for pulse activation and recording the CAP along with changes to the rat's heart rate and blood pressure. The inset shows the proximal ligature, the CorTec bipolar cuff electrode used to initiate the CAP and recording needle electrode on the rat's left cervical vagus nerve. Adapted from Ref. [12]. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 2

Example of CAP recordings. Shown using conduction delay to the 5th peak for ease of visualization. This study used the thresholds of peaks 1 & 2 which correlate to myelinated A‐fibers with calibers of ~11.5 and ~6 μm were analyzed. Artifacts have been removed by shape subtraction and the waveforms are vertically shifted for viewing purposes, with dark gray showing the remaining stimulus artifact and lighter gray indicating the section shown to the right. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 3

(A) Shows the in‐silico geometry of the cuff about the nerve in a two‐dimensional axial symmetric model. (B) Shows the steady state results of the simulation and the iso‐potential lines. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 4

In‐air, low frequency resistivity (ρ ₀), measured relative permittivity (ɛ _r,measured) 400 μs pulse width weight and activating functions. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 5

In‐saline environment (A) and in‐air environment (B) 400 μs pulse width weight and activating functions. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 6

A log–log plot of the in‐vivo thresholds for peaks 1–3 and in‐silico strength‐duration curves using the 11.5 μm fiber diameter MRG model for different parameters. The peaks of the in‐vivo measurements were estimated to have conduction velocities as follows: Peak 1 (34.3 m/s), Peak 2 (8.1 m/s), and Peak 3 (1.4 m/s). The fiber diameter estimated for Peak 1 is ~9.5 μm. However, peaks 2 and 3 are outside the conduction velocity range of simulated MRG fibers. There is clear separation between the saline based simulations and the in‐air based simulation. Although all curves are within a factor of 10 to those measured in‐vivo the best fit is with the in‐air case for the SD curves for peaks 1 and 2. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 7

A log–log plot of the in‐vivo and in‐silico in‐air environment only SD curves using the 11.5 and 5.7 μm fiber diameter MRG model. [Color figure can be viewed at wileyonlinelibrary.com]

FIGURE 8

|Z| as a function of pulse width for in‐saline environment (top) and in‐air environment (bottom) weight functions. [Color figure can be viewed at wileyonlinelibrary.com]