Search tips
Search criteria 


Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
Brain Stimul. Author manuscript; available in PMC 2011 April 1.
Published in final edited form as:
Brain Stimul. 2010 April 1; 3(2): 65–67.
doi:  10.1016/j.brs.2010.01.003
PMCID: PMC2895675

Patient-specific models of deep brain stimulation: Influence of field model complexity on neural activation predictions


Deep brain stimulation (DBS) of the subthalamic nucleus (STN) has become the surgical therapy of choice for medically intractable Parkinson’s disease. However, quantitative understanding of the interaction between the electric field generated by DBS and the underlying neural tissue is limited. Recently, computational models of varying levels of complexity have been used to study the neural response to DBS. The goal of this study was to evaluate the quantitative impact of incrementally incorporating increasing levels of complexity into computer models of STN DBS. Our analysis focused on the direct activation of experimentally measureable fiber pathways within the internal capsule (IC). Our model system was customized to an STN DBS patient and stimulation thresholds for activation of IC axons were calculated with electric field models that ranged from an electrostatic, homogenous, isotropic model to one that explicitly incorporated the voltage-drop and capacitance of the electrode-electrolyte interface, tissue encapsulation of the electrode, and diffusion-tensor based 3D tissue anisotropy and inhomogeneity. The model predictions were compared to experimental IC activation defined from electromyographic (EMG) recordings from eight different muscle groups in the contralateral arm and leg of the STN DBS patient. Coupled evaluation of the model and experimental data showed that the most realistic predictions of axonal thresholds were achieved with the most detailed model. Furthermore, the more simplistic neurostimulation models substantially overestimated the spatial extent of neural activation.

Keywords: deep brain stimulation, computational modeling, neural activation, Parkinson's disease


The subthalamic nucleus (STN) is an integral component of the basal ganglia and is known to play an important role in the pathophysiology of Parkinson’s disease (PD) (1, 2). Chronic high frequency deep brain stimulation (DBS) of the STN and its surrounding structures has become an established therapy for the treatment of PD (3, 4). However, scientific understanding of the cellular effects and physiological mechanisms of DBS remains largely incomplete (5, 6).

The fundamental purpose of DBS is to modulate neural activity with applied electric fields. Unfortunately, quantitative understanding of the effects of manipulating various stimulation parameters (frequency, pulse-width, and amplitude) on the neural response to DBS is lacking. In turn, numerous investigators have worked to develop a wide range of computational models to predict the electric field and stimulating influence generated by DBS (727). Recently, these DBS modeling efforts have focused on the integration of patient imaging data and finite element models of the DBS electric field (28). However, an important issue that has not been explicitly addressed is the degree of model complexity that is needed to make accurate predictions on the neural response to DBS.

When building a model of a complex system such as DBS, it is important to know the impact of various simplifications and assumptions on the predictive capabilities of the model. Neurostimulation electric field models have traditionally assumed electrostatic conditions with perfect voltage coupling between the electrode and tissue medium, as well as simplifying the conductivity of the surrounding tissue to be homogeneous and isotropic. However, clinical DBS electrodes are placed in an anisotropic and inhomogenous tissue medium (12, 18, 21), capacitive components of the electrode-tissue interface limit the applicability of the electrostatic assumption (14, 26), and a substantial voltage drop occurs at the electrode interface during charge transduction from the electrode surface to the ionic medium (19, 29). Therefore, the focus of this study was to evaluate the quantitative role of these factors on stimulation predictions in DBS models.

This study concentrated on DBS of the corticospinal tract (CST) with electrodes implanted within the subthalamic region. The CST is a major fiber pathway within the internal capsule, which defines the lateral border of the STN. Consequently, motor evoked responses from activation of larger-diameter CST fibers can be elicited with relatively low thresholds during STN DBS (3033). Clinically, CST activation is an unwanted side effect of DBS (34). However, the generation of muscle contractions via stimulation of the CST represents a direct link between DBS, known neural substrates, and clinically measurable behaviors. In turn, we used experimental measurements of CST activation during DBS to address the level of model complexity required to accurately predict stimulation induced neural activation. Preliminary portions of this study were previously presented as a conference paper (35).

Materials and Methods

This study used a series of DBS computational models, customized to an individual human patient, following and expanding upon methodology previously described in Butson et al. (12). Our fundamental goal was to evaluate the quantitative importance of four components of voltage-controlled DBS electric field models: 1) electrode interface voltage drop, 2) electrode interface capacitance, 3) tissue encapsulation of the electrode, and 4) tissue anisotropy/inhomogeneity. The modeling system combined both anatomical and diffusion tensor magnetic resonance imaging (MRI) data. Pre- and post-operative T1 images were used to position and align 3D surfaces representing anatomical nuclei of interest (i.e. thalamus and STN), as well as to determine the position of the DBS electrode within the patient’s brain. Diffusion tensor imaging (DTI) data (36) was used to both define axonal trajectories of the internal capsule (IC), and to estimate 3D tissue anisotropy and inhomogeneity in the tissue region surrounding the DBS electrode (12).

The electric field generated by DBS was calculated with finite element models (FEMs). We created five variants of the DBS FEM to address differences associated with the degree of model complexity. Model I was the most simplistic, an electrostatic model that ignored the interface voltage drop, electrode capacitance, encapsulation, and tissues anisotropy/inhomogeneity. Models II-V incrementally added explicit representations of the electrode interface with the brain, and the tissue anisotropy/inhomogeneity (Table 1).

