Skip to main content
Advertisement
  • Loading metrics

Kilohertz waveforms optimized to produce closed-state Na+ channel inactivation eliminate onset response in nerve conduction block

  • Guosheng Yi,

    Roles Data curation, Formal analysis, Funding acquisition, Methodology, Visualization, Writing – original draft, Writing – review & editing

    Affiliations Department of Biomedical Engineering, Pratt School of Engineering, Duke University, Durham, North Carolina, United States of America, School of Electrical and Information Engineering, Tianjin University, Tianjin, People’s Republic of China

  • Warren M. Grill

    Roles Conceptualization, Funding acquisition, Methodology, Project administration, Supervision, Visualization, Writing – review & editing

    warren.grill@duke.edu

    Affiliations Department of Biomedical Engineering, Pratt School of Engineering, Duke University, Durham, North Carolina, United States of America, Department of Electrical and Computer Engineering, Pratt School of Engineering, Duke University, Durham, North Carolina, United States of America, Department of Neurobiology, School of Medicine, Duke University, Durham, North Carolina, United States of America, Department of Neurosurgery, School of Medicine, Duke University, Durham, North Carolina, United States of America

Abstract

The delivery of kilohertz frequency alternating current (KHFAC) generates rapid, controlled, and reversible conduction block in motor, sensory, and autonomic nerves, but causes transient activation of action potentials at the onset of the blocking current. We implemented a novel engineering optimization approach to design blocking waveforms that eliminated the onset response by moving voltage-gated Na+ channels (VGSCs) to closed-state inactivation (CSI) without first opening. We used computational models and particle swarm optimization (PSO) to design a charge-balanced 10 kHz biphasic current waveform that produced conduction block without onset firing in peripheral axons at specific locations and with specific diameters. The results indicate that it is possible to achieve onset-free KHFAC nerve block by causing CSI of VGSCs. Our novel approach for designing blocking waveforms and the resulting waveform may have utility in clinical applications of conduction block of peripheral nerve hyperactivity, for example in pain and spasticity.

Author summary

Many neurological disorders, including pain and spasticity, are characterized by undesirable increases in sensory, motor, or autonomic nerve activity. Local application of kilohertz frequency alternating currents (KHFAC) can effectively and completely block the conduction of undesired hyperactivity through peripheral nerves and could be a therapeutic approach for alleviating disease symptoms. However, KHFAC nerve block produces an undesirable initial burst of action potentials prior to achieving block. This onset firing may result in muscle contraction and pain and is a significant impediment to potential clinical applications of KHFAC nerve block. We present a novel engineering optimization approach for designing a blocking waveform that completely eliminated the onset firing in peripheral axons by moving voltage-gated Na+ channels to closed-state inactivation. Our results suggest that the resulting KHFAC waveform can generate electric nerve block without an onset response. Our approach for optimizing blocking waveforms represents a novel engineering design methodology with myriad potential applications and has relevance for the conduction block of peripheral nerve hyperactivity.

Introduction

Multiple neurological disorders are characterized by undesirable increases in sensory, motor, or autonomic nerve activity, for example pain and spasticity. Blocking the conduction of activity in peripheral nerves has the potential to alleviate disease symptoms [14]. Kilohertz frequency alternating current (KHFAC) stimulation is an effective method for producing a rapid, controlled, and locally-acting conduction block [48]. KHFAC nerve block is quickly reversible, and thus enables both temporal and spatial control of the target tissue with minimal side effects. However, an important drawback of KHFAC block is the transient burst of action potentials (APs) activated at the onset of blocking currents [3, 6, 812]. This onset firing is likely to lead to intense muscle contraction and or painful sensations and is a significant impediment to potential clinical applications of KHFAC nerve block. Several approaches are proposed to reduce or eliminate the onset response, including transitioning KHFAC waveform from a high amplitude and high frequency to a low amplitude and low frequency [13], slowly increasing the amplitude of the KHFAC waveform [14, 15], introducing direct current (DC) nerve block [12, 16], and modifying electrode geometries [17]. However, these approaches are either unable to eliminate completely the onset firing or are challenging for clinical implementation [1217]. We employed engineering optimization to design a stimulus waveform to achieve conduction block while preventing the onset responses.

Nerve excitation is regulated by voltage-gated Na+ channels (VGSCs). During an AP, the VGSCs first activate from non-conducting closed states to open states when membrane depolarization exceeds a threshold. The open channels allow influx of Na+ ions, which further depolarizes the cell and drives the upstroke of the AP. Membrane depolarization also causes inactivation of VGSCs, which prevents Na+ influx and facilitates repolarization to the resting potential. Na+ inactivation is coupled to activation [1820], and typically occurs from the open state at strongly depolarized membrane potentials [21, 22], also called open-state inactivation (OSI). However, recordings in cardiac cells [20, 23, 24], skeletal muscle [25], neuroblastoma cells [19, 26], squid giant axon [21], medullary raphe neurons [27], and dorsal root ganglion neurons [28, 29] confirm that VGSCs can also inactivate from pre-open closed states without first opening, called closed-state inactivation (CSI). This path occurs at hyperpolarized and modestly depolarized potentials without generation of an AP [18, 22, 24, 26]. Based on these findings, we hypothesized that moving VGSCs into CSI without first causing activation would eliminate the initiation of APs and thus prevent the onset response.

The objective of this study was to design a current waveform to produce KHFAC nerve block without onset firing. We used Markov-type kinetic models [30] to simulate Na+ channel CSI, and the models were previously validated with electrophysiological data from all VGSC isoforms. We first modified the McIntyre-Richardson-Grill (MRG) model [31] of a mammalian peripheral myelinated axon to incorporate Markov-type Nav 1.1 and Nav 1.6 channels. We then used particle swarm optimization (PSO) [32] first to design the profile of membrane voltage required to move the two types of VGSCs into CSI, and, subsequently, to determine the stimulation waveform shape required to drive that voltage profile. PSO is an effective search method for dealing with non-linear optimization problems by moving the particles to a “preferred” position in their constrained domains [32]. It is effective for determining sequences of real numbers in non-linear systems [33, 34], and was well suited to designing our voltage and current waveforms. The outcome of our engineering optimization was a KHFAC current waveform that produced nerve conduction block without generating an onset response. Our results revealed that moving VGSCs into the CSI without first opening completely suppressed the onset firing and suggest that using this novel approach for designing blocking waveforms may have relevance for the conduction block of peripheral nerve hyperactivity.

Results

Design of transmembrane voltage trajectory to generate closed-state inactivation of VGSCs

Our first step was to design a voltage trajectory that generated CSI of VGSCs. We used a 21-node 10 μm diameter MRG myelinated axon model of mammalian peripheral nerve fiber (Fig 1A). The fast Na+ (INaf) and persistent Na+ (INap) currents in each node were originally modeled using Hodgkin-Huxley (HH)-type channels [31], which represented a VGSC as an assembly of several independent gating variables and were unable to represent the dependence of Na+ inactivation on activation. Since KHFAC waveforms generate a locally-acting conduction block in peripheral nerves, the original INaf and INap were replaced respectively by Markov-type Nav 1.1 and Nav 1.6 channels [30] in the middle nine nodes (6 to 14) to create a hybrid model that enabled transitions of VGSCs between states. The CSI path of Nav 1.1 and Nav 1.6 channels was C1 → I1 → I2, and the OSI path was C1 → C2 → O1 → I1 → I2 (Fig 1B). The original conductance density of INaf was too small for INav11 to propagate faithfully an AP from node 0 to node 20, and we increased the conductance density of INav11 to enable high fidelity propagation along the hybrid axon model at rates of up to 400 spikes per second (S1 Fig).

