1. Introduction

The localization of the generators of electric or magnetic fields recorded at/near the scalp have been for decades a basic goal of brain researchers. Not even the quick development of functional imaging techniques (PET, SPECT, fMRI) has reduced the importance of this problem. So far, no functional technique can compete with neuroelectromagnetic measurements in terms of temporal resolution. This is an indisputable advantage of the EEG/MEG over other functional techniques, since many normal mental processes are known to occur within an interval of 50-500 milliseconds after stimulus presentation. In this fast processing of information many different brain structures are known to be activated serially or in parallel. The major problem with scalp recorded electromagnetic fields is that precise anatomic localization of the activity requires the solution of a rather complex mathematical problem: the neuroelectromagnetic inverse problem (NIP). The NIP lacks a unique solution and thus anatomical, physical, neurophysiological, mathematical or functional constraints have to be incorporated. The reliability of the reconstruction of the neuronal generators obtained will be dictated by the veracity of the constraints incorporated to the solution. Also, the spatial resolution of the reconstructed activity will be mainly driven by the resolution of the scalp measurements and consequently rather gross when compared with single neuron studies. 

In summary, neuroelectromagnetic brain imaging has excellent temporal resolution while lacks spatial resolution. The contrary applies to functional neuroimaging, in particular, to functional magnetic resonance (FMRI). It is therefore extremely appealing to essay to combine these techniques in order to take advantage of their relative strengths. Such integration has been facilitated by the recent development of a technique denominated “event related fMRI“ (see Rosen et al, 1998, for a review). In event related fMRI, regional responses to single sensory or cognitive events can be detected using experimental paradigms identical to the ones employed in neurophysiological or cognitive studies. Accordingly, it is nowadays possible to design experimental setups in which functional and electric/magnetic brain responses to identical stimuli are available. This paper discusses theoretical alternatives to formalize this integration as well as some of its pitfalls and caveats.

The simplest alternative to extract information common to both modalities is to compare the maps of the anatomically constrained solution to the NIP with the statistical result of the analysis of the fMRI. Using the individual subject MRI a linear or a nonlinear solution can be applied to the scalp recorded data to estimate the generators in a solution space restricted to the gray matter. Analysis of the temporal curves provided by the inverse solution for activation blobs, which are common in both modalities, can lead to information about the temporal development of the ongoing processes. A second variant is to directly incorporate the functional images as part of the a priori information needed to solve the NIP. The basic assumption underlying latter approach is the existence of a coupling between neuronal activation and BOLD responses, a point still largely subject to debate. 

The first part of the paper discusses some alternatives to incorporate functional constraints into the solution of the inverse problem. The second part discusses the main practical aspects that limit the integration of modalities. Some examples of the divergence in the anatomical localization found using functional techniques and electrocortical stimulation are presented. The examples are used to motivate a discussion about the potential dangers of attributing excessive weights to the functional constraints. Finally, a comparison is presented between the results of an inverse solution without functional constraints in an epileptic patient and the EEG triggered fMRI. The example shows that the localization accuracy of the inverse solution is comparable to that of the fMRI while surpassing it in terms of temporal resolution. Some techniques to improve the spatial resolution of the inverse solution are also discussed. 

2. Methods

2.1 Some alternatives to constraint the NIP by functional information

The discrete electromagnetic inverse problem can be stated as: 

(1)

where column vector d represents the data , i.e., the measurements obtained over the N_s electric or magnetic sensors at a fixed time and vector e, represents the noise contribution. Matrix L, is usually termed the lead field matrix (Hämäläinen et al. 1993), and describes the physical relationship existing between the source model and the measurements in the selected head shape model. Vector j stands for the N_p unknown parameters that determine the sources. 

Since the number of unknowns is generally bigger than the number of sensors, the solution of (1) requires some additional information. One alternative to obtain the additional information is the use of functional images (PET, SPECT, fMRI). A first limitation for the integration of these images into the solution is that while functional images provide information about a scalar field, the inverse problem (1) is commonly stated in terms of a vector field (e.g. current density vector). Consequently, we prefer here to restate problem (1) in terms of the estimation of a scalar field for both the electrical and magnetic field cases. 

1) Electrical case. 