Table 1
Characteristics of the DBS electric field models.

The electric field generated by each DBS FEM was applied to detailed multi-compartment cable models of myelinated axons which had trajectories defined by DTI tractography (Fig. 1). Stimulation thresholds were calculated for each IC axon by applying the extracellular voltage distribution generated along the axon trajectory for each of the five variants of the FEM. These model results were then compared to experimentally defined CST thresholds acquired using electromyogram (EMG) recordings from the patient.

Fig. 1
Patient-specific DBS model

Clinical Data

The study received prior approval from the Cleveland Clinic Institutional Review Board, and the patient provided informed written consent. The subject was a 63-year old male patient with Parkinson’s disease, previously implanted with a Soletra pulse generator and 3387 DBS electrode (Medtronic Inc, Minneapolis, MN) in the STN region, who exhibited good therapeutic benefit from the device based on the Unified Parkinson’s Disease Rating Scale. Our models used the patient’s pre-operative and post-operative high-resolution T1-weighted MRIs, acquired on a Siemens Symphony 1.5 T scanner and on a Siemens 1.5 T Magnetom Vision, respectively. Both images were acquired with a 256 mm × 256 mm field of view and were interpolated to have a 1 mm3 isotropic voxel resolution. The post-operative MRI was performed with imaging parameters previously defined as being safe by extensive phantom testing in the specific scanner used to acquire the images (37).

The clinical experiments were conducted at a time point greater than one year post-surgery. Differential EMG recordings were made with electrode pairs placed over the biceps, triceps, flexor carpi ulnaris, extensor carpi radialis, quadriceps, tibialis anterior, and lateral gastrocnemius. During these experiments, 20-second recording epochs were gathered while the patient experienced unilateral low frequency, monopolar stimulation (5 Hz, 60 µsec pulse width, 0 to −10 V in −1 V increments). These recordings were individually performed with stimulation applied through each of the four DBS electrode contacts. The results reported in this paper are for the right-side DBS electrode (measurements and recordings were made on the left arm and leg).

EMG activity was recorded with a Biotop 6R12 amplifier with the following settings: low frequency filter at 5 Hz and a high frequency filter at 1500 Hz, with a 1 mV/full-scale, where the full-scale was 6 V. Signals were subsequently filtered with a ninth-order Butterworth high-pass filter with a cutoff frequency of 50 Hz to remove any low frequency baseline drift. The stimulus artifact was recorded using a surface electrode from the connecting lead on the patient’s neck. Using this as the trigger event, time-triggered average EMG signals were computed from the other channels. These signals were analyzed for a threshold response indicating a muscle twitch, which we interpreted as stimulation spillover into the IC and activation of CST fibers. To simplify the data presentation, EMG thresholds are reported with muscles divided into two distinct groups: arm (bicep, tricep, flexor carpi ulnaris, extensor carpi radialis) and leg (quadriceps, tibialis anterior, lateral gastrocnemius).

Image Registration and Anatomical Nuclei

MRI data formed the basis for the patient-specific DBS computer models. The patient’s pre- and post-operative MRI datasets were co-registered with the Wakana et al. (36) diffusion tensor atlas brain using Analyze 6.0 (Lenexa, KS). The 3D co-registration algorithm was used within the Insight Toolkit feature of Analyze, and it implemented an intensity-based stochastic approach (38). The DTI atlas brain was acquired with a 2-mm3 isotropic voxel size with a diffusion gradient weighting of 700 mm2/s (36).

The general structure of the patient-specific DBS computer model was created following the methodology of Butson et al (12) (Fig. 1). Graphical representations of relevant anatomical nuclei (STN and thalamus) were defined by warping 3D surfaces to fit the patient’s pre-operative MRI data using a non-linear algorithm (39), originally developed by Surgical Navigation Technologies (now Medtronic Navigation, Louisville, CO). The electrode tip location and insertion trajectory were determined by segmenting the electrode from the post-operative MRI. This procedure used an image thresholding method to extract out the dark, hollowed cavity artifact created in the MR by the electrode. A virtual replica of the Medtronic 3387 DBS electrode, linked to a multi-resolution finite element mesh, was then placed at that location within the image volume. The anatomical and diffusion tensor MRI data, along with the anatomical nuclei and patient-specific electrode position were all loaded into a common 3D visualization and simulation environment SCIRun/BioPSE (Scientific Computing and Imaging Institute, University of Utah, Salt Lake City, UT). (Fig. 1).

Internal Capsule Tractography

Individual tensors from the DTI atlas brain, along with a fiber tractography algorithm within SCIRun/BioPSE (40), were used to extract individual model axon trajectories within the internal capsule. Two hundred forty seed points were placed equidistantly from one another within a 6 × 0.5 × 1 mm rectangular region just lateral to the STN. The resulting 240 trajectories were used to represent a population of IC axons in the patient-specific DBS model (Fig. 1E).

Multi-compartment cable models of myelinated axons were created for each of the 240 IC fiber trajectories. These cable models, 10 µm in diameter, included detailed representations of the nodes of Ranvier, paranodal, and internodal sections of the individual axons (41). Each axon had 51 nodes of Ranvier and 551 total compartments along its 50 mm path length.