thumbnail
Fig 1. Hybrid model of a peripheral myelinated axon incorporating Markov-type VGSCs.

(A) Schematic of the model, which included 21 nodes of Ranvier and 20 internodal segments. In the original MRG model, each node included fast Na+ (INaf), persistent Na+ (INap), slow K+ (IKS), and linear leakage (IL) currents, and membrane capacitance (Cn), and the ionic currents were described using HH-type kinetics. We substituted INaf and INap in node 6 to node 14 with Markov-type Nav 1.1 (INav11) and Nav 1.6 (INav16) channels, respectively. The internodal segments used a double-cable concentric structure to represent the myelin attachment segment (MYSA), paranode main segment (FLUT), and internode segment (STIN). The internodal segments were modeled with linear conductances (Gm, Gi) and membrane capacitances (Cm, Ci). Adjacent segments were connected by axoplasmic (Ga) and periaxonal conductances (Gp). Resting membrane potential Vrest was –80 mV, and the fiber diameter was 10 μm. The axon-specific electrical and geometrical parameters followed McIntyre et al [31] and are provided in S1 and S2 Tables. (B) Schematic of the states and transitions of Markov-type VGSCs [30]. C1 and C2 were closed states, O1 and O2 were open states, I1 and I2 were inactivation states, and Aw1w2 was the transition rate from state w1 to state w2. Pink lines represented a CSI path, and blue lines represented an OSI path. The parameters of the transition rate for Nav 1.1 and Nav 1.6 channels are provided in S3 Table.

https://doi.org/10.1371/journal.pcbi.1007766.g001

We first examined the responses in Nav 1.1 and Nav 1.6 channels when a standard voltage step (from –80 mV to –10 mV) was applied to the central node as a voltage clamp (Fig 2A). At the resting potential of –80 mV, the two VGSCs were in closed states C1. After the onset of the voltage step to –10 mV, both Nav 1.1 and Nav 1.6 channels first entered open states O1 and then inactivated (Fig 2B), i.e., OSI occurred. The magnitude of INav11 and INav16 reached –81.85 mA/cm2 and –0.10 mA/cm2, respectively. The conductance density of Nav 1.6 channels (i.e., 0.01 mS/cm2) was much smaller than that of Nav 1.1 channels (i.e., 11.9 mS/cm2), and INav11 was much larger than INav16. An AP was initiated and propagated in the hybrid model axon at the onset of voltage step (Fig 2C).

thumbnail
Fig 2. Trajectory of transmembrane voltage to cause Na+ CSI.

(A) Top panel: Voltage step profile Vstep. The transmembrane voltage was increased from the resting potential of -80 mV to -10 mV at 25 ms. Center panel: Nav 1.1 current during Vstep. Bottom: Nav 1.6 current during Vstep. (B) Top panel: fractions of Nav 1.1 channels in states C1 (C11), I1 (I11), I2 (I12), and O1 (O11) during Vstep. Bottom panel: fractions of Nav 1.6 channels in states C1 (C61), I1 (I61), I2 (I62), and O1 (O61) during Vstep. (C) Voltage response recorded in each node when applying Vstep at node 10. (D) Top panel: PSO-generated optimized voltage profile VPSO. Center panel: Nav 1.1 current during VPSO. Bottom panel: Nav 1.6 current during VPSO. (E) Top panel: C11, I11, I12, and O11 of Nav 1.1 channels during VPSO. Bottom panel: C61, I61, I62, and O61 of Nav 1.6 channels during VPSO. (F) Voltage response recorded in each node when applying VPSO as a series of voltage clamps at node 10.

https://doi.org/10.1371/journal.pcbi.1007766.g002

We subsequently applied PSO to design an optimized trajectory of transmembrane voltage to drive transitions along the CSI path. PSO is a search algorithm inspired by the social behavior of swarms. It finds the optimal or near-optimal solutions to numerical and qualitative problems by adjusting the trajectories of individual particles in the defined search space of an objective function [32], and is particularly applicable to solving our waveform optimization. We defined a population of 50 particles as the membrane voltages that were used as the voltage profile to control a series of voltage clamps at node 10. The objective of the PSO was to maximize the fractions of Nav 1.1 and Nav 1.6 channels in the inactivated states while minimizing the fractions in the open states. The PSO-generated voltage profile VPSO increased nonlinearly from the resting potential (–80 mV) to a suprathreshold level (–10 mV) over a duration of 46 ms (Fig 2D, top) and drove transitions of Nav 1.1 and Nav 1.6 channels to CSI.

The transitions from the closed states to the inactivated states caused by voltage profile VPSO were similar for Nav 1.1 and Nav 1.6 channels. At the resting potential, VGSCs were in C1 state. After depolarization, the fraction of either channel in C1 was gradually reduced (Fig 2E), the fractions in I2 increased, and there was little increase in the fraction in O1. Thus, more Nav 1.1 and Nav 1.6 channels were driven to inactivation directly from the resting state, which reduced the availability of each channel to conduct and participate in AP generation, and when the membrane voltage reached suprathreshold levels, the VGSCs were unable to enter the open states. The Na+ conductance was determined by the fraction of VGSCs in the open states, and since few Nav 1.1 channels entered the O1 state during VPSO, INav11 reached only –0.62 mA/cm2 (Fig 2D, bottom). In contrast to the voltage clamp step, the PSO-generated VPSO successfully moved both Nav 1.1 and Nav 1.6 channels to CSI, and VPSO delivered as a series of voltage clamps to the central node did not evoke an AP (Fig 2F).

PSO-generated DC waveform drives VPSO and generates onset-free conduction block

Having generated a transmembrane voltage trajectory that caused CSI of Nav 1.1 and Nav 1.6 channels, our next step was to design a stimulus waveform to drive this voltage profile in the axon. We again applied PSO, with the particles defined as the electrode currents at each time step, to design an extracellular stimulus to drive VPSO in node 10 using an extracellular monopolar electrode placed 1 mm above the central node (Fig 3A). The goal of the PSO was to minimize the difference between the membrane voltage recorded in node 10 and the optimized voltage profile VPSO. The resulting current IPSO increased nonlinearly from 0 mA to a plateau of –0.69 mA (Fig 3B) and drove the transmembrane potential to track VPSO. No AP was generated in the axon by IPSO (Fig 3C) and it successfully blocked nerve conduction. When the time delay (Δτ) between the onset of an intracellular test pulse and the onset of IPSO was ≤ 35 ms, the single spike propagated through the axon (Fig 3D). With Δτ ≥ 36 ms, the transmission of single spike was blocked (Fig 3D), as were spike trains at 50 Hz, 100 Hz or 200 Hz (Fig 3E–3G). Thus, the PSO-generated waveform IPSO blocked nerve conduction without activating an onset response in the hybrid model axon.