If the analysis is confined to electrical data, there is no reason to expect solenoidal (silent) sources. Thus problem (1) can be stated in terms of a scalar field j which represents the discretization of the electrical potential in depth. L is the product of the discrete lead field operator times the gradient operator. This model expresses that only irrotational currents can be at the origin of scalp measured electric fields and has been termed ELECTRA. More details are given in Grave et al (2000). 

2) Electrical and /or magnetic case. 

A vector field v(r) can be always decomposed as: 

(2)

Where d(r) is a direction vector field with norm (modulus) one and m(r) is a scalar field with absolute value equal to the modulus of v(r). Here L stands for the discretization of the product of the discrete lead field and the direction vector field d(r) while j denotes the discrete scalar field m(r). The direction vector d(r) can be selected on anatomical basis, i.e., normal to the brain surface or could be estimated from an inverse solution, e.g., the minimum norm. Note that if the estimated  m(r) is negative then the direction in d(r) should be reverted.
 

After this transformation, it is easier to combine the scalar field provided by the functional image, f(r), and the scalar field associated to the electromagnetic inverse problem j(r). Let’s assume that electrophysiological events and functional events are coupled, i.e., there is a correlation (linear or not ) between j(r) and f(r). Then the following strategies could be applied: 

a) Construction of a parametric model and identification of the parameters. 

Let’s assume the existence of a correlation between both modalities, i.e., a certain spatio-temporal model can be assumed for both j(r) and f(r) in the form, j(r)=j(r,Q) and f(r)=f(r,Q).  In such case, we could identify the parameters Q of the model  from the functional image and afterwards estimate j using equation (I). One example of the parametric models that can be considered is the Markov Random field model (Kinderman and Snell 1980, Geman and Geman 1984)  that consider local spatial models and that allow the incorporation of temporal information in a simple manner.  Other possibilities are autoregressive spatial models (Ripley, 1988), general statistical (e.g., Bayessian non-gaussian) models with (a priori) distributions estimates based on f(r), etc. Note that the availability of image f(r) allows for the construction of more detailed and restrictive models (linear and non-linear) that will act as binding conditions in the estimation of j(r). Note also that while the model for j(r)=j(r,Q) can be non-linear with respect to the model parameters Q,  the estimation of j can still remain linear. 

b) Direct linear inverse estimation. 

The general solution of equation I, can be written for the linear case as (Grave and Gonzalez, 1998),: 

(3)

Where W_d and W_j are symmetric positive definite matrices representing the metrics associated with the measurement space and the source space respectively. Vector j_p denotes any a priori value of the unknown function and denotes the regularization parameter. 

The functional information in f(r) can be directly incorporated into Equation 3 to produce an estimation of j(r). For example, if the metric is selected using the  linear probabilistic approach where W_j is interpreted as the covariance matrix of the sources (Dale and Sereno 1993), such matrix can be computed from the fMRI data as in Babiloni et al. 1999a. 

The matrix W_j can be also directly interpreted in terms of a metric which incorporates specific features on the source space. One of the simplest examples of this approach is the lead field column scaling that leads to the weighted minimum norm solution. Still, one could be interested in combining this column scaling W_j1 with an additional diagonal matrix W_j2 that weights the points in the solution space according to the functional image. Assuming that both matrices are diagonal and taking into account that when =0 an scaling factor in W_j does not influence the solution, then different scaling strategies can be combined normalizing the weighting matrices in the following way: 

Wj = {  Wj1 /max (Wj1) } * {  Wj2 / max (Wj2) }

Where max (W) denotes the maximum element in matrix W. The resulting weighting matrix with elements scaled to be no bigger than one can be used in (III) to get a weighted minimum norm estimate of j.  Note that this approach can include any arbitrary combination of diagonal weighting matrices. Other alternatives to constraint the inverse solution based on the fMRI data have been reported in Dale et al. 2000 and Babiloni et al 1999b among others. 

Another appealing strategy to incorporate the functional image f into the solution j is the use of space varying regularization methods. In this approach, a different regularization strategy can be locally designed on the basis of some local spatial properties of the image f(r) (e.g., the local variance). 