Electrical Model

A multi-resolution finite element mesh of the DBS electrode and surrounding tissue medium was constructed using FEMLAB 3.1 (Comsol Inc., Burlington, MA). This 3D mesh consisted of over 4.2 million nodes, most of which were located circumferentially around the electrode to provide for greater resolution near the stimulating contacts. The same mesh was used for all variants of the DBS FEM (Table 1). The Poisson equation was solved in 3D to determine the voltage distribution generated in the tissue medium by the DBS electrode (Fig. 1D). The FEM solutions were performed on an 8-processor 32 GB shared-memory SGI Prism (Silicon Graphics Inc., Mountain View, CA).

The voltage solutions from each variant of the DBS FEM were linearly interpolated onto the center of every compartment of each axon trajectory (Fig. 1F). Simulations of the neural response to the applied field were preformed in NEURON 6.1.2 (42). The stimulus waveform (60 µsec pulse width) applied to the axon models mimicked the output of the Soletra pulse generator implanted in the patient (43). Each of the 240 model axons had an activation threshold for each model variant that was defined as the minimum stimulus voltage necessary to generate a propagating action potential.

Model Evaluations

Internal capsule axon activation during monopolar DBS was evaluated at the patient’s clinically defined therapeutic electrode contact (contact 3), as well as each of the other three contacts. Five models with increasing levels of complexity were examined for each of the four contacts (Table 1). The simplest model (Model I) consisted of a homogeneous and isotropic tissue medium (0.3 S/m) with no electrode encapsulation, and stimulation was applied under electrostatic conditions with no voltage drop at the electrode interface (Fig. 2A). Model II was identical to Model I, but included the 42% voltage drop at the electrode interface (see Appendix) (Fig. A1) (Fig. 2B). A slightly more complex model (Model III) integrated electrode capacitance (3.3 µF), producing a more realistic simulation waveform in the tissue medium (14) (Fig. 2C). The fourth model (Model IV) incorporated a 0.5 mm tissue encapsulation layer (0.18 S/m) around the electrode, to account for the chronic electrode impedance (~900 Ω) estimated by the patient’s implanted pulse generator (13) (Fig. 2D). Finally, the most complex model (Model V) added the diffusion tensor based tissue conductivities to represent the anisotropic and inhomogenous bulk tissue medium (12, 18, 19, 44) (Fig. 2E). A simple linear transform (0.8 (S-s)/mm2 scaling factor) was used to convert the diffusion tensors into conductivity tensors (19, 44).

Fig. 2
DBS FEM comparison
Fig. A1
Characterizing the electrode-electrolyte interface voltage drop


Voltage Distribution Generated by DBS

Figure 2 describes how the various model characteristics affected the voltage distribution in the tissue medium when a −1 V stimulus (as programmed into the pulse generator) was applied. Incorporation of the 42% voltage drop at the electrode-tissue interface produced the greatest attenuation of stimulus spread (Model II), while adding the electrode encapsulation had the next largest effect (Model IV). Including electrode capacitance modified the shape of the stimulus pulse. Finally, representing the tissue medium with DTI-based conductivity tensors produced an asymmetric, non-spherical voltage spread (Model V).

Action Potential Initiation

Extracellular electrical stimulation of myelinated axons with a monopolar cathode generates both depolarization and hyperpolarization along the axon. Depolarization occurs in the nodes of Ranvier closest to the active electrode contact; whereas, hyperpolarization occurs in flanking regions of the axon (45). Action potential initiation occurs in the node of Ranvier with the greatest stimulation induced depolarization, and then propagates in both directions (antidromically and orthodromically). For each model variant and each of the four electrode contacts, we determined the threshold stimulation voltage required for action potential initiation and propagation in each of the 240 modeled IC axons. We also calculated the shortest distances from the respective axon to the center of the stimulating contact on the DBS electrode. The voltage-distance relationships for each variant of the FEM are illustrated in Figure 3. These results show that axonal threshold stimulation voltages were nonlinearly related to electrode-to-axon distance. In addition, the different FEMs generated very different activation predictions, with the more simple models predicting lower thresholds than the more detailed models.

Fig. 3
Voltage-distance relationship

Voltage-Distance Relationship

The results presented in Figure 3 provide the opportunity to evaluate the potential utility of analytical equations to predict the spread of stimulation. The voltage-distance equation:


where Vth is the estimated threshold voltage necessary to activate an IC axon, V0 is an offset from baseline, r is the distance from the center of the stimulating contact to the axon, and the constant k is the slope, has been used in many applications to predict the radius of activation from extracellular stimulation. Therefore, we performed least-squares fits of the voltage-distance equation to the voltage-distance relationships of the models (i.e. axon thresholds from each contact for each FEM variant) (Table 2). These fits were completed using MATLAB, and V0 was forced to be ≥ 0, while k was left unconstrained. Given those parameter bounds, the best fits for all models were achieved with V0 = 0. We noted substantial variability in the slope (k) of the voltage-distance equation when comparing models of different complexity. Further, when using the most detailed model (Model V), the fitted parameters were not consistent across the different electrode contacts in the same patient, due to the complex electrical properties of the 3D tissue medium. These results bring into question the accuracy of using the voltage-distance equation, especially with an arbitrarily defined k value, to quantify stimulation spread on a patient-specific basis.