thumbnail
Fig 3. PSO-generated DC waveform drives VPSO and generates onset-free conduction block.

(A) Simulation setup. A monopolar block electrode was placed 1 mm over node 10, and suprathreshold intracellular test pulses were delivered at node 0 to generate propagating APs. (B) Top panel: Voltage response recorded at node 10 (red solid line) and the PSO-generated voltage profile VPSO (black dotted line). Bottom panel: PSO-generated DC waveform IPSO for driving VPSO in node 10. (C) Voltage response recorded in each node after applying IPSO. No test pulse was delivered at node 0. (D) Propagation or block of a single test AP along the axon initiated at a time delay Δτ after the onset of IPSO. Top panel: Δτ = 35 ms, bottom panel: Δτ = 37 ms. (E) Block of train of 50 Hz test APs by IPSO. (F) Block of train of 100 Hz test APs by IPSO. (G) Block of train of 200 Hz test APs by IPSO.

https://doi.org/10.1371/journal.pcbi.1007766.g003

KHFAC waveform for onset-free conduction block

The DC IPSO waveform was not charge-balanced, and thus can cause irreversible faradaic reactions at the electrode—electrolyte interface [35], and these may result in electrode and or nerve damage when applied for long periods of time [12, 16, 36]. We therefore used a scaled version of the DC IPSO waveform as the envelope of symmetric rectangular charge-balanced biphasic pulses [12, 37], applied at 10 kHz to block nerve conduction (Fig 4A, bottom). However, there was a time delay after the onset of resulting KHFAC waveform IBI and complete conduction block (Fig 5A), and this delay depended on the amplitude scale factor (Fig 5B). The minimum scale factor for producing conduction block was 1.0, but we used a scale factor of 1.5, which enabled the KHFAC waveform to generate conduction block after a minimum time delay.

thumbnail
Fig 4. PSO-based charge-balance biphasic KHFAC waveform for onset-free conduction block.

(A) Transmembrane voltage (top) and Nav 1.1 current INav11 (center) recorded in node 10 in response to KHFAC current IBI (bottom). The frequency of IBI was 10 kHz, and each pulse was rectangular, symmetric, and charge-balanced. The envelope of IBI was determined by multiplying IPSO (red line) by a scale factor of 1.5. (B) Fractions of Nav 1.1 channels in the states C1 (C11), I2 (I12), and O1 (O11). (C) Voltage response recorded in each node after applying IBI with a monopolar electrode placed 1 mm above the central node. No test pulse was injected at node 0. The red lines in (A) and (B) were the responses associated with DC waveform IPSO.

https://doi.org/10.1371/journal.pcbi.1007766.g004

thumbnail
Fig 5. Onset-free conduction block by PSO-based KHFAC waveform.

(A) Propagation or block of a single test AP along the axon initiated at a time delay Δτ after the onset of IBI. Top panel: Δτ = 47 ms, bottom panel: Δτ = 49 ms. A scale factor of 1.5 was used to design the envelope of IBI. (B) Minimum time delay Δτmin for producing conduction block as a function of scale factor. (C) Block of train of 50 Hz test APs by IBI. (D) Block of train of 100 Hz test APs by IBI. (E) Block of train of 200 Hz test APs by IBI. PSO-based KHFAC waveform IBI was applied using a monopolar block electrode placed over the central node at an electrode-fiber distance of 1 mm.

https://doi.org/10.1371/journal.pcbi.1007766.g005

The KHFAC waveform IBI, delivered using a monopolar electrode placed 1 mm over node 10, generated passive oscillations in transmembrane voltage at node 10, and the upper envelope of the transmembrane voltage response followed the optimized voltage profile VPSO (Fig 4A, top). The intrinsic low-pass filtering and rectification of AC electrical signals by the cell membrane prevented neural firing from following the KHFAC stimulus, while the envelop of the PSO-based waveform IBI drove the Nav 1.1 channels to CSI (Fig 4B). The fraction C11 of Nav 1.1 channels in C1 state was higher than that in response to the DC waveform IPSO (Fig 4B), which enabled more channels to enter O1 state instead of I2 and resulted in more Na+ influx through the Nav 1.1 channels, causing a depolarized peak in the transmembrane voltage of node 10. Importantly, this depolarized peak did not generate a propagating AP (Fig 4C), and the KHFAC waveform enabled onset-free conduction block. When the time delay (Δτ) between the onset of an intracellular test pulse and the onset of IBI was ≤ 47 ms, the single spike propagated through the axon (Fig 5A). With Δτ ≥ 48 ms, the transmission of the single spike was blocked (Fig 5A), and the minimum delay Δτmin for producing conduction block was dependent on the scale factor (Fig 5B), reaching a minimum when the scale factor was 1.5 ~ 1.9. Further, IBI completely blocked conduction of spike trains at 50 Hz, 100 Hz, or 200 Hz (Fig 5C–5E). Thus, the PSO-based charged-balanced KHFAC waveform blocked nerve conduction without generating an onset response.

Discussion

We implemented a novel design process to develop a KHFAC stimulus waveform to enable nerve conduction block without the onset response that occurs with conventional KHFAC waveforms. We implemented Markov-type kinetic models to simulate CSI of Nav 1.1 and Nav 1.6 channels in a novel hybrid model of a mammalian myelinated nerve fiber. We used engineering optimization with a PSO algorithm first to design an optimized transmembrane voltage trajectory to drive the two VGSCs to CSI without first causing activation, and subsequently, to generate an extracellular DC stimulation waveform that drove the voltage profile in the central node of the model to follow the optimized trajectory. Finally, we designed a charge-balanced 10 kHz biphasic waveform based on the PSO-generated DC waveform, and the PSO-based KHFAC waveform produced robust conduction block without an onset response.

In in vivo experiments, a KHFAC waveform is commonly delivered using an extracellular electrode and can affect a number of axons with different electrode-fiber distances. We calculated the PSO-based biphasic waveform IBI in a 21-node model at distance of 1 mm, and this waveform did not produce onset firing. The block threshold increases with increasing electrode-fiber distance [17, 38, 39], and we examined whether our PSO-based KHFAC waveform optimized at a specific location produced onset-free conduction block in the axon at other locations. Applying IBI with same plateau amplitude generated an onset response in the model at an electrode-fiber distance of 0.3 mm (S2 Fig). These closer nerve fibers had a block threshold lower than the plateau amplitude of IBI, and as the KHFAC amplitude approached the plateau, the node underneath the blocking electrode was strongly depolarized and APs were initiated in adjacent nodes. While the PSO-based biphasic waveform eliminated the onset response in nerve fibers at certain locations, it may still activate onset firing at the other locations.

