Removal of FMRI environment artifacts from EEG data using optimal basis sets

Abstract

The combination of functional magnetic resonance imaging (FMRI) and electroencephalography (EEG) has received much recent attention, since it potentially offers a new tool for neuroscientists that makes simultaneous use of the strengths of the two modalities. However, EEG data collected in such experiments suffer from two kinds of artifact. First, gradient artifacts are caused by the switching of magnetic gradients during FMRI. Second, ballistocardiographic (BCG) artifacts related to cardiac activities further contaminate the EEG data. Here we present new methods to remove both kinds of artifact. The methods are based primarily on the idea that temporal variations in the artifacts can be captured by performing temporal principal component analysis (PCA), which leads to the identification of a set of basis functions which describe the temporal variations in the artifacts. These basis functions are then fitted to, and subtracted from, EEG data to produce artifact-free results. In addition, we also describe a robust algorithm for the accurate detection of heart beat peaks from poor quality electrocardiographic (ECG) data that are collected for the purpose of BCG artifact removal. The methods are tested and are shown to give superior results to existing methods. The methods also demonstrate the feasibility of simultaneous EEG/FMRI experiments using the relatively low EEG sampling frequency of 2048 Hz.

Introduction

Functional neuroimaging techniques such as positron emission tomography (PET) and functional magnetic resonance imaging (FMRI) have made it possible to map specific areas of the brain that are involved in carrying out tasks of differing complexities. However, such techniques measure brain activity indirectly, i.e., they detect a secondary effect of brain activation rather than the neural activity itself. Furthermore, although such imaging modalities provide good spatial localization, they suffer from relatively low temporal sampling frequency. Electrophysiological mapping methods such as electroencephalography (EEG) and magnetoencephalography (MEG) measure brain electrical activity directly, and in real time, but suffer from spatial blurring of activation. This has encouraged the neuroimaging community to investigate multimodal imaging in order to combine the strengths of the individual techniques. In particular, the simultaneous combination of FMRI and EEG has received recent attention. Initial studies aimed to demonstrate the safety, potential and data quality possible using this technique (Ives et al., 1993, Lemieux et al., 1997, Krakow et al., 2000, Goldman et al., 2000). Subsequently, the use of simultaneous EEG and FMRI has been used to study the generators of the alpha rhythm (Goldman et al., 2002, Laufs et al., 2003, Moosmann et al., 2003, Niazy et al., 2004), event-related brain responses (Bonmassar et al., 1999, Kruggel et al., 2000, Liebenthal et al., 2003), brain activation during different sleep stages (Czisch et al., 2002, Liebenthal et al., 2003) and epileptic activities (Seeck et al., 2001, Krakow et al., 2001, Lemieux et al., 2001, Iannetti et al., 2002, Bénar et al., 2003).

Common to all simultaneous EEG and FMRI experiments, however, are the MRI environment artifacts that contaminate EEG data. There are two kinds of MRI environment artifact. The first kind is the gradient (or imaging) artifact caused by the switching of the magnetic field gradients (Felblinger et al., 1999, Allen et al., 2000). The second kind is the ballistocardiographic (BCG) artifact, which is caused by heart-related blood and electrode movements inside the static magnetic field of the MRI scanner (Allen et al., 1998, Bonmassar et al., 2002), regardless of whether or not MR scanning is being performed. In this paper, we briefly discuss these artifacts, review some of the methods previously proposed to remove them and then propose new methods for their removal. Additionally, a companion paper by Iannetti et al. (this issue) further validates our methods by applying them to simultaneous laser-evoked potentials (LEPs) and FMRI experiments and demonstrating the quality of the acquired data. In the remainder of this section, more detailed information is given about the different kinds of artifact, current methods for their removal, as well as a summary of our proposed methods. Following the Introduction, the Methods section details our proposed algorithms. This is followed by the Validation section, which will describe how the algorithms were tested. Results and Discussion will then follow.

During MR imaging, the magnetic field inside the MRI scanner continuously changes as a result of the switching of the magnetic field gradients. The gradients change according to the imaging sequence being used. In an echo planar imaging (EPI) sequence typically used in FMRI, gradient switching is repeated each time a new slice is collected, resulting in artifacts that repeat with the collection of each new FMRI slice. The amplitude of such artifact can be 100 times greater than the EEG signal and its frequency content overlaps that of the EEG, thus gradient artifacts cannot be simply filtered out. The artifact shape and amplitude varies from one EEG channel to another depending on the location of the electrodes and the wire connections. For more details, see Hoffmann et al., 2000, Anami et al., 2003, Garreffa et al., 2003.

Different approaches have been proposed to remove gradient artifacts from biological signals collected during MRI scanning. Hoffmann et al. (2000) proposed a frequency domain method, where the amplitude and phase of the data were set to zero at frequencies matching an artifact power spectrum template. However, this approach suffers from the typical ‘ringing’ effect common to such frequency domain filters (Bénar et al., 2003). The most used method is average artifact subtraction (Allen et al., 2000, Bénar et al., 2003). This utilizes the repetitiveness of the artifact to form an average artifact template, which is then subtracted from the EEG data. The efficacy of this approach has been demonstrated in the literature (Allen et al., 2000, Bénar et al., 2003, Salek-Haddadi et al., 2003), though a number of quality and practicality issues still remain. Firstly, some residual artifacts remain on some channels. Allen et al. (2000) proposed the use of adaptive noise cancellation (ANC) to remove these residuals; however, this approach does not remove all residual artifacts. Secondly, in order to minimize the residuals, a high sampling frequency is needed. From our experience, some unsatisfactory results are obtained from commercial implementation of this algorithm even at sampling rates of 10 kHz. However, even if better quality data were to be achieved at such high sampling rates, the amount of generated data (especially in high electrode density experiments) limits the length of the experiments and causes practical problems when the data need to be analyzed using third party software such as Matlab® (The MathWorks, Inc., MA, USA).

Fig. 1 shows an example of the gradient artifacts in EEG data and their origin. After subtraction of the average artifact, some residuals remain, which result in sharp deflection in the data after low-pass filtering. This is mainly due to the inaccuracy of the artifact template being subtracted. This inaccuracy is partly caused by the fact that the MRI machine and the EEG system are typically driven by separate clocks, which means that the artifact is not always sampled at exactly the same location. This introduces a slight variation in the shape of the artifact from one slice to another. Also, the average template calculation and subtraction processes are usually dependent on the triggers received from the MRI machine that indicate the start of each artifact segment. The location of these triggers is often slightly inconsistent from one segment to another, which causes a temporal jitter in the onset of the different artifact segments and thus degrades the accuracy of the calculated template. These problems are exacerbated if the sampling frequency is decreased.

In this paper, we propose a new method for the removal of gradient artifacts; FMRI artifact slice template removal (FASTR). In FASTR, a unique artifact template for each slice artifact in each EEG channel is constructed and then subtracted. Each slice template is constructed as the local moving average plus a linear combination of basis functions that describe the variation of residuals. The basis functions are derived by performing temporal principal component analysis (PCA) on the artifact residuals and selecting the dominant components to serve as a basis set. This technique is demonstrated to be superior to imaging artifact reduction (IAR) (Allen et al., 2000) and applicable at a sampling rate as low as 2048 Hz.

Recent independent work by Negishi et al. (2004) has proposed a similar approach to ours. In contrast to our approach, Negishi's work requires the collection of extra EEG data without scanning to serves as a reference. This could well introduce a number of problems, as it assumes the two data sets differ only in the introduction of the gradient artifacts. However, their method utilizes all estimated principal components (rather than the strongest few) and has the advantage of automatically weighting each component depending on its projection on both the contaminated and clean data; Negishi's use of PCA is fairly different and, to an extent, complementary to what we have developed. More discussion of the differences between the two approaches is given later.

The BCG artifact is a distortion in the EEG data caused by cardiac-related activities. In a normal, ideal EEG environment it is usually caused by an electrode being directly above a pulsating scalp vessel, and the problem can be avoided by changing the electrode position. Inside the MRI magnet, this problem is greatly magnified. Causes and characteristics of the BCG artifact have been described in the literature (Allen et al., 1998, Bonmassar et al., 2002), and in general it is caused by electrode movement due to pulsatile scalp and blood movement related to the cardiac cycle. This movement of electrodes and conductive blood inside the magnetic field induces the artifacts (Allen et al., 1998). The magnitude of the BCG artifact may be as much as 200 μV at 1.5 T (3–4 times that of the EEG) (Allen et al., 1998), it is spread throughout the heart beat period and it can be observed across the scalp (Allen et al., 1998), although its magnitude and shape can vary considerably from one EEG channel to another. In contrast to the imaging artifact, although the basic shape of the BCG artifact is similar from one occurrence to the next in any single EEG channel, there exists considerable variation in the artifact shape, amplitude and scale over time.

Several approaches have been proposed to remove the BCG artifact. Adaptive filtering has been proposed by Bonmassar et al. (2002): a piezoelectric sensor was used to generate a reference BCG signal, which was then used to filter out BCG contributions from the EEG. This method is computationally expensive, requires the use of an extra sensor and assumes that no EEG correlated information is present in the sensor signal. Spatial PCA and independent component analysis (ICA) filters have also been proposed (Bénar et al., 2003, Srivastava et al., 2005). One problem with these approaches is that they necessitate the presence of a large number of sensors. Also, the identification of artifact components can be subjective and is usually done manually. Most importantly, spatial filters assume that all the sensors are contaminated by common sources, which is not the case. The BCG artifact derives from sources that are rotating/moving, which contaminate different sensors at different points during the cardiac cycle with different effects. The most commonly used method for removing the BCG artifact is the average artifact subtraction (AAS) (Allen et al., 1998), in which a moving average artifact template is computed from successive artifact occurrences, then subtracted from the data. This assumes that the BCG artifact is a slowly changing signal that can be accurately captured by a moving average. This can result in residual artifacts in the data. A variation of this approach was used by Goldman et al. (2000), where instead of computing a simple moving average, a weighted average is used such that artifacts that lie further from the one being processed are less emphasized. Sijbers et al. (2000) and Ellingson et al. (2004) used a median filter to construct a template which is then scaled in time and amplitude to fit each artifact instance. Again, all these approaches assume a temporal relationship between the different occurrences of the artifact. Another central issue to such subtraction-based methods is the accurate detection of heartbeat locations. EEG systems often provide limited ECG recording facilities , e.g., one single bipolar channel. In addition, the ECG is usually distorted inside the MR machine due to blood conductivity (Wendt et al., 1988). These factors can lead to inaccurate detection of QRS peaks in the ECG, especially when simple thresholding detection methods are used (Allen et al., 1998).

Similar to the removal of the residual gradient artifacts, we propose a method where a basis set is constructed by performing temporal PCA on each EEG channel data. The basis set is then fitted to, and subtracted from, each artifact occurrence. This approach has the advantage of not assuming any temporal relation between the different occurrences of the BCG artifact in a given EEG channel. Rather, the assumption is that over a sufficient period of EEG recording from any single EEG channel, the different BCG artifact occurrences in that channel are all sampled from a constant pool of possible shapes, amplitudes and scales. The principal components of all the occurrences can then describe most of the variations of the BCG artifact in that channel. This method is shown to be superior to AAS (Allen et al., 1998). In addition, an accurate, robust procedure based on the work of Christov (2004) followed by a correction algorithm is proposed for the accurate detection of QRS complexes in ECG data collected inside the magnet. The work of Negishi et al. (2004) also sheds some light on the usability of temporal PCA to remove BCG artifacts. However, they affirm that their results were somewhat unsatisfactory and more work was needed. More discussion about this is given later.