Table 2
Regression fit parameters for the slope k (V/mm2) of the voltage-distance equation for all five models and four contacts.

Clinical Thresholds and Model Comparison

EMG defined threshold voltages for stimulation-induced muscle twitches were calculated for each recorded muscle. To simplify the data presentation, the clinical thresholds were combined into one of two generalized muscle groups (arm or leg). The clinical activation thresholds for the arm and leg muscle groups were then compared to the model predictions of IC axonal activation (Table 3).

Table 3
EMG threshold results for the arm and leg muscle groups, and their respective model predictions for the recruitment of IC axons.

Physiologically it may be possible for activation of a single cortico-motoneuronal fiber to generate a muscle twitch (46). However, it is likely that the muscle responses we recorded were the result of simultaneous activation of multiple axons, especially given the resolution of the experimental testing (1 V increments). Therefore, our expectation was that the DBS model should predict ~5–15% IC activation at the experimentally defined threshold voltage. For example, the EMG threshold for the triceps muscle for contact 3 in this patient was at −5 V. Model V predicted no IC activation at −4 V (sub-threshold), 15% activation at −5 V (threshold), and 36% activation at −6 V (super-threshold) (Fig. 4). As the stimulation voltage increased, additional fibers were recruited in a non-linear fashion (Fig. 3).

Fig. 4
Comparison of model and experimental results

The more simplistic models (Models I, II, and III) all generated predictions of excessive axonal activation at the experimentally defined thresholds. These models excluded some or all of the major components defining the electrode interface with the brain. All three models also lacked a sheath of resistive tissue that typically encompasses the electrode after chronic implantation (13, 47). Consequently, these models suffered from an underestimation the electrode impedance. For example, at the clinically defined therapeutic contact (contact 3) the impedance measured by the Medtronic IPG was 956 Ω. Models I–V exhibited a corresponding impedance of 409, 636, 636, 1129, and 960 Ω, respectively.

Model IV exhibited reduced axonal activation when compared to Models I-III, which was attributed to the inclusion of electrode encapsulation. However, Model IV also showed substantial variability in the percentage of axons activated through the range of experimental thresholds measured at the various electrode contacts. This variability was credited to failure to account for the 3D tissue conductivity differences surrounding the different contacts. Model V represented the most detailed model and generated predictions that most consistently corresponded to the level of axonal activation expected at the experimental thresholds defined for each electrode contact.

Sensitivity Analysis of Model Activation Predictions

Based on our history of developing patient-specific computational models for DBS applications, we believe that given an appropriately parameterized electrical model, limitations in defining the actual electrode position in the brain represent the next most important source of error. Post-operative imaging artifacts limit certainty in defining the actual electrode location (48), and the potential for brain shift during surgery (49, 50) limits confidence in using frame-based stereotactic coordinates relative to a pre-operative image. Therefore, we performed a sensitivity analysis on electrode position using the two most relevant DBS FEMs (Models IV and V). We calculated the axonal activation induced by −5 V stimulation through contact 3 (corresponding to Figure 4) after moving the electrode a total of ±1.0 mm in 0.25 mm increments within the transverse (axial) xy-plane (Figure 5). Not surprisingly, Model IV showed that as the electrode placement got closer to the fiber pathway the percentage of activated axons increased (i.e. electrode movement in the anteriolateral direction). However, incorporation of DTI-based tissue conductivities (Model V) minimized the impact of electrode placement uncertainty on the activation predictions. Over the evaluated space, the activation ranges for Model IV and Model V were 16–44% and 7–17%, respectively. These results further reinforce the role of 3D tissue electrical properties in DBS FEMs, as stimulation induced activation is dictated by the second spatial derivative of the voltage distribution along the axon. Voltage spread is enhanced parallel to (hindered perpendicular to) the preferred direction of anisotropy (i.e. the IC axon trajectory) (Fig. 1C,D,E). Hence, axons within highly anisotropic fiber pathways are somewhat shielded from activation, relative to axons within more isotropic brain regions.

Fig. 5
Sensitivity analysis


The clinical success of DBS has prompted substantial scientific interest in characterizing its underlying effects on the nervous system. However, experimental analyses of DBS often encounter substantial difficulties in controlling all of the relevant variables, and interpretation of the results can be ambiguous. In turn, computational analyses have been employed to provide quantitative guidance on the response of neurons to DBS electric fields (28). Recently, numerous DBS models have been developed, spanning a wide spectrum of detail and complexity. The goal of this study was to quantitatively address the impact of DBS electric field model complexity on the spread of stimulation in a clinically relevant context. The results show that each of the four features examined in this study (electrode interface voltage drop, electrode capacitance, electrode encapsulation, and bulk tissue anisotropy/inhomogeneity) impacted the model predictions, and should be considered when creating and/or using DBS models. Future computational studies intending to draw correlations between patient-specific DBS parameter settings and clinical outcomes should make every effort to use the most accurate model possible, and at a minimum avoid the gross simplifications of Models I-III or the voltage-distance equation. However, it should also be noted that while Model V was the most detailed, it still has substantial room for improvement assuming parallel advances in medical imaging technology.

Data Integration and Study Limitations