In in vivo experiments, a KHFAC waveform delivered using an extracellular electrode can also affect axons with different fiber diameters. The block threshold decreases with increasing fiber diameter [17, 38, 39], and we examined whether the PSO-based KHFAC waveform produced onset-free conduction block in smaller diameter (2 μm and 8.7 μm) model axons. The internodal lengths of 2 μm or 8.7 μm axons are much shorter than those of the 10 μm fiber (S1 Table), which significantly reduced the total length of the 21-node model axon. Since the electrode-fiber distance was a fixed value (1 mm), we increased the number of nodes to 251 to examine conduction block using the PSO-based KHFAC waveform and to avoid potential end effects, as the ends of the 21-node model axons with shorter internodal lengths were much closer to the electrode. When scaling the plateau amplitude of IBI to the block threshold of a 251-node model at a diameter of 2 μm, two APs were generated as the KHFAC amplitude approached the block threshold (S3 Fig). Applying IBI with same plateau amplitude to a 251-node model at a diameter of 8.7 μm also produced an onset response. While the PSO-based KHFAC waveform eliminated onset response in the fibers for which it was designed, it may generate onset responses in other fiber diameters.

The original MRG models used the HH-type formalism to describe fast Na+ and persistent Na+ channels in each node [31], and we designed the waveform using the hybrid model incorporating Markov-type channels. Applying our PSO-based biphasic waveform to a 21-node 10 μm HH-type MRG axon, i.e., without the Markov-type VGSCs in the central nodes, generated eight APs as KHFAC amplitude approached the block threshold (S4 Fig). The initial ramp of biphasic stimulus drove the central node under the block electrode to suprathreshold depolarization, which enabled the initiation of onset firing. The VGSC gating variables m and h operate independently of each other, and no transition occurred between them as KHFAC amplitude increased. In contrast, the Markov-type Nav 1.1 channels moved from the closed state to inactivation without fully opening in response to the PSO-based biphasic waveform, and thus no onset response was generated. To reduce overall computational demands, we implemented Nav 1.1 and Nav 1.6 channels only in the central nine nodes of a 21-node MRG model, but the specific number of Markov-type nodes did not alter the block threshold or the lack of an onset response in the hybrid model (S5 Fig).

An earlier study of somatic voltage clamp [40] reported that applying a voltage pre-pulse inactivated axonal channels after evoking an axonal spike, and the voltage steps following the pulse no longer triggered axonal spikes. Our study showed that directly applying voltage clamp of the PSO-generated voltage profile VPSO at the central node generated CSI of VGSCs, and no APs were escaped from the node. Further, we used a monopolar point source electrode to design PSO-based DC and KHFAC waveforms in the hybrid MRG model, while bipolar and tripolar cuff electrodes are commonly used in vivo [6, 14, 15, 17]. Linearly ramping KHFAC amplitude from zero to block threshold in rat sciatic nerve did not eliminate the onset response with a tripolar electrode [14], while linearly ramping KHFAC waveforms from non-zero amplitudes reduced the onset response with a bipolar electrode [15]. We re-calculated the PSO-generated DC and biphasic waveforms using representations of bipolar and tripolar cuff electrodes, and the resulting KHFAC waveform was similar for each of three electrodes. Using a bipolar cuff, there was onset firing activated in the 21-node hybrid model as the KHFAC amplitude approached the block threshold (S6 Fig), but using the tripolar cuff, no onset response was generated by the initial ramp of the KHFAC waveform, as observed with the original monopolar electrode.

Other waveforms were proposed to eliminate the onset response. Tai et al. [41] predicted that ramping KHFAC amplitude from 0 mA to a level sufficient to produce block eliminated the onset response in a HH model of unmyelinated giant squid axons. Miles et al. [14] showed that the KHFAC ramp that started from zero amplitude did not eliminate the onset response in either the HH-type MRG model or in experiments on rat sciatic nerve. Vrabec et al. [15] recently reported that ramping KHFAC amplitudes up to 125% of block threshold from non-zero amplitudes were successful in reducing the onset response in rat sciatic nerve. In vivo recordings demonstrated that the onset response was minimized with high amplitudes and high frequencies [9, 17, 42], and Gerges et al [13] showed that a transition from a high amplitude and high frequency to a low amplitude and low frequency during KHFAC nerve block minimized or eliminated the onset response with longer transition times. A DC block was also introduced to eliminate the onset response. Ackermann et al. [43] reported that a brief DC block in combination with KHFAC prevented the onset firing, but such DC was found to produce nerve damage and prolonged conduction failure. Subsequently, Franke et al. [12] suggested that combining KHFAC and charge-balanced DC block significantly reduced or completely prevented the onset response without a change in nerve conduction by using a block electrode of high capacitance materials. The PSO-generated DC waveform also eliminated the onset response in the hybrid MRG axon, but our DC block was not charge-balanced.

We presented an engineering optimization method to design a charge-balanced biphasic current waveform that eliminated the onset response during KHFAC nerve block. Our PSO-based KHFAC waveform completely prevented the onset firing in axons at specific locations and with specific diameters, which suggested that it is possible to achieve KHFAC nerve block without an onset response by moving VGSCs to CSI.

Methods

Axon model

We used the MRG double-cable myelinated axon model to simulate KHFAC conduction block, and this model was originally developed to reproduce a wide range of excitation properties of mammalian peripheral nerve fibers [31]. The model was composed of equivalent electrical circuits for both nodal and internodal segments (Fig 1A). The nodal circuit included fast Na+ (INaf), persistent Na+ (INap), slow K+ (IKS), and linear leakage (IL) currents, and the membrane capacitance Cn, and the ionic currents were originally modeled by using HH-type kinetics. INaf was responsible for AP initiation in each node (maximum conductance: 3 S/cm2), and INap was responsible for generating a non-inactivating current at the afterpotentials (maximum conductance: 0.01 S/cm2). The resting potential in each node was –80 mV. The axon-specific parameters were implemented following McIntyre et al [31, 44] (S1 and S2 Tables) and were previously validated to be consistent with experimental observations of KHFAC conduction block [38, 39, 45]. The axon models were solved in the NEURON simulation environment [46] (version 7.7) using a time step of 0.005 ms.

Markov-type models of VGSCs

There are two types of models commonly used to reproduce the electrophysiological behaviors of VGSCs, i.e., HH- and Markov-type models. The HH-type models represent an ionic channel as an assembly of several independent gating particles, which are not sufficient to capture specific kinetics of VGSCs [47, 48]. Unlike the HH formalism, Markov-type models represent an ion channel as a series of conformational states of the channel protein and a collection of transitions between them [18, 22], which allowed us to reproduce the dependence of Na+ inactivation on activation. We used Markov-type kinetic models of VGSCs to design the transmembrane voltage trajectory for causing Na+ CSI.

The Markov-type models used in our simulations were developed by Balbi et al [30] and were validated against electrophysiological data from all VGSC isoforms. The model framework was detailed, unifying, and computationally efficient, and can account for different features of each human VGSC isoform with a minimal set of states and transitions. The model included six states, and the state diagram and transitions of the Markov-type channels are shown in Fig 1B. C1 and C2 were two closed states, O1 and O2 were two open states, and I1 and I2 were two inactivation states. The transition from O1 to I1 was irreversible, and the transitions between the other consecutive states were reversible. The dynamics of fraction of VGSCs in each state were governed by following Markovian equations [30] (1) (2) (3) (4) (5) (6)

