Deep brain stimulation modeling for several anatomical and electrical considerations

Modelos de estimulación cerebral profunda para diferentes consideraciones anatómicas y eléctricas

Cristian Alejandro Torres-Valencia¹, Genaro Daza-Santacoloma², Mauricio Alexander Álvarez-López³, Álvaro Ángel Orozco-Gutiérrez⁴

1 M.Sc. Universidad Tecnológica de Pereira. Pereira, Colombia. cristian.torres@utp.edu.co.
2 Ph.D Instituto de Epilepsia y Parkinson del Eje Cafetero - NEUROCENTRO. Pereira, Colombia. research@neurocentro.com.co.
3 Ph.D Universidad Tecnológica de Pereira. Pereira, Colombia. malvarez@utp.edu.co.
4 Ph.D Universidad Tecnológica de Pereira. Pereira, Colombia. aaog@utp.edu.co.

ABSTRACT

Deep Brain Stimulation (DBS) is a clinical treatment for Parkinson disease symptoms. DBS consists in the implantation of a stimulation electrode into the Subthalamic nucleus (STN) for the excitation of specific regions inside the STN. The stimulation potential has a few parameters that should be adjusted in order to achieve the desired treatment effect. The adjust is performed by the neurologist in several sessions with the patients and is not an exact procedure. In recent years there have been several works on the construction of propagation models of DBS, including head geometries and medium properties in order to visualize the possible effects of DBS while the stimulation parameters are adjusted. This work presents the construction of propagation models using the Finite Element Method (FEM) for the solution of Laplace or Poisson equations that govern the propagation phenomena. By the construction of these models, the shape and magnitude of the electric propagation inside the objective structures can be obtained.

KEYWORDS: Deep Brain Stimulation, electric propagation, finite element method, Laplace equation, Parkinson disease, Poisson equation, subthalamic nucleus.

RESUMEN

La estimulación cerebral profunda (DBS) es una terapia quirúrgica validada para el tratamiento de los síntomas asociados con la enfermedad de Parkinson. Consiste en la implantación de un electrodo de estimulación generalmente en la región del núcleo subtalámico (STN), con el cual se excitan regiones específicas a partir de un potencial eléctrico con ciertos parámetros específicos. El ajuste de los parámetros de estimulación es un proceso realizado por parte del neurólogo y puede tardar varios meses hasta alcanzar los resultados deseados. Es por esto que en años recientes se ha estudiado la construcción de modelos de propagación eléctrica de las estructuras objetivo de la DBS con el fin de visualizar los posibles resultados de la distribución de campo eléctrico y la activación del tejido cerebral que sirven como guía para el ajuste de los parámetros de estimulación, optimizando el procedimiento de configuración. En este trabajo se presenta la comparación de modelos de simulación que incluyen la definición de geometrías complejas representando diferentes estructuras cerebrales con propiedades de diferentes tejidos, con los cuales se obtienen los patrones de propagación eléctrica cerebral por medio del método de elementos finitos (FEM) aplicado a la solución de las ecuaciones de Laplace y Poisson.

PALABRAS CLAVE: Ecuación de Laplace, ecuación de Poisson, enfermedad de Parkinson, estimulación cerebral profunda, método de elementos finitos, núcleo subtalámico, propagación eléctrica.

1. INTRODUCTION

Deep Brain Stimulation (DBS) was introduced in 1995 as a clinical treatment for Parkinson disease symptoms [1]. DBS consists on the placement of an electrode into a specific brain structure (e.g. subthalamic nucleus - STN) in order to deliver a continuous electrical stimulation [2]. The effects of DBS were subject of investigation in several works in order to assess the real improvement and side effects for Parkinson patient treatment [2][3]. From the analysis developed over DBS in several works, it was concluded that benefits could cover brain diseases as psychiatric conditions, refractory epilepsy and dystonia [3].

Since the DBS effects are directly related to the electric stimulation of a specific brain structure, the stimulation pulse parameters have to be configured by the neurologist in several sessions with the patients [4]. The magnitude, frequency, pulse width and the stimulation mode (monopolar or bipolar) are the parameters that define the shape and behavior of the potential into the STN [5]. The electric propagation into non desired areas in or outside the STN could derive in adverse effects for the patient [2], [6].