The methodology implemented in this study required the integration of multiple forms of both computational and clinical data. Although several different software tools and computer algorithms were used in various phases of the project, SCIRun enabled us to place everything into a common coordinate system and visualization platform. Nevertheless, an inherent issue of such a study is the accuracy in which different data sets are co-registered. Our MRI data and 3D nuclei were co-registered with state-of-the-art algorithms, but such registrations still have associated errors on the order of 1 mm (38, 39). Further, DTI tractography only provides a rough estimate of the fiber tract and/or individual axon trajectories (51). It should also be noted that estimation of the voltage drop at the electrode-electrolyte interface as a simple linear percentage fails to account for the complex interactions of charge transduction between metal electrodes and the ionic medium (19). Given these caveats, we were able to integrate detailed anatomical information on DBS electrode location, relative to a known fiber tract, and simulate direct activation of that fiber tract. This allowed for quantitative comparisons between model predictions of axonal activation and clinically measured EMG responses in a DBS patient. Our results show that when an appropriate model is used, accurate stimulation predictions can be made, but inappropriate models (i.e. Models I-III) provide very poor predictions.

The methodology presented in this study concentrated on direct stimulation of the IC near the STN. We chose this neural population as the focus of our analysis because measurements of the activation of these axons could be performed using simple EMG recordings in awake, permanently implanted human patients. While other neural populations may be more relevant to the desired clinical effects of DBS (e.g. STN projection neurons, GPi fibers of passage, SNc fibers of passage, cortical afferent inputs) (20, 52, 53), none of them represent a neural entity which can be easily reconstructed via DTI tractography, or have a simple and direct behavioral effect from stimulation that can be measured non-invasively. It should also be noted that the basic neural response to extracellular electrical stimulation is dictated by the axon (54), and the basic biophysics of how axonal activation occurs is independent of neuron type or fiber diameter (45, 55). In turn, the IC represents the largest and most easily accessible population of axons near DBS electrodes; thereby representing an excellent medium for studying the neural response to DBS.

Our simulated IC fiber bundle was comprised of 240 uniformly distributed fibers, each 10 µm in diameter. The real IC actually contains a wide range of fiber diameters and other associated axonal properties which affect the threshold for action potential generation (41, 55). In general, larger diameter fibers have lower thresholds than smaller diameter fibers in response to extracellular stimulation. While most fibers within the human IC have diameters less than 4 µm, a substantial number of fibers have diameters of ~10 µm (56). These large diameter fibers represent some of the most excitable neural elements in the STN region. The compound action potential and subsequent motor responses commonly measured from IC stimulation with STN DBS electrodes are associated with activation of these larger diameter fibers (3034), and consequently were the focus of our analysis.

The fundamental goal of this study was to demonstrate the relative impact of each DBS electric field model characteristic (i.e. tissue encapsulation, electrode capacitance, voltage-drop at the interface, and tissue inhomogeneities) on a functionally relevant outcome, specifically IC axonal activation. For simplicity sake, we intentionally avoided detailed parameter sensitivity analyses of each parameter of each model variant, relying on the fact that the parameter values we used were good estimates based on the available experimental data in the literature. However, it should be noted that given a priori knowledge of the desired axonal activation thresholds, it would be possible to optimally fit each model variant to the experimental data, albeit with model parameter values that are perhaps outside of realistic experimental ranges. Unfortunately, this circular exercise does not validate the model, and such a practice could lead to extraneous and/or inaccurate predictions when the model is used to interpret new experimental/clinical data.

Clinical Significance

We examined stimulation of the internal capsule in this study for two basic reasons. First, the generation of muscle contractions via stimulation of the IC represents a direct link between STN DBS, known neural substrates, and clinically measurable behaviors. In turn, we were able to make a connection between our patient-specific DBS models and experimental data recorded from that patient. Second, because IC activation is a relatively common unwanted side effect of STN DBS it is important to understand the stimulation conditions that control it. In turn, one possible application of the modeling techniques presented in this study would be to provide visual feedback to the clinician and help them identify techniques to avoid IC activation with DBS (12, 57). Software technology employing such models could be used intra-operatively to assist in optimizing DBS electrode placement (58), and post-operatively to assist in the definition of therapeutic stimulation parameter settings (59).

Our simulation results suggest that DBS induced axonal activation depends on a long list of factors, many of which can be accounted for with an appropriate model. We have previously come to similar conclusions using a range of coupled simulations and experiments in both humans and monkeys, but DBS models still have numerous limitations. Nonetheless, the DBS modeling community is beginning to develop the computational infrastructure and scientific methodology required to account for most of the relevant factors impacting model prediction accuracy. In turn, we foresee many new opportunities to utilize the coupled analysis of clinical data and computational models to evaluate the effects of DBS (60). However, the results of this study show that when performing such analysis it is important to use the right model for the task at hand. Our results indicate that many of the commonly employed simplifying assumptions in neurostimulation modeling generate a substantial overestimation of stimulation spread. In turn, the standard for scientific studies attempting to relate DBS FEMs to clinical data should be to use models at the level of at least Model IV or V, and new efforts are warranted to further improve the predictive capabilities of these models.


This work was supported by grants from the National Institutes of Health (R01 NS059736, R21 NS050449, F32 NS052042). The authors would also like to thank Jaimie Henderson for providing the 3D nuclei surfaces, Susumu Mori for providing the diffusion tensor image brain atlas, Christopher Maks and Svjetlana Miocinovic for assistance with the model simulations, and Barbara Wolgamuth for assistance with clinical threshold data collection.