Here Aw1w2 was the transition rate from state w1 to state w2, which was computed by (7)

The Na+ conductance was determined by the fraction of open channels, and the Na+ current was computed by (8) where was the maximum conductance, V was local membrane voltage, and ENa = 50 mV was the Na+ reversal potential.

Implementation of Markov-type VGSCs

KHFAC waveforms are applied in practice using comparatively small electrodes and generate a locally-acting conduction block in peripheral nerves. In our simulations, there were 21 nodes in the model axon and the simulated block electrode was placed above the central node (i.e., node 10). To reduce overall computational demands, we implemented Markov-type VGSCs in the central nine nodes, which were closest to the block electrode. The other nodes far away from the site of block were modeled with HH-type VGSCs. Note that the specific number of Markov-type nodes did not alter the block threshold or lack of onset response in the hybrid model (S5 Fig).

Experimental data indicated that the Nav 1.1 and Nav 1.6 channels are preferentially expressed in the nodes of mammalian peripheral nerve fibers [4951]. We substituted the HH-type models of INaf and INap in node 6 to node 14 with the Markov-type models of Nav 1.1 and Nav 1.6 channels, respectively. The parameters of each transition rate for two VGSCs are provided in S3 Table. We modified only the maximum conductances of INav11 and INav16 and did not change any parameter of their transition rates between consecutive states. The maximum conductance of Nav 1.6 channels was set to = 0.01 S/cm2, which was identical to that of INap. For Nav 1.1 channels, if we set to the minimal value (4.7 S/cm2) for faithfully propagating a single spike from node 0 to node 20, there was substantial attenuation in the AP as it propagated from HH-type nodes to Markov-type nodes (S1 Fig). In this case, the was too small to propagate APs at high frequencies. For example, when firing rate was 200 spikes per second, APs were lost as they propagated along the axon. We set to the minimal value (11.9 S/cm2) to propagate faithfully spike trains at rates of up to 400 spikes per second. At = 11.9 S/cm2, there was no attenuation in the APs as they propagated along the axon (S1 Fig). Although was larger than , INav11 and INaf had similar magnitude underlying an AP.

Particle swarm optimization

We applied PSO algorithms to design the transmembrane voltage trajectory for causing Na+ CSI and to generate current waveform to drive this voltage trajectroy, which were implemented in Matlab (version R2016a). PSO is a meta-heuristic optimization method inspired by the information circulation and social behavior of swarms [32]. It uses a population-based search stratergy to solve nonlinear optimization problems, based on cooperation and competition among the particles of a swarm. In PSO, each particle is a point in the search space, which represents a candidate solution of the optimization problem. The movement of a particle in the search space is driven by its own fitness as well as the highest fitness in the swarm. The algorithm iteratively updates the velocity of each particle towards the position of the highest fitness. Through communication between particles over multiple iterations, the particles explore the problem space to identify the optimal solution of the problem.

Mathematically, the position of particle i (1 ≤ i ≤ M) at the kth iteration was described by a vector , where M was the number of particles and D was the dimension of the search space. Each dimension of the particle was bounded between xmin and xmax. The personal best position prevoiusly found by particle i was , and the best position expericed by all particles in the whole swarm was indicated by . The velocity of particle i in the kth iteration was , which determined the movement direction of the particle through the search space. At the kth iteration, the d-dimension (1 ≤ d ≤ D) of particle i was updated based on following equations (9) (10) where c = 1.5 was a learning coefficient, and r1 and r2 were uniform random numbers within [0, 1] used to help promote exploration.

Implementation of PSO

We applied PSO to generate an optimized voltage trajectory VPSO for causing the CSI of Nav 1.1 and Nav 1.6 channels in the 21-node model (Fig 6A). We defined a population of 50 particles as the membrane voltages used to control a time series of voltage clamps at node 10. We recorded the Na+ currents through Nav 1.1 and Nav 1.6 channels as well as their fractions in states C1, I1, I2, and O1 during voltage clamp. We defined a multi-objective optimization problem: find a voltage trajectory to maximize the fractions of Nav 1.1 and Nav 1.6 channels in inactivated states while minimizing the fractions of the two VGSCs in open states. Mathematically, the following cost function f1 was defined to evaluate each solution, (11) where O11, O12, and I12 were the fractions of Nav 1.1 channels in states O1, O2 and I2, O61, O62, and I62 were the fractions of Nav 1.6 channels in states O1, O2 and I2, j corresponded to the time point, and B = 20000 was a weight used to reduce the fractions in open states. The lower bound of each particle was xmin = −80 mV (i.e., resting potential), and the upper bound was xmax = −10 mV. The time constant of inactivation state I1 was smaller than I2, and the optimization results were sensitive to the changes in I1. Thus, we excluded the fraction of state I1 from the fitness function f1. The initial positions of 50 particles were defined to be randomly distributed around a straight line, which linearly increased from xmin to xmax.

thumbnail
Fig 6.

Flowchart of the PSO algorithm for (A) designing a transmembrane voltage trajectory for causing Na+ CSI and (B) generating DC waveform for driving that trajectory.

https://doi.org/10.1371/journal.pcbi.1007766.g006

Subsequently, we used PSO to generate a current waveform for driving the optimized voltage trajectory in the central node of the 21-node model (Fig 6B). We placed a point source electrode over node 10 at an electrode-fiber distance of 1 mm. We defined 50 particles as a time series of electrode currents. The extracellular medium was assumed to be infinite homogenous and isotropic with a conductivity of σ = 0.303 S/m [52, 53]. For a segment at location (x, y, z), the extracellular potential Vext was calculated as [53, 54] (12) where Iext was the amplitude of electrode current, and (x0, y0, z0) was the location of the point source electrode. We recorded the resulting transmembrane voltage in node 10 at each iteration. We set VPSO as the desired voltage trajectory, and defined an one-objective optimization problem: find electrode currents to minimize the difference between transmembrane voltage Vnode10 recorded in node 10 and desired voltage trajectory VPSO. Mathematically, the following cost function f2 was defined to evaluate the quality of this problem solution (13) where j corresponded to the time point. By using the point source electrode, cathodic stimulus depolarized node10 and anodic stimulus hyperpolarized node10. Since VPSO started from the resting potential, the upper bound of each particle was xmax = 0 mA. The lower bound was defined to xmin = −0.7 mA, which was sufficient to drive transmembrane voltage in node 10 to a plateau level higher than –10 mV. The initial positions of particles were defined to be randomly distributed around a straight line from the xmin to xmax.

When solving the optimal solutions for f1 and f2, the time step of PSO was 2 ms and the search space dimension was D = 50. The termination criterion was that the algorithm generation reached a limit of 1000, which was sufficient for convergence. When designing the optimized voltage trajectory VPSO for causing Na+ CSI, the minimum f1 converged to within 1% of the final iteration by 330 iterations (S7 Fig). When generating the DC waveform IPSO for driving the optimized trajectory VPSO, the minimum fitness f2 after 499 iterations converged to within 1% of the final iteration (S8 Fig).

Simulation protocols for KHFAC conduction block