2.2 Some studies about the coupling between functional and electrophysiological images.

 The basic rationale behind the idea of integrating functional information into the inverse solution is that functional information can help to overcome the lack of uniqueness inherent to this problem. This means that functional images will orient the inverse algorithm to select one of the infinite possible solutions by providing a reliable independent a priori information about some features of the sources. It is self evident that the reliability of the inverse reconstruction will depend upon the credibility of the a priori information. Also, if different sources of a priori information are combined, the specific weight given to each source should be proportional to the level of confidence we have on this information. 

In what follows we describe some results of studies validating the localization of functional regions provided by fMRI against the gold standard of invasive electrophysiological studies. This cannot be considered an overview on this topic and we would like to remark that there are many studies considering this topic with contradictory results. Our purpose is to show that in spite of the recognized high spatial resolution of the fMRI, its functional localization results are not always coincident with the ones obtained with electrophysiological methods. This is a key aspect to consider since in the solution of the NIP, the purpose is to search for electrophysiological generators that behave according to electrodynamics laws and which reflect the electrochemical processes which are at the origin of the electromagnetic fields measured on the scalp. 

Blanke at. Al., (2000): These authors compared the localization of the frontal eye field (FEF) obtained using electrical cortical stimulation in six epileptic patients with that reported in the literature using fMRI. After normalization to the Tailarach atlas these authors report differences between both modalities of up to 3 centimeters. Electrophysiological responses inducing motor or sensory effects were found within the area functionally defined as the FEF. 

Stippich et al. (1999): These authors compared the localization provided by (fMRI) with the one obtained with dipole localization techniques applied to MEG data in six subjects during self-paced finger movement performance, tactile somatosensory stimulation and binaural auditory stimulation using identical stimulation paradigms. The mean distances found in this study between fMRI activity and the corresponding MEG dipoles were 10.1 mm (motor), 10.7 mm (somatosensory), 13.5 mm (auditory right hemisphere) and 14.3 mm (auditory left hemisphere). They concluded that the differences found may reflect the different underlying substrates of neurophysiology measured by fMRI and MEG. It should be noted that this study was carried out using an spherical head model for the source localization. The reconstructed sources were a posteriori matched onto the individual anatomical MRI. It is difficult to decide in such case which part of the differences could be attributed to this matching since a simple co-registration between a real brain and a sphere can not be considered as an anatomical constraint. 

Castellano-Smith (2000): This Ph.D. thesis carefully reviews and dissects the literature that compares functional localization with electrocortical stimulation of the sensorimotor cortex in epileptic patients. The main conclusion in this study is that for this group of epilepsy surgery patients with lesions in or near the sensorimotor cortex the fMRI cannot be considered a reliable tool for localizing the eloquent regions of the cortex. They partially attribute their results to the difficulty experienced by epileptic patients in remaining still in 

the MR scanner during the fMRI acquisitions and also to deformation of the brain surface during surgery to implant the electrode mats used for electrophysiological mapping. 

Disbrow et al. (2000): These authors developed an animal model appropriate for the study of the relationship between bold responses and electrophysiological events. They present a study of cortical maps generated using both fMRI and electrophysiological methods in the same animals under identical stimulus conditions for the topography of somatosensory areas 3a, 3b, 1 and 2 located in  anterior parietal cortex. fMRI and electrophysiologically defined maps were considered concordant if the centroids of the fMRI volume of activation fell within the electrophysiologically defined map. With this definition they found a concordance rate between the fMRI and electrophysiological maps of  55%. In the other 45% variability was     highest in the anterior-posterior plane, perpendicular to large local vessels. In this A-P plane, the centroid of activated pixels,  defined by fMRI, was an average of 8.6mm (SD=2.8mm) from the center of the electrophysiologically defined map. In the slice plane (superior-inferior) and the medial-lateral plane the centroids fell within the electrophysiologically defined maps. The area of fMRI activation was larger than the electrophysiological map, and was inversely related to anesthetic concentration. They conclude that the disparity between maps may be attributed to the hemodynamic source of the fMRI signal, which is only an indirect correlate of neuronal activity. They also point out that “fMRI, as typically performed, should be correlated with neurophysiology with caution”.

These four studies report discrepancies between fMRI and electrophysiological results for different brain areas, which are in the order of a few millimeters up to several centimeters. Let’s briefly see which is, in average, the spatial accuracy reported for inverse solutions that do not use functional constraints. 