In order to achieve an accurate description of the propagation of the electric potential into the STN, simulated models of the STN, brain or the entire head has been developed [7]. These simulation models allow the configuration of the stimulation parameters, by defining the medium properties and the visualization of the propagated potential in a defined geometry. Several models are proposed in the literature that include different brain structures and idealized 2D and 3D head geometries including medium properties [7], [8], [9]. The geometry of the DBS lead is also included in the models for monopolar or bipolar stimulation modes using the active contacts of the lead. Even with the inclusion of a realistic head model, different tissues, structures and medium properties, the modeling of the DBS problem is still an open research field.

After the model is built, a definition of the boundary conditions and the stimulation parameters configuration is performed [7]. Since the electric propagation in any medium is governed by the Maxwell equations [10], Laplace and Poisson equations are used to compute the solution in the DBS simulation models [7][4]. Finding the electric potential propagated into the model is not a straightforward process: an analytic solution is hard to compute due the geometry and the different mediums considered for each model. Numerical approximations for the solution of the propagation equation for the DBS model should be implemented [4]. The Finite Differences Method (FDM) and the Finite Element Method (FEM) are numerical methods used in DBS simulations. These methods are based on the definition of a mesh in the geometry of the model, boundary and initial conditions, the value of the variable to be computed is depicted from each node in the mesh following arithmetic operations [11], [12]. FEM is preferred over FDM since the mesh in FDM must be uniform, giving less resolution in areas where there are smaller structures inside the global geometry [11].

The aim of this work is to compare several DBS simulation models with different medium properties, including brain structures and head layers, DBS lead geometry, Dirichlet and Newmman boundary conditions related to the solution of Laplace and Poisson equations. These models allows a comparative analysis of the differences between the electric potential propagation by the inclusion of different features in the model. A realistic head model including the DBS generator is assumed as the more accurate model for the neurologist to assess the stimulation parameters that allows the effective application of DBS therapy.

2. PREVIOUS WORK

The study of the effects of DBS has become an important research field in recent years. Several works have developed studies of the effects of DBS and simulation models in order to obtain an approximation of the electric propagation inside the STN. In 1995, Patricia Limousin et al. published a study of three patients that were operated in order to treat the akinetic-rigid Parkinson's disease with STN-DBS. A posterior evaluation of the symptoms of these patients shows an improving around 42% to 84% following the unified Parkinson's disease rating scale [1]. In 2003, Benabid et al. presented a work about the clinical results of DBS when the frequency parameter changes in the adjustment of the stimulator [2]. From the study of the effects, it was concluded that the benefits derived from DBS cover different motor symptoms and the quality of life [2].In another study, a simulation of three layers of neurons was developed by Tarparelli et al. in order to analyze the activation of the surrounding tissue due DBS stimulation [3]. The simulation process was developed using a commercial software for neurologic modelling and the results allow an initial calibration of the parameters configuration of the DBS pulse [3].

From the analysis developed over the effects of DBS, several works have been proposed in order to achieve a specific description of the electric potential propagated into the STN and adjacent structures. McIntyre et al. in [5] built a simulation model of the neurons near to the electrode in order to examine the influence of the stimulus waveform and frequency. In an additional study [7], McIntyre developed a study of the axonal tissue activated by DBS in the STN. 3D cubic geometries are built following tissue considerations for a homogeneous and an inhomogeneous medium of propagation [7]. The Laplace equation is used for finding the potential propagation into the cubic model with conductivity assumptions for the medium. Dirichlet boundary conditions are defined for all the surface of the cube and a commercial software is used for the implementation of the Finite Element Method that allows the numerical computing of the solution of the propagation model [7]. Following different configurations for the stimulation electrode into the models, the results shows that electric propagation is influenced by the medium properties [7]. In a similar study, Butson in 2007, studied the variations in the electric potential propagation in a simulation model of the STN and the DBS lead, incorporating the conductivity information from Magnetic Resonance Imaging [13]. The results were computed using FEM and the volume of tissue activated (VTA) was used for measuring the effects of a specific configuration of the stimulation parameters [13].

Liberti presented a different approach for the simulation modelling in 2007 [9]. Several 2D models were built in order to test the influence of the size of the model and the boundary conditions. Box models of 5x5, 15x15, 50x50 and 100x100 mm with the definition of the whole boundary or the base side of the box as the ground were developed [9]. The results show that dimension and ground definition of the model influences the electric potential propagation shape and magnitude and the activating function [9]. A more realistic approach was developed by Grant et al. in 2009 [8]. Indeed, four different 3D-spherical models were built, including the geometry for the entire head with definitions of different structures and tissues [8]. Each subdomain of the model was configured with the conductivity value associated to a specific head structure. The DBS lead geometry was also included, the electric potential propagation was computed using FEM. Isotropic, and anisotropic mediums were considered in the models. For three out of the four models the ground definition was settled on a cylindrical extension of the sphere, designated as the reference surface. The results show differences between the electric potential computed using FEM due the inclusion of different brain structures and head layers as well as the ground placement [8].

