Edited by: Mark Holmes, University of Washington, USA
Reviewed by: Jose F. Tellez-Zenteno, University of Saskatchewan, Canada; Don Tucker, Electrical Geodesics, Inc. and the University of Oregon, USA
*Correspondence: Patrícia Figueiredo, Department of Bioengineering, Institute for Systems and Robotics, Instituto Superior Técnico, Av. Rovisco Pais, 1, 1049-001 Lisboa, Portugal. e-mail:
This article was submitted to Frontiers in Epilepsy, a specialty of Frontiers in Neurology.
This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in other forums, provided the original authors and source are credited and subject to any copyright notices concerning any third-party graphics etc.
Simultaneous electroencephalogram (EEG)-functional Magnetic Resonance Imaging (fMRI) recordings have seen growing application in the evaluation of epilepsy, namely in the characterization of brain networks related to epileptic activity. In EEG-correlated fMRI studies, epileptic events are usually described as boxcar signals based on the timing information retrieved from the EEG, and subsequently convolved with a hemodynamic response function to model the associated Blood Oxygen Level Dependent (BOLD) changes. Although more flexible approaches may allow a higher degree of complexity for the hemodynamics, the issue of how to model these dynamics based on the EEG remains an open question. In this work, a new methodology for the integration of simultaneous EEG-fMRI data in epilepsy is proposed, which incorporates a transfer function from the EEG to the BOLD signal. Independent component analysis of the EEG is performed, and a number of metrics expressing different models of the EEG-BOLD transfer function are extracted from the resulting time courses. These metrics are then used to predict the fMRI data and to identify brain areas associated with the EEG epileptic activity. The methodology was tested on both ictal and interictal EEG-fMRI recordings from one patient with a hypothalamic hamartoma. When compared to the conventional analysis approach, plausible, consistent, and more significant activations were obtained. Importantly, frequency-weighted EEG metrics yielded superior results than those weighted solely on the EEG power, which comes in agreement with previous literature. Reproducibility, specificity, and sensitivity should be addressed in an extended group of patients in order to further validate the proposed methodology and generalize the presented proof of concept.
Over the years, the electroencephalogram (EEG) has been the tool of choice for the diagnosis and characterization of epilepsy. With the possibility to acquire the EEG simultaneously with functional Magnetic Resonance Imaging (fMRI), studies of Blood Oxygen Level Dependent (BOLD) signals correlated with epileptic activity proliferated (Ives et al.,
Although a great amount of both experimental and theoretical work has been dedicated to the clarification of the relationship between neural activity and associated hemodynamic changes, neurovascular coupling mechanisms remain an active area of research (Rosa et al.,
In most epilepsy EEG-fMRI studies, the goal is to find brain networks exhibiting hemodynamic changes associated with interictal and/or ictal activity, which are expected to be correlated with the epileptogenic areas (Hoffmann et al.,
In this paper, we address the issue of modeling the BOLD dynamics associated with epileptic EEG activity by extracting metrics from the EEG spectrum expressing different models of the transfer function between neuronal and hemodynamic signals. A methodology is proposed for the analysis of EEG-fMRI data in epilepsy, consisting of ICA decomposition of the EEG followed by component selection based the reproducibility across different acquisition runs, Morlet wavelet spectral analysis, and EEG metric extraction. The resulting time courses are convolved with a canonical HRF and used as regressors of interest in a GLM analysis of the fMRI data. The proposed methodology is applied to the study of a patient with epilepsy associated with a hypothalamic hamartoma. The different EEG-fMRI transfer functions are compared with each other, as well as with a conventional GLM methodology based on the identification of ictal and interictal events on the EEG by the neurophysiologist, and also with an fMRI-ICA approach.
We focused on the simultaneous EEG-fMRI data recorded from a 2-year-old patient with a giant hypothalamic hamartoma suffering from gelastic epilepsy, as part of the pre-surgical evaluation under the Program for Epilepsy Surgery of the Hospital Center of West Lisbon. This case study has been previously described (Leal et al.,
Seizures occurred more than 50 times per day and typically lasted for 20–30 s, involving almost exclusively the left hemisphere. The clinical manifestations of the seizures consisted of slowing of motor activity, variable interruption of consciousness, eyelid rhythmic movements with bilateral nystagmus to the right, and occasionally gelastic laughter. After the acquisition of the EEG-fMRI data, the patient underwent a two-stage hamartoma disconnection surgery, 1 year after which the seizures were reduced to 1–3 per day.
The EEG interictal activity demonstrated a persistent slow-wave abnormality over the left temporal-occipital area, associated with abundant spike activity with occasional contralateral propagation. Abundant spike activity also occurred over the left hemisphere frontal lobe. Topological analysis of the interictal spikes presented a spatial stationarity for the posterior spikes, whereas the frontal ones changed significantly in configuration from an occipital dipolar potential at spike onset to a dipolar frontal potential at spike peak. The ictal EEG pattern was very monotonous and consisted of early diffuse desynchronization, followed by the build-up of spike activity over the left occipital and temporal areas and, in the later stages of the seizure, over the frontal area. Occasionally secondary propagation of spike activity to the right temporal areas occurred.
The EEG was recorded using an MR-compatible 37-channel system (Maglink, Neuroscan, Charlotte, NC, USA) with two ECG channels, using a sampling rate of 1000 Hz with a bandwidth DC-250 Hz and reference at electrode FCz. The imaging was performed on a 1.5 T MRI scanner (GE Cvi/NVi). Six fMRI runs were collected using a gradient-echo echoplanar imaging (EPI) sequence, with TR = 2.275 s, 3.75 mm × 3.75 mm × 5.00 mm voxel size and a total of 154 volumes (the first four were subsequently rejected in order to discard T1 relaxation unstationarities). A T1-weighted structural image was also acquired using a 3D spoiled gradient recovery (SPGR) sequence, with 0.94 mm × 0.94 mm in-plane resolution and 0.6 mm slice thickness. During the scanning session, the patient was administered light anesthesia with Sevoflurane at 1% (Abbott Laboratories, Abbot Park, IL, USA), through mask, as established by the MRI protocol for small children and uncooperative patients. A 5 min EEG recording was also performed outside the scanner and before anesthesia, in order to allow for cross-validation of the simultaneous EEG-fMRI recordings in terms of the presence of MR related artifacts and the effects of anesthesia. Periods of interictal and ictal activity were identified on the artifact-corrected EEG by the neurophysiologist. For runs 2 and 3, two and five ictal events occurred, respectively. For the remainder runs, only interictal spikes were detected on the EEG. The visual inspection of the EEG inside and outside the scanner revealed a slight increase in both beta and theta background activity under anesthesia, but there was no significant change in the morphology of interictal spikes.
The EEG pre-processing was executed using the EEGLAB toolbox (Delorme and Makeig,
In an attempt to separate out the activity of interest in the EEG, the pre-processed data were decomposed by ICA using the
A reproducibility analysis of the ICs was performed in order to identify the associated topographies that were consistent across the six acquisition runs.
The spectral profiles of the selected ICs were obtained by time-frequency analysis through convolution of the respective time courses with Morlet wavelets. The subset of metrics found to be more relevant in Rosa et al. (
The difference between frequency averaging measures (RMSF or MF) was found to have a negligible effect on BOLD signal prediction; hence the results emerging from the MF and uMF metrics will be omitted, as they were not significantly different from those of the RMSF and uRMSF metrics, respectively.
The fMRI data were analyzed using FSL
Pre-processing consisted on: motion correction using MCFLIRT (Jenkinson et al.,
Two GLM analyses were specified in order to identify BOLD signal changes associated with: (1) the
An ICA of each fMRI run was also performed using Probabilistic ICA as implemented in MELODIC (Multivariate Exploratory Linear Decomposition into Independent Components) Version 3.09, part of FSL (FMRIB’s Software Library (see text footnote 1); Beckmann and Smith,
The comparison between EEG metrics was performed through GLM analysis followed by inferences based on
The comparison between consistent IC groups of interest was performed in a similar way, for each acquisition run, applying the selected EEG metric. Here, we considered a single GLM comprising the three consistent IC groups of interest (I, II, and V, as will be shown in the Results) with three regressors of interest each (canonical HRF, time, and dispersion derivatives), totaling nine regressors of interest, in one three-way comparison.
In order to assess the plausibility and consistency of the EEG metric-derived BOLD maps, the EEG metric/IC group combination yielding the largest maps (and hence best at predicting BOLD signal changes) was selected for subsequent comparison with the GLM analysis using the EEG EA regressors and also with the ICA analyses. The consistency between each two regressors was measured as their temporal correlation. The consistency between two maps was measured as their spatial overlap, i.e., the ratio of the volume of the map intersection with the volume of the map union.
In this section, the results obtained through the proposed methodology will be presented.
An exerpt of the EEG trace obtained after pre-processing is presented in Figure
The consistent IC groups, obtained as a function of the reproducibility of the associated IC topographies across runs, are presented in Figure
The spectral profiles obtained with Morlet wavelet decomposition for one IC (IC10) in run 3, as well as the EEG channel with highest absolute weight for this component (P3), are shown in Figure
The BOLD statistical maps obtained using each EEG metric and consistent IC group were generally consistent with each other, across runs and also with the patient’s seizure semiology, but differed considerably in terms of the number of voxels showing statistically significant EEG metric-related BOLD changes. Firstly, the results of the comparison between the fMRI analysis using the different EEG metrics and consistent IC groups will be presented. A single metric and IC group will then be selected for comparison with the EA GLM and fMRI-ICA analyses.
The three-way analysis did not show meaningful activations to allow for the comparison of the individual EEG metrics; however, the two-way comparisons yielded a clear superiority of the RMSF metric when compared to the TP and uRMSF metrics, for every IC group and run, in terms of the number of voxels showing statistically significant EEG metric-related BOLD changes, as shown in Figure
The three-way comparison regarding the IC group showed a significant superiority of IC group II, in terms of the number of voxels exhibiting statistically significant EEG metric-related BOLD changes, as shown in Figure
The BOLD statistical maps obtained using the RMSF metric for IC groups I and V are shown in Figures
The comparison between the fMRI results obtained using the EEG RMSF metric and the EA GLM analysis, for IC group II, are summarized in Figure
In general, the maps obtained using the EEG metric approach were consistent across runs and also with the patient’s seizure semiology. In fact, clusters were found on the left occipital/parietal and left frontal lobes, consistently with the observation on the EEG of spike buildups in left occipital/parietal channels followed by the left frontal channels, after a generalized desynchronization. In run 3, the EA approach was also capable of identifying these brain areas, but the corresponding statistical maps were more significant and extensive for the EEG metric approach. Furthermore, with the EEG metric approach, clusters involving left thalamic, left hippocampal, and left frontal ventral areas, as well as the hamartoma itself, were also found, which are generally consistent with the brain network associated with seizures originating in a hypothalamic hamartoma. For the remaining runs, in contrast with the EEG metric approach, the fMRI statistical maps obtained with the EA approach were in general not consistent across runs nor with the expected epileptic network.
The results regarding the EEG metric GLM analysis and the fMRI-ICA decomposition of all the runs, for IC group II, are presented in Figure
We proposed a new methodology for BOLD signal prediction in EEG-correlated fMRI studies in epilepsy, by incorporating a model of the EEG-BOLD transfer function. Specifically, independent components of the EEG associated with consistent topographies were translated into BOLD signal predictions by a set of model-based metrics. Interestingly, we found that increases in the MF of the EEG were better than its power at predicting BOLD increases, in support of the heuristic proposed in (Kilner et al.,
The heuristic proposed by Kilner and colleagues puts forward the MF (in an RMS sense) of the EEG spectrum as a surrogate of the neuronal activity eliciting the BOLD signal, following a broad physiological inspiration: the BOLD eliciting signal is assumed to be proportional to the electric work dissipated by the ionic currents across the cell membranes; this can in turn be shown to be proportional to the RMSF of the LFP/EEG spectrum if the covariance of the membrane potentials of the cells are assumed constant (Kilner et al.,
The EEG metric-related BOLD change maps were consistent with the ones obtained using the epileptic activity regressors defined by the neurophysiologist, whenever ictal activity was identified on the EEG. For the runs on which no ictal events were detected on the EEG, the proposed methodology was able to identify the same brain network that was involved in the ictal runs. This was achieved by regressors based on the same IC scalp topography as those from the ictal runs, suggesting the presence of underlying epileptic activity in the identified epileptic network, which only occasionally manifests itself with an ictal character. In contrast with the proposed EEG metric based approach, however, the conventional analysis of the interictal runs yielded inconsistent or no results for the same statistical significance threshold. These findings suggest that the interictal events detected on the EEG may not fully reflect the activity of the underlying brain network, while the selected EEG metrics may be more powerful in depicting it. The same network was also identified by fMRI-ICA on both ictal and interictal runs, which further supports this idea. However, the fully data-driven method of fMRI-ICA lacks an explanatory model for the data, in contrast with the proposed methodology, which is based on modeling the EEG-BOLD transfer function.
For the purpose of verifying the specific epileptic character of the identified brain networks, their consistency across runs and their plausibility as the epileptic brain network underlying the patient’s seizure semiology were considered. The overlap of the corresponding BOLD change maps across runs was very high, both for the EEG metrics GLM and the fMRI-ICA approaches. Moreover, these maps were in good agreement with the brain network known to be involved in the epileptic activity of this patient, which includes the hamartoma, as well as left hemisphere hypothalamus, hippocampus, parietal–occipital lobe, cingulate gyrus, and dorsal–lateral frontal lobe (Leal et al.,
An ICA decomposition of the EEG was used here with the purpose of separating out the activity of interest into a univariate timecourse, which was expected to exhibit a consistent and meaningful topography. ICA is a popular technique for the removal of muscular activity or eye movement artifacts in EEG data processing (Vigario,
Other approaches for the extraction of a univariate EEG time course, representative of the epileptic activity, have been proposed in the literature, namely continuous Electrical Source Imaging (Vulliemoz et al.,
A limitation of the proposed approach is the bias of the EEG measures toward superficial cortical activity, which possibly precludes the identification of hemodynamic changes associated with deep brain activity. This limitation is however common to all EEG-fMRI “integration through prediction” approaches. Nonetheless, the particular syndrome of epilepsy associated with hypothalamic hamartomas is rather well described in the literature, with clear evidence for the hypothalamic hamartoma being the seizure focus, and this region was in fact found in our work in five out of six EEG-fMRI datasets. Although the proposed EEG-BOLD transfer functions are not specific to epileptic activity, this specificity is achieved in the presented methodology by selecting the EEG topography used to extract the metric based on an ICA procedure. Nevertheless, a possible limitation of the proposed transfer functions is that they are not specific to the start of the seizure, and hence to the epileptogenic focus. However, our aim was the description of the full epileptic network, leaving the specific identification of the seizure focus and seizure propagation dynamics for other, related lines of research (Murta et al.,
The study presented here focused on data from a single patient with the aim of providing a proof of concept for the potential usefulness of the proposed methodology for the identification of epileptic networks. The choice of this patient was based on the quality of the EEG data that could be achieved due to the absence of movement associated with the beginning of the seizures. Moreover, epilepsy cases associated with hypothalamic hamartomas are relatively stereotypical in terms of the electrophysiological patterns of seizure propagation, which makes the interpretation of the results relatively more robust in comparison to other types of epilepsy. In general, the clinical utility of EEG-fMRI is yet to be established. Nevertheless, in this case study, the results of the EEG-fMRI investigation indicate possible alternative treatment approaches involving the surgical interruption of the seizure propagation pathways (Leal et al.,
In conclusion, we presented a new approach for EEG-fMRI integration in the field of epilepsy, which incorporates and tests different models of the transfer function between EEG and BOLD signals, hence allowing better predictions of the hemodynamic changes associated with epileptic activity. This work therefore provides a contribution to our understanding of the link between EEG and BOLD signals as well as for improving the yield of EEG-fMRI studies in epilepsy.
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
We acknowledge the Portuguese Science Foundation (FCT) for financial support through Project PTDC/SAU-BEB/65977/2006, Project PTDC/SAU-ENB/112294/2009, and Project PEst-OE/EEI/LA0009/2011.
1