After generating PSO-based DC and KHFAC waveforms for conduction block, we first examined whether they activated onset responses. We used an extracellular point source electrode to deliver blocking currents to the central node (10) of the model, and the electrode-fiber distance was 1 mm. Subsequently, we determined whether the PSO-generated waveform blocked nerve conduction. We applied an intracellular monopolar electrode to deliver suprathreshold test pulses (width: 0.1 ms and amplitude: 2.5Ith) at one end of the axon (node 0), which was used to stimulate APs in node 0 after the onset of blocking currents and test whether these APs propagated through the site of the block electrode. Ith = 0.43 mA was the threshold of intracellular test pulse for activation of an AP in node 0.

A 251-node model was used to examine the conduction block at the fiber diameters (2 μm and 8.7 μm) less than 10 μm using PSO-based KHFAC waveform. Each fiber diameter had specific lengths for all segments (S1 Table), which were identified based on experimental data. The Markov-type Nav 1.1 and Nav 1.6 channels were implemented in node 121 to node 129, and the block electrode was placed 1 mm above node 125.

Supporting information

S1 Fig. Setting maximum conductance of Nav 1.1 channels.

Propagation of a test AP initiated in node 0 along the 21-node and 10 μm model at (A) = 4.7 S/cm2 and (B) 11.9 S/cm2. = 4.7 S/cm2 was the minimum conductance for faithfully propagating a single AP along the axon, and = 11.9 S/cm2 was the minimum conductance for faithfully propagating spike trains at rates of up to 400 Hz along the axon. INaf was recorded in node 3, and INav11 was recorded in node 10. Left panels: a single test pulse (width: 0.1 ms and amplitude: 2.5Ith) was delivered at the node 0. Right panels: 200 Hz test pulses were delivered at the node 0. Maximum conductance of HH-type INaf was = 3.0 S/cm2, maximum conductance of HH-type INap was = 0.01 S/cm2, and maximum conductance of Markov-type INav16 was = 0.01 S/cm2.

https://doi.org/10.1371/journal.pcbi.1007766.s001

(TIF)

S2 Fig. Conduction block at electrode-fiber distance of 0.3 mm by PSO-based KHFAC waveform.

(A) Simulation setup. A monopolar block electrode was placed 0.3 mm over the central node of 21-node 10 μm diameter model nerve fiber. An intracellular test pulse (width: 0.1 ms and amplitude: 2.5Ith) was delivered at node 0 to generate a propagating AP at t = 120 ms. (B) Transmembrane voltages (top) recorded in node 0 and node 20 in response to KHFAC waveform IBI (bottom). A scale factor of 1.5 was used to design the envelope of IBI, and the green dotted line was the block threshold at electrode-fiber distance of 0.3 mm. (C) Voltage responses recorded in node 10 to node 14.

https://doi.org/10.1371/journal.pcbi.1007766.s002

(TIF)

S3 Fig. Effects of fiber diameter on conduction block by PSO-based KHFAC waveform.

(A) Simulation setup. A monopolar block electrode was placed 1 mm over node 125 and delivered KHFAC waveform IBI to a 251-node model. An intracellular test pulse (width: 0.1 ms and amplitude: 2.5Ith) was delivered at node 0 to generate a propagating AP at t = 120 ms. Markov-type Nav 1.1 and Nav 1.6 channels were implemented in node 121 to node 129. (B) Block of a test AP along the axon by IBI with a fiber diameter of 2.0 μm. (C) Block of a test AP along the axon by IBI with a fiber diameter of 8.7 μm. In (B) and (C), the green dotted lines were the block threshold at each fiber diameter.

https://doi.org/10.1371/journal.pcbi.1007766.s003

(TIF)

S4 Fig. Conduction block in HH-type model by PSO-based KHFAC waveform.

(A) Block of a test AP along the 21-node 10 μm diameter model nerve fiber by PSO-based waveform IBI. No Markov-type VGSCs were implemented in the central nodes. The block electrode was placed 1 mm above node 10, and a single test pulse (width: 0.1 ms and amplitude: 2.5Ith) was delivered at node 0 to generate a propagating AP at t = 70 ms. (B) Top panel: transmembrane voltage recorded in node 10 in response to KHFAC waveform IBI. Center panel: activation gating variable m and inactivation gating variable h of fast Na+ current. Bottom panel: PSO-based KHFAC waveform IBI with a plateau amplitude scaled to the block threshold.

https://doi.org/10.1371/journal.pcbi.1007766.s004

(TIF)

S5 Fig. Effects of number of Markov-type nodes on conduction block by PSO-based KHFAC waveform.

(A) Simulation setup. We implemented Nav 1.1 and Nav 1.6 channels in node (10−n) to node (10+n), where 0 ≤ n ≤ 10. A monopolar block electrode was placed 1 mm over node 10, and an intracellular test pulse (width: 0.1 ms and amplitude: 2.5Ith) was delivered at node 0 to generate a propagating AP at t = 120 ms. (B) Block threshold as a function of the number of Markov-type nodes. At 3 ≤ n ≤ 10, PSO-based KHFAC waveform IBI produced onset-free conduction block. At 0 ≤ n ≤ 2, no conduction block occurred and onset firing was activated by IBI. (C) Propagation of a test AP along the axon by IBI (n = 1). A scale factor of 1.0 was used to design the envelope of IBI, and the number of Markov-type nodes was three.

https://doi.org/10.1371/journal.pcbi.1007766.s005

(TIF)

S6 Fig. PSO-based KHFAC waveform designed using bipolar and tripolar cuff electrodes.

(A) We interpolated the extracellular potential at each segment of 21-node hybrid model from a finite element model (FEM) of the rat tibial nerve, which had a bipolar nerve cuff and a model nerve fiber that was 0.398 mm away from the electrode. Since the bipolar cuff resulted in the largest depolarization in node 11, we applied PSO algorithm to generate a DC waveform IPSO to drive the optimized voltage trajectory VPSO in node 11. The envelop of KHFAC waveform IBI was determined by multiplying IPSO by a scale factor of 0.73, and the plateau amplitude of IBI was the block threshold. The resulting KHFAC waveform blocked nerve conduction with activation of an onset response. (B) We used a FEM of a tripolar cuff on the rat tibial nerve to interpolate the extracellular potentials of 21-node hybrid model, and the electrode-fiber distance was 0.398 mm. The DC waveform IPSO was designed to drive VPSO in node 10. A scale factor of 0.93 was applied to scale the plateau amplitude of IBI to the block threshold. The resulting KHFAC waveform produced onset-free conduction block. A single test pulse (width: 0.1 ms and amplitude: 2.5Ith) was injected in node 0 at t = 70 ms. Fiber diameter was 10 μm.

https://doi.org/10.1371/journal.pcbi.1007766.s006

(TIF)

S7 Fig. Progression of PSO for generating transmembrane voltage trajectory to cause Na+ CSI.

(A) Changes of transmembrane voltage trajectory across iterations. The sequence of plots show the voltage profile with the minimum cost f1 at each indicated iteration. (B) Minimum cost f1 of 50 particles at each iteration. Red dotted line is 1.01 times the minimum cost at the final generation. (C) Voltage response recorded in each node when applying the voltage profile at each indicated iteration as a series of voltage clamps at node 10. Fiber diameter was 10 μm.

