ERD/ERS when there is only one event

đź’¬ Support & Discussions
,
1

Hello,

for a project, I need to build a machine learning algorithm that can predict MI-BCI performance based on mu suppression. So first, I need to calculate the ERD/ERS of the mu band during user’s BCI tasks. The dataset that I’m using contains data from participants during a MI-BCI task where they had to imagine left or right-hand movement. Each trial took 8 seconds where the first 4 seconds are the rest period and the last 4 seconds are the imagery period. Each period is indicated by a TimeStamp from -3 to 5, to at 0 the imagery period starts. I want to calculate the average power of the resting period and of the imagery period and then use the formula for ERD/ERS. To do so, I used code I found on

So far it looks like this:

# import data set 
df = pd.read_csv('TestSubject.csv')

#tmin, tmax 
tmin, tmax = min(df['TimeStamp']),max(df['TimeStamp'])

#transpose data set 
data_trans = df.transpose()

# information for raw file 
ch_names = ['C3','C4','Trigger','class']
sfreq = 250 

#create raw mne file 
info = mne.create_info(ch_names = ch_names,sfreq = sfreq)
raw = mne.io.RawArray(data_trans,info)
raw.plot()

#get picks from raw data 
picks = mne.pick_channels(raw.info["ch_names"], ["C3", "C4"])

events = mne.make_fixed_length_events(raw, id=1, duration=8) #get epochs of 8 seconds

#filter 
raw.filter(l_freq = 8.0 ,h_freq = 13.0,picks=picks,fir_design='firwin')

#create epochs 
epochs = mne.Epochs(raw, events,None, tmin, tmax,
                    picks=picks, baseline=None, preload=True)

event_id, tmin, tmax = 1,-3, 5  # define id and epochs around events (in s)

freqs = np.arange(8, 13, 1)
n_cycles = freqs  # use constant t/f resolution
vmin, vmax = -1, 1.5  # set min and max ERDS values in plot
baseline = [0, 1]  # baseline interval (in s)

cmap = center_cmap(plt.cm.RdBu, vmin, vmax)  # zero maps to white
kwargs = dict(n_permutations=100, step_down_p=0.05, seed=1,
              buffer_size=None, out_type='mask')  # for cluster test

tfr = tfr_multitaper(epochs, freqs=freqs, n_cycles=n_cycles,
                      use_fft=True, return_itc=False, average=False,
                      decim=2,picks=picks)

power = tfr_multitaper(epochs[event_id], freqs=freqs, n_cycles=n_cycles,
                           use_fft=True, return_itc=False, decim=2,picks=picks)
power.crop(tmin, tmax)

How do calculate the ERD/ERS for each epoch?

Thank you in advance.

2

That’s already a good start. A few notes regarding your project:

  • ERD/ERS is defined as an average activity over epochs. If you want to use that as a feature, you will probably have to use band power instead. Otherwise, you will not be able to apply your pipeline to a BCI, where data is continually streaming in. I’m assuming that you are probably computing ERD/ERS maps because you want to select individual frequency bands for later use in band power features. Is this correct? Otherwise, you can try standard bands such as alpha and beta for all subjects.
  • You mention that there are 4s of rest and 4s of MI, but later you are using 3s and 5s. What are the correct durations?
  • If you only have one event ID, you don’t need the for event in event_id loop at all.
  • Currently you compute ERD/ERS twice for C3.
  • The code lacks some important function calls and definitions such as tfr; start and stop are used before they are defined. Maybe it got lost when you copied it here.
3

Thank you so much for your fast response!

  • I will not actually apply my code to a BCI, I only want to use the current dataset to see whether there is a correlation between personality traits and mu suppression (and therefore BCI performance). So, I need to calculate the average mu suppression for each participant during their MI task, that is what I tried to use ERS/ERS for. Is that the correct approach?

  • There are 4 seconds of rest and 4 seconds of MI, but the timestamp starts at -3 and stops at 5.

  • It’s true, I didn’t copy the whole code actually. I implemented start and stop at the beginning.

4

I edited my post, so you can see my whole code for now. I am only using one subject for now, I actually have 54.

5

Can you share a data set so that I can run the code?

If you are not actually going to implement a BCI, ERD/ERS is perfectly OK.

7

It looks like your data is already epoched. Here is some code to create an mne.Epochs object more efficiently:

import pandas as pd
import mne


df = pd.read_csv("TestSubject.csv", index_col=0)
data = df.iloc[:, 1:17].values.reshape((40, 2000, 16)).transpose(0, 2, 1)
fs = 250
info = mne.create_info(list(df.columns[1:17]), fs, "eeg")

epochs = mne.EpochsArray(data, info, tmin=-3 + 1/fs)

Note that I’m using times from -3 to 5, this should be adapted if you really want it to be -4 to 4 (meaning that the timestamps in the file are incorrect).

Now you can go ahead and compute ERD/ERS maps, but before that I was wondering if you really want that. You are filtering the data between 8-13 Hz, so maybe you just want ERD/ERS in this alpha band? This would be much simpler than creating a time/frequency representation.

Also, I’m not sure what the columns Trigger, t, t2, and class actually contain - can you explain?

8

Yes, I only want to create the ERD/ERS of the mu waves (8 – 13 Hz) to calculate the mu suppression for each trial. I thought that I can do that by first filtering the data and then applying a CSP. What is the easier method you are talking about?

These were values I used to re-create the timestamps and trialnumbers in the data.

Trigger indicates the runnumber: -0.1 is calibration, -0.2 is first feedback run, -0.3 second feedback run, -0.4 the third feedback run. t defines stimulus onset by changing when the classes change. t2 is 1 before stimulus onset and 2 after.

Actually, these are not important anymore as I already added the TimeStamp and divided the data into trials.

Thank you!

2809B852AAE04DF6B86B85B9A895EF91.png
2809B852AAE04DF6B86B85B9A895EF91.png702Ă—1 83 Bytes

9

Instead of computing ERD/ERS over many frequencies, you can filter the data in the desired band, square it, and then average over epochs to compute one ERD/ERS time course.

However, for filtering it would be better if you could do that on continuous data, because you will get edge artifacts when you filter epoched data. Do you have access to the continuous data?

10

No, I don’t have the continuous data unfortunately. I extracted the mu waves using that mne.filter function. So my plan was to continue like this: I divided the data into two kinds of epochs, namely imagery and rest. Then I want to calculate the average power of the mu during rest and during imagery and then use this formula: (mu - gamma) / gamma, where gamma represents the average activity during rest and mu the average power during motor imagery. But I am stuck with this step at the moment.
Why did you say that I should square the data first?

Thank you.

11

Hello again,

I tried to calculate the average power of the mu band during rest and imagery with the tfr_morlet module like this (I also implemented the code that you sent me, but I didn’t know how to use it when I want to have two kinds of epochs?)

def calculate_ERDS(raw_filter):
raw_filter, tmin, tmax, start4, picks = create_raw()
#find events based one rest/imagery
events = make_fixed_length_events(raw_filter, duration=8)
event_id = 1

#define active periods 
epochs_active = Epochs(raw_filter, events, event_id, start4, tmax, proj=False,
            picks=picks, baseline=None, preload=True,verbose=False)

# define rest epochs 
epochs_rest = Epochs(raw_filter, events, event_id, tmin, start4, proj=False,
                picks=picks, baseline=None, preload=True,verbose=False)

freqs = np.arange(8, 13, 1) #mu waves 
n_cycles = freqs  # use constant t/f resolution

power_active = epochs_active.average(picks=picks)
tfr = tfr_morlet(epochs_active, freqs=freqs, n_cycles=n_cycles, return_itc=False,picks=picks)
power_active_average = np.mean(tfr.data)
# print(power_active_average)

power_rest = epochs_rest.average(picks=picks)
tfr = tfr_morlet(epochs_rest, freqs=freqs, n_cycles=n_cycles, return_itc=False,picks=picks)
power_rest_average = np.mean(tfr.data)
# print(power_rest_average)

ERDS = (power_active_average - power_rest_average) / power_rest_average

But sometimes I get positive values which is not correct because mu suppression should always be negative… Do you see what I did wrong?

Thank you and kind regards,
Laura

12

You can compute power in a specific frequency band by squaring the (filtered) EEG signals. You never did this, and this is why you’re seeing positive and negative values.

I created just one Epochs object, which contains time points in both the rest (times < 0) and activity periods (times > 0). By selecting only the relevant time periods you can compute both averages from this object. Of course, you can also stick with two separate Epochs objects, but this seems more complicated.