A statistically robust EEG re-referencing procedure to mitigate reference effect

Abstract

Background: The electroencephalogram (EEG) remains the primary tool for diagnosis of abnormal brain activity in clinical neurology and for in vivo recordings of human neurophysiology in neuroscience research. In EEG data acquisition, voltage is measured at positions on the scalp with respect to a reference electrode. When this reference electrode responds to electrical activity or artifact all electrodes are affected. Successful analysis of EEG data often involves re-referencing procedures that modify the recorded traces and seek to minimize the impact of reference electrode activity upon functions of the original EEG recordings.

New method: We provide a novel, statistically robust procedure that adapts a robust maximum-likelihood type estimator to the problem of reference estimation, reduces the influence of neural activity from the re-referencing operation, and maintains good performance in a wide variety of empirical scenarios.

Results: The performance of the proposed and existing re-referencing procedures are validated in simulation and with examples of EEG recordings. To facilitate this comparison, channel-to-channel correlations are investigated theoretically and in simulation.

Comparison with existing methods: The proposed procedure avoids using data contaminated by neural signal and remains unbiased in recording scenarios where physical references, the common average reference (CAR) and the reference estimation standardization technique (REST) are not optimal.

Conclusion: The proposed procedure is simple, fast, and avoids the potential for substantial bias when analyzing low-density EEG data.

Fig. 1

Example synthetic re-referenced data illustrate the utility of rCAR. Simulation results are summarized in Section 3.2 and described in detail in Section 3.1. Left: simulated α rhythms in posterior electrodes. Right: simulated K-complex (a brief, oscillatory event). Top row: ideal voltage, s_t defined in Section 2.1. Remaining rows: competing re-referencing procedures. The proposed re-reference procedure (rCAR) (Section 2.4) out-performs all of the re-reference procedures in the posterior alpha rhythm and K-complex cases (see also Fig. 3). The proposed re-reference procedure (rCAR) does not “bleed” channel signal across channels. This feature of the rCAR procedure restricts large signals on some channels from affecting other channels and is, heuristically speaking, the essence of its success. It is demonstrated more clearly in Fig. 2, where the ideal EEG recordings are removed from the re-referenced data and the resulting difference is displayed.

Fig. 2

Example synthetic re-referenced data resulting from the simulations described in Section 3.1 reduced by the ideal channel recordings to highlight their differences. The simulations are identical to those used to produce Fig. 1. Ideally, these images would be zero for all times and channels. A departure from zero indicates an error in the re-referencing procedure. The across-channel “bleeding” effect described in the text and in Fig. 1 is apparent. The proposed re-referencing procedure (rCAR) is the most accurate. The vertical stripes result from the fact that the reference effects all of the channels equally in the data model, Eq. (1). When the estimate of the reference deviates from the actual reference the result is an equal error on all of the channels.

Fig. 3

The proposed robust referencing procedure (rCAR, red) out-performs the competing re-referencing procedures by avoiding badly biased estimates of the reference, r_t. Both plots: color indicates re-referencing procedure: C2 (green), CAR (blue), REST-C2 (magenta), REST-CAR (black), and the proposed re-referencing procedure, rCAR (red). Symbol indicates the simulation scenario: K-complex (x), alpha rhythm (right-oriented triangle), active half-cortex (upward-oriented triangle), and random cortical source activity over the whole cortex labeled background only, (downward-oriented triangle). Together these scenarios investigate the performance of the re-referencing procedures on increasingly non-focal cortical activity. Each mark indicates the result of a single simulation for a single re-referencing procedure. Forty simulations are performed for each simulation scenario. For a complete description of the simulations, including a description of the simulation scenarios, see Section 3.1. Left: actual sum of the absolute channel-to-channel correlation (computed from the synthetic silent-reference recordings) vs. the sum of the absolute channel-to-channel correlation estimated from the re-referenced synthetic data produced in simulation. The proposed robust reference procedure (rCAR, red) avoids badly biased correlation estimates computed using the competing references. Right: the sum of the absolute channel-to-channel correlation (vertical axis) tends to increase with increasing reference estimator error. The standard deviation normalized reference estimator error is estimated (horizontal axis) from synthetic data as the scaled sum of the squared deviations of the reference estimate,r̂_t, from the actual reference, r_t, divided by the product of the channel standard deviations, σ_kσ_k′. This quantity directly estimates $\frac{MS{E}_e}{{\sigma }_k{\sigma }_{k^{\prime }}}$ appearing in Eq. (C.11). As expected from Eq. (C.11), the sum of the absolute channel-to-channel correlations tends to increase with the scaled reference estimator error. This result, established in simulation for a known r_t, is used to establish a performance metric for the competing re-referencing procedures when re-referencing real data where r_t is not known. The sum of the absolute channel-to-channel correlations for actual data are plotted in Fig. 9. The proposed re-referencing procedure, rCAR, is always fairly accurate and does not exhibit the large bias demonstrated by the competing re-referencing procedures in at least one simulation scenario. (For interpretation of the references to color in this legend, the reader is referred to the web version of the article.)