In 2009, Walckiers presented his doctoral thesis with a detailed analysis of all the stages of DBS procedure and modelling. The definition of a complex geometry that represents the head, neck and the top of the chest were the DBS generator is implanted was proposed. Several Dirichlet and Newmman boundary conditions were defined to bring an accurate description of the propagation problem into the head [11]. The results from this study show differences between models with boundary conditions that affect the impedance of the model. It was concluded that the complete closed loop of the DBS excitation should be considered into monopolar and bipolar simulations [4].

3. MATERIAL AND METHODS

In this section, the methods used in the development of this work are depicted. First, the DBS modelling problem is well explained and the laws that govern the propagation phenomena inside the brain are introduced. A brief introduction to the Finite Element Method (FEM) is included in this section and the implementation details of the models are included as well.

3.1. Deep Brain Stimulation Modeling

From the knowledge of the DBS therapy and the procedure for configuring the stimulation parameters, several models have been proposed in the state of the art, as it was presented in section 2. The fundamentals of the DBS are based on the construction of a geometrical representation of the brain structures involved or geometrical representations of the entire head. Initially the models were developed only for smaller geometries that could represent the subthalamic nucleus (STN), that is the objective structure in DBS [6]. In recent years, the inclusion of the different brain tissues characteristics and head structures as the scalp and the skull has been considered in the model building. The inclusion of these additional components will give more realistic considerations into the propagation problem to be solved [8], [10].

From the geometric model and the structures considered in each model, other assumptions related to the conductivity of the tissue and the dimensions should be taken account. The equation that governs the electric phenomena inside the geometry is also another parameter that is related to the solution of the model. Laplace and Poisson equations are the two mathematical models derived from the Maxwell equations used in DBS modeling. It can be demonstrated from the analysis of these equations that with a quasistatic approximation, the conductivity values of the tissues when the Laplace equation is used is not influencing the solution. The mathematic formulation for the Laplace equation is presented in (1).

where σ, is the medium conductivity. Alternatively, the Poisson equation allows the definition of conductivity values of each tissue due the function term that is included in the mathematical model, as (2) shows.

3.2. Finite Element Method (FEM)

The Finite Element Method FEM is a numerical method for approximating the solution of differential equations associated to a problem with complexes geometries i.e. DBS modelling [11]. FEM allows the computing of the approximate solution of a problem over a domain, which is govern by differential equations. The computation is performed by first dividing the original domain into a high number of subdomains known as finite elements. Inside each finite element can be distinguished a set of representative points called nodes [12], [14]. Two nodes are adjacent if they are contained by one finite element and one node on the boundary of an element can belongs to several elements. The compilation of nodes with all the defined relationships between them is called mesh [12]. The FEM method has two basic divisions, the first one uses finite elements to obtain the joint displacements and member forces of a structural framework or mesh [12]. From the adjacency conditions of the mesh, for each node can be proposed a set of variables denominated freedom degrees. The relationship between the values of each variable in the nodes can be expressed as a linear system of equations. The number of equations is proportional to the nodes quantity [12], [14]. Then FEM combines several mathematical concepts to produce a large system of linear or nonlinear equations. Numerical methods should be used in order to solve these type of problems [12].

Fig. 1 shows an example of a 2D domain with a field variable, Ø (x,y) that has to be represented in each point P(x,y) [12]. The development of an algorithm of finite elements could be expressed in four stages:

• The problem must be formulated as a variational problem.
• The domain of independent variables is divided into finite elements. With all the finite elements a vectorial space is built, being the solution a lineal combination of elements in such vectorial space.
• The variational problem is projected over the finite element space. This derive into a finite system of equations. The number of unknowns is equal to the vectorial space dimension and as high is the dimension of the space as accurate will be the numerical approximation.
• Finally the equation system is solved.

Due the number of different subdomains that can be defined into the DBS modelling, as well as the geometry of the DBS lead and different head structures, FEM has proved to be an effective tool for computing the electric potential.

3.3. Implementation Details

For model building, a spherical geometry is selected in order to reproduce a head like model. The dimensions of the sphere are a radius of 150 mm, since it reproduces the average head size [15], [16]. The DBS lead geometry is included in the models, based in the Medtronic 3387 configuration of electrodes and dimensions. The geometry of the model is presented in Fig. 2.

