Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
Trends in EEG source localization Zoltan J. Koles* Department of Biomedical Engineering, University of Alberta, Edmonton, Alberta, Canada
Abstract The concepts underlying the quantitative localization of the sources of the EEG inside the brain are reviewed along with the current and emerging approaches to the problem. The concepts mentioned include monopolar and dipolar source models and head models ranging from the spherical to the more realistic based on boundary and finite elements. The forward and inverse problems in electroencephalography are discussed, including the non-uniqueness of the inverse problem. The approaches to the solution of the inverse problem described include single and multiple time-slice localization, equivalent dipole localization and the weighted minimum norm. The multiple time-slice localization approach is highlighted as probably the best available at this time and is discussed in terms of the spatiotemporal model of the EEG. The effect of noise corruption, artifacts and the number of recording electrodes on the accuracy of source localization is also mentioned. It is suggested that the main appeal of the minimum norm is that it does not assume a model for the sources and provides an estimate of the current density everywhere in the three dimensional volume of the head. 1998 Elsevier Science Ireland Ltd. Keywords: Current dipole; Electroencephalography; Head models; Minimum norm; Source localization; Spatiotemporal modeling; The equivalent dipole
1. Introduction Interpretation of the clinical EEG almost always involves speculation as to the possible locations of the sources inside the brain that are responsible for the observed activity on the scalp. The earliest efforts to quantitatively locate the sources of the EEG occurred more than 40 years ago when researchers began to correlate the existing body of electrophysiologic knowledge about the brain to the basic physical principles governing volume currents in conductive media (Brazier, 1949; Shaw and Roth, 1955; Plonsey, 1969; Schneider, 1972; Henderson et al., 1975; Nunez, 1981). The basic principle that seemed to apply was that an active source of current inside a finite conductive medium would produce volume currents throughout the medium and lead to potential differences on its surface. Given the columnar arrangement of the pyramidal cells in the cerebral cortex, it was suggested that, if enough of these cells in localized regions were synchronized in their oscillations between states of hyper- and depolarization, then volume currents large enough to produce measurable potential differences on the scalp might well be generated. * Tel.: +1 403 4926302; fax: +1 403 4928259; e-mail:
[email protected] 0013-4694/98/$19.00 1998 Elsevier Science Ireland Ltd. All rights reserved PII S0013-4694 (97 )0 0115-6
The process of predicting scalp potentials from current sources inside the brain is generally referred to as the forward problem in electroencephalography. If the configuration and distribution of sources inside the brain are known at every instant in time and the conductive properties of the tissues are known everywhere within the volume of the head, the potentials everywhere on the scalp can be calculated from basic physical principles. Conversely, the process of predicting the locations of the sources of the EEG from measurements of the scalp potentials is called the inverse problem. This is now a widely researched topic and there are many documented attempts to quantitatively predict the locations inside the brain of both event-related potentials and the background activity. At first sight, it might seem that the forward and inverse problems are similar. However, it turns out that, given a finite number of sites on the scalp at which the scalp potentials are measured at some instant in time, it is possible for an infinite number of source configurations to account for those measurements (Plonsey, 1963). Theoretically, only an infinite number of measurement sites on the surface of the brain would enable a unique determination of the locations of the responsible sources inside. In view of the non-uniqueness of the solution to the
EEG 33
128
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
inverse problem, it would seem that attempting to localize the sources of the EEG is a futile exercise. However, if the problem is approached by recognizing that, in the collection of possible source configurations inside the brain that can account for the measurements on the scalp, there is only one that is correct, then a fruitful way to begin might be to eliminate at least those regions in the brain were sources are known not to exist. For example, sources cannot be located in the skull or in the ventricles of the brain. These constraints alone will eliminate many of the source configurations that can account for the measurements. The conclusion therefore is that the solution to the inverse problem involves the collection of sufficient knowledge, some of which might not even come from the EEG, that will enable the exclusion of all but the correct source locations. Anatomic information obtained using one of the imaging modalities, for example, could be used to constrain sources to only certain locations around cortical folds (Scherg and Berg, 1991). Also, methods such as principal component analysis (Golub and van Loan, 1983), applied to the temporal aspect of the EEG, can be used to estimate the number of active sources and to isolate those that are relevant to the disease process of interest (Koles et al., 1995). With the necessary constraints in hand, source localization proceeds with an analysis of the spatial patterning of the EEG. In this analysis, the number of electrodes that should be used to obtain the scalp recording becomes a question. Indeed, it has been suggested that perhaps 200 electrodes on the scalp are necessary to capture all of the spatial information in the EEG (Gevins, 1984). Moreover, simulations have shown that variations in the thickness of the skull will have a significant effect upon the shape of the potential surface on the scalp. Since skull tissue is two orders of magnitude less conductive than brain tissue, the potential surface on the scalp is very much distorted from what it is on the surface of the brain. The effect of this low conductivity is to spread the volume currents out over the scalp and to attenuate the magnitude of the potential differences there. Therefore it would seem evident that knowledge of the variations in skull thickness is important for the accurate prediction of the source locations. Although Gevins (1997) does not attempt to localize the sources of cognitive evoked potentials, he is clearly concerned about the distortion of the potential surfaces on the skull. He employs a method that utilizes finite elements to model the variations in skull thickness. The method leads to topographic maps of potential surfaces as they might appear on the surface of the brain. The simplest realistic model of a source in the brain is the current dipole. It can easily be shown, for example, that a monopolar source cannot produce potential differences on the scalp. The current dipole consists of a source and a sink of electrical charge, with the source part delivering the charge to the brain volume and the sink part accepting the same amount of charge back from the volume. Therefore, with dipolar sources, the brain as a whole always remains electrically neutral. With elementary current dipoles as
source models, the localization process requires the iterative modification of 3 location parameters, two orientation parameters and one so-called moment parameter, for each active dipole. The moment parameter can be thought of as the strength of the dipole. Fig. 1 shows the potentials on the surface of a coronal slice through the head at the ears resulting from radially and tangentially oriented current dipole sources. The head model is a 3 shell sphere of unit radius with the scalp represented by the thick circular line. The asterisks on the scalp show the locations of the recording electrodes. The fine lines passing outside and inside the circle represent the potential profiles due to the two sources on the scalp. When the line is outside the scalp, its magnitude represents positive potentials; when it is inside the scalp, its magnitude represents negative potentials. The asterisks on the potential profiles show the potentials at the electrode sites. The figure indicates that the peak positive potential resulting from a radially oriented source is at the top of the head (Cz) while the positive potentials resulting from the tangentially oriented source occur on the left side of the head. It is interesting to note that, even though both the radially and tangentially oriented sources are at the same location inside the brain, the potential distributions that are produced on the scalp are very different. It would seem evident then that quantitative localization of the sources of the EEG requires accurate models of the generators of the electrical activity inside the brain, knowl-
Fig. 1. The potential on the circumference of a slice along the x-z plane of a spherical model of the head resulting from radially and tangentially oriented dipole current sources. The sphere is of unit radius and consists of three concentric shells. The origin of the coordinate system is at the center of the sphere with the +x axis pointing through the right ear, the +y axis pointing through the nose and the +z axis pointing through the vertex of the head. The thick circular line represents the scalp while the finer lines outside and inside the scalp represent positive and negative potentials respectively. The asterisks mark the electrode positions on the scalp, while on the potential profiles they mark the values of the potentials at the electrode sites.
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
129
head model that is not realistic in terms of shape, skull thickness and tissue types will certainly have an adverse affect on the accuracy of the predicted source locations.
2. Single time-slice source localization
Fig. 2. The potential profiles on the surface of the coronal slice of Fig. 1 resulting from a single tangential dipole source (dark fill, dark line) and from eight dipole sources located in a cluster nearer the center of the slice (gray fill, gray line). The two potential profiles match exactly at the electrode sites.
edge of the number of active generators and an accurate model for the shape and conductive properties of the head. Given these, the forward solution can be used, in a computerized algorithm, to iteratively modify the location parameters of the source models, subject to the available constraints, until the potentials predicted at the measurement sites on the surface of the head model match that observed at the same sites on the scalp. Finally, it should be re-emphasized that there is not a unique solution to the unconstrained inverse problem in electroencephalography. This non-uniqueness is illustrated in Figs. 2 and 3. In Fig. 2, the potential profile on the coronal slice produced by a single tangential dipole (black fill) is shown with the black line. The potential profile produced on the slice by the 8 dipoles (gray fill) is shown with the gray line. The profiles match exactly only at the electrode sites. This situation is illustrated more dramatically in Fig. 3 where the 8 dipoles (gray fill) have been moved to different locations. The disagreement between the profiles is now highly exaggerated; however, the match is still exact at the measurement sites. These figures clearly indicate that two combinations of sources in totally different regions of the brain can account for an observed EEG. The reality is that knowledge of the number of active sources responsible for the brain process of interest and whether these are point or distributed sources is necessary to obtain the correct solution to the problem. Specification of a source as a current dipole when it is actually an extended dipolar sheet may not lead to useful location information. Mis-specification of the number of active sources will usually lead to inaccurate estimates of the locations of all of the sources. Using a
This approach to the localization of the sources of the EEG was first applied to event-related potentials and more recently to epileptiform spikes (Ebersole and Wade, 1990; Wong, 1998). The main assumption is that time points in the recording corresponding to important landmarks such as, for example, at the P300 or the peak of the epileptiform spike, are largely the result of the activity of a known number of dipole current sources. These landmarks are usually points in time at which the waveform amplitude is maximum or minimum at one electrode site in the recording montage. This site is referred to as the location of maximum field power. As shown by Fig. 1, a single maximum in the field power distribution on the scalp part of the head is an indication of a radially oriented source while positive-negative pair indicates a tangential orientation of the source. The shape of the potential surface around the point of maximum field power determines the depth of the source. A rapid decline around a large maximum indicates a shallow source while a slow decline of the field power from a smaller maximum indicates a deeper source. With this insight from the recording as the starting point, a head model, an estimate of the number dipole sources involved and an adequate number of recording electrodes, the inverse problem can be solved for single time-slices of the EEG.
Fig. 3. The potential profiles on the surface of the coronal slice of Fig. 1 resulting from a single tangential dipole source (dark fill, dark line) and from 8 dipole sources located in a cluster in the left region of the slice (gray fill, gray line). The two potential profiles match exactly at the electrode sites.
130
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
To localize dipole sources, given that the number of parameters which must be estimated to localize a single current dipole is 6, the EEG must provide at least six independent views of the source at each time-slice. This means that the recording montage must contain at least 6 electrodes for each source to be localized and that each of these electrodes must provide some new information about the location of the source. Therefore, the electrodes must be well separated on the potential surface generated by each source on the scalp. Having only a few electrodes closely placed on this surface does not provide these views well and results in what is termed an ill-conditioned problem. Solutions to ill-conditioned problems are known to be unstable and easily biased by, for example, artifacts in the recording. In practice, where artifacts and extraneous noise contaminate the recordings, electrodes additional to the minimum 6 for each source particularly on regions of the scalp where the potential surface is changing most rapidly is desirable for accurate source localization. In other words, to deal with aspects of the recording which are not related to sources in the brain, the images of the sources on the scalp should be oversampled as much as possible. Noise and artifacts in recordings are particularly problematic in single time-slice source localization since there is no way of averaging out their effects other than by utilizing additional electrodes. In mathematical terms, the single time-slice source localization problem is formulated as follows. Given the potentials at the single time-slice, these can be collected into a vector v which has one column with N rows, each row corresponding to an electrode site. The problem then is to find a vector m, the collection of potentials at the same electrode sites but generated by the assumed sources at known locations inside the brain, which is the same as v. In practice, the procedure is first to estimate a starting point (from the scalp field as suggested above) and then to iteratively move the sources around in the allowed locations inside the brain (using, for example, anatomic constraints) in the attempt to produce the match. This involves solving the forward problem repetitively and calculating a measure of the difference between the observed and the predicted potential vectors at each step. The measure most often used to do this is the mean-squared distance between the two vectors. This is written as: J = kv − mk2
3. Equivalent dipole source localization The assumption that the sources of the EEG are current dipoles will probably seldom be correct in practice. Given the anatomy of the brain, particularly the cerebral cortex, a more likely model for the source is the current dipole sheet. The approach in this case could be one of attempting to match the measurement v on the scalp with that produced by a current dipole sheet and modifying the parameters of this model source until the error function in Eq. (1) is minimized. This would seem like a reasonable approach since the image of a dipole sheet on the surface of a model head can be calculated by simply summing the contributions of the individual dipoles in the sheet. However, each dipole added to the sheet brings with it six new parameters which must be iteratively estimated in the localization process. Supposing that the dipole sheet is made of up say 10 elementary dipoles, 60 parameters would be required to specify the location of the sheet. While the previous discussion would suggest that 60 electrode sites are required for measurements that are well enough conditioned to do this, it may be that far fewer would be sufficient. This is because the 60 parameters are constrained to form a dipole sheet and, therefore, they are not totally independent of one another in terms of their location and moment parameters. Therefore, attempts to localize dipole sheets may not be as ill conceived as might first seem. An alternative approach to the above is to attempt to minimize Eq. (1) using a single current dipole irrespective of the actual configuration of the source. The results obtained with this procedure yield the location of what is called the equivalent dipole for the extended source. It is assumed that the equivalent dipole is located at the center of activity of the actual source. While this may have some clinical value, it could be somewhat removed from point
(1)
where J is the error function which is to be minimized. This measure amounts to the sum of the squared differences between the observed and predicted potentials at all the electrode sites. When this function is minimized, subject to whatever constraints are available, the process of moving the sources stops. It should be noted that, in practice, J will not be zero at its minimum. This could happen because the number of sources was not estimated correctly, the sources are not dipoles, the head model was not accurate and there was environmental/instrumentation noise and/or artifacts present in the measurements on the scalp.
Fig. 4. The coronal slice in the x-z plane of the head model with an extended source consisting of 13 current dipoles. Together the dipoles span an arc length of 0.46 units in the sphere of unit radius. The location of the equivalent dipole is indicated by the black dot.
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
131
Fig. 5. The topography of the transformed error function − log(J) in the coronal slice obtained during the localization of the equivalent dipole shown in Fig. 4. The meshes show the values of the error function for candidate sources on a 0.05 unit grid in the slice. Side A of the figure shows the topography with no noise added to the potential measurements on the scalp while side B shows the topography with a scalp signal to noise variance ratio of 0.05.
of activity relevant to the disease process of interest. Fig. 4 shows a dipole sheet used as a model source in a spherical head model in some of our work. The sheet is made up of 13 unit dipoles in an arc arrangement and is placed in the same coronal slice shown in Figs. 1–3. Fig. 4 indicates that the location of the equivalent dipole, found as described above, is as expected but away from any of regions of the source that are actually active. The main difficulty in locating the equivalent dipole of an extended source is related to finding the minimum in the error function. In this case, since the vector v is the image on the scalp of the extended source and m is the image of a dipole source, it is almost certain that the minimum for J will never be zero. A bigger problem in practice is that, while the surface representing J does in fact come to a minimum at the location of the equivalent dipole, the region around the minimum can be very shallow. Because of this, the minimum is sometimes difficult to find accurately. Moreover, this difficulty can be exacerbated by noise (or extraneous sources) in the recorded scalp potential. To illustrate this problem, we computed the potential on the surface of a spherical head model produced by the source of Fig. 4 and computed values for J by placing a candidate dipole source at several hundred grid points in the coronal slice. Instead of plotting the values of J directly we plotted −log(J). This transformation results in a better visual presentation of the topography of J in the slice because it emphasizes any values of J near zero and results in the minimum in J appearing as a maximum. Fig. 5a shows the topography of J in the slice with a clear maximum at the location of the equivalent dipole source. Fig. 5b shows this surface after noise was added to the scalp potential produced by the extended source. The flatness of the surface and the poor definition of the peak are now very evident. In Fig. 5, the x and z directions are the same as those used in Figs. 1–4.
4. Multiple time-slice source localization The multiple time-slice approach to the localization of the sources of the EEG utilizes both the spatial and the temporal components of the EEG (Scherg and Ebersole, 1993). The basic premise is that the EEG can be modeled as a number of current sources whose locations, during the time interval of observation, are fixed (stationary) inside the brain and that the variations in scalp potential are due only to variations in the strengths of these sources. The coupling of the sources to the scalp is by volume conduction and the assumption is that the contributions of each of the sources to the scalp potential are additive. In other words, the conductivity of the head tissue is assumed to be linear. Given the above, the EEG can be separated into a spatial component, consisting of the images of unit strength sources on the scalp, and a temporal component, indicating how these strengths change through time. The image of a unit source at the electrode sites is also called the lead field of the source. In multiple time-slice source localization, the sources are positioned iteratively in the head model until they account for the maximum portion of the temporal variations present in the EEG. In terms of the spatiotemporal model, the EEG is described mathematically as: Q
vi (t) = ∑ mij sj (t) j=1
(2)
where vi(t) is the potential at electrode i at time t, mij is the lead field of unit source j at electrode i and sj(t) is the moment waveform of source j. There are Q sources in the model. Therefore, the EEG at each electrode site is the sum of the weighted lead fields from each source. The weighting of the lead fields changes with time and is reflected in the source waveforms. If the EEG is obtained in digital form, Eq. (2) can be written in matrix notation as:
132
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
Fig. 6. The topography of the transformed error function − log(J) for two dipole sources located at x = ±0.05 and z = 0.60 in the coronal slice. The meshes show the values of the error function for candidate sources on a 0.05 unit grid. Side A of the figure was obtained with 5 electrodes between the ears on the circumference of the slice while side B was obtained using 19 electrodes. Note that the x and z axes in this figure are reversed from those in Fig. 5.
V=M S
(3)
where V is the N electrodes (rows) by T time-slices (columns) EEG data matrix, M is the N electrodes by Q sources lead field matrix and S is the Q sources by T time-slices source waveform matrix. Matrix multiplication of M and S yields the EEG data matrix. Multiple time-slice source localization involves attempting to match the right side of Eq. (3) to the measured data matrix by, once again, iterating the location parameters of the sources. This iteration is continued until the error function: J = kV − M Sk2
(4)
is minimized. In practice, assuming Q sources and their location parameters at some starting point inside the head, the M matrix is calculated from the source and head models. The source waveform matrix S which can best account for the data matrix V in the least-squares sense is then computed and the value of J in Eq. (4) determined. The iteration of source parameters is continued until a minimum for J is found. The appealing aspects of the multiple time-slice approach to the localization of the sources of the EEG are related firstly to the reduced number of source parameters that must be estimated by iteration and secondly to its immunity to noise and artifacts in the EEG data matrix. The source parameter reduction property comes from the fact that the lead field matrix M contains the images of unit strength sources. Therefore, in the case of current dipole sources, only 5 and not 6 parameters need to be iterated in the localization process for each source assumed to be present. Given M for the unit sources, the moment matrix S can be calculated deterministically to fit the measured data. In other words, iteration is not necessary to find S. In addition, since the S matrix at each step is related to what the sources at particular locations inside the brain can account for maximally in V, environmental sources and muscle artifacts with their sources presumably
outside the head will tend to be excluded. Therefore, the least-squares determination of S at each iterative step can be viewed as a filtering process that rejects sources of noncerebral origin from the recording of the EEG. The source waveform matrix S which best accounts for the data matrix V in the least-squares sense is obtained from: p
S=M V
(5)
where * is the so-called pseudoinverse of the estimated lead field matrix M at each iteration. If Eq. (5) is substituted into Eq. (4) and a small amount of algebra is performed, it can be shown that: p
J = k(I − M M )Vk2
(6)
and the explicit determination of S is removed from the process. The pseudoinverse is a well-defined matrix operation derived from what is called the singular value decomposition of M. In terms of Eq. (6), it is said that the projection of V (the measured EEG) into the left null space of the lead field matrix M should be minimized in the source localization process. We note here that, if the number of assumed sources Q is equal to N the number of electrodes, the pseudoinverse becomes the proper inverse of M and J = 0. The proper inverse will exist as long as the N sources assumed to be active are distinct in the sense that none of their images can be expressed as the weighted sum of the images of the others. Therefore, if Q = N, N sources placed anywhere in the head will be able to account fully for the observed EEG! Evidently, for the multiple timeslice approach to work, the number of independent sources assumed to account for the observed EEG must be less than the number of electrodes used for its measurement. This statement is not inconsistent with the earlier statement that, for each dipole to be localized, there should be at least 6 electrodes in the recording montage. However, in the sense that the number of sources must be less than the number of
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
measurement sites, the multiple time-slice source localization approach is said to be an overdetermined problem. An examination of Eq. (6) indicates the major shortcoming of the multiple time-slice source localization procedure as described to this point. In order to compute the pseudoinverse of M, the location parameters of all the sources must be known. Therefore, in order to implement this approach, these parameters must all be iterated at once and the mislocation of one source will probably lead to the mislocation of all of them. To overcome this problem, a modified form of the multiple time-slice source localization procedure has recently been developed based on the multiple signal classification method (Koles et al., 1995; Koles and Soong, 1998). This method is usually referred to as MUSIC (Mosher et al., 1992). The modified MUSIC method requires that the EEG be first numerically decomposed into spatial and temporal components as: V=X Y
(7)
where in this case X is a set of arbitrary spatial patterns and Y is the corresponding set of temporal patterns. The spatial patterns are obtained from an analysis of the covariance of the EEG between recording sites. These can be obtained using, for example, the method of principal component analysis (Achim et al., 1988) or from that of common spatial patterns (Koles, 1991) or from others (Harman, 1976). The objective is to isolate the brain process of interest into a subset of the spatial patterns in X. This is done by examining the temporal waveforms in Y and selecting those that, for example, contain epileptiform features. We have done this by finding the common spatial patterns between an EEG from an epilepsy patient and from a normal control and, after rotating these using the Varimax procedure, we were able to recognize interictal spikes in only a few of the temporal waveforms (Koles et al., 1995). To localize the sources of these spikes, the error function used was: p
J = k(I − Xs Xs )Mk2
(8)
where Xs was the column subset of X whose corresponding temporal waveforms contained the spikes. The main advantage of Eq. (8) for source localization is that the lead field matrix is now projected into the left null signal space containing the epileptiform activity of interest. Given the rules of matrix algebra, it is evident that this projection of M is the sum of the projections of its columns. Therefore, the minimum for J will occur when the projection of each of the columns (the lead field of each source) is also minimum. In practice, then, the null signal space can be scanned with just a single source and the first Q minima will correspond to the sources sought. In other words, when formulated this way, the multiple time-slice source location method enables sources to be localized one at a time and this number does not need to be estimated beforehand. In a study designed to validate the method of multiple time-slice source localization described by Eq. (8),
133
we placed two current dipoles in the coronal plane of the spherical head model and simulated source waveforms for each of them (Koles and Soong, 1998). Using Eq. (3) we computed the EEG on the circumference of the slice at the electrode locations. By analyzing the spatial covariance of this simulated EEG, we were able to isolate the signal space Xs containing the sources. A candidate dipole source was once again placed at several hundred grid points in the slice and the topography of J determined according to Eq. (8). Fig. 6a shows the topography of −log(J) for 5 electrodes on the slice, while Fig. 6b shows it for 19 electrodes. Both figures show two sharp peaks at the locations at which the sources were actually placed. Fig. 6a, however, shows several additional peaks apparently unrelated to the locations of the actual sources. These spurious peaks are thought to be related to the fact that at least 8 electrodes must be used to uniquely locate two dipole sources in a slice through the head (two location, one orientation and one moment parameter for each source). A practical difficulty with the above method is that multiple minima in J must be found. This can be a problem if some of the minima are not well defined. An approach that would appear to overcome this problem consists of finding only the global minimum for J and adding the lead field vector of the located source to Xs as a new column. The process is then repeated over and over at each step finding a new global minimum and adding a new column to Xs . This is similar to the recursively applied MUSIC (RAP-MUSIC) algorithm recently described by Mosher (1996). The major difficulty with RAP-MUSIC is that, if one of the sources is mislocated, all those following it will probably also be mislocated. In order to gain some insight into the precision of the predicted location of an equivalent dipole in practical situations using the multiple time-slice approach, we have recently completed a simulation study in which the extended source shown in Fig. 4 was placed in a 3 concentric shell spherical head model (Koles and Soong, 1998). The lead field produced by this source at 43 sites on the surface of the model was calculated and a source waveform from an actual EEG used with Eq. (2) to create the data matrix. Various amounts of random noise were then added to surfaces potentials and the location of the equivalent dipole source estimated by minimizing Eq. (8). The results of this study are summarized in Fig. 7. This figure indicates, for example, that at a signal-tonoise variance ratio of 0.07 and using 3 s of the temporal component of the simulated EEG, an error of about 3 mm in the location of the equivalent dipole could be achieved.
5. Minimum norm localization The source localization approaches described above are
134
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
it now contains 3 columns for every grid point in the head model. At a single time-slice, the moments of the sources at the grid points are related to the potentials measured at the electrode sites by the vector-matrix expression: v=M i
Fig. 7. The results of a Monte Carlo study designed to determine the precision of the predicted location of an equivalent dipole source using the multiple time-slice source localization approach. The source of the simulated EEG was the extended dipole sheet shown in Fig. 4. The source waveform for the simulation was derived from an actual EEG. Various amounts of random noise were added to the scalp potentials and the location of the equivalent dipole source estimated by minimizing Fig. 8. In this study, the signal space Xs was taken to be the first principal component of the simulated EEG data matrix. The epoch length used to obtain this principal component is that plotted in this figure.
dependent on the model chosen for the source. They involve scanning the head volume with the model source until an error function is minimized. Because in most applications the exact form of the source is not known, extended sources are usually located with an equivalent current dipole. In the multiple time-slice approach, the number of active sources can be estimated from a covariance analysis of the recorded data. This cannot be done with single time-slice data. When the number of sources is not correctly estimated, it is likely that none of the predicted source locations will be correct. These problems have led to the development of what are called minimum norm approaches to the localization of the sources of the EEG (Pascual-Marqui et al., 1994). These approaches lead to estimates of the current density everywhere in the 3 dimensional volume of the head. To localize the sources of the EEG, it is then necessary to scan the head volume to locate the regions of highest current density. These regions are then taken to be the locations of the sources. In the minimum norm approach to source localization, the head model is first mapped onto a 3 dimensional grid. The lead field at each electrode site is calculated from 3 mutually perpendicular dipole current sources of unit strength placed at each of the grid points. The 3 model sources are oriented with respect to one another along the axes of the measurement system (usually denoted x, y and z) at each grid point. The lead fields from all of the grid points are collected into N electrodes by 3 times the number of grid points matrix. In practical situations, the number of columns in this matrix could be 100 000 or more. This matrix is the same as the source image matrix M defined above in the multiple time-slice source localization approach except that
(9)
where v is the scalp potential vector (one column length equal to the number of electrodes N) and i the current density vector (one column of length equal to 3 times the number of grid points). In this equation, the vector i can be viewed as that weighting of M, the lead fields due to unit sources, required to produce the actual potentials measured on the scalp. The objective in this localization approach is to determine i which leads to the measurement v. We note here the similarity between Eq. (3) and Eq. (9). In the minimum-norm approach to source localization, the objective is therefore to solve for perhaps 100 000 unknown currents given a relatively small number of observation points on the scalp. This is a vastly underdetermined problem and one therefore for which there is an infinite number of solutions. However, it is argued that of all the solutions possible, a reasonable one to choose is the one that contains the least energy. Minimum energy here means that the overall current density inside the head is minimum. Stated mathematically: J = kik2
(10)
should be minimized. In minimum norm source localization, Eq. (10) is minimized subject to the constraint that the solution i satisfies Eq. (9). In other words, the minimum norm solution to Eq. (9) is chosen from amongst the infinitely many that are possible. It can be shown that the minimum norm solution to Eq. (9) is given by: p
i=M v
(11)
where, once again, * indicates the pseudoinverse of the matrix M. Simulation studies using dipole current sources in spherical head models indicate that minimization of Eq. (10) leads to predicted current densities that are larger near the surface of the head model than is actually the case. It is said, therefore, that the minimum norm solution leads to current densities that are biased towards the measurement electrodes. In an attempt to eliminate this bias, the weighted minimum norm source localization approach has been suggested. In this approach, the function minimized instead is: J = kW ik2
(12)
where W is a square diagonal weighting matrix with dimensions equal to the number of rows in i. Each diagonal element of W is chosen to downweight the grid points in the head model that are closer to the surface. In this way, the bias toward large current densities near the surface is reduced.
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
Results obtained using the weighted minimum norm approach applied to single dipole source models indicate that these sources can indeed be localized but that the 3 dimensional current density predicted inside the head is very blurred. This makes it difficult to accurately locate the position of the sources. In the attempt to rectify this situation, a recursive application of the weighted minimum norm approach has recently been proposed (Gorodnitsky et al., 1995). This algorithm, termed FOCUSS (FOCal Underdetermined System Solution), involves repetitively minimizing Eq. (11) but at each step reducing the weighting on those grid points that have the lowest current densities. This process is continued until the current densities at most of the grid points are zero and no change occurs in the 3 dimensional current density after each step. It is shown that remarkably accurate localization of dipole current sources including dipole sheets can be obtained with this approach.
6. Head models The forward problem in electroencephalography requires a model of the head to connect the activity of the sources to the measurements on the scalp. In the head model, the electrical properties of the various tissues involved, namely the brain, the skull and the scalp, are mathematically described. As indicated earlier, these properties are assumed to be purely resistive. It is evident that the accuracy with which the head is modeled will determine the accuracy with which the source locations can be estimated. The head models that have been used for source localization have been categorized as either spherical or non-spherical (Rosenfeld et al., 1996). The spherical models used have been both homogeneous and non-homogeneous; however, the non-homogeneities in these models have generally been confined to concentric shells representing the skull and the scalp. Within the volume in each shell, the tissue is assumed to be homogeneous and isotropic. The reason for the choice of the spherical model is obviously simplicity; however, there is also the fact that the potential image on the scalp due a dipole current source located anywhere inside can be predicted analytically (Ary et al., 1981). That is, given the location parameters of the dipole, a mathematical formula can be used to calculate the potential at any point on the surface of the model. Given that the model is a volume conductor and that it is linear, the images on the scalp of more complex extended sources such as dipole sheets can be calculated by simply summing the contributions of the individual dipoles. Generally, the scalp images of sources in more realistic head models which incorporate actual head shapes and volumes of variable conductivity cannot be computed analytically. They must instead be computed using a more complex numerical algorithm. This usually means a significant
135
increase in computational burden and of storage memory. The numerical approach requires that the volume occupied by the head be discretized. This discretization ranges from the relatively coarse, such as in what are called the boundary element methods (Cuffin, 1990), to the fine, such as in the finite element methods (Yan et al., 1991). In both the boundary element and finite element models, analytic methods are applied in each region or volume, but numerical methods are required to deal with the interaction between the boundary or volume elements. In these numerical methods, the potential on all the boundary elements and the volume current in all of the finite elements as a result of the activity of the source are involved simultaneously in the computation of the potential on the outside surface of the model. The boundary element methods are based on head models which are essentially distortions of the concentric spherical shell models. These distortions are introduced to better fit the shape of the actual head involved as indicated, for example, by magnetic resonance or computed tomographic imagery. These methods allow a non-uniform skull thickness to be modeled. Boundaries between the tissue types, for example the brain and the skull, are represented by tiled surfaces, with a point on each tile forming a single grid point in the model. Therefore, in the boundary element methods, it is the shell surfaces that are discretized. In most applications of the boundary element method, about 1000 triangular tiles are used for each of 3 boundaries in the model: the brain/ skull, the skull/scalp and scalp/air boundaries. The main disadvantage of the boundary element methods is that each bounded volume must be electrically homogeneous and isotropic. In the finite element methods, the head is divided into 3 dimensional elemental volumes. In these methods, each of the elemental volumes forms a grid point in the head model. This enables both the actual head shape and conductivity variations between the gray and white matter and the cerebrospinal fluid and of the variations in skull thickness to be included in the model. In addition, anisotropy in the conductivity of any of these tissues can be accounted for. In applications of the finite elements method, as many as 340 000 grid points have been used to model the head. Solving the forward problem in electroencephalography using the finite elements method will therefore place considerable demands on the computing resources available in terms of both the memory requirement and the computational burden. The main problem with the finite elements method is that the density of grid points in the region of the source must be high to account for rapidly changing current densities. It is sometimes not practical to maintain this density throughout the entire head volume. Therefore, when attempting to solve the inverse problem by moving sources iteratively around inside the head, regridding of the head volume may be required at each step. This situation will increase computational demands and slow the localization process.
136
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137
7. Discussion There has been significant progress over the last 15 years in the development of new methodologies for localizing the sources of the EEG. In the 1970s single dipole sources were localized in the attempt to account for single time-slices usually in event-related potentials. Since that time, based on a careful consideration of the basic physical principles underlying the EEG and upon mathematical models developed on the basis of the new understanding obtained in this way, researchers are now attempting to account for the full spatial and temporal characteristics of the EEG and for the 3 dimensional distribution of current densities inside the head. The latter is a highly underdetermined problem and the need for more electrodes in the collection of EEG data would seem to be obvious. The basic difficulty will, however, remain and that is that there is no unique solution to the inverse problem in electroencephalography. The conclusion is that, to localize the sources of the EEG, it is necessary to choose the source configuration that can satisfy all of the known constraints on the problem. In multiple time-slice source localization, an analysis of the spatial covariance of the EEG enables the estimation of the number of active sources responsible for the recording. This is the main advantage of the multiple time-slice approach compared to the single time-slice approach where there is no quantitative way of estimating the multiplicity of the sources. The main advantage of the minimum norm approaches is that knowledge of the source multiplicity is not required. These approaches, however, produce very blurred images of the source locations and are, at present, only applicable to single time-slices from the EEG. The challenge that remains then is to extend the minimum norm approach to multiple time-slices and to apply it with realistic models of the head. Methods of deblurring the current distributions predicted by the minimum norm approaches are beginning to emerge and it may be that one of the ways in which this will happen is through a combination of this approach with the dipole-based localization approaches. In other words, the minimum norm might well serve as the starting point for the iterative localization of current dipole sheets. Probably the most exciting aspect of the trends in EEG source localization is that there has been progress toward being able to predict the underlying current densities everywhere in the 3 dimensional volume of the head. As these approaches come to utilize more realistic models of the head, the day can be foreseen when tomographic images of the current densities in the brain that are correlated to epilepsy, for example, can be coregistered with structural and functional images obtained with computed tomography, magnetic resonance imaging, single photon emission computed tomography and positron emission tomography. Indeed, these modalities will become an integral part of the source localization process in that they will provide some of the constraints necessary
to provide clinically meaningful solutions to the inverse problem. In conclusion, it can be said that the accuracy of the localization of the sources of the EEG is strongly influenced by the accuracy of the head model. It should also be said, however, that the use of the most realistic finite elements head model may not, in every case, be justified, given the cost in computing hardware involved. In pychopathology, for example, localization of equivalent sources to even one of the lobes of the cerebrum may be a significant step forward in mental health research. Therefore, in this case, the 3 shell spherical model is probably sufficient.
References Achim, A., Richer, F. and Saint-Hilaire, J.M. Methods for separating temporally overlapping sources of neuroelectric data. Brain Topogr., 1988, 1: 22–28. Ary, J.P., Klein, S.A. and Fender, D.H. Location of sources of evoked scalp potentials: corrections for skull and scalp thicknesses. IEEE Trans. Biomed. Eng., 1981, 6: 447–452. Brazier, M.A.B. A study of the electrical fields at the surface of the head. Electroenceph. clin. Neurophysiol. Suppl., 1949, 2: 38–52. Cuffin, B.N. Effects of head shape on EEGs and MEGs. IEEE Trans. Biomed. Eng., 1990, 37: 44–52. Ebersole, J.S. and Wade, M.D. Spike voltage topography and equivalent dipole localization in complex partial epilepsy. Brain Topogr., 1990, 3: 21–34. Gevins, A.S. Analysis of the electromagnetic signals of the human brain: milestones, obstacles and goals. IEEE Trans. Biomed. Eng., 1984, 31: 833–850. Gevins. Electroenceph. clin. Neurophysiol., 1997. Gorodnitsky, I.F., George, J.S. and Rao, B.D. Neuromagnetic source imaging with FOCUSS: a recursive weighted minimum norm algorithm. Electroenceph. clin. Neurophysiol., 1995, 95: 231–251. Golub, G.H. and van Loan, C.F. Matrix Computations. The Johns Hopkins University Press, Baltimore, MD, 1983. Harman, H.H. Modern Factor Analysis, 3rd edn. University of Chicago Press, Chicago, 1976. Henderson, C.J., Butler, S.R. and Glass, A. The localization of the equivalent dipoles of EEG sources by the application of electric field theory. Electroenceph. clin. Neurophysiol., 1975, 39: 117–130. Koles, Z.J. The quantitative extraction and topographic mapping of abnormal components in the clinical EEG. Electroenceph. clin. Neurophysiol., 1991, 79: 440–447. Koles, Z.J. and Soong, A.C.K. EEG source localization: implementing the spatio-temporal decomposition approach. Electroenceph. clin. Neurophysiol., 1998, in press. Koles, Z.J., Lind, J.C. and Soong, A.C.K. Spatio-temporal decomposition of the EEG: a general approach to the isolation and localization of sources. Electroenceph. clin. Neurophysiol., 1995, 95: 219–230. Mosher, J.C. Recursively applied MUSIC (RAP-MUSIC), MEG and EEG for human brain mapping, fundamentals through advanced issues. In: Proceedings 2nd International Conference on Functional Mapping of the Human Brain, Boston, MA, 1996. Mosher, J.C., Lewis, P.S. and Leachy, R.M. Multiple dipole modeling and localization from spatio-temporal MEG data. IEEE Trans. Biomed. Eng., 1992, 39: 541–557. Nunez, P.L. Electric Fields of the Brain: The Neurophysics of EEG. Oxford University Press, New York, 1981. Pascual-Marqui, R.D., Michel, C.M. and Lehmann, D. Low resolution electromagnetic tomography: a new method for localizing electrical activity in the brain. Int. J. Psychophysiol., 1994, 18: 49–65.
Z.J. Koles / Electroencephalography and clinical Neurophysiology 106 (1998) 127–137 Plonsey, R. Reciprocity applied to volume conductors and the EEG. IEEE Trans. Biomed. Eng., 1963, 10: 9–12. Plonsey, R. (Ed.). Bioelectric Phenomena. McGraw-Hill, New York, 1969, pp. 304–308. Rosenfeld, M., Tanami, R. and Abboud, S. Numerical solution of the potential due to dipole sources in volume conductors with arbitrary geometry and conductivity. IEEE Trans. Biomed. Eng., 1996, 43: 679–689. Scherg, M. and Berg, P. Use of prior knowledge in brain electromagnetic source analysis. Brain Topogr., 1991, 4: 143–150. Scherg, M. and Ebersole, J.S. Models of brain sources. Brain Topogr., 1993, 5: 419–423.
137
Schneider, M. A multistage process for computing virtual dipolar sources of EEG discharges from surface information. IEEE Trans. Biomed. Eng., 1972, 19: 1–12. Shaw, J.C. and Roth, M. Potential distribution analysis II: A theoretical consideration of its significance in terms of electrical field theory. Electroenceph. clin. Neurophysiol., 1955, 7: 285–292. Wong, P.K. Potential fields, EEG maps, and cortical spike generators. Electroenceph. clin. Neurophysiol., 1998, 106: 138–141. Yan, Y., Nunez, P.L. and Hart, R.T. Finite-element model of the head: scalp potentials due to dipole sources. Med. Biol. Eng. Comput., 1991, 29: 475–481.