3. About the Accuracy of the Solutions to the Neuroelectromagnetic Inverse Problem

Electroencephalography (EEG) and Magnetoencephalography (MEG) have been validated as non-invasive methods to study the functional principles of the human brain. Information on the strength and the localization of the neural activity can be derived from the distribution of these fields on the scalp surface by applying source localization algorithms. Many different source localization methods have been developed or applied in the last years (see e.g., Fuchs et al., 1999, Mosher et al., 1999; Grave de Peralta et al., 1997; Grave de Peralta et al., 2000; Goronidtsky and Rao, 1997; Sekihara and Scholtz, 1998). Also, source and head models have considerably evolved by incorporating more detailed anatomical and physiological information (Yan et al., 1991; Yvert et al., 1995).  At the present stage of development and using reasonable models, the localization accuracy reported for dipolar inverse solutions is smaller than 13 millimeters (Leahy et al., 1999). Unfortunately, similar phantom studies dedicated to evaluate the spatial accuracy that can be reached with distributed solutions are scarce. 

A reasonable alternative to experimentally evaluate distributed source models is to compare the outcome of the localization procedure with the localization established by invasive electrophysiological techniques in epileptic patients. The major difficulty for such studies is that simultaneous recordings of intracranial and extracranial potentials is a technically difficult problem. Thus, most of the studies have to be confined to compare the localization results obtained for presurgically recorded data with the position of the resected area in patients that are totally or partially seizure free after surgery. A major problem is that the size of the area removed during surgery is usually larger than the extent of the activated area as detected by the inverse solution. 

Only for illustration purposes we present here an example of the localization results that can be obtained for epileptic patients using a distributed solution (ELECTRA, Grave de Peralta et al., 2000) that restricts the source model to the type of currents physically capable of generate the measured maps. This example also constraints the feasible solution space to the gray matter detected by a semiautomatic segmentation procedure applied to the high resolution anatomical MRI of the patient. The patient is a 18 years old female which had a seizure frequency of 5-12/week, often followed by secondary generalization. The habitual seizures began with an impairment of consciousness and were followed by manual automatisms and a rightward deviation of the eyes and the head. Presurgical evaluation, including continuos video-EEG recording, nuclear imaging, and neuropsychological testing indicated left frontal epilepsy. Invasive monitoring was demanded in order to precisely localize the epileptic focus and to differentiate it from eloquent cortex. However, functional mapping resulted in speech arrest in close proximity to the epileptogenic focus and not complete resection of the left frontal lobe could be carried out. During a follow-up period of 6 months, five habitual seizures occurred, corresponding to a marked reduction of her preoperative seizure frequency. More than 12 seizures were recorded in this patient using 28 electrodes (sampling rate 128 Hz). Spikes selected by a specialist were aligned and averaged as shown in Figure 1. 
 

Figure 1. Spikes selected by a specialist, aligned and averaged

Figure 2. Localization results for the times marked as green vertical lines. Maxima are encircled in red and minima in blue.

 

Figure 3. The two slices with more significant activity found after the analysis of the EEG spike triggered fMRI.

There are several conclusions that can be drawn from this simple visual comparison. First of all, the distributed solution produces quite focal results suggesting the left frontal lobe as the epileptogenic site. This result is concordant with the invasive electrocorticography studies. Also the propagation of the epileptogenic activity observed in the solution coincides with the electrocorticography findings. In contrast to this focalized sequence of maps, the functional image shows highly focal activity at many spatially separated brain areas. There are no clues in the functional images to prefer one of these sites as the main epileptogenic focus.  It is important to mention, that this is not an isolated example but a consistent finding over a larger population of patients (Gorantz et al., submitted). In examples like this one, we see no reason to expect a substantial enhancement of the inverse solution by introducing the functional constraint. 

4. Discussion