The solid blue line observed in Fig. 2 correspond to the geometry of the DBS lead. An addition to the model geometry presented in Fig 2, is the inclusion of a cylinder in the inferior portion of the ellipse in order to simulate the neck of a human being and to place a different ground condition that emulates the DBS generator. This additional feature is presented in Fig. 3. The active contact of the electrode is settled at 1V for simulation pruposes.

The definitions of the subdomains inside the model, correspond to different brain structures and other head components that will bring a major realism to the simulation process. From the paper presented by Grant in 2009 [8], the definition of the subdomains include the tissues and structures documented in Table 1.

Several model are proposed by the inclusion of the different head layers described in Table 1 for comparison purposes. For each included layer, the solution is computed using FEM and stored in an xml file. The solution of the model is computed using the Laplace formulation following Dirichlet boundary conditions for the ground and the stimulation potential. When the solution is computed using the Poisson equation, the source is included as a constant function in the region of the active contact of the electrode.

When the mesh of the geometry is generated, a previous refinement stage before computation is performed. The selection of the nodes to be refined corresponds to the regions that are closer to the electrode in which a higher resolution is needed. The mesh is stored for further computation since the general geometry of the model does not change.

4. RESULTS

This section presents the results obtained from the solution of the proposed geometry when DBS is performed using FEM for computing the electric potential. Since the results are hard to compare visually, a brief description of the potential curves is presented first.

The presented Fig. 4 corresponds to the solution of the DBS model using FEM and the Laplace equation, including different structures. Fig 4a corresponds to the model with the definition of the domain as the cerebrospinal and brain tissue with the active contact of the electrode. As can be seen, the potential propagation has a uniform shape from the center of the electrode to the exterior of the sphere. Fig. 4b, 4c and 4d, show the solution for the electric potential of the models including the rest of the DBS lead and the encapsulation tissue (Fig. 4b), the inner cortical skull, cancellous bone and outer cortical skull (Fig. 4c), the fat and the skin (Fig. 4d). For the last models, it can be seen that the potential propagation has a non uniform shape from the definition of the different subdomains with conductivity values. The propagation seems to reach further regions of the model for the case with no definition of the additional head simulation layers.

In Fig. 5 some of the results obtained for the same models than in Fig. 4 are presented. For this solution of the DBS problem, FEM is computed using the Poisson equation. Fig. 5a. presents the case of the model with only the active contact of the electrode. Fig 5b. presents the complete model with the definition of the other brain structures and head layers.

It can be seen from Fig 5 that a similar results are obtained compared against the results from the Laplace equation. From each potential level on the figures presented before, the minimum distance from the active contact center of the electrode is computed. These minimum distances are presented in the following figures for the Laplace and Poisson equations.

In Fig. 6 can be observed that higher potential levels were reached in further regions of the model when only the subdomain for the active contact of the lead was modelled. For the other three models the minimum distances are comparable between each other, with an initial potential around 0.9V in regions closer than 10mm and a final electric potential around 0.05V in regions up to 40mm.

From Fig. 7, the same behavior presented for the models when Laplace was used is presented for the Poisson equation. In this case the difference between the models that include more subdomains against the model with only the active contact of the electrode can be noticed. Although a significant difference is observed for further distances from the electrode, in closer regions the potential propagation can be comparable. For a detailed analysis, the percentage difference between the different models is presented in Table II. This percentage difference corresponds to the computation of the relative percentage error between two measures of the same quantity.

From the data presented in Table II, it can be seen the difference of the potential propagation due the inclusion of the properties of the medium in the defined subdomains. For the Laplace equation, the results from the different models are similar since the conductivity has no influence mathematically when the Laplace equation is solved. More evident differences are noticed when the Poisson equation is solved for the models with the inclusion of the conductivity for each subdomain.

Table II

5. CONCLUSIONS

Since the results presented in section 3 clearly demonstrated the difference in the potential propagation between the models with different tissues and structures configuration, it can be concluded that a more realistic head model will derive in simulations that can be used effectively by the neurologist in the DBS parameter configuration. The comparative analysis of the potential curves shows that the boundary conditions definition and subdomains properties plays an important role in the simulation and have to be considered when DBS modeling is developed, as Fig 6, Fig 7 and Table II shows.

The results from the Laplace and Poisson equations were comparable and it can be concluded that those mathematical models has a reliable application to the DBS problem. The conductivity values from the different subdomains have a clear influence in the shape of the propagated potential when both Laplace and Poisson equations are used.

