- Python version: 3.10.2
- MNE version: 0.24.1
- OS: Linux Manjaro 21.2.4
Hi everyone,
I looked for a similar question but couldn’t find it in the forum (apologies if I missed it). I am brand new to Python and MNE, so I’m possibly biting more than I can chew… any help would be greatly appreciated!
I have a preprocessed EEG dataset with 64 electrodes and two conditions, and I would like to identify ERP components in a data-driven way. Specifically, I want to average the stimulus-locked ERPs of both conditions and run permutation t-test on this grand-averaged dataset, in order to identify statistically significant activity at each electrode and time point. This would result in an unbiased localization of brain activity that should correspond to well-known ERP components (e.g., N1), based on the significant time window and scalp topography. In a second phase, I plan to use the electrode ROI and time window to test condition differences.
I can already do this in MATLAB with the Mass Univariate ERP Toolbox, i.e., calculate one-sample t-tests vs. 0 at each time point and electrode and control the familywise error rate using a t_max approach. According to this tutorial, the permutation_t_test
function in mne
should accomplish a similar goal.
Unfortunately, I am not able to retain the electrode dimension: the permutation_t_test
function seems to automatically average activity across electrodes. Below a (hopefully) reproducible example:
###########################
# import libraries
###########################
import random
import numpy as np
import glob
import re
import mne
from mne.stats import permutation_t_test
###########################
# set parameters
###########################
project_seed = 999 # RNG seed
random.seed(project_seed) # set seed to ensure computational reproducibility
preproc_path = '/path/to/data' # directory with preprocessed files
filenames = glob.glob(preproc_path + '/*.fif') # list of .fif files in directory
conditions = ['cond1', 'cond2'] # condition names
###########################
# load data & compute evoked responses
# (this example only with dataset from one participant)
###########################
# load data, select scalp electrodes only
epochs = mne.read_epochs(filenames[0], preload = True).pick_types(eeg = True, exclude = ['M1', 'M2', 'VEOG', 'HEOG'])
# equalize epochs across conditions (drop epochs from condition with more data)
epochs_equal, dropped_epochs = epochs.equalize_event_counts()
###########################
# permutation t-test
###########################
# select scalp electrodes only
picks = mne.pick_types(epochs_equal.info, eeg = True, exclude = ['M1', 'M2', 'VEOG', 'HEOG'])
# separate condition-specific epochs
epochs_equal_cond1 = epochs_equal['cond1']
epochs_equal_cond2 = epochs_equal['cond2']
# extract all data and put them in one array
all_epochs = [epochs_equal_cond1.get_data(),
epochs_equal_cond2.get_data()]
# define time window of interest
# in this example I don't need to run t-tests on all time points, because I'm only interested in the N1 (this also lowers the chance of Type II error)
t_min = 0
t_max = 0.3
# extract only time points of interest
temporal_mask = np.logical_and(t_min <= epochs_equal.times, epochs_equal.times <= t_max)
# average across conditions, then across epochs, then extract values only from time points of interest
all_epochs_perm = np.mean(np.mean(all_epochs, axis = 0), axis = 0)[:, temporal_mask]
# verify that the array dimensions are as expected (64 channels x 154 time points)
np.shape(all_epochs_perm)
# number of permutations for t-test
n_perm = 5000
# t-test (at last!)
T0, p_values, H0 = permutation_t_test(
all_epochs_perm,
n_permutations = n_perm,
tail = 0,
n_jobs = 1,
seed = project_seed)
The expected output is (for both T0
and p_values
) an array of dimensions 64 x 154 (channels x time points), one value per electrode x time point.
The actual output is (for both T0
and p_values
) an array of dimensions 154 x 1 (time points), one value per time point
Am I using the wrong function? If not, is there a reason why permutation_t_test
does not retain the electrode dimension? Is there a hack that would allow me to do that?
Thanks!