In Vitro Characterization of the DBS Voltage Distribution

Following and expanding upon methodology described in Miocinovic et al. (19), we characterized the voltage drop at the electrode-electrolyte interface with in vitro experiments on a Medtronic (Minneapolis, MN) 3387 human DBS electrode (Fig. A1). The lead was suspended inside a glass beaker filled with saline and placed on an electric hot plate within an electrically shielded Faraday cage. The Faraday cage allowed for recording at a higher signal-to-noise ratio, and prevented additional noise potentially caused by external electric fields. The 600 mL glass beaker was 8 cm in diameter and was filled with a solution of 0.9% NaCl, heated to 37 °C. A stainless steel coil was used for the return electrode, and was loosely wound around the inner wall of the beaker. Both, an Ag|AgCl wire and a tungsten microelectrode (FHC, Bowdoin, ME) were suspended within the saline solution to serve as the reference and recording electrodes, respectively.

Stimulus pulses were delivered through the DBS electrode with a 0.7 V, 5 Hz, 60 µs pulse-width waveform, generated by a Medtronic Itrel II implantable pulse generator (IPG). The low stimulation voltage was used to prevent saturation in the recording amplifier. Voltages were recorded at specific points along seven different microelectrode recording tracks parallel to the DBS electrode (Fig. A1–A). Each point in each track was acquired sequentially by moving the recording microelectrode relative to the DBS electrode using a micromanipulator (World Precision Instruments, Sarasota, FL). High-resolution photographs were taken to verify the microelectrode location relative to the DBS electrode. The recorded signals were band-pass filtered between 1 Hz and 20 kHz using a differential amplifier (A-M Systems, Model 3000, Sequim, WA), digitized at a sampling rate of 100 kHz, and stored for offline analysis (Cambridge Electronic Design, Power 1401 and Spike2 software, Cambridge, UK). The analysis involved averaging the peak voltages during each 20-second acquisition.

An FEM was created to mimic the in vitro experiments. The in vitro FEM relied on a mesh consisting of nearly 3.7 million elements. The mesh density was highest at the electrode contact and element size increased further away from the electrode. A cylindrical boundary was defined 8 cm from the DBS electrode, mimicking the dimensions of the glass beaker, and it was set to ground. Nodes on the active electrode surface (contact 1) were used as voltage sources, which was consistent with the voltage-controlled stimulation employed by the Itrel II pulse generator. The Poisson equation was solved to determine voltage as a function of space within the saline medium which was assumed to be homogenous and isotropic (2 S/m). The FEM was iteratively solved to identify the voltage drop at the electrode-electrolyte interface that minimized the error between the experimentally recorded voltages and the model solutions, as previously described (19) (Fig. A1). The model predicted a 42% voltage drop at the electrode interface for the Medtronic 3387 human DBS electrode at 37 °C.


Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Conflicts of interest

CCM and CRB authored intellectual properties related to the project methodology, and are shareholders in Intelect Medical Inc. CCM, CRB, and AC are paid consultants for Intelect Medical Inc.