Future works could be focused on the inclusion of the anisotropic configuration of the brain tissue. Information from Magnetic Resonance Imaging (MRI) could derive in the computation of the tensor information that could be included in the model.

AGRADECIMIENTOS

The authors would like to thank the Universidad Tecnológica de Pereira for providing the resources needed to develop this work. This work was developed under the project 111056934461 funded by COLCIENCIAS. Author GDS was partially supported by "Patrimonio Autónomo Fondo Nacional de Financiamiento para la Ciencia, la Tecnología y la Innovación, Francisco José de Caldas", project 499153530997.

REFERENCIAS

[1] P. Limousin, P. Pollak, A. Benazzouz, D. Hoffmann, J. F. Le Bas, E. Broussolle, J. E. Perret, and A. L. Benabid, "Effect of parkinsonian signs and symptoms of bilateral subthalamic nucleus stimulation." Lancet, vol. 345, no. 8942, pp. 91-95, 1995.         [ Links ]

[2] Benabid, A. L. Deep brain stimulation for Parkinson's disease. Current Opinion in Neurobiology, no, 13, pp. 696-706, 2003.         [ Links ]

[3] Tarparelli, R.; Iovine, R.; La Spada, L.; Vegni, L., "Electromagnetic analysis of deep brain stimulation," Electromagnetics in Advanced Applications (ICEAA), International Conference, pp. 1153-1155, 9-13, 2013.         [ Links ]

[4] G. Walckiers, B. Fuchs, J.-P. Thiran, J. R. Mosig, and C. Pollo, "Influence of the implanted pulse generator as reference electrode in finite element model of monopolar deep brain stimulation," Journal of Neuroscience Methods, vol. 186, pp. 90-96, 2010.         [ Links ]

[5] C. C. McIntyre and W. M. Grill, "Extracellular stimulation of central neurons: Influence of stimulus waveform and frequency on neuronal output," Journal of Neurophysiology, vol. 88, 2002, pp. 1592-1604.         [ Links ]

[6] A.L. Benabid, S. Chabardes, J. Mitrofanis, P. Pollak, "Deep brain stimulation of the subthalamic nucleus for the treatment of Parkinson's disease", Lancet Neurol, vol. 8, no. 1, 2009, pp. 67-81.         [ Links ]

[7] C. C. McIntyre, S. Morib, D. L. Shermanc, N. V. Thakorc, and J. L. Vitek, "Electric field and stimulating influence generated by deep brain stimulation of the subthalamic nucleus," Clinical neurophysiology, vol. 115, pp. 589-595, 2004.         [ Links ]

[8] P. F. Grant and M. M. Lowery, "Electric field distribution in a finite-volume head model of deep brain stimulation," Medical engineering and physics, vol. 31, pp. 1095-1103, 2009.         [ Links ]

[9] M. Liberti, F. Apollonio, A. Paffi, M. Parazzini, F. Maggio, T. Novellino, P. Ravazzani, and G. D'Inzeo, "Fundamental electrical quantities in deep brain stimulation: Influence of domain dimensions and boundary conditions." Lyon, France: 29th annual international conference of the IEEE EMBS, pp. 6668-6671, 2007.         [ Links ]

[10] M. N. O. Sadiku, "Elements of electromagnetics." USA: Oxford University Press, 2002.         [ Links ]

[11] G. Walckiers, "Bio-electromagnetic model of deep brain stimulation." IEL, 2009.         [ Links ]

[12] D. Hutton. "Fundamentals of Finite Element Analysis," United States: McGraw-Hill, 2004.         [ Links ]

[13] C. Butson, S. E. Cooper, J. M. Henderson and C. C. McIntyre, "Patient-specific analysis of the volume of tissue activated during deep brain stimulation", NeuroImage, 2007.         [ Links ]

[14] L.J. Segerind, "Applied Finite Element Analysis", New York: John Wiley and Sons, 1984.         [ Links ]

[15] J. W. Kim and P. A. Robinson, "Compact dynamical model of brain activity", American Physical Society, Phys. Rev. E. 031907, Vol. 75, i. 3, p. 031907-031917, 2007.         [ Links ]

[16] A. Chaturvedi, C. R. Butson, F. Scott, F. Lempka, E. Cooper, and C. Cameron, "Patient-specific models of deep brain stimulation: Influence of field model complexity on neural activation predictions", Brain Stimulation, vol. 3, no. 2, pp. 65-77, 2010.         [ Links ]