In kusumikakd/EEG: A collection of utilities for EEG analysis

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  eval = nzchar(Sys.getenv("EVAL_EEG_VIGNETTE"))
)

eegUtils is a package for performing basic EEG preprocessing and plotting of EEG data. Many of these functions are wrappers around existing R functions to make them behave in consistent ways and produce output that is more amenable to many of the subsequent steps in EEG analysis.

The package implements custom objects to contain EEG data and associated metadata. Some of its functions depend on data being stored in this format, but part of the philosophy of the package is that any object stored in the custom eeg_data and associated formats will always be convertible to a standard data.frame or tibble for subsequent use in whatever way the user desires. Plotting functions will typically work on both eeg_data objects and standard formats, while more complex processing functions will require an eeg_data or related custom object (such as eeg_tfr).

Basic EEG processing

There is currently suport for loading raw data in the .BDF (typically BioSemi), .CNT (32-bit; associated with Neuroscan), and .vhdr/.vmrk/.dat Brain Vision Analyzer 2.0 file formats using the import_raw() command. Loading data in these formats results in an eeg_data object - a structure that contains the raw data and a variety of metadata.

library(eegUtils)
eeg_example <- import_raw("F:\\Dropbox\\EEGData\\RawEEGData\\BasicCat\\S2B1.bdf")
eeg_example

The data itself is stored in the signals data frame within the object; the timing of events is stored in the events data frame, and timings for each sample are stored in the timings data frame. Times are stored in samples and seconds.

A common first step would be to rereference the data, which can be done using the reref_eeg() command. By default, if no electrodes are specified, the data will be referenced to a common average. The median rather than the mean can be used by specifiying robust == TRUE. Electrodes can also be excluded from the average reference calculation; they will nevertheless still be referenced to the same reference as the other channels.

eeg_example <- eeg_reference(eeg_example, ref_chans = "average")
eeg_example <- eeg_reference(eeg_example, ref_chans = "average", robust = TRUE)
eeg_example <- eeg_reference(eeg_example, ref_chans = "average",
                             exclude = c("EXG7", "EXG8"),
                             robust = TRUE)

EEGLAB can be loaded using the import_set() command. The command detects whether the data is epoched or not and, by default, arranges the data in the eeg_data or eeg_epochs format.

Filtering can be performed using the eeg_filter() command. This uses IIR or FIR filters to modify the frequency response of the signal, remove low or high frequency fluctuations as requested.

eeg_example <- eeg_filter(eeg_example,
                          method = "iir",
                          low_freq = .1,
                          high_freq = 40,
                          filter_order = 3) #specify a bandpass filter

Creating epochs

Data can be epoched around events/triggers using epoch_data(), which outputs an eeg_epochs object. A list of the event triggers found in the data can be retrieved using list_events(). In this case we’ll epoch around event code 201.

list_events(eeg_example)
eeg_example <- epoch_data(eeg_example,
                          events = 201,
                          time_lim = c(-.3, .75))

Plotting

eeg_epochs can then be plotted using plot_butterfly() or plot_timecourse(). Both plot_butterfly() and plot_timecourse() average over epochs. plot_timecourse() will also average over electrodes - all electrodes if none are specified, or over any specified electrodes. Baseline correction can also be applied at this point if desired. Baseline correction can also be applied using rm_baseline(), or for plotting only using the baseline parameter in the plotting call.

plot_butterfly(eeg_example,
               time_lim = c(-.2, .5),
               legend = FALSE)
plot_butterfly(eeg_example,
               time_lim = c(-.2, .5),
               legend = FALSE,
               baseline = c(-.1, 0))
plot_timecourse(eeg_example,
                time_lim = c(-.2, .5),
                baseline = c(-.1, 0),
                electrode = "A30") # Plot POz
plot_timecourse(eeg_example,
                time_lim = c(-.2, .5),
                baseline = c(-.1, 0),
                electrode = "B6") # Plot Fz

Standard channel locations can be added using the electrode_locations() command. This function supplies default locations for over 300 typical locations accroding to the 10-05 system. There are several specific montages provided, and can be specified using the montage parameter. In this case, the data was recorded using Biosemi with alpha-numeric channel labels (e.g. A1 rather than Fp1), which are built-in with biosemi64alpha. topoplot() can then be used to plot topographical representation of selected data. Note that it is not compulsory to use locations from electrode_locations(); if the data has x and y columns when it is a data frame, or added to chan_info element of the eeg_data/eeg_epochs object, then those will be used.

eeg_example <- electrode_locations(eeg_example,
                                   montage = "biosemi64alpha")
topoplot(eeg_example,
         time_lim = c(.1, .2))

Note that one electrode (Fp1/A1) was malunctioning during the the recording. We can interpolate that using spherical spline interpolation (interp_elecs()).

eeg_example <- interp_elecs(eeg_example,
                            bad_elecs = "A1")
topoplot(eeg_example,
         time_lim = c(.1, .2))

At any point, eegUtils objects can be transformed into data frames for use with functions that don't natively support them.

library(ggplot2)
library(dplyr)
eeg_example %>%
  select_epochs(epoch_no = 1:5) %>%
  select_elecs(c("A29", "B6")) %>%
  as.data.frame(long = TRUE) %>%
  ggplot(aes(x = time, y = amplitude)) +
  geom_line(aes(group = epoch), alpha = 0.2) + 
  stat_summary(fun.y = mean,
               geom = "line",
               size = 2,
               aes(colour = electrode)) + 
  facet_wrap(~electrode) + 
  theme_classic()

Tidyverse functions

In addition, there are overloaded versions of some dplyr functions that operate on the signals element of eeg_data and eeg_epochs objects. For example, select() can be used to choose particular electrodes, and filter() can be used to filter out epochs or timepoints. mutate() can be used to add new columns (e.g. creating ROIs from multiple electrodes).

eeg_example %>%
  mutate(occipital = (A29 + A27 + B32) / 3) %>%
  select(A29, B6, occipital) %>%
  filter(epoch <= 60, time < .5) %>%
  as.data.frame(long = TRUE) %>%
  ggplot(aes(x = time, y = amplitude)) +
  geom_line(aes(group = epoch), alpha = 0.2) + 
  stat_summary(fun.y = mean,
               geom = "line",
               size = 2,
               aes(colour = electrode)) +
  facet_wrap(~electrode) + 
  scale_colour_viridis_d() +
  theme_classic()