The use of functional images as constraints for the solution of the electroencephalographic inverse problem relies upon the assumption that both imaging techniques are closely linked. Otherwise, using functional images as hard constraints for the inverse solution could mislead the inversion procedure. The term hard constraints, refers here to inverse solutions that force the electromagnetic solution to agree strictly with the functional data. The studies described in this paper suggest that the substantiation of a tight link between neuronal processes of interest and BOLD responses remains on shaky foundations. The fMRI bold signal is sensitive to parameters reflecting energy consumption, in particular to the cerebral rate of oxygen metabolism and the cerebral rate of blood flow (Ogawa et al., 1998) While the neuronal activities of interest are those involved in the communication of information between neurons, the brain consumes energy for many more processes which are not directly linked to it (Rothmann et al., 1999). Neurotransmitter release and uptake, vescicular recycling and maintenance of membrane potentials are examples of processes which consume energy (Shepherd, 1994). Usually neglected, glia also require energy which might explain the existence of activation often found in the fMRI in the white matter, which is usually excluded from the feasible inverse solution space in anatomically constrained head models. 

There are additional practical limitations to the integration fMRI inverse solutions. The first aspect is that the volume of the head covered by the fMRI and the inverse solution space used by realistic head models, which is selected from anatomical images, is not necessarily the same. Quite often and due to technical limitations, the fMRI is confined to preselected slices while the inverse solutions intends to cover as much as possible the gray matter detected in high resolution MRI images. Although this technical limitation of the fMRI can be circumvented in the nearby future the differences in the temporal resolution attainable by both techniques will certainly subsist. Time resolution of fMRI will remain low even when the technique becomes more advanced because changes in cerebral blood oxygenation occur at a slower time scale compared with relevant neuronal events that may only take milliseconds. Thus, the same static functional image has to serve as a constraint for a rather large set of scalp maps which can be very dissimilar in topography. 

The last years have not only lead to important developments in the field of functional imaging but also in the field of inverse solutions. The source and head models currently in use are more sophisticated and accurate. The basic limitations of these models have been described (Grave and Gonzalez, 1998) and there is today a more clear picture of what are the limits of these techniques.  The introduction of reasonable constraints in the spatial and temporal features of the generators have allowed in the last years to obtain interesting results in the analysis of event related data and epileptic activity. Also there are some recent proposals of techniques to improve and assess the reliability of the distributed solution maps. One of this approaches relies on the automatically isolation of scalp potential maps which are simple enough to expect reasonable results after applying a distributed source localization procedure. The isolation technique is based on the time frequency decomposition of the scalp measured data by means of the Short Time Fourier Transform (STFT). The basic rationale behind the approach is that neural generators synchronize during short time periods over given frequency bands for the codification of information and its transmission. Consequently potential patterns specific for certain time frequency pairs should be simpler than those appearing at single times but for all frequencies. The method considers the general case of distributed source models with non-stationary time behavior. 

At this stage, the essential question to be answered is whether the introduction of functional constraints can help to further increase the accuracy and reliability of inverse solutions. Considering the above mentioned difficulties and particularly the problems to totally correlate energy consumption and neuronal activity we believe that no improvements can be obtained by using the functional image as a hard constraint, i.e., as constraint that needs to be fulfilled. In practical terms this means, that a model that assigns a zero weight to a point non activated in the fMRI is more restrictive and dangerous than one assigning a lower probability of activation at the same point. While in the first case the activity at the point is banned in latter one is simply penalized. 

5. Conclusions

In this paper we considered the problem of introducing functional images as constraints into the solution of the neuroelectromagnetic inverse problem. In the first part we discussed some of the alternative manners to state this integration. Some studies about the relationship between functional and electrical activation were included. An example of the accuracy that can be obtained in the localization of focal sources with anatomically and physically constrained inverse solutions was presented.  In our opinion these results do not limit the recognized value of fMRI to study brain function, but suggest that the integration fMRI-inverse solutions needs to be applied and rigorously tested before concrete conclusions can be drawn as to its utility. To summarize we consider as the two more reasonable alternatives at this moment: 1) to combine EEG and fMRI and perform the same investigation with both methods, to cross-validate assumptions in either the time or spatial domain of each method (Ahlfors, et al., 1999; George et al, 1996) and 2) to introduce functional images as soft constraints in the inverse solutions by using a parametric approach that relates both modalities in terms of certain parameters. Latter approach may be more feasible in future once the issues of concern have been addressed and more specific models to explain the coupling are developed.

Acknowledgments

Work supported by the Programme commun de recherche en genie biomedical 1999-2002 and the Swiss National Foundation Grants 2053-059341.99/1 and 20-59341.99.