https://doi.org/10.1371/journal.pcbi.1007766.s007

(TIF)

S8 Fig. Minimum cost function f2 of 50 particles at each iteration.

Red dotted line was 1.01 times the minimum cost at the final generation.

https://doi.org/10.1371/journal.pcbi.1007766.s008

(TIF)

S1 Table. Geometric parameters of MRG models

https://doi.org/10.1371/journal.pcbi.1007766.s009

(DOCX)

S2 Table. Electrical parameters of MRG models

https://doi.org/10.1371/journal.pcbi.1007766.s010

(DOCX)

S3 Table. Parameters of the transition rates for Nav 1.1 and Nav 1.6 channels

https://doi.org/10.1371/journal.pcbi.1007766.s011

(DOCX)

Acknowledgments

We would like to thank Edgar Peña, Ph.D., for comments on designing PSO-based KHFAC waveforms and Nathan Titus for implementing the PSO algorithm on the Duke Shared Compute Cluster.

References

  1. 1. Amatya B, Khan F, La Mantia L, Demetrios M, Wade DT. Non pharmacological interventions for spasticity in multiple sclerosis. Cochrane Database Syst Rev. 2013 Feb 28; (2): CD009974. pmid:23450612
  2. 2. Tiede J, Brown L, Gekht G, Vallejo R, Yearwood T, Morgan D. Novel spinal cord stimulation parameters in patients with predominant back pain. Neuromodulation. 2013 Jul–Aug; 16(4): 370–5. pmid:23433237
  3. 3. Boger A, Bhadra N, Gustafson KJ. Bladder voiding by combined high frequency electrical pudendal nerve block and sacral root stimulation. Neurourol Urodyn. 2008; 27(5): 435–9. pmid:18041769
  4. 4. Avendano-Coy J, Serrano-Munoz D, Taylor J, Goicoechea-Garcia C, Gomez-Soriano J. Peripheral nerve conduction block by high-frequency alternating currents: a systematic review. IEEE Trans Neural Syst Rehabil Eng. 2018 Jun; 26(6): 1131–1140. pmid:29877837
  5. 5. Tanner JA. Reversible blocking of nerve conduction by alternating-current excitation. Nature. 1962 Aug 18; 195: 712–3. pmid:13919574
  6. 6. Kilgore KL, Bhadra N. Reversible nerve conduction block using kilohertz frequency alternating currents. Neuromodulation. 2014 Apr; 17(3): 242–54. pmid:23924075
  7. 7. Bhadra N, Vrabec TL, Bhadra N, Kilgore KL. Reversible conduction block in peripheral nerve using electrical waveforms. Bioelectron Med (Lond). 2018 Jan; 1(1): 39–54. pmid:29480897
  8. 8. Pelot NA, Grill WM. In vivo quantification of excitation and kilohertz frequency block of the rat vagus nerve. Journal of Neural Engineering 2020 Mar 9; 17(2):026005. pmid:31945746
  9. 9. Patel YA, Butera RJ. Challenges associated with nerve conduction block using Kilohertz electrical stimulation. J Neural Eng. 2018 Jun; 15(3): 031002. pmid:29415877
  10. 10. Bhadra N, Kilgore KL. High-frequency electrical conduction block of mammalian peripheral motor nerve. Muscle Nerve. 2005 Dec; 32(6): 782–90. pmid:16124008
  11. 11. Ackermann DM Jr, Foldes EL, Bhadra N, Kilgore KL. Conduction block of peripheral nerve using high frequency alternating currents delivered through an intrafascicular electrode. Muscle Nerve. 2010 Jan; 41(1): 117–9. pmid:19813186
  12. 12. Franke M, Vrabec T, Wainright J, Bhadra N, Bhadra N, Kilgore K. Combined KHFAC+DC nerve block without onset or reduced nerve conductivity after block. J Neural Eng. 2014 Oct; 11(5): 056012. pmid:25115572
  13. 13. Gerges M, Foldes EL, Ackermann DM, Bhadra N, Bhadra N, Kilgore KL. Frequency and amplitude transitioned waveforms mitigate the onset response in high frequency nerve block. J Neural Eng. 2010 Dec; 7(6): 066003. pmid:20966536
  14. 14. Miles JD, Kilgore KL, Bhadra N, Lahowetz EA. Effects of ramped amplitude waveforms on the onset response of high-frequency mammalian nerve block. J Neural Eng. 2007 Dec; 4(4): 390–8. pmid:18057506
  15. 15. Vrabec TL, Eggers TE, Foldes EL, Ackermann DM, Kilgore KL, Bhadra N. Reduction of the onset response in kilohertz frequency alternating current nerve block with amplitude ramps from non-zero amplitudes. J Neuroeng Rehabil. 2019 Jun; 16(1): 80. pmid:31253152
  16. 16. Bhadra N, Kilgore KL. Direct current electrical conduction block of peripheral nerve. IEEE Trans Neural Syst Rehabil Eng. 2004 Sep; 12(3): 313–24. pmid:15473193
  17. 17. Ackermann D. M. Jr., Bhadra N., Foldes E. L., Wang X. F. & Kilgore K. L. Effect of nerve cuff electrode geometry on onset response firing in high-frequency nerve conduction block. IEEE Trans Neural Syst Rehabil Eng. 2010 Dec; 18(6): 658–65. pmid:20813650
  18. 18. Bähring R, Covarrubias M. Mechanisms of closed-state inactivation in voltage-gated ion channels. J Physiol. 2011 Feb 1; 589(Pt 3): 461–79. pmid:21098008
  19. 19. Aldrich RW, Corey DP, Stevens CF. A reinterpretation of mammalian sodium channel gating based on single channel recording. Nature. 1983 Dec 1–7; 306(5942): 436–41. pmid:6316158
  20. 20. Scanley BE, Hanck DA, Chay T, Fozzard HA. Kinetic analysis of single sodium channels from canine cardiac Purkinje cells. J Gen Physiol. 1990 Mar; 95(3): 411–37. pmid:2157791
  21. 21. Armstrong CM. Na channel inactivation from open and closed states. Proc. Proc Natl Acad Sci U S A. 2006 Nov 21; 103(47): 17991–6. pmid:17101981
  22. 22. Pal K, Gangopadhyay G. Dynamical characterization of inactivation path in voltage-gated Na(+) ion channel by non-equilibrium response spectroscopy. Channels. 2016 Nov; 10(6): 478–97. pmid:27367642
  23. 23. Dharmawan T, Nakajima T, Iizuka T, Tamura S, Matsui H, Kaneko Y, Kurabayashi M. Enhanced closed-state inactivation of mutant cardiac sodium channels (SCN5A N1541D and R1632C) through different mechanisms. J Mol Cell Cardiol. 2019 May; 130: 88–95. pmid:30935997
  24. 24. Lawrence JH, Yue DT, Rose WC, Marban E. Sodium channel inactivation from resting states guinea-pig ventricular myocytes. J Physiol. 1991 Nov; 443: 629–50. pmid:1668345
  25. 25. Groome J, Lehmann-Horn F, Holzherr B. Open- and closed-state fast inactivation in sodium channels: differential effects of a site-3 anemone toxin. Channels 2011 Jan-Feb; 5(1): 65–78. pmid:21099342
  26. 26. Goldman L. Sodium channel inactivation from closed states: evidence for an intrinsic voltage dependency. Biophys J. 1995 Dec; 69(6): 2369–77. pmid:8599643
  27. 27. Milescu LS, Yamanishi T, Ptak K, Smith JC. Kinetic properties and functional dynamics of sodium channels during repetitive spiking in a slow pacemaker neuron. J Neurosci. 2010 Sep 8; 30(36): 12113–27. pmid:20826674
  28. 28. Cummins TR, Howe JR, Waxman SG. Slow closed-state inactivation: a novel mechanism underlying ramp currents in cells expressing the hNE/PN1 sodium channel. J Neurosci. 1998 Dec 1; 18(23): 9607–19. pmid:9822722
  29. 29. Herzog RI, Cummins TR, Ghassemi F, Dib-Hajj SD, Waxman SG. Distinct repriming and closed-state inactivation kinetics of Nav1.6 and Nav1.7 sodium channels in mouse spinal sensory neurons. J Physiol. 2003 Sep; 551(Pt 3): 741–50. pmid:12843211
  30. 30. Balbi P, Massobrio P, Hellgren Kotaleski J. A single Markov-type kinetic model accounting for the macroscopic currents of all human voltage-gated sodium channel isoforms. PLoS Comput Biol. 2017 Sep 1; 13(9): e1005737. pmid:28863150
  31. 31. McIntyre CC, Richardson AG, Grill WM. Modeling the excitability of mammalian nerve fibers: influence of afterpotentials on the recovery cycle. J Neurophysiol. 2002 Feb; 87(2): 995–1006. pmid:11826063
  32. 32. Kennedy J, Eberhart R. Particle swarm optimization. IEEE Int Conf Neural Networks Proc. 1995; 4: 1942–8.
  33. 33. Peña E, Zhang S, Patriat R, Aman JE, Vitek JL, Harel N, Johnson MD. Multi-objective particle swarm optimization for postoperative deep brain stimulation targeting of subthalamic nucleus pathways. J Neural Eng. 2018 Dec; 15(6): 066020. pmid:30211697
  34. 34. Peña E, Zhang S, Deyo S, Xiao Y, Johnson MD. Particle swarm optimization for programming deep brain stimulation arrays. J Neural Eng. 2017 Feb; 14(1): 016014. pmid:28068291
  35. 35. McCreery DB, Agnew WF, Yuen TG, Bullara LA. Comparison of neural damage induced by electrical stimulation with faradaic and capacitor electrodes. Ann Biomed Eng. 1988; 16(5): 463–81. pmid:3189974
  36. 36. Cogan SF. Neural stimulation and recording electrodes. Annu Rev Biomed Eng. 2008; 10: 275–309. pmid:18429704
  37. 37. Gwilliam JC, Horch K. A charge-balanced pulse generator for nerve stimulation applications. J Neurosci Methods. 2008 Feb 15; 168(1): 146–50. pmid:17950907
  38. 38. Pelot NA, Behrend CE, Grill WM. Modeling the response of small myelinated axons in a compound nerve to kilohertz frequency signals. J Neural Eng. 2017 Aug; 14(4): 046022. pmid:28361793
  39. 39. Bhadra N, Lahowetz EA, Foldes ST, Kilgore KL. Simulation of high-frequency sinusoidal electrical block of mammalian myelinated axons. J Comput Neurosci. 2007 Jun; 22(3): 313–26. pmid:17200886
  40. 40. Milescu LS, Bean BP, Smith JC. Isolation of somatic Na+ currents by selective inactivation of axonal channels with a voltage prepulse. J Neurosci. 2010 Jun; 30(22): 7740–8. pmid:20519549
  41. 41. Tai C, de Groat WC, Roppolo JR. Simulation analysis of conduction block in unmyelinated axons induced by high-frequency biphasic electrical currents. IEEE Trans Biomed Eng. 2005 Jul; 52(7): 1323–32. pmid:16041996
  42. 42. Bhadra N, Bhadra N, Kilgore K, Gustafson KJ. High frequency electrical conduction block of the pudendal nerve. J Neural Eng. 2006 Jun; 3(2): 180–7. pmid:16705274
  43. 43. Ackermann DM Jr, Bhadra N, Foldes EL, Kilgore KL. Conduction block of whole nerve without onset firing using combined high frequency and direct current. Med Biol Eng Comput. 2011 Feb; 49(2): 241–51. pmid:20890673
  44. 44. McIntyre CC, Grill WM, Sherman DL, Thakor NV. Cellular effects of deep brain stimulation: model-based analysis of activation and inhibition. J Neurophysiol. 2004 Apr; 91(4): 1457–69. pmid:14668299
  45. 45. Ackermann DM Jr, Foldes EL, Bhadra N, Kilgore KL. Effect of bipolar cuff electrode design on block thresholds in high-frequency electrical neural conduction block. IEEE Trans Neural Syst Rehabil Eng. 2009 Oct; 17(5): 469–77. pmid:19840914
  46. 46. Hines ML, Carnevale NT. Neuron: a tool for neuroscientists. Neuroscientist. 2001 Apr; 7(2): 123–35. pmid:11496923
  47. 47. Armstrong CM. Sodium channels and gating currents. Physiol Rev. 1981 Jul; 61(3): 644–83. pmid:6265962
  48. 48. Bezanilla F. Gating of sodium and potassium channels. J Membr Biol. 1985; 88(2): 97–111. pmid:2419568
  49. 49. Yang F, Anderson M, He S, Stephens K, Zheng Y, Chen Z, et al. Differential expression of voltage-gated sodium channels in afferent neurons renders selective neural block by ionic direct current. Sci Adv. 2018 Apr 11; 4(4): eaaq1438. pmid:29651458
  50. 50. Black JA, Dib-Hajj S, McNabola K, Jeste S, Rizzo MA, Kocsis JD, et al. Spinal sensory neurons express multiple sodium channel alpha-subunit mRNAs. Brain Res Mol Brain Res. 1996 Dec 31; 43(1–2): 117–31. pmid:9037525
  51. 51. Caldwell JH, Schaller KL, Lasher RS, Peles E, Levinson SR. Sodium channel Nav1.6 is localized at nodes of ranvier, dendrites, and synapses. Proc Natl Acad Sci U S A. 2000 May 9; 97(10): 5616–20. pmid:10779552
  52. 52. McIntyre CC, Grill WM. Excitation of central nervous system neurons by nonuniform electric fields. Biophys J. 1999 Feb; 76(2): 878–88. pmid:9929489
  53. 53. Yi GS, Grill WM. Frequency-dependent antidromic activation in thalamocortical relay neurons: effects of synaptic inputs. J Neural Eng. 2018 Oct; 15(5): 056001. pmid:29893711
  54. 54. Zhang TC, Grill WM. Modeling deep brain stimulation: point source approximation versus realistic representation of the electrode. J Neural Eng. 2010 Dec; 7(6): 066009. pmid:21084730