1. DeLong MR, Wichmann T. Circuits and circuit disorders of the basal ganglia. Arch Neurol. 2007;64:20–24. [PubMed]
2. Rivlin-Etzion M, Marmor O, Heimer G, et al. Basal ganglia oscillations and pathophysiology of movement disorders. Curr Opin Neurobiol. 2006;16:629–637. [PubMed]
3. Deuschl G, Schade-Brittinger C, Krack P, et al. A randomized trial of deep-brain stimulation for Parkinson's disease. N Engl J Med. 2006;355:896–908. [PubMed]
4. Obeso JA, Olanow CW, Rodriguez-Oroz MC, et al. Deep-brain stimulation of the subthalamic nucleus or the pars interna of the globus pallidus in Parkinson's disease. N Engl J Med. 2001;345:956–963. [PubMed]
5. Lozano AM, Dostrovsky J, Chen R, Ashby P. Deep brain stimulation for Parkinson's disease: disrupting the disruption. Lancet Neurol. 2002;1:225–231. [PubMed]
6. McIntyre CC, Savasta M, Kerkerian-Le Goff L, Vitek JL. Uncovering the mechanism(s) of action of deep brain stimulation: activation, inhibition, or both. Clin Neurophysiol. 2004;115:1239–1248. [PubMed]
7. Astrom M, Johansson JD, Hariz MI, Eriksson O, Wardell K. The effect of cystic cavities on deep brain stimulation in the basal ganglia: a simulation-based study. J Neural Eng. 2006;3:132–138. [PubMed]
8. Arle JE, Mei LZ, Shils JL. Modeling parkinsonian circuitry and the DBS electrode. I. Biophysical background and software. Stereotact Funct Neurosurg. 2008;86:1–15. [PubMed]
9. Kuncel AM, Cooper SE, Grill WM. A method to estimate the spatial extent of activation in thalamic deep brain stimulation. Clin Neurophysiol. 2008;119:2148–2158. [PMC free article] [PubMed]
10. Kuncel AM, Grill WM. Selection of stimulus parameters for deep brain stimulation. Clin Neurophysiol. 2004;115:2431–2441. [PubMed]
11. Astrom M, Zrinzo LU, Tisch S, et al. Method for patient-specific finite element modeling and simulation of deep brain stimulation. Med Biol Eng Comput. 2009;47:21–28. [PubMed]
12. Butson CR, Cooper SE, Henderson JM, McIntyre CC. Patient-specific analysis of the volume of tissue activated during deep brain stimulation. Neuroimage. 2007;34:661–670. [PMC free article] [PubMed]
13. Butson CR, Maks CB, McIntyre CC. Sources and effects of electrode impedance during deep brain stimulation. Clin Neurophysiol. 2006;117:447–454. [PubMed]
14. Butson CR, McIntyre CC. Tissue and electrode capacitance reduce neural activation volumes during deep brain stimulation. Clin Neurophysiol. 2005;116:2490–2500. [PubMed]
15. Hemm S, Mennessier G, Vayssiere N, et al. Deep brain stimulation in movement disorders: stereotactic coregistration of two-dimensional electrical field modeling and magnetic resonance imaging. J Neurosurg. 2005;103:949–955. [PubMed]
16. Johnson MD, McIntyre CC. Quantifying the neural elements activated and inhibited by globus pallidus deep brain stimulation. J Neurophysiol. 2008;100:2549–2563. [PubMed]
17. McIntyre CC, Grill WM, Sherman DL, Thakor NV. Cellular effects of deep brain stimulation: model-based analysis of activation and inhibition. J Neurophysiol. 2004;91:1457–1469. [PubMed]
18. McIntyre CC, Mori S, Sherman DL, Thakor NV, Vitek JL. Electric field and stimulating influence generated by deep brain stimulation of the subthalamic nucleus. Clin Neurophysiol. 2004;115:589–595. [PubMed]
19. Miocinovic S, Lempka SF, Russo GS, et al. Experimental and theoretical characterization of the voltage distribution generated by deep brain stimulation. Exp Neurol. 2009;216:166–176. [PMC free article] [PubMed]
20. Miocinovic S, Parent M, Butson CR, et al. Computational analysis of subthalamic nucleus and lenticular fasciculus activation during therapeutic deep brain stimulation. J Neurophysiol. 2006;96:1569–1580. [PubMed]
21. Sotiropoulos SN, Steinmetz PN. Assessing the direct effects of deep brain stimulation using embedded axon models. J Neural Eng. 2007;4:107–119. [PubMed]
22. Vasques X, Cif L, Hess O, et al. Stereotactic model of the electrical distribution within the internal globus pallidus during deep brain stimulation. J Comput Neurosci. 2009;26:109–118. [PubMed]
23. Walckiers G, Fuchs B, Thiran JP, Mosig JR, Pollo C. Influence of the implanted pulse generator as reference electrode in finite element model of monopolar deep brain stimulation. J Neurosci Methods. 2009 [PubMed]
24. Wei XF, Grill WM. Current density distributions, field distributions and impedance analysis of segmented deep brain stimulation electrodes. J Neural Eng. 2005;2:139–147. [PubMed]
25. Yousif N, Bayford R, Wang S, Liu X. Quantifying the effects of the electrode-brain interface on the crossing electric currents in deep brain recording and stimulation. Neuroscience. 2008;152:683–691. [PMC free article] [PubMed]
26. Yousif N, Bayford R, Liu X. The influence of reactivity of the electrode-brain interface on the crossing electric current in therapeutic deep brain stimulation. Neuroscience. 2008;156:597–606. [PMC free article] [PubMed]
27. Yousif N, Liu X. Investigating the depth electrode-brain interface in deep brain stimulation using finite element models with graded complexity in structure and solution. J Neurosci Methods. 2009;184:142–151. [PMC free article] [PubMed]
28. McIntyre CC, Miocinovic S, Butson CR. Computational analysis of deep brain stimulation. Expert Rev Med Devices. 2007;4:615–622. [PubMed]
29. Mayer S, Geddes LA, Bourland JD, Ogborn L. Faradic resistance of the electrode/electrolyte interface. Med Biol Eng Comput. 1992;30:538–542. [PubMed]
30. Ashby P, Kim YJ, Kumar R, Lang AE, Lozano AM. Neurophysiological effects of stimulation through electrodes in the human subthalamic nucleus. Brain. 1999;122(Pt 10):1919–1931. [PubMed]
31. Baker KB, Montgomery EB, Jr, Rezai AR, Burgess R, Luders HO. Subthalamic nucleus deep brain stimulus evoked potentials: physiological and therapeutic implications. Mov Disord. 2002;17:969–983. [PubMed]
32. Costa J, Valls-Sole J, Valldeoriola F, Rumia J, Tolosa E. Motor responses of muscles supplied by cranial nerves to subthalamic nucleus deep brain stimuli. Brain. 2007;130:245–255. [PubMed]
33. MacKinnon CD, Webb RM, Silberstein P, et al. Stimulation through electrodes implanted near the subthalamic nucleus activates projections to motor areas of cerebral cortex in patients with Parkinson's disease. Eur J Neurosci. 2005;21:1394–1402. [PubMed]
34. Tommasi G, Krack P, Fraix V, et al. Pyramidal tract side effects induced by deep brain stimulation of the subthalamic nucleus. J Neurol Neurosurg Psychiatry. 2008;79:813–819. [PubMed]
35. Chaturvedi A, Butson CR, Cooper SE, McIntyre CC. Subthalamic nucleus deep brain stimulation: accurate axonal threshold prediction with diffusion tensor based electric field models. Conf Proc IEEE Eng Med Biol Soc. 2006;1:1240–1243. [PubMed]
36. Wakana S, Jiang H, Nagae-Poetscher LM, van Zijl PC, Mori S. Fiber tract-based atlas of human white matter anatomy. Radiology. 2004;230:77–87. [PubMed]
37. Baker KB, Tkach JA, Phillips MD, Rezai AR. Variability in RF-induced heating of a deep brain stimulation implant across MR systems. J Magn Reson Imaging. 2006;24:1236–1242. [PubMed]
38. Viola P, Wells WM., III Alignment by maximization of mutual information; Proceedings of the Fifth International Conference on Computer Vision: IEEE Computer Society; 1995.
39. Christensen GE, Joshi SC, Miller MI. Volumetric transformation of brain anatomy. IEEE Trans Med Imaging. 1997;16:864–877. [PubMed]
40. Weinstein D, Kindlmann G, Lundberg E. Tensorlines: advection-diffusion based propagation through diffusion tensor fields. Visualization '99: IEEE. 1999:249–530.
41. McIntyre CC, Richardson AG, Grill WM. Modeling the excitability of mammalian nerve fibers: influence of afterpotentials on the recovery cycle. J Neurophysiol. 2002;87:995–1006. [PubMed]
42. Hines ML, Carnevale NT. The NEURON simulation environment. Neural Comput. 1997;9:1179–1209. [PubMed]
43. Butson CR, McIntyre CC. Differences among implanted pulse generator waveforms cause variations in the neural response to deep brain stimulation. Clin Neurophysiol. 2007;118:1889–1894. [PMC free article] [PubMed]
44. Tuch DS, Wedeen VJ, Dale AM, George JS, Belliveau JW. Conductivity tensor mapping of the human brain using diffusion tensor MRI. Proc Natl Acad Sci U S A. 2001;98:11697–11701. [PubMed]
45. McNeal DR. Analysis of a model for excitation of myelinated nerve. IEEE Trans Biomed Eng. 1976;23:329–337. [PubMed]
46. Rosen I, Asanuma H. Peripheral afferent inputs to the forelimb area of the monkey motor cortex: input-output relations. Exp Brain Res. 1972;14:257–273. [PubMed]
47. Haberler C, Alesch F, Mazal PR, et al. No tissue damage by chronic deep brain stimulation in Parkinson's disease. Ann Neurol. 2000;48:372–376. [PubMed]
48. Yelnik J, Damier P, Demeret S, et al. Localization of stimulating electrodes in patients with Parkinson disease by using a three-dimensional atlas-magnetic resonance imaging coregistration method. J Neurosurg. 2003;99:89–99. [PubMed]
49. Halpern CH, Danish SF, Baltuch GH, Jaggi JL. Brain shift during deep brain stimulation surgery for Parkinson's disease. Stereotact Funct Neurosurg. 2008;86:37–43. [PubMed]
50. Khan MF, Mewes K, Gross RE, Skrinjar O. Assessment of brain shift related to deep brain stimulation surgery. Stereotact Funct Neurosurg. 2008;86:44–53. [PubMed]
51. Mori S, Kaufmann WE, Davatzikos C, et al. Imaging cortical association tracts in the human brain using diffusion-tensor-based axonal tracking. Magn Reson Med. 2002;47:215–223. [PubMed]
52. Gradinaru V, Mogri M, Thompson KR, Henderson JM, Deisseroth K. Optical deconstruction of parkinsonian neural circuitry. Science. 2009;324:354–359. [PubMed]
53. Lee KH, Blaha CD, Harris BT, et al. Dopamine efflux in the rat striatum evoked by electrical stimulation of the subthalamic nucleus: potential mechanism of action in Parkinson's disease. Eur J Neurosci. 2006;23:1005–1014. [PubMed]
54. McIntyre CC, Grill WM. Excitation of central nervous system neurons by nonuniform electric fields. Biophys J. 1999;76:878–888. [PubMed]
55. Rattay F. Analysis of models for external stimulation of axons. IEEE Trans Biomed Eng. 1986;33:974–977. [PubMed]
56. Graf von Keyserlingk D, Schramm U. Diameter of axons and thickness of myelin sheaths of the pyramidal tract fibres in the adult human medullary pyramid. Anat Anz. 1984;157:97–111. [PubMed]
57. Kamada K, Todo T, Masutani Y, et al. Combined use of tractography-integrated functional neuronavigation and direct fiber stimulation. J Neurosurg. 2005;102:664–672. [PubMed]
58. Lujan JL, Noecker AM, Butson CR, et al. Automated 3-dimensional brain atlas fitting to microelectrode recordings from deep brain stimulation surgeries. Stereotact Funct Neurosurg. 2009;87:229–240. [PMC free article] [PubMed]
59. Frankemolle AMM, Wu J, Noecker AM, et al. Reversing cognitive-motor impairments in Parkinson’s disease patients using a computational modelling approach to deep brain stimulation programming. Brain. 2010 (in press). [PMC free article] [PubMed]
60. Maks CB, Butson CR, Walter BL, Vitek JL, McIntyre CC. Deep brain stimulation activation volumes and their association with neurophysiological mapping and therapeutic outcomes. J Neurol Neurosurg Psychiatry. 2009;80:659–666. [PMC free article] [PubMed]