Fig. 4

The simulation results for the proposed robust re-reference procedure (rCAR) presented in Fig. 3, reproduced for clarity. The proposed robust referencing procedure yields unbiased estimates of the channel-to-channel correlation in the four simulation scenarios.

Fig. 5

The periodogram of the CAR and rCAR re-referenced data. The periodogram of rCAR is more similar to the periodogram of the ideal silent-reference recordings than the periodogram of the CAR. In particular, the low-frequency streaking across the channels present in the CAR re-referenced data is absent in the rCAR re-referenced data.

Fig. 6

The phase of the discrete Fourier transformed CAR and rCAR re-referenced data. The phase of rCAR is more similar to the phase of the ideal silent-reference recordings than the phase of the CAR recordings.

Fig. 7

The 19-channel EEG data recorded with respect to the C2 physical reference (upper-left), CAR (upper-right), and rCAR (bottom). For each panel, clockwise from upper left: a left centrotemporal spike (t = 0.2 s), a right centrotemporal spike (t = 0.21 s), a posterior dominant alpha rhythm, an example of EKG activity, a run of generalized spike and wave activity, and a K-complex (t = .65 s). The data has been bandpass filtered prior to re-referencing to the 1–50 Hz band of frequencies. Vertical streaking prominent when using the CAR re-referenced data for the right centrotemporal spike, the K-complex and the alpha rhythm (t = .41 s), is reduced or absent when re-referencing with rCAR, while the good behavior of CAR is maintained when re-referencing the EKG artifact.

Fig. 8

The difference in the proposed (rCAR) estimate of the unknown reference time-series r_t and the estimate of the r_t time-series associated with CAR. Substantial differences exist between${\check{r}}_t^{\left(rCAR\right)}$ and${\check{r}}_t^{\left(CAR\right)}$ for the six different 19-channel EEG data examples.

Fig. 9

The sum of the absolute channel-to-channel correlations with the proposed re-referenced data (rCAR) is lower than the absolute channel-to-channel correlations estimated in nearly all scenarios. The two panels illustrate the same results with symbols indicating different references (left) and symbols indicating different EEG activity patterns (right).

Fig. 10

The proposed sliding window version of ${\widehat{r}}_t^{\left(rCAR\right)}$ successfully re-references dramatically non-stationary data. The CAR (upper-left) and rCAR (upper-right) re-referenced onset of the 3 Hz bifrontal spike wave complex. This analysis is identical to that used to re-reference the spike wave complex depicted for seconds 10–12.5 in Figs. 7 and 8, except that ${\widehat{r}}_t^{\left(rCAR\right)}$ is computed using the sliding window (50 ms duration) procedure discussed in Section 2.4, for seconds 6–10 as opposed to ${\widehat{r}}_t^{\left(rCAR\right)}$ computed using all of the data from second 10 to second 12.5. Upper left: CAR re-referenced data. Upper right: rCAR re-referenced data. The normalized sum of absolute correlations (NSAC), plotted for the synthetic data in Fig. 3 and for the real data in Fig. 9, is again larger for CAR than for the sliding version of rCAR and is similar to values found in both the previous simulations and in the previous data analysis. Middle left: the difference between the rCAR and CAR reference estimates. The variance of the differences is consistent with the change in the across-channel variance (bottom left), as measured with the median of the across-channel absolute deviations from the across-channel median (MAD). The change in MAD at around second 8 signifies a dramatic change in the across-channel variation. Note that the MAD measure of variability is not affected by dramatic variability present on less than half of the channels. Lower right: Whisker-box plots of the 171 channel-to-channel pairwise correlations computed from the CAR re-referenced data 1st column), computed from the sliding rCAR re-referenced data (2nd column) and the difference between the pairwise correlations: rCAR - CAR (3rd column). Note that the very large negative correlations present with CAR are absent with sliding rCAR, and the median of the pairwise correlations computed using the CAR re-referenced data is lower than the median of the pairwise correlations computed using sliding rCAR. This manifests in the pairwise correlation difference (rCAR-CAR) exhibiting a positive bias – as expected due to the ability of rCAR to avoid the spatially distributed negative correlation exhibited by data re-referenced with CAR.

Fig. B.11

Tukey's bisquare function plotted for c = 2. The re-descending tails down-weight data values that deviate from the theoretical location parameter,–R_f.