Effects of head models and dipole source parameters on EEG fields.

Head model and an efficient method for computing the forward EEG (electroencephalography)problem are essential to dipole source localization(DSL). In this paper, we use less expensive ovoid geometry to approximate human head, aiming at investigating the effects of head shape and dipole source parameters on EEG fields. The application of point least squares (PLS) based on meshless method was introduced for solving EEG forward problem and numerical simulation is implemented in three kinds of ovoid head models. We present the performances of the surface potential in the face of varying dipole source parameters in detail. The results show that the potential patterns are similar for different dipole position in different head shapes, but the peak value of potential is significantly influenced by the head shape. Dipole position induces a great effect on the peak value of potential and shift of peak potential. The degree of variation between sphere head model and non-sphere head models is seen at the same time. We also show that PLS method with the trigonometric basis is superior to the constant basis, linear basis, and quadratic basis functions in accuracy and efficiency.

INTRODUCTION

Electroencephalogrphy (EEG) is a non-invasive method of measuring the electrical activity of the brain. Dipole source localization (DSL) has been widely used in analyzing electroencephalograms (EEG’s), which is a valuable tool in the pre-surgical evaluation of patients suffering from epilepsy [1] and basic brain research. In order to solve this inverse problem, one must solve the EEG forward problem to determine the surface potentials from the dipole current sources in the head volume. For a quantitative interpretation of electroencephalograms(EEG’s), a mathematical geometry is needed to approximate the human head. Spherical head models[2-4] have been commonly used in EEG problems. However, some key drawbacks are shown using spherical head models in [5, 6]. Nonspherical head models, such as ellipsoid models [6-9] and realistic head shape [10-13] based on segmented 3D magnetic resonance images are also applied. These realistic head models are more adequate to describe the human head than spherical models. However, excessive computing time is required. Therefore, such models are not convenient for the theoretical analysis.

For accurate localization of dipole current source, the information about the effects of the head shape and dipole source parameters on EEG field is necessary. The related work [6, 13-15] developed some information from different aspects when using spherical or realistically shaped models. Boundary element methods (BEM’s) have been widely studied in EEG forward problem [16-20]. The performance of BEM methods is dependent on the selection of the basis and weighting functions. The comparison of accuracy was made between the BEM methods with constant, linear and second order basis function for the potential using spherical models in [16,17]. Currently, the efficiency of EEG forward solutions based on BEM for non-spherical head models is limited [5, 21]. In this paper, we use three analytical models (Figs. (-)), lengthwise ovoid, oblate ovoid model, and sphere to approximate human head aiming at investigating the effects of head shape and dipole source parameters on EEG fields. Point least squares(PLS) based on meshless method for solving EEG forward problem is developed and numerically implemented. Four kinds of basis functions; constant basis, linear basis, quadratic basis and trigonometric basis functions with compact support are applied. We include numerical comparisons of four kinds of basis functions in accuracy and efficiency. We organize the paper as follows. Section 2 provides the EEG field equations and our head models. In Section 3, we derive formula of Point Least Square (PLS) for computing the scalp electric potential. The description of numerical experiments is presented in section 4. In section 5, simulation results and discussion are presented. In section 6, conclusion is given.

EEG FIELD EQUATION AND HEAD MODELS

The typical head model assumes that it is made up of a set of compartments with homogeneous and isotropic conductivity. Let  be the space region occupied by the head and  be the boundary of  Assume that  is divided by surfaces  (from the beginning of the scalp) into some subregions  with constant isotropic conductivity  in each .  and  denote the conductivities inside  and outside , respectively.  and  are the coordinate vectors in .

The physiology basis for EEG is that the electrical activity of active nerve cells in the brain produces small currents which reach the scalp surface and result into voltage difference. EEG forward problem aims at computing the scalp potential that is produced by assumed primary current in the head volume. The relationship between the scalp potential and dipole source can be derived from quasistatic Maxwell’s equations [22], that is,

(1)

where ,  is outward directed unit vector normal to the jth surface  and  is the potential in an infinite homogeneous medium with conductivity , that is,

(2)

The primary current is generally modeled as a current dipole with moment  located at , then we have

(3)

The head model we will use is a new kind of analytical geometry, ovoid geometry. The surface equation of ovoid geometry can be expressed as

(4)

where  is the length of longitudinal axis, and  is called the index of ovoid. Note that

1.   When ,  is lengthwise ovoid (Fig. ()).

2.   When ,  is oblate ovoid (Fig. ()).

3.   When ,  is a sphere with the radius  and center on z-axis (see Fig. ()).

4. The oval cross-section with different index is shown in Fig. ().

Three-layer models can be seen in Fig. (), representing scalp, skull and brain.

POINT LEAST SQUARES BASED ON MESHLESS METHOD WITH COMPACTLY SUPPORTED FUNCTIONS

The techniques widely used to solve the EEG forward problem are the finite element method (FEM) and boundary element method (BEM). FEM for the anisotropic conducting compartments have an increased computational demand [23, 24]. In contrast, BEM for isotropic conducting compartments can reduce much work. But the drawback of BEM is still the computational cost and the singularity in the vicinity of a surface separating two compartments. The meshless method with compactly supported functions has been widely applied in computing mechanics and fluid dynamics communities [25, 26]. Meshless method we use, a finite point method, is to use Weighted Residual method with compactly supported trial functions, which is constructed by discrete points. According to the formula (1), we define the operate

(5)

Then the formula (1) is equivalent to:

Let the trial function of the potential be:

(6)

where  are basis functions, and  is a set of corresponding unknown coefficients.

We need to set the residual error:

Then the general form of the method of weighted residual can be described as follows:

(7)

where  is weighting function. Generally, the weighting function is taken as the following form:

where  are coefficients and  are basis functions. Thus, inserting each basis function into the formula (8), we have unified form for the method of weighted residuals as follows:

(8)

If we choose the weighting functions as:

then general least squares takes the form:

Point of least squares is that we replace the integral by sum as follows:

where  are node points,  is the number of basis functions.

Equivalently, we write the matrix form as:

(9)

where  and  are  vectors, and  and  are as follows:

(10)

where,

(11)

and

(12)

It is seen that  is a symmetric, positive definite and dense matrix.

Noting that the surface equation of the ovoid in spherical coordinate system is

and the parametric form of the surface is the following

(13)

we transform the domain of boundary surface integration in Eq.(11) into the rectangular region in spherical coordinate parameters , and , that is,

(14)

According to the node point we divide the surface  into subsurfaces . The corresponding domain is sub rectangular . Except for the collocation points which are used to construct the trial function, a number of auxiliary points are also adopted (Fig. ). Here we choose the constant basis, linear basis, quadratic basis and trigonometric basis functions with compact support as follows, respectively.

(15)

(16)

(17)

(18)

where  and  are the parameters of surface in spherical coordinates. According to the parametric form of surface Eq. (13), the collocation points are written as

(19)

Noting that the elements of matrix in Eq.(11) need to handle boundary integrals, which show a strong singularity of the order  because of the highly varying kernel  in Eq.(12), we apply boundary regularization method which means that the singular points are replaced by the corresponding points on imaginary boundary surface [27]. Therefore, the formula (11) can be rewritten as:

(20)

is the element of area ,  is the index of ovoid,  and  are defined as the form (13), and  are collocation points defined by Eq.(19)in which  is replaced by . Here  is called boundary regularization parameter.

DESCRIPTION OF THE NUMERICAL EXPERIMENTS

To demonstrate the performance of potentials in different head models versus the varying dipole source parameters, we perform simulations in three head models with indexes ,  and . In simulations, we not only use the collocation points but also adopt a number of auxiliary points to construct the trial function. See Fig. ().

The parameters are taken as the following: the conductivities: , , the length of longitude axis , and the boundary regularization parameter . The collocation points are distributed over the part of surface with . We compute potentials in the following dipole source parameters. Dipole direction: . The dipole  was located on  axis moving from  to .

The forward calculation for EEG is divided into two steps. In the first step, we formulate the matrix  in Eq. (11), which is independent of the parameters of dipole source. We need to compute boundary regularized double integrals in Eq. (20). We implement quadl in dbquadl based on Matlab which is a high order method using an adaptive Gauss/Lobatto quadrature rule. We solve linear systems by LU decomposition. In the second step, the potential in infinite homogeneous space related to the parameter of dipole source is computed. In order to test the accuracy of our method, we use the relative error as follows

(21)

where  is the value based on the approximation of higher order and  is the value based on the approximation of lower order.

RESULTS AND DISCUSSION

The simulations are conducted in three head shapes. Under the same dipole direction(), the maps of the potentials versus the sensor positions parameter () for different dipole position on z-axis are shown in Figs. (-). It is shown that their patterns are similar but there is a difference in magnitude. The peak potential value shifts from the right to the left as the dipole moves from the bottom toward the scalp in each head shape. The maps of potentials versus dipole position on Z axis for three head shapes are seen in Fig. (), which shows the difference of the magnitude of potentials between three head models in the case of the same dipole position and direction. We also understand the degree of variation between the different head shapes. The difference between the oblate ovoid shape and the sphere shape is less than that between the lengthwise ovoid shape and sphere model. In Fig. (), the peak potentials versus the dipole position in three head models are compared. We can find that the trends of peak potentials as dipole moving toward the surfaces are similar in different head models. The relative differences of potentials computed by four basis functions are presented in Fig. (). It is seen that the relative difference of potential based on trigonometric basis is much better than that based on other three kinds of basis functions. The comparison of the computing time based on four kinds of basis functions can be seen in Table . It is seen that trigonometric basis has the advantage over other three kinds of basis functions in computing time. This result is true of the three head models is true.

CONCLUSION

The analytic head model we introduced is less expensive and convenient to make theoretical analysis, which can be used to make further studies of other conductor properties of the head. Point of Least Squares (PLS) based on meshless method are applied. Simulations are implemented in four kinds of basis functions for three head models. We obtained the detail information about the effects of head models and dipole source parameters on EEG. The results indicate that potential patterns in different head shapes are similar, but there is a huge difference in the magnitude of potentials in different head shapes. The peak value of potential is significantly influenced by the head shape. Dipole location has great effect on the peak value of potential and shift of peak potential. We not only obtained the detailed information about the effect of head shapes and dipole source parameters on EEG fields, but also understand the degree of variation between sphere head model and non-sphere head models. This provides necessary information for accurate localization of dipole current source. In comparison to PLS based on four kinds of basis functions, constant basis, linear basis, quadratic basis and trigonometric basis functions, we find that the trigonometric basis is superior to the other three basis functions in accuracy and efficiency. It is a significant reference for further computation of EEG forward solutions with anisotropic conductivity head model.

Table 1.

Running time (s) for four kinds of basis functions.

Method	Const.	Linear	Quadratic	Trig.
n0 = 1/3	149	159	161	90
n0 = 1	148	158	160	89
n0 = 3	147	158	159	88

time: timing benchmarks were determined using 2.00 GHz Pentium Duo CPU PC with 2.00 GB RAM running double-precision MATLAB based programs.