From d79a199a40c3711420148735be2374b26213cd6a Mon Sep 17 00:00:00 2001 From: Okba Bekhelifi Date: Fri, 29 Nov 2019 20:21:06 +0100 Subject: [PATCH] - issue #5 : - auto lint with autopep8 - issue #3 : - Added visualization module: - Learning curve: train_val - ROC curve - Confustion Matrix On branch dev Changes to be committed: modified: aawedha/analysis/preprocess.py modified: aawedha/evaluation/base.py modified: aawedha/evaluation/cross_subject.py modified: aawedha/evaluation/single_subject.py modified: aawedha/io/base.py modified: aawedha/io/bci_comp_iv_2a.py modified: aawedha/io/exoskeleton.py modified: aawedha/io/freiburg_online.py modified: aawedha/io/inria_ern.py modified: aawedha/io/laresi.py modified: aawedha/io/physionet_mi.py modified: aawedha/io/san_diego.py modified: aawedha/io/tsinghua.py modified: aawedha/models/EEGModels.py modified: aawedha/paradigms/base.py modified: aawedha/paradigms/erp.py modified: aawedha/paradigms/hybrid.py modified: aawedha/paradigms/motor_imagery.py modified: aawedha/paradigms/ssvep.py modified: aawedha/paradigms/subject.py modified: aawedha/process/base.py modified: aawedha/scripts/cross_subj.py modified: aawedha/scripts/single_subj.py new file: aawedha/visualization/__init__.py new file: aawedha/visualization/viz.py --- aawedha/analysis/preprocess.py | 21 +- aawedha/evaluation/base.py | 122 ++++--- aawedha/evaluation/cross_subject.py | 155 ++++----- aawedha/evaluation/single_subject.py | 152 +++++---- aawedha/io/base.py | 21 +- aawedha/io/bci_comp_iv_2a.py | 85 ++--- aawedha/io/exoskeleton.py | 121 +++---- aawedha/io/freiburg_online.py | 54 ++-- aawedha/io/inria_ern.py | 100 +++--- aawedha/io/laresi.py | 276 ++++++++-------- aawedha/io/physionet_mi.py | 53 ++-- aawedha/io/san_diego.py | 83 ++--- aawedha/io/tsinghua.py | 119 +++---- aawedha/models/EEGModels.py | 459 +++++++++++++-------------- aawedha/paradigms/base.py | 9 +- aawedha/paradigms/erp.py | 14 +- aawedha/paradigms/hybrid.py | 29 +- aawedha/paradigms/motor_imagery.py | 29 +- aawedha/paradigms/ssvep.py | 16 +- aawedha/paradigms/subject.py | 8 +- aawedha/process/base.py | 13 +- aawedha/scripts/cross_subj.py | 34 +- aawedha/scripts/single_subj.py | 26 +- aawedha/visualization/__init__.py | 0 aawedha/visualization/viz.py | 76 +++++ 25 files changed, 1088 insertions(+), 987 deletions(-) create mode 100644 aawedha/visualization/__init__.py create mode 100644 aawedha/visualization/viz.py diff --git a/aawedha/analysis/preprocess.py b/aawedha/analysis/preprocess.py index 3379e20..bb9c7b5 100644 --- a/aawedha/analysis/preprocess.py +++ b/aawedha/analysis/preprocess.py @@ -3,24 +3,27 @@ def bandpass(eeg, band, fs, order=2): - B,A = sig.butter(order, np.array(band)/(fs/2), btype='bandpass') - return sig.filtfilt(B, A, eeg, axis=0) + B, A = sig.butter(order, np.array(band) / (fs / 2), btype='bandpass') + return sig.filtfilt(B, A, eeg, axis=0) def eeg_epoch(eeg, epoch_length, markers): ''' - input + input : eeg : continuous eeg : 2D numpy array samples x channels - : epoch_length : 1D numpy array (size 2 )epoch start and stop in samples + : epoch_length : 1D numpy array (size 2 )epoch start and stop + in samples : markers : event markers onset in samples 1D array 1,n_markers - - output + + output : eeg_epochs : 3D array of epoched EEG : samples x channels x trials (Fortran ordering aka MATLAB format) ''' channels = int(eeg.shape[1]) epoch_length = np.around(epoch_length) - dur = np.arange(epoch_length[0], epoch_length[1]).reshape((np.diff(epoch_length)[0],1)) * np.ones( (1, len(markers)),dtype=int) + dur = np.arange(epoch_length[0], epoch_length[1]).reshape( + (np.diff(epoch_length)[0], 1)) * np.ones((1, len(markers)), dtype=int) samples = len(dur) epoch_idx = dur + markers - eeg_epochs = np.array(eeg[epoch_idx,:]).reshape((samples, len(markers), channels), order='F').transpose((0,2,1)) - return eeg_epochs \ No newline at end of file + eeg_epochs = np.array(eeg[epoch_idx, :]).reshape( + (samples, len(markers), channels), order='F').transpose((0, 2, 1)) + return eeg_epochs diff --git a/aawedha/evaluation/base.py b/aawedha/evaluation/base.py index 3c499c6..8bebd87 100644 --- a/aawedha/evaluation/base.py +++ b/aawedha/evaluation/base.py @@ -9,23 +9,22 @@ import abc import os + class Evaluation(object): - def __init__(self, n_subjects=0, partition=[], - folds=[], dataset=None, - model=None, verbose=2 - ): - + def __init__(self, n_subjects=0, partition=[], folds=[], dataset=None, + model=None, verbose=2): + ''' + ''' self.n_subjects = n_subjects self.partition = partition self.folds = folds self.dataset = dataset self.model = model - self.cm = [] # confusion matrix per fold - self.results = {} # dict + self.cm = [] # confusion matrix per fold + self.results = {} # dict self.model_history = {} self.verbose = verbose - @abc.abstractmethod def generate_split(self, nfolds): @@ -33,37 +32,34 @@ def generate_split(self, nfolds): ''' pass - @abc.abstractmethod + @abc.abstractmethod def run_evaluation(self): ''' ''' pass - def measure_performance(self, Y_test, probs): ''' ''' preds = probs.argmax(axis=-1) - y_true = Y_test.argmax(axis=-1) + y_true = Y_test.argmax(axis=-1) classes = Y_test.shape[1] - acc = np.mean(preds == y_true) - + acc = np.mean(preds == y_true) + self.cm.append(confusion_matrix(Y_test.argmax(axis=-1), preds)) - if classes == 2: - fp_rate, tp_rate, thresholds = roc_curve(y_true, probs[:,1]) - auc_score = auc(fp_rate, tp_rate) + if classes == 2: + fp_rate, tp_rate, thresholds = roc_curve(y_true, probs[:, 1]) + auc_score = auc(fp_rate, tp_rate) return acc.item(), auc_score.item(), fp_rate, tp_rate else: return acc.item() - - def results_reports(self, res, tfpr={}): ''' ''' results = {} - + # if tfpr: # res : (metric, subjects, folds) # means = res.mean(axis=-1) # mean across folds @@ -72,105 +68,104 @@ def results_reports(self, res, tfpr={}): results['acc_mean_per_fold'] = r1.mean(axis=0) results['acc_mean_per_subj'] = r1.mean(axis=1) results['acc_mean'] = r1.mean() - + # r2 = np.array(res['auc']) results['auc'] = r2 results['auc_mean_per_fold'] = r2.mean(axis=0) results['auc_mean_per_subj'] = r2.mean(axis=1) results['auc_mean'] = r2.mean() - + # results['fpr'] = tfpr['fp'] results['tpr'] = tfpr['tp'] else: # res : (subjects, folds) results['acc'] = res - results['acc_mean_per_fold'] = res.mean(axis=0) # mean across folds - results['acc_mean_per_subject'] = res.mean(axis=1) # mean across subjects and folds + # mean across folds + results['acc_mean_per_fold'] = res.mean(axis=0) + # mean across subjects and folds + results['acc_mean_per_subject'] = res.mean(axis=1) results['acc_mean'] = res.mean() - return results + return results def get_folds(self, nfolds, population, tr, vl, ts, exclude_subj=True): ''' ''' folds = [] - if hasattr(self.dataset, 'test_epochs'): + if hasattr(self.dataset, 'test_epochs'): if self.__class__.__name__ == 'CrossSubject': # independent test set - # list : nfolds : [nsubjects_train] [nsubjects_val] + # list : nfolds : [nsubjects_train] [nsubjects_val] for subj in range(self.n_subjects): selection = np.arange(0, self.n_subjects) if exclude_subj: # fully cross-subject, no subject train data in fold selection = np.delete(selection, subj) - for fold in range(nfolds): - np.random.shuffle(selection) - folds.append([np.array(selection[:tr]), np.array(selection[tr:])]) + for fold in range(nfolds): + np.random.shuffle(selection) + folds.append([np.array(selection[:tr]), + np.array(selection[tr:])]) elif self.__class__.__name__ == 'SingleSubject': - # generate folds for test set from one set + # generate folds for test set from one set pop = population t = tr v = vl - s = ts + s = ts for f in range(nfolds): - if type(population) is np.ndarray: - # inconsistent numbers of trials among subjects + if isinstance(population, np.ndarray): + # inconsistent numbers of trials among subjects sbj = [] - for subj in range(self.n_subjects): + for subj in range(self.n_subjects): pop = population[subj] t = tr[subj] v = vl[subj] s = ts[subj] - tmp = np.array(random.sample(range(pop), pop)) - sbj.append([tmp[:t], tmp[t:t+v], tmp[-s:]]) + tmp = np.array(random.sample(range(pop), pop)) + sbj.append([tmp[:t], tmp[t:t + v], tmp[-s:]]) folds.append(sbj) - else: - # same numbers of trials for all subjects + else: + # same numbers of trials for all subjects tmp = np.array(random.sample(range(pop), pop)) - folds.append([tmp[:t], tmp[t:t+v], tmp[-s:]]) - - #for _ in range(nfolds): - # tmp = np.array(random.sample(range(population), population)) - # folds.append([tmp[:tr], tmp[tr:tr+vl], tmp[-ts:]]) + folds.append([tmp[:t], tmp[t:t + v], tmp[-s:]]) else: # generate folds for test set from one set for _ in range(nfolds): tmp = np.array(random.sample(range(population), population)) - folds.append([tmp[:tr], tmp[tr:tr+vl], tmp[-ts:]]) - - return folds - + folds.append([tmp[:tr], tmp[tr:tr + vl], tmp[-ts:]]) + # + return folds + def fit_scale(self, X): mu = X.mean(axis=0) sigma = X.std(axis=0) - X = np.subtract(X, mu[None,:,:]) - X = np.divide(X, sigma[None,:,:]) + X = np.subtract(X, mu[None, :, :]) + X = np.divide(X, sigma[None, :, :]) return X, mu, sigma - def transform_scale(self, X, mu, sigma): - X = np.subtract(X, mu[None,:,:]) - X = np.divide(X, sigma[None,:,:]) + X = np.subtract(X, mu[None, :, :]) + X = np.divide(X, sigma[None, :, :]) return X - def class_weights(self, y): ''' ''' cl_weights = {} classes = np.unique(y) - n_perclass = [np.sum(y==cl) for cl in classes] + n_perclass = [np.sum(y == cl) for cl in classes] n_samples = np.sum(n_perclass) - ws = np.array([np.ceil(n_samples / cl).astype(int) for cl in n_perclass]) + ws = np.array([np.ceil(n_samples / cl).astype(int) + for cl in n_perclass]) if np.unique(ws).size == 1: # balanced classes - cl_weights = {cl:1 for cl in classes} + cl_weights = {cl: 1 for cl in classes} else: # unbalanced classes if classes.size == 2: - cl_weights = {classes[ws == ws.max()].item():ws.max(), classes[ws < ws.max()].item():1} + cl_weights = {classes[ws == ws.max()].item( + ): ws.max(), classes[ws < ws.max()].item(): 1} else: - cl_weights = {cl:ws[idx] for idx,cl in enumerate(classes)} + cl_weights = {cl: ws[idx] for idx, cl in enumerate(classes)} return cl_weights def labels_to_categorical(self, y): @@ -180,18 +175,17 @@ def labels_to_categorical(self, y): if np.isin(0, classes): y = to_categorical(y) else: - y = to_categorical(y-1) + y = to_categorical(y - 1) return y - + def save_model(self, folderpath=None): ''' ''' if not os.path.isdir('trained'): os.mkdir('trained') - if not folderpath: - folderpath = 'trained' + if not folderpath: + folderpath = 'trained' prdg = self.dataset.paradigm.title dt = self.dataset.title - filepath = folderpath + '/' + '_'.join(['model',prdg,dt,'.h5']) + filepath = folderpath + '/' + '_'.join(['model', prdg, dt, '.h5']) self.model.save(filepath) - diff --git a/aawedha/evaluation/cross_subject.py b/aawedha/evaluation/cross_subject.py index 95ab4d2..c241e6d 100644 --- a/aawedha/evaluation/cross_subject.py +++ b/aawedha/evaluation/cross_subject.py @@ -1,10 +1,8 @@ from aawedha.evaluation.base import Evaluation -from sklearn.metrics import roc_auc_score import numpy as np -import random -class CrossSubject(Evaluation): +class CrossSubject(Evaluation): def generate_split(self, nfolds=30, excl=True): ''' @@ -12,19 +10,20 @@ def generate_split(self, nfolds=30, excl=True): # folds = [] n_phase = len(self.partition) train_phase, val_phase = self.partition[0], self.partition[1] - + if n_phase == 2: - # independent test set available + # independent test set available test_phase = 0 elif n_phase == 3: # generate a set set from the dataset - test_phase = self.partition[2] + test_phase = self.partition[2] else: # error : wrong partition raise AssertionError('Wrong partition scheme', self.partition) - self.folds = self.get_folds(nfolds, self.n_subjects, train_phase, val_phase, test_phase, exclude_subj=excl) - + self.folds = self.get_folds( + nfolds, self.n_subjects, train_phase, val_phase, + test_phase, exclude_subj=excl) def run_evaluation(self, clbs=[]): ''' @@ -32,38 +31,35 @@ def run_evaluation(self, clbs=[]): # generate folds if folds are empty if not self.folds: self.folds = self.generate_split(nfolds=30) - # - res_acc , res_auc = [], [] + # + res_acc, res_auc = [], [] res_tp, res_fp = [], [] for fold in range(len(self.folds)): - rets = self._cross_subject(fold,clbs) - if type(rets) is tuple: + rets = self._cross_subject(fold, clbs) + if isinstance(rets, tuple): res_acc.append(rets[0]) - res_auc.append(rets[1]) + res_auc.append(rets[1]) res_fp.append(rets[2]) res_tp.append(rets[3]) else: - res_acc.append(rets) - + res_acc.append(rets) + # Aggregate results if res_auc: res = np.array([res_acc, res_auc]) tfpr = np.array([res_fp, res_tp]) else: res = np.array(res_acc) - tfpr = [] - - self.results = self.results_reports(res, tfpr) + tfpr = [] + + self.results = self.results_reports(res, tfpr) def _cross_subject(self, fold, clbs=[]): ''' ''' - # - res_acc = [] - res_auc = [] - split = self._split_set(fold) - classes = split['classes'] + # + split = self._split_set(fold) # normalize data X_train, mu, sigma = self.fit_scale(split['X_train']) X_val = self.transform_scale(split['X_val'], mu, sigma) @@ -74,14 +70,15 @@ def _cross_subject(self, fold, clbs=[]): # cws = self.class_weights(np.argmax(Y_train, axis=1)) # evaluate model on subj on all folds - - self.model_history = self.model.fit(X_train, Y_train, batch_size = 96, epochs = 500, - verbose = self.verbose, validation_data=(X_val, Y_val), - class_weight = cws, callbacks=clbs) - # train/val + + self.model_history = self.model.fit(X_train, Y_train, batch_size=96, + epochs=500, verbose=self.verbose, + validation_data=(X_val, Y_val), + class_weight=cws, callbacks=clbs) + # train/val probs = self.model.predict(X_test) - rets = self.measure_performance(Y_test, probs) - + rets = self.measure_performance(Y_test, probs) + return rets def _split_set(self, fold): @@ -89,69 +86,81 @@ def _split_set(self, fold): ''' split = {} kernels = 1 - if type(self.dataset.epochs) is list: + if isinstance(self.dataset.epochs, list): # TODO - X_train, Y_train = self._cat_lists(fold,0) - X_val, Y_val = self._cat_lists(fold,1) - X_test, Y_test = self._cat_lists(fold,2) + X_train, Y_train = self._cat_lists(fold, 0) + X_val, Y_val = self._cat_lists(fold, 1) + X_test, Y_test = self._cat_lists(fold, 2) classes = np.unique(Y_train) - samples, channels, _ = X_train.shape - X_train = X_train.transpose((2,1,0)).reshape((X_train.shape[2], kernels, channels, samples)) - X_val = X_val.transpose((2,1,0)).reshape((X_val.shape[2], kernels, channels, samples)) - X_test = X_test.transpose((2,1,0)).reshape((X_test.shape[2], kernels, channels, samples)) + samples, channels, _ = X_train.shape + X_train = X_train.transpose((2, 1, 0)).reshape( + (X_train.shape[2], kernels, channels, samples)) + X_val = X_val.transpose((2, 1, 0)).reshape( + (X_val.shape[2], kernels, channels, samples)) + X_test = X_test.transpose((2, 1, 0)).reshape( + (X_test.shape[2], kernels, channels, samples)) Y_train = self.labels_to_categorical(Y_train) - Y_val = self.labels_to_categorical(Y_val) - Y_test = self.labels_to_categorical(Y_test) - + Y_val = self.labels_to_categorical(Y_val) + Y_test = self.labels_to_categorical(Y_test) + else: - - x = self.dataset.epochs - subjects, samples, channels, trials = x.shape + + x = self.dataset.epochs + subjects, samples, channels, trials = x.shape y = self.dataset.y - x = x.transpose((0,3,2,1)) + x = x.transpose((0, 3, 2, 1)) # classes = np.unique(y) - y = self.labels_to_categorical(y) - tr, val = self.partition[0], self.partition[1] # n_subjects per train/val - if len(self.partition) == 3: - ts = self.partition[2] - - X_train = x[self.folds[fold][0],:,:,:].reshape((tr*trials,kernels,channels,samples)) - X_val = x[self.folds[fold][1],:,:,:].reshape((val*trials,kernels,channels,samples)) + y = self.labels_to_categorical(y) + # n_subjects per train/val + tr, val = self.partition[0], self.partition[1] + if len(self.partition) == 3: + ts = self.partition[2] + + X_train = x[self.folds[fold][0], :, :, :].reshape( + (tr * trials, kernels, channels, samples)) + X_val = x[self.folds[fold][1], :, :, :].reshape( + (val * trials, kernels, channels, samples)) ctg_dim = y.shape[2] - Y_train = y[self.folds[fold][0],:].reshape((tr*trials, ctg_dim)) - Y_val = y[self.folds[fold][1],:].reshape((val*trials, ctg_dim)) - + Y_train = y[self.folds[fold][0], :].reshape((tr * trials, ctg_dim)) + Y_val = y[self.folds[fold][1], :].reshape((val * trials, ctg_dim)) + if hasattr(self.dataset, 'test_epochs'): trs = self.dataset.test_epochs.shape[3] - X_test = self.dataset.test_epochs.transpose((0,3,2,1)) - X_test = X_test.reshape((trs*self.n_subjects , kernels , channels , samples)) + X_test = self.dataset.test_epochs.transpose((0, 3, 2, 1)) + X_test = X_test.reshape( + (trs * self.n_subjects, kernels, channels, samples)) # Y_test = tf_utils.to_categorical(self.dataset.test_y) - Y_test = self.labels_to_categorical(self.dataset.test_y.reshape((self.n_subjects*trs))) + Y_test = self.labels_to_categorical( + self.dataset.test_y.reshape((self.n_subjects * trs))) else: - X_test = x[self.folds[fold][2],:,:,:].reshape((ts*trials,kernels,channels,samples)) - Y_test = y[self.folds[fold][2],:].reshape((ts*trials, ctg_dim)) + X_test = x[self.folds[fold][2], :, :, :].reshape( + (ts * trials, kernels, channels, samples)) + Y_test = y[self.folds[fold][2], :].reshape( + (ts * trials, ctg_dim)) split['X_train'] = X_train - split['X_val'] = X_val - split['X_test'] = X_test + split['X_val'] = X_val + split['X_test'] = X_test split['Y_train'] = Y_train - split['Y_val'] = Y_val - split['Y_test'] = Y_test + split['Y_val'] = Y_val + split['Y_test'] = Y_test split['classes'] = classes return split def _cat_lists(self, fold=0, phase=0): ''' ''' - if phase == 2: # test phase + if phase == 2: # test phase if hasattr(self.dataset, 'test_epochs'): - X = np.concatenate([self.dataset.test_epochs[idx] for idx in range(self.n_subjects)], axis=-1) - Y = np.concatenate([self.dataset.test_y[idx] for idx in range(self.n_subjects)], axis=-1) - return X,Y - X = np.concatenate([self.dataset.epochs[idx] for idx in self.folds[fold][phase]], axis=-1) - Y = np.concatenate([self.dataset.y[idx] for idx in self.folds[fold][phase]], axis=-1) - return X,Y - - + X = np.concatenate([self.dataset.test_epochs[idx] + for idx in range(self.n_subjects)], axis=-1) + Y = np.concatenate([self.dataset.test_y[idx] + for idx in range(self.n_subjects)], axis=-1) + return X, Y + X = np.concatenate([self.dataset.epochs[idx] + for idx in self.folds[fold][phase]], axis=-1) + Y = np.concatenate([self.dataset.y[idx] + for idx in self.folds[fold][phase]], axis=-1) + return X, Y diff --git a/aawedha/evaluation/single_subject.py b/aawedha/evaluation/single_subject.py index 76efb7d..88c4186 100644 --- a/aawedha/evaluation/single_subject.py +++ b/aawedha/evaluation/single_subject.py @@ -1,56 +1,55 @@ from aawedha.evaluation.base import Evaluation from sklearn.model_selection import KFold -# from sklearn.metrics import roc_auc_score -from sklearn.metrics import roc_curve, auc -from tensorflow.keras.callbacks import EarlyStopping import numpy as np -import random class SingleSubject(Evaluation): - - + # def generate_split(self, nfolds=30, strategy='Kfold'): ''' ''' # folds = [] n_phase = len(self.partition) train_phase, val_phase = self.partition[0], self.partition[1] - - if type(self.dataset.y) is list: - # - n_trials = np.array([self.dataset.y[i].shape[0] for i in range(self.n_subjects)]) + # + if isinstance(self.dataset.y, list): + # + n_trials = np.array([self.dataset.y[i].shape[0] + for i in range(self.n_subjects)]) else: n_trials = self.dataset.y.shape[1] - + # if n_phase == 2: - # independent test set available + # independent test set available test_phase = 0 elif n_phase == 3: # generate a set set from the dataset - test_phase = self.partition[2] + test_phase = self.partition[2] else: # error : wrong partition - raise AssertionError('Wrong partition scheme', self.partition) - + raise AssertionError('Wrong partition scheme', self.partition) + # part = np.round(n_trials / np.sum(self.partition)).astype(int) - + # train_phase = train_phase * part val_phase = val_phase * part test_phase = test_phase * part if strategy == 'Kfold': - self.folds = self._get_folds(nfolds, n_trials, train_phase, val_phase, test_phase) + self.folds = self._get_folds(nfolds, n_trials, train_phase, + val_phase, test_phase) elif strategy == 'Shuffle': - self.folds = self.get_folds(nfolds, n_trials, train_phase, val_phase, test_phase, exclude_subj=False) + self.folds = self.get_folds(nfolds, n_trials, train_phase, + val_phase, test_phase, + exclude_subj=False) def run_evaluation(self, subject=None, clbs=[]): ''' - ''' + ''' # generate folds if folds are empty if not self.folds: self.folds = self.generate_split(nfolds=30) - # + # res_acc = [] res_auc, res_tp, res_fp = [], [], [] @@ -63,41 +62,40 @@ def run_evaluation(self, subject=None, clbs=[]): independent_test = True else: # concatenate train & test data - # test data are different subjects - n_subj = self._fuse_data() + # test data are different subjects + n_subj = self._fuse_data() self.n_subjects = n_subj if subject: operations = [subject] - else: + else: operations = range(self.n_subjects) for subj in operations: - # res_per_subject, avg_res = self._single_subject(subj, independent_test) - # acc_per_subject, auc_per_subject, fp, tp = self._single_subject(subj, independent_test, clbs) + # res_per_subject, avg_res = self._single_subject(subj, + # independent_test) + # acc_per_subject, auc_per_subject, fp, tp = self._single_subject + # (subj, independent_test, clbs) rets = self._single_subject(subj, independent_test, clbs) - if type(rets[0]) is tuple: + if isinstance(rets[0], tuple): res_acc.append([elm[0] for elm in rets]) res_auc.append([elm[1] for elm in rets]) res_fp.append([elm[2] for elm in rets]) - res_tp.append([elm[3] for elm in rets]) + res_tp.append([elm[3] for elm in rets]) else: - res_acc.append(rets) - + res_acc.append(rets) # Aggregate results - tfpr = {} - if res_auc: + if res_auc: res = {} res['acc'] = res_acc - res['auc'] = res_auc + res['auc'] = res_auc tfpr['fp'] = res_fp tfpr['tp'] = res_tp else: - res = np.array(res_acc) - + res = np.array(res_acc) + # self.results = self.results_reports(res, tfpr) - def _single_subject(self, subj, indie=False, clbs=[]): ''' @@ -106,14 +104,11 @@ def _single_subject(self, subj, indie=False, clbs=[]): x = self.dataset.epochs[subj][:, :, :] y = self.dataset.y[subj][:] samples, channels, trials = x.shape - x = x.transpose((2,1,0)) - kernels = 1 # + x = x.transpose((2, 1, 0)) + kernels = 1 # x = x.reshape((trials, kernels, channels, samples)) - # - classes = np.unique(y) + # y = self.labels_to_categorical(y) - # res_acc = [] - # res_auc = [] rets = [] # get in the fold!!! for fold in range(len(self.folds)): @@ -125,43 +120,45 @@ def _single_subject(self, subj, indie=False, clbs=[]): X_test = self.transform_scale(split['X_test'], mu, sigma) # Y_train = split['Y_train'] - Y_test = split['Y_test'] - Y_val = split['Y_val'] + Y_test = split['Y_test'] + Y_val = split['Y_val'] # class_weights = self.class_weights(np.argmax(Y_train, axis=1)) # evaluate model on subj on all folds # rename the fit method - self.model_history = self.model.fit(X_train, Y_train, batch_size = 64, epochs = 500, - verbose = self.verbose, validation_data=(X_val, Y_val), - class_weight = class_weights, callbacks=clbs) - + self.model_history = self.model.fit(X_train, Y_train, + batch_size=64, epochs=500, + verbose=self.verbose, + validation_data=(X_val, Y_val), + class_weight=class_weights, + callbacks=clbs) + # probs = self.model.predict(X_test) - rets.append(self.measure_performance(Y_test, probs)) + rets.append(self.measure_performance(Y_test, probs)) - return rets - # res = [] # shape: (n_folds, 2) - + return rets def _fuse_data(self): ''' ''' # TODO : when epochs/y/test_epochs/y_test are lists??? # make lists of lists - if type(self.dataset.epochs) is list: + if isinstance(self.dataset.epochs, list): # TODO pass else: - self.dataset.epochs = np.vstack((self.dataset.epochs, self.dataset.test_epochs)) + self.dataset.epochs = np.vstack( + (self.dataset.epochs, self.dataset.test_epochs)) self.dataset.y = np.vstack((self.dataset.y, self.dataset.test_y)) - return self.dataset.epochs.shape[0] # n_subject + return self.dataset.epochs.shape[0] # n_subject def _equale_subjects(self): ''' ''' ts = 0 tr = len(self.dataset.epochs) - if hasattr(self.dataset, 'test_epochs'): - ts = len(self.dataset.test_epochs) + if hasattr(self.dataset, 'test_epochs'): + ts = len(self.dataset.test_epochs) return tr == ts def _split_set(self, x=None, y=None, subj=0, fold=0, indie=False): @@ -171,47 +168,46 @@ def _split_set(self, x=None, y=None, subj=0, fold=0, indie=False): # folds[0][0] : same trials numbers for all subjects split = {} trials, kernels, channels, samples = x.shape - - if type(self.dataset.epochs) is list: + # + if isinstance(self.dataset.epochs, list): f = self.folds[fold][subj][:] else: f = self.folds[fold][:] - - X_train = x[f[0],:,:,:] - X_val = x[f[1],:,:,:] + # + X_train = x[f[0], :, :, :] + X_val = x[f[1], :, :, :] Y_train = y[f[0]] - Y_val = y[f[1]] + Y_val = y[f[1]] if indie: trs = self.dataset.test_epochs[0].shape[2] - X_test = self.dataset.test_epochs[subj][:,:,:].transpose((2,1,0)) - X_test = X_test.reshape((trs, kernels, channels, samples)) + X_test = self.dataset.test_epochs[subj][:, :, :].transpose( + (2, 1, 0)) + X_test = X_test.reshape((trs, kernels, channels, samples)) Y_test = self.labels_to_categorical(self.dataset.test_y[subj][:]) else: - X_test = x[f[2],:,:,:] - Y_test = y[f[2]] + X_test = x[f[2], :, :, :] + Y_test = y[f[2]] split['X_train'] = X_train - split['X_val'] = X_val - split['X_test'] = X_test split['Y_train'] = Y_train - split['Y_val'] = Y_val - split['Y_test'] = Y_test - + split['X_test'] = X_test + split['Y_test'] = Y_test + split['X_val'] = X_val + split['Y_val'] = Y_val return split def _get_folds(self, nfolds=4, n_trials=0, tr=0, vl=0, ts=0): ''' ''' - if type(n_trials) is np.ndarray: + if isinstance(n_trials, np.ndarray): t = [np.arange(n) for n in n_trials] sbj = [] - folds = [] - sbj = [self._get_split(nfolds, n, tr, vl) for n in t] + folds = [] + sbj = [self._get_split(nfolds, n, tr, vl) for n in t] folds.append(sbj) else: t = np.arange(n_trials) folds = self._get_split(nfolds, t, tr, vl) - return folds def _get_split(self, nfolds, t, tr, vl): @@ -223,9 +219,7 @@ def _get_split(self, nfolds, t, tr, vl): if len(self.partition) == 2: # independent test set folds.append([train, test]) - elif len(self.partition) ==3: + elif len(self.partition) == 3: # generate test set from the entire set - folds.append([train[:tr], train[tr:tr+vl], test]) - + folds.append([train[:tr], train[tr:tr + vl], test]) return folds - \ No newline at end of file diff --git a/aawedha/io/base.py b/aawedha/io/base.py index 68921cb..f1e9ac6 100644 --- a/aawedha/io/base.py +++ b/aawedha/io/base.py @@ -3,13 +3,12 @@ """ from abc import ABCMeta, abstractmethod import os -import sys import pickle class DataSet(metaclass=ABCMeta): """DataSet - + Attributes ---------- title : str @@ -49,6 +48,7 @@ class DataSet(metaclass=ABCMeta): generate_set() """ + def __init__(self, title='', ch_names=[], fs=None, doi=''): self.epochs = [] self.y = [] @@ -63,7 +63,7 @@ def __init__(self, title='', ch_names=[], fs=None, doi=''): @abstractmethod def load_raw(self): pass - + @abstractmethod def generate_set(self): pass @@ -82,10 +82,10 @@ def save_set(self, save_folder=None): if not self.title: self.title = 'unnamed_set' - fname = save_folder + '/' + self.title +'.pkl' - print(f'Saving dataset {self.title} to destination: {fname}') + fname = save_folder + '/' + self.title + '.pkl' + print(f'Saving dataset {self.title} to destination: {fname}') f = open(fname, 'wb') - pickle.dump(self, f, protocol=pickle.HIGHEST_PROTOCOL) + pickle.dump(self, f, protocol=pickle.HIGHEST_PROTOCOL) f.close() # log if verbose @@ -98,11 +98,4 @@ def load_set(self, file_name=None): else: raise FileNotFoundError f.close() - return data - - - - - - - + return data diff --git a/aawedha/io/bci_comp_iv_2a.py b/aawedha/io/bci_comp_iv_2a.py index 0d72edd..c90df26 100644 --- a/aawedha/io/bci_comp_iv_2a.py +++ b/aawedha/io/bci_comp_iv_2a.py @@ -23,12 +23,13 @@ class Comp_IV_2a(DataSet): 2012. Available: http://journal.frontiersin.org/article/10.3389/fnins.2012.00055 """ + def __init__(self): super().__init__(title='Comp_IV_2a', - ch_names=['Fz','FC3','FC1','FCz', - 'FC2','FC4','C5','C3','C1','Cz', - 'C2','C4','C6','CP3','CP1','CPz', - 'CP2','CP4','P1','Pz','P2','POz'], + ch_names=['Fz', 'FC3', 'FC1', 'FCz', + 'FC2', 'FC4', 'C5', 'C3', 'C1', 'Cz', + 'C2', 'C4', 'C6', 'CP3', 'CP1', 'CPz', + 'CP2', 'CP4', 'P1', 'Pz', 'P2', 'POz'], fs=250, doi='https://doi.org/10.3389/fnins.2012.00055' ) @@ -36,14 +37,15 @@ def __init__(self): self.test_y = [] self.test_events = [] - def load_raw(self, path=None, mode='', - epoch_duration=2, - band=[4.0, 40.0], - order=3): + def load_raw(self, path=None, mode='', + epoch_duration=2, + band=[4.0, 40.0], + order=3): ''' ''' set_log_level(verbose=False) - epoch_duration = np.round(np.array(epoch_duration) * self.fs).astype(int) + epoch_duration = np.round( + np.array(epoch_duration) * self.fs).astype(int) labels_folder = path + '/true_labels' if mode == 'train': @@ -52,7 +54,7 @@ def load_raw(self, path=None, mode='', elif mode == 'test': data_files = glob.glob(path + '/*E.gdf') labels_files = glob.glob(labels_folder + '/*E.mat') - + data_files.sort() labels_files.sort() @@ -61,67 +63,71 @@ def load_raw(self, path=None, mode='', Y = [] for subj in subjects: - x, y = self._get_epoched(data_files[subj], + x, y = self._get_epoched(data_files[subj], labels_files[subj], epoch_duration, band, order) X.append(x) - Y.append(y) - + Y.append(y) + samples, channels, trials = X[0].shape X = np.array(X) Y = np.array(Y) return X, Y def generate_set(self, load_path=None, epoch=2, - band=[4.,40.], - order=3, - save_folder=None): + band=[4., 40.], + order=3, + save_folder=None): + ''' ''' - ''' self.epochs, self.y = self.load_raw(load_path, 'train', - epoch, band, - order - ) - - self.test_epochs, self.test_y = self.load_raw(load_path, - 'test', epoch, - band, order - ) - + epoch, band, + order + ) + + self.test_epochs, self.test_y = self.load_raw(load_path, + 'test', epoch, + band, order + ) + self.paradigm = self._get_paradigm() self.save_set(save_folder) - + def get_path(self): NotImplementedError - def _get_epoched(self, data_file, label_file, + def _get_epoched(self, data_file, label_file, dur, band, order): ''' ''' Left = 769 Right = 770 - Foot = 771 + Foot = 771 Tongue = 772 Unkown = 783 raw = read_raw_edf(data_file) lb = loadmat(label_file) y = lb['classlabel'].ravel() - signal = np.divide(raw.get_data()[:22,:].T, raw._raw_extras[0]['units'][:22]).T # keep EEG channels only, multiple by gain + # keep EEG channels only, multiple by gain + signal = np.divide( + raw.get_data()[:22, :].T, raw._raw_extras[0]['units'][:22]).T # get events events_raw = raw._raw_extras[0]['events'] events_pos = events_raw[1] events_desc = events_raw[2] - ev_idx = (events_desc == Left) | (events_desc == Right) |(events_desc == Foot) | (events_desc == Tongue) | (events_desc == Unkown) + ev_idx = (events_desc == Left) | (events_desc == Right) | ( + events_desc == Foot) | (events_desc == Tongue) | (events_desc == Unkown) ev_desc = events_desc[ev_idx] ev_pos = events_pos[ev_idx] # filter signal = bandpass(signal, band, self.fs, order=order) # epoch - epochs = eeg_epoch(signal.T, dur, ev_pos) - self.subjects.append(self._get_subjects(raw._raw_extras[0]['subject_info'])) + epochs = eeg_epoch(signal.T, dur, ev_pos) + self.subjects.append(self._get_subjects( + raw._raw_extras[0]['subject_info'])) return epochs, y def _get_labels(self): @@ -131,21 +137,20 @@ def _get_labels(self): def _get_subjects(self, raw_subject_info): ''' - ''' + ''' return Subject(id=raw_subject_info['id'], gender=raw_subject_info['sex'], age=raw_subject_info['age'], handedness=raw_subject_info['handedness'], - condition=raw_subject_info['medication'] + condition=raw_subject_info['medication'] ) def _get_paradigm(self): ''' ''' - return MotorImagery(title='Comp_IV_2a', - stimulation=3000, - break_duration=2000, - repetition=72, stimuli=4, + return MotorImagery(title='Comp_IV_2a', + stimulation=3000, + break_duration=2000, + repetition=72, stimuli=4, phrase='' ) - \ No newline at end of file diff --git a/aawedha/io/exoskeleton.py b/aawedha/io/exoskeleton.py index ff9b2ec..05604fc 100644 --- a/aawedha/io/exoskeleton.py +++ b/aawedha/io/exoskeleton.py @@ -10,46 +10,47 @@ class Exoskeleton(DataSet): """ - SSVEP-based BCI of subjects operating an upper limb exoskeleton during a shared control - - References : - [1] Emmanuel K. Kalunga, Sylvain Chevallier, Olivier Rabreau, Eric Monacelli. Hybrid interface: - Integrating BCI in Multimodal Human-Machine Interfaces. IEEE/ASME International Conference + SSVEP-based BCI of subjects operating an upper limb exoskeleton during a shared control + + References : + [1] Emmanuel K. Kalunga, Sylvain Chevallier, Olivier Rabreau, Eric Monacelli. Hybrid interface: + Integrating BCI in Multimodal Human-Machine Interfaces. IEEE/ASME International Conference on Advanced Intelligent Mechatronics (AIM), 2014, Besancon, France. - [2] Emmanuel Kalunga, Sylvain Chevallier, Quentin Barthelemy. Data augmentation in - Riemannian space for Brain-Computer Interfaces, STAMLINS (ICML workshop), + [2] Emmanuel Kalunga, Sylvain Chevallier, Quentin Barthelemy. Data augmentation in + Riemannian space for Brain-Computer Interfaces, STAMLINS (ICML workshop), 2015, Lille, France. - [3] Emmanuel K. Kalunga, Sylvain Chevallier, Quentin Barthelemy. Online SSVEP-based BCI + [3] Emmanuel K. Kalunga, Sylvain Chevallier, Quentin Barthelemy. Online SSVEP-based BCI using Riemannian Geometry. Neurocomputing, 2016. arXiv research report on arXiv:1501.03227. """ + def __init__(self): super().__init__(title='Exoskeleton_SSVEP', - ch_names=['Oz','O1','O2','PO3','POz','PO7','PO8','PO4'], - fs=256, - doi='http://dx.doi.org/10.1109/AIM.2014.6878132') + ch_names=['Oz', 'O1', 'O2', 'PO3', + 'POz', 'PO7', 'PO8', 'PO4'], + fs=256, + doi='http://dx.doi.org/10.1109/AIM.2014.6878132') self.test_epochs = [] self.test_y = [] self.test_events = [] - - def load_raw(self, path=None, mode='', epoch_duration=[2,5], - band=[5.0, 45.0], order=6, augment=False): + def load_raw(self, path=None, mode='', epoch_duration=[2, 5], + band=[5.0, 45.0], order=6, augment=False): ''' ''' - files_list = sorted(glob.glob(path + '/s*')) + files_list = sorted(glob.glob(path + '/s*')) n_subjects = 12 X = [] Y = [] for subj in range(n_subjects): - session = glob.glob(files_list[subj]+'/*.pz') + session = glob.glob(files_list[subj] + '/*.pz') if mode == 'train': - records = np.arange(0,len(session)-1).tolist() + records = np.arange(0, len(session) - 1).tolist() elif mode == 'test': - records = [len(session)-1] + records = [len(session) - 1] x, y = self._get_epoched(session, records, - epoch_duration, band, - order, augment) + epoch_duration, band, + order, augment) X.append(x) Y.append(y) @@ -59,35 +60,35 @@ def load_raw(self, path=None, mode='', epoch_duration=[2,5], # Y = np.array(Y) return X, Y - def generate_set(self, load_path=None, epoch=[2,5], - band=[5.,45.], - order=6, - save_folder=None, - augment=False): + def generate_set(self, load_path=None, epoch=[2, 5], + band=[5., 45.], + order=6, + save_folder=None, + augment=False): ''' ''' self.epochs, self.y = self.load_raw(load_path, 'train', - epoch, band, - order, augment - ) + epoch, band, + order, augment + ) - self.test_epochs, self.test_y = self.load_raw(load_path, - 'test', epoch, - band, order, augment - ) + self.test_epochs, self.test_y = self.load_raw(load_path, + 'test', epoch, + band, order, augment + ) self.subjects = self._get_subjects(n_subjects=12) self.paradigm = self._get_paradigm() self.save_set(save_folder) - + def _get_epoched(self, files=[], records=[], - epoch=[2,5], band=[5.,45.], - order=6, augment=False): + epoch=[2, 5], band=[5., 45.], + order=6, augment=False): ''' ''' OVTK_StimulationId_VisualStimulationStart = 32779 labels = [33024, 33025, 33026, 33027] - if type(epoch) is list: + if isinstance(epoch, list): epoch = np.array(epoch).astype(int) * self.fs else: epoch = np.array([0, epoch]).astype(int) * self.fs @@ -96,52 +97,52 @@ def _get_epoched(self, files=[], records=[], for sess in records: f = gzip.open(files[sess], 'rb') df = pickle.load(f, encoding='latin1') - raw_signal = df['raw_signal'] * 1000 # samples, channels - event_pos = df['event_pos'].reshape((df['event_pos'].shape[0])) + raw_signal = df['raw_signal'] * 1000 # samples, channels + event_pos = df['event_pos'].reshape((df['event_pos'].shape[0])) event_type = df['event_type'].reshape((df['event_type'].shape[0])) - desc_idx = np.logical_or.reduce([event_type == lb for lb in labels]) + desc_idx = np.logical_or.reduce( + [event_type == lb for lb in labels]) desc = event_type[desc_idx] # y.append(desc - 33023) - pos = event_pos[event_type == OVTK_StimulationId_VisualStimulationStart] + pos = event_pos[event_type == + OVTK_StimulationId_VisualStimulationStart] raw_signal = bandpass(raw_signal, band, self.fs, order) if augment: stimulation = 5 * self.fs - augmented = np.floor(stimulation /np.diff(epoch))[0].astype(int) - v = [eeg_epoch(raw_signal, epoch + np.diff(epoch)*i, pos) for i in range(augmented)] + augmented = np.floor( + stimulation / np.diff(epoch))[0].astype(int) + v = [eeg_epoch(raw_signal, epoch + np.diff(epoch) * i, pos) + for i in range(augmented)] epchs = np.concatenate(v, axis=2) - y.append(np.tile(desc-33023, augmented)) + y.append(np.tile(desc - 33023, augmented)) else: y.append(desc - 33023) epchs = eeg_epoch(raw_signal, epoch, pos) - x.append(epchs) - + x.append(epchs) + if len(x) > 1: - x = np.concatenate(x, axis=-1) + x = np.concatenate(x, axis=-1) else: x = np.array(x).squeeze() - + y = np.array(self.flatten(y)) - - return x, y #epochs and labels - def _get_subjects(self, n_subjects=0): - return [Subject(id='S'+str(s),gender='M',age=0, handedness='') - for s in range(1, n_subjects+1)] + return x, y # epochs and labels + def _get_subjects(self, n_subjects=0): + return [Subject(id='S' + str(s), gender='M', age=0, handedness='') + for s in range(1, n_subjects + 1)] def _get_paradigm(self): return SSVEP(title='SSVEP_LED', control='Async', - stimulation=5000, - break_duration=3000, repetition=8, - stimuli=3, phrase='', - stim_type='ON_OFF', frequencies=['idle','13','21','17'], - ) + stimulation=5000, + break_duration=3000, repetition=8, + stimuli=3, phrase='', + stim_type='ON_OFF', frequencies=['idle', '13', '21', '17'], + ) def flatten(self, list_of_lists=[]): return [item for sublist in list_of_lists for item in sublist] def get_path(self): NotImplementedError - - - \ No newline at end of file diff --git a/aawedha/io/freiburg_online.py b/aawedha/io/freiburg_online.py index b35c7bc..c8e53bc 100644 --- a/aawedha/io/freiburg_online.py +++ b/aawedha/io/freiburg_online.py @@ -20,41 +20,39 @@ class FreiburgOnline(DataSet): def __init__(self): super().__init__(title='Freiburg_ERP_online', - ch_names=['C3','C4', 'CP1','CP2', - 'CP5','CP6','Cz','F10', - 'F3','F4','F7','F8', - 'F9','FC1','FC2','FC5', - 'FC6','Fp1', 'Fp2', 'Fz', - 'O1','O2','P10','P3', - 'P4','P7','P8','P9', - 'Pz','T7', 'T8'], - fs=100, - doi='https://doi.org/10.1371/journal.pone.0175856' - ) - - + ch_names=['C3', 'C4', 'CP1', 'CP2', + 'CP5', 'CP6', 'Cz', 'F10', + 'F3', 'F4', 'F7', 'F8', + 'F9', 'FC1', 'FC2', 'FC5', + 'FC6', 'Fp1', 'Fp2', 'Fz', + 'O1', 'O2', 'P10', 'P3', + 'P4', 'P7', 'P8', 'P9', + 'Pz', 'T7', 'T8'], + fs=100, + doi='https://doi.org/10.1371/journal.pone.0175856' + ) def load_raw(self, path=None): ''' ''' - files_list = glob.glob(path+'/S*.mat') - files_list.sort() + files_list = sorted(glob.glob(path + '/S*.mat')) n_subjects = 13 X = [] Y = [] for subj in range(n_subjects): data = loadmat(files_list[subj]) - X.append(data['epo']['x'][0][0][20:,:,:]) + X.append(data['epo']['x'][0][0][20:, :, :]) Y.append(data['epo']['y'][0][0].argmin(axis=0)) del data - + samples, channels, trials = X[0].shape - X = np.array(X).reshape((n_subjects, samples, channels, trials), order='F') + X = np.array(X).reshape( + (n_subjects, samples, channels, trials), order='F') Y = np.array(Y).reshape((n_subjects, trials), order='F') - + return X, Y - - def generate_set(self, load_path=None, + + def generate_set(self, load_path=None, save_folder=None): ''' ''' @@ -66,23 +64,17 @@ def generate_set(self, load_path=None, def _get_subjects(self, n_subjects=0): ''' ''' - return [Subject(id='S'+str(s),gender='M',age=0, handedness='') - for s in range(1, n_subjects+1)] - + return [Subject(id='S' + str(s), gender='M', age=0, handedness='') + for s in range(1, n_subjects + 1)] def _get_paradigm(self): ''' ''' - return ERP(title='ERP_FRIEBURG', stimulation=100, + return ERP(title='ERP_FRIEBURG', stimulation=100, break_duration=150, repetition=68, - phrase='Franzy jagt im komplett verwahrlosten Taxi quer durch Freiburg', + phrase='Franzy jagt im komplett verwahrlosten Taxi quer durch Freiburg', flashing_mode='SC', speller=[]) def get_path(self): NotImplementedError - - - - - diff --git a/aawedha/io/inria_ern.py b/aawedha/io/inria_ern.py index f9e701d..fce67f9 100644 --- a/aawedha/io/inria_ern.py +++ b/aawedha/io/inria_ern.py @@ -21,87 +21,91 @@ class Inria_ERN(DataSet): def __init__(self): super().__init__(title='Inria_ERN', - ch_names=['Fp1', 'Fp2', 'AF7','F3', 'AF4','F8','F7','F5','F3','F1', - 'Fz', 'F2', 'F4', 'F6','F8','FT7','FC5','FC3','FC1', - 'FCz', 'FC2', 'FC4', 'FC6', 'FT8','T7','C5','C3','C1', - 'Cz','C2','C4','C6','T8','TP7','CP5','CP3','CP1','CPz', - 'CP2','CP4','CP6','TP8','P7', 'P5','P3','P1','Pz', - 'P2','P4','P6','P8','PO7','POz','P08','O1','O2'], - fs=200, - doi = 'http://dx.doi.org/10.1155/2012/578295' - ) + ch_names=['Fp1', 'Fp2', 'AF7', 'F3', 'AF4', 'F8', 'F7', 'F5', 'F3', 'F1', + 'Fz', 'F2', 'F4', 'F6', 'F8', 'FT7', 'FC5', 'FC3', 'FC1', + 'FCz', 'FC2', 'FC4', 'FC6', 'FT8', 'T7', 'C5', 'C3', 'C1', + 'Cz', 'C2', 'C4', 'C6', 'T8', 'TP7', 'CP5', 'CP3', 'CP1', 'CPz', + 'CP2', 'CP4', 'CP6', 'TP8', 'P7', 'P5', 'P3', 'P1', 'Pz', + 'P2', 'P4', 'P6', 'P8', 'PO7', 'POz', 'P08', 'O1', 'O2'], + fs=200, + doi='http://dx.doi.org/10.1155/2012/578295' + ) self.test_epochs = [] self.test_y = [] self.test_events = [] - - def load_raw(self, path=None, epoch_duration=1, - band=[1.0, 40.0], order=5): + def load_raw(self, path=None, epoch_duration=1, + band=[1.0, 40.0], order=5): """ """ - list_of_tr_files = glob.glob(path + '/train/Data_*.csv') - list_of_tr_files.sort() + list_of_tr_files = sorted(glob.glob(path + '/train/Data_*.csv')) list_of_ts_files = glob.glob(path + '/test/Data_*.csv') - list_of_ts_files.sort() - + list_of_ts_files.sort() n_subjects = 26 - epoch_duration = round(epoch_duration * self.fs) + epoch_duration = round(epoch_duration * self.fs) channels = len(self.ch_names) epochs = 340 - + X = self._get_epoched(list_of_tr_files, epoch_duration, band, order) - X_test = self._get_epoched(list_of_ts_files, epoch_duration, band, order) - + X_test = self._get_epoched( + list_of_ts_files, epoch_duration, band, order) + train_subjects = 16 test_subjects = 10 - + samples = X[0].shape[0] # X = np.array(X).transpose((2,1,0)).reshape((n_subjects, epoch_duration, channels, epochs)) # X = np.array(X).transpose((2,1,0)).reshape((n_subjects, epoch_duration, channels, epochs), order='F') - X = np.array(X).reshape((train_subjects, epochs, samples, channels), order='F').transpose((0,2,3,1)) - X_test = np.array(X_test).reshape((test_subjects, epochs, samples, channels), order='F').transpose((0,2,3,1)) + X = np.array(X).reshape((train_subjects, epochs, samples, + channels), order='F').transpose((0, 2, 3, 1)) + X_test = np.array(X_test).reshape( + (test_subjects, epochs, samples, channels), order='F').transpose((0, 2, 3, 1)) - labels = np.genfromtxt(path + '/train//TrainLabels.csv',delimiter=',',skip_header=1)[:,1] + labels = np.genfromtxt( + path + '/train//TrainLabels.csv', delimiter=',', skip_header=1)[:, 1] Y = labels.reshape((train_subjects, epochs)) - labels_test = np.genfromtxt(path+'/test//true_labels.csv',delimiter=',') + labels_test = np.genfromtxt( + path + '/test//true_labels.csv', delimiter=',') Y_test = labels_test.reshape((test_subjects, epochs)) - - return X, Y, X_test, Y_test # : 4-D array : (subject, samples, channels, epoch) + + # : 4-D array : (subject, samples, channels, epoch) + return X, Y, X_test, Y_test def _get_epoched(self, files_list=None, - epoch=1, band=[1.,40.], order=5): + epoch=1, band=[1., 40.], order=5): ''' ''' if not files_list: return None - + X = [] for f in files_list: sig = np.array(pd.io.parsers.read_csv(f)) - eeg = sig[:,1:-2] - trigger = sig[:,-1] + eeg = sig[:, 1:-2] + trigger = sig[:, -1] signal = bandpass(eeg, band, self.fs, order) - idxFeedback = np.where(trigger==1)[0] + idxFeedback = np.where(trigger == 1)[0] for idx in idxFeedback: - X.append(signal[idx:idx+epoch, :]) - - del sig # saving some RAM + X.append(signal[idx:idx + epoch, :]) + + del sig # saving some RAM return X - def generate_set(self, load_path=None, epoch=1, - band=[1.0, 40.0], - order=5, save_folder=None): + def generate_set(self, load_path=None, epoch=1, + band=[1.0, 40.0], + order=5, save_folder=None): """ """ - self.epochs, self.y, self.test_epochs, self.test_y = self.load_raw(load_path, epoch, band, order) + self.epochs, self.y, self.test_epochs, self.test_y = self.load_raw( + load_path, epoch, band, order) self.subjects = self._get_subjects(n_subjects=16) self.paradigm = self._get_paradigm() self.save_set(save_folder) - + def get_path(self): NotImplementedError @@ -109,19 +113,9 @@ def _get_events(self): NotImplementedError def _get_subjects(self, n_subjects=0): - return [Subject(id='S'+str(s),gender='M',age=0, handedness='') - for s in range(1, n_subjects+1)] + return [Subject(id='S' + str(s), gender='M', age=0, handedness='') + for s in range(1, n_subjects + 1)] def _get_paradigm(self): return ERP(title='ERP_ERN', stimulation=60, break_duration=50, repetition=12, - phrase='', flashing_mode='RC') - - - - - - - - - - + phrase='', flashing_mode='RC') diff --git a/aawedha/io/laresi.py b/aawedha/io/laresi.py index 3a6241c..dedfee5 100644 --- a/aawedha/io/laresi.py +++ b/aawedha/io/laresi.py @@ -9,14 +9,15 @@ import glob import pickle import os - + OV_ExperimentStart = 32769 -OV_ExperimentStop = 32770 -Base_Stimulations = 33024 +OV_ExperimentStop = 32770 +Base_Stimulations = 33024 OV_Target = 33285 -OV_NonTarget = 33286 +OV_NonTarget = 33286 OV_TrialStart = 32773 -OV_TrialStop = 32774 +OV_TrialStop = 32774 + class LaresiHybrid: """ @@ -24,68 +25,69 @@ class LaresiHybrid: Reference: [] coming soon... """ - def __init__(self): self.title = 'Laresi_hybrid' - self.ch_names = ['Pz','PO5','PO3','POz','PO4', - 'PO6','O1','Oz','O2','Fz', - 'Cz','P3','P4','FCz'] + self.ch_names = ['Pz', 'PO5', 'PO3', 'POz', 'PO4', + 'PO6', 'O1', 'Oz', 'O2', 'Fz', + 'Cz', 'P3', 'P4', 'FCz'] self.fs = 512 - self.doi='' - self.erp_set = LaresiERP(title='Hybrid_ERP', ch_names=self.ch_names, fs=self.fs) - self.ssvep_set = LaresiSSVEP(title='Hybrid_SSVEP', ch_names=self.ch_names, fs=self.fs) - self.paradigm = HybridLARESI() + self.doi = '' + self.erp_set = LaresiERP( + title='Hybrid_ERP', ch_names=self.ch_names, fs=self.fs) + self.ssvep_set = LaresiSSVEP( + title='Hybrid_SSVEP', ch_names=self.ch_names, fs=self.fs) + self.paradigm = HybridLARESI() def load_raw(self, path=None): ''' ''' - files_list = sorted(glob.glob(path + '/s*')) + files_list = sorted(glob.glob(path + '/s*')) n_subjects = 1 cnts, infos = [], [] for subj in range(n_subjects): - train_session = glob.glob(files_list[subj]+'/*calib*.csv') - test_session = glob.glob(files_list[subj]+'/*online*.csv') + train_session = glob.glob(files_list[subj] + '/*calib*.csv') + test_session = glob.glob(files_list[subj] + '/*online*.csv') cnt, info = self._read_raw(train_session[0]) cnt_test, info_test = self._read_raw_online(test_session[0]) cnts = [cnt, cnt_test] - infos = [info, info_test] - return cnts, infos # cnt: EEG, info : dict - + infos = [info, info_test] + return cnts, infos # cnt: EEG, info : dict + def generate_set(self, load_path=None, save_folder=None, - epoch=[0.,0.7, 0.,1.0], band=[1,10,5,45], - order=[2,6]): + epoch=[0., 0.7, 0., 1.0], band=[1, 10, 5, 45], + order=[2, 6]): ''' ''' - n_subjects = 1 # 1 for now + n_subjects = 1 # 1 for now for subj in range(n_subjects): cnts, infos = self.load_raw(load_path) - self.erp_set.generate_set(cnts, [infos[0]['ERP'],infos[1]], - epoch=[epoch[0],epoch[1]], - band=[band[0],band[1]], order=order[0]) - self.ssvep_set.generate_set(cnts, [infos[0]['SSVEP'],infos[1]], - epoch=[epoch[2],epoch[3]], - band=[band[2],band[3]], order=order[1]) - + self.erp_set.generate_set(cnts, [infos[0]['ERP'], infos[1]], + epoch=[epoch[0], epoch[1]], + band=[band[0], band[1]], order=order[0]) + self.ssvep_set.generate_set(cnts, [infos[0]['SSVEP'], infos[1]], + epoch=[epoch[2], epoch[3]], + band=[band[2], band[3]], order=order[1]) + self.erp_set._cat_lists() self.ssvep_set._cat_lists() self.save_set(save_folder) - + def save_set(self, save_folder=None): ''' ''' - # save dataset + # save dataset if not os.path.isdir(save_folder): os.mkdir(save_folder) if not self.title: self.title = 'unnamed_set' - fname = save_folder + '/' + self.title +'.pkl' - print(f'Saving dataset {self.title} to destination: {fname}') + fname = save_folder + '/' + self.title + '.pkl' + print(f'Saving dataset {self.title} to destination: {fname}') f = open(fname, 'wb') - pickle.dump(self, f, protocol=pickle.HIGHEST_PROTOCOL) + pickle.dump(self, f, protocol=pickle.HIGHEST_PROTOCOL) f.close() def load_set(self, file_name=None): @@ -98,24 +100,24 @@ def load_set(self, file_name=None): raise FileNotFoundError f.close() return data - + def _read_raw(self, path=None): ''' - ''' + ''' cnt = [] info = {} if path: cnt, raw_desc, raw_pos = self._convert_raw(path) - start_interv = raw_pos[raw_desc==OV_ExperimentStart] - end_interv = raw_pos[raw_desc==OV_ExperimentStop] + start_interv = raw_pos[raw_desc == OV_ExperimentStart] + end_interv = raw_pos[raw_desc == OV_ExperimentStop] if end_interv.size > 1: # Hybrid paradigm - start2 = start_interv[start_interv> end_interv[0]][0] - starts = np.array([start_interv[0], start2]) - else: + start2 = start_interv[start_interv > end_interv[0]][0] + starts = np.array([start_interv[0], start2]) + else: starts = start_interv[0] - info = self._get_info(raw_desc, raw_pos, starts, end_interv) - + info = self._get_info(raw_desc, raw_pos, starts, end_interv) + return cnt, info def _read_raw_online(self, path=None): @@ -127,56 +129,61 @@ def _read_raw_online(self, path=None): cnt_test, raw_desc, raw_pos = self._convert_raw(path) info_test['raw_desc'] = raw_desc info_test['raw_pos'] = raw_pos - + return cnt_test, info_test - def _convert_raw(self, path=None): ''' ''' markers = {} raw = pd.read_csv(path) - cnt = raw[self.ch_names].to_numpy() # samples x channels + cnt = raw[self.ch_names].to_numpy() # samples x channels raw_desc = self._get_events(raw, 'Event Id') raw_pos = self._get_events(raw, 'Event Date') - return cnt, raw_desc, raw_pos - + return cnt, raw_desc, raw_pos + def _get_info(self, raw_desc, raw_pos, starts, ends): ''' ''' - + info = {} for i in range(len(starts)): - - idx = np.logical_and(raw_pos >=starts[i], raw_pos <=ends[i]) + + idx = np.logical_and(raw_pos >= starts[i], raw_pos <= ends[i]) desc = raw_desc[idx] pos = raw_pos[idx] - evs = np.unique(desc) - stimuli = np.sum(np.logical_and(evs > Base_Stimulations, evs < OV_Target)) - id_stim = np.logical_and(desc > Base_Stimulations, desc <= Base_Stimulations + stimuli) + evs = np.unique(desc) + stimuli = np.sum(np.logical_and( + evs > Base_Stimulations, evs < OV_Target)) + id_stim = np.logical_and( + desc > Base_Stimulations, desc <= Base_Stimulations + stimuli) desc = desc[id_stim] - Base_Stimulations pos = np.floor(pos[id_stim] * self.fs).astype(int) - - if np.any(evs == OV_Target): # erp + + if np.any(evs == OV_Target): # erp info['ERP'] = {} - y = raw_desc[np.logical_or(raw_desc==OV_Target, raw_desc==OV_NonTarget)] - y[y==OV_Target] = 1 - y[y==OV_NonTarget] = 0 - info['ERP']['desc'] = desc + y = raw_desc[np.logical_or( + raw_desc == OV_Target, raw_desc == OV_NonTarget)] + y[y == OV_Target] = 1 + y[y == OV_NonTarget] = 0 + info['ERP']['desc'] = desc info['ERP']['pos'] = pos info['ERP']['y'] = y - info['ERP']['session_interval'] = (np.floor( [starts[i], ends[i]]) * self.fs).astype(int) - else: #ssvep + info['ERP']['session_interval'] = ( + np.floor([starts[i], ends[i]]) * self.fs).astype(int) + else: # ssvep y = desc - info['SSVEP'] = {} - info['SSVEP']['desc'] = desc + info['SSVEP'] = {} + info['SSVEP']['desc'] = desc info['SSVEP']['pos'] = pos info['SSVEP']['y'] = y - info['SSVEP']['session_interval'] = (np.floor( [starts[i], ends[i]]) * self.fs).astype(int) - - # session_interval = np.floor( [starts, ends] ) * self.fs # begin, end of session + info['SSVEP']['session_interval'] = ( + np.floor([starts[i], ends[i]]) * self.fs).astype(int) + + # session_interval = np.floor( [starts, ends] ) * self.fs # begin, + # end of session return info - + def _get_events(self, dataframe, key): events_id = dataframe[key].notna() events = dataframe[key].loc[events_id] @@ -186,8 +193,8 @@ def _get_events(self, dataframe, key): return ev def _get_subjects(self, n_subjects=0): - return [Subject(id='S'+str(s),gender='M',age=0, handedness='') - for s in range(1, n_subjects+1)] + return [Subject(id='S' + str(s), gender='M', age=0, handedness='') + for s in range(1, n_subjects + 1)] # single paradigm class @@ -198,32 +205,34 @@ def __init__(self, title='', ch_names=[], fs=None, doi=''): self.test_epochs = [] self.test_y = [] self.test_events = [] - - def load_raw(self, cnts=[], infos=[], epoch_duration=[0,.7], - band=[1., 10.], order=2, - augment=False): + + def load_raw(self, cnts=[], infos=[], epoch_duration=[0, .7], + band=[1., 10.], order=2, + augment=False): ''' ''' - epochs = self._get_epochs(cnts[0],infos[0], epoch_duration, band, order) - test_epochs, test_y = self._get_epochs_test(cnts[1], infos[1], epoch_duration, band, order) + epochs = self._get_epochs( + cnts[0], infos[0], epoch_duration, band, order) + test_epochs, test_y = self._get_epochs_test( + cnts[1], infos[1], epoch_duration, band, order) return epochs, test_epochs, test_y def generate_set(self, - cnts=[], - infos=[], - epoch=[0,.7], - band=[1., 10], - order=2, - augment=False): + cnts=[], + infos=[], + epoch=[0, .7], + band=[1., 10], + order=2, + augment=False): ''' ''' epochs, test_epochs, y_test = self.load_raw(cnts, infos, - epoch, band, - order, augment - ) + epoch, band, + order, augment + ) self.epochs.append(epochs) - self.y.append(infos[0]['y']) - self.test_epochs.append(test_epochs) + self.y.append(infos[0]['y']) + self.test_epochs.append(test_epochs) self.test_y.append(y_test) self.subjects = [] self.paradigm = [] @@ -234,16 +243,17 @@ def get_path(self): def _get_epochs(self, cnt, info, epoch, band, order): ''' ''' - epoch = np.round(np.array(epoch)* self.fs).astype(int) - signal = cnt[info['session_interval'][0]:info['session_interval'][1] , :] - signal = bandpass(signal, band, self.fs, order) - cnt[info['session_interval'][0]:info['session_interval'][1] ,:] = signal + epoch = np.round(np.array(epoch) * self.fs).astype(int) + signal = cnt[info['session_interval'] + [0]:info['session_interval'][1], :] + signal = bandpass(signal, band, self.fs, order) + cnt[info['session_interval'][0]:info['session_interval'][1], :] = signal eps = eeg_epoch(cnt, epoch, info['pos']) return eps @abstractmethod def _get_epochs_test(self): - pass + pass def _cat_lists(self): ''' @@ -255,79 +265,87 @@ def _cat_lists(self): class LaresiERP(LaresiEEG): - + def _get_epochs_test(self, cnt=None, markers={}, - epoch=[0.,0.7], band=[1,10], order=2): + epoch=[0., 0.7], band=[1, 10], order=2): ''' ''' raw_desc = markers['raw_desc'] raw_pos = markers['raw_pos'] - - start = raw_pos[raw_desc==OV_TrialStart] - end = raw_pos[raw_desc==OV_TrialStop] + + start = raw_pos[raw_desc == OV_TrialStart] + end = raw_pos[raw_desc == OV_TrialStop] erp_start = np.round(start[::2] * self.fs).astype(int) erp_end = np.round(end[::2] * self.fs).astype(int) - erp_tr = np.logical_or(raw_desc==OV_Target, raw_desc==OV_NonTarget) - erp_mrk = np.round(raw_pos[erp_tr]* self.fs).astype(int) - + erp_tr = np.logical_or(raw_desc == OV_Target, raw_desc == OV_NonTarget) + erp_mrk = np.round(raw_pos[erp_tr] * self.fs).astype(int) + y_erp = raw_desc[erp_tr] - y_erp[y_erp==OV_Target] = 1 - y_erp[y_erp==OV_NonTarget] = 0 + y_erp[y_erp == OV_Target] = 1 + y_erp[y_erp == OV_NonTarget] = 0 erp_epochs = [] - epoch = np.round(np.array(epoch)* self.fs).astype(int) + epoch = np.round(np.array(epoch) * self.fs).astype(int) for tr in range(len(erp_start)): # filter, epoch, append - idx = np.logical_and(erp_mrk >= erp_start[tr], erp_mrk <= erp_end[tr]) - mrk = erp_mrk[idx] - erp_start[tr] - erp_signal = bandpass(cnt[erp_start[tr]:erp_end[tr],:], band, self.fs, order) - erp_epochs.append(eeg_epoch(erp_signal, epoch, mrk) ) - - samples, channels, trials = erp_epochs[0].shape + idx = np.logical_and( + erp_mrk >= erp_start[tr], erp_mrk <= erp_end[tr]) + mrk = erp_mrk[idx] - erp_start[tr] + erp_signal = bandpass( + cnt[erp_start[tr]:erp_end[tr], :], band, self.fs, order) + erp_epochs.append(eeg_epoch(erp_signal, epoch, mrk)) + + samples, channels, trials = erp_epochs[0].shape blocks = len(erp_epochs) - eps = np.array(erp_epochs).transpose((1,2,3,0)).reshape((samples, channels, trials*blocks)) + eps = np.array(erp_epochs).transpose((1, 2, 3, 0)).reshape( + (samples, channels, trials * blocks)) return eps, y_erp + class LaresiSSVEP(LaresiEEG): - + def _get_epochs_test(self, cnt=None, markers={}, - epoch=[0.,1.0], band=[5,45], order=6): + epoch=[0., 1.0], band=[5, 45], order=6): ''' ''' raw_desc = markers['raw_desc'] raw_pos = markers['raw_pos'] - start = raw_pos[raw_desc==OV_TrialStart] - end = raw_pos[raw_desc==OV_TrialStop] + start = raw_pos[raw_desc == OV_TrialStart] + end = raw_pos[raw_desc == OV_TrialStop] ssvep_start = start[1::2] ssvep_end = end[1::2] ev_desc = [] ev_pos = [] base = Base_Stimulations - for i in range(len(ssvep_start)): - idx1 = np.logical_and(raw_pos>ssvep_start[i],raw_posbase, raw_desc<=base+10) # 10 is an arbitary number - idx = np.logical_and(idx1,idx2) - ev_desc.append(raw_desc[idx][0]-base) - ev_pos.append(np.round(raw_pos[idx][0]*self.fs).astype(int)) - - ev_desc = np.array(ev_desc) # y + for i in range(len(ssvep_start)): + idx1 = np.logical_and( + raw_pos > ssvep_start[i], raw_pos < ssvep_end[i]) + # 10 is an arbitary number + idx2 = np.logical_and(raw_desc > base, raw_desc <= base + 10) + idx = np.logical_and(idx1, idx2) + ev_desc.append(raw_desc[idx][0] - base) + ev_pos.append(np.round(raw_pos[idx][0] * self.fs).astype(int)) + + ev_desc = np.array(ev_desc) # y ev_pos = np.array(ev_pos) ssvep_epochs = [] - epoch = np.round(np.array(epoch)* self.fs).astype(int) + epoch = np.round(np.array(epoch) * self.fs).astype(int) ssvep_start = np.round(ssvep_start * self.fs).astype(int) ssvep_end = np.round(ssvep_end * self.fs).astype(int) for tr in range(len(ev_desc)): # filter, epoch, append - mrk = ev_pos[tr] - ssvep_start[tr] - ssvep_signal = bandpass(cnt[ssvep_start[tr]:ssvep_end[tr],:], band, self.fs, order) - ssvep_epochs.append(eeg_epoch(ssvep_signal, epoch, np.array([mrk])).squeeze()) - - ssvep_epochs = np.array(ssvep_epochs).transpose((1,2,0)) - return ssvep_epochs, ev_desc \ No newline at end of file + mrk = ev_pos[tr] - ssvep_start[tr] + ssvep_signal = bandpass( + cnt[ssvep_start[tr]:ssvep_end[tr], :], band, self.fs, order) + ssvep_epochs.append( + eeg_epoch(ssvep_signal, epoch, np.array([mrk])).squeeze()) + + ssvep_epochs = np.array(ssvep_epochs).transpose((1, 2, 0)) + return ssvep_epochs, ev_desc diff --git a/aawedha/io/physionet_mi.py b/aawedha/io/physionet_mi.py index 8d51adf..ef19484 100644 --- a/aawedha/io/physionet_mi.py +++ b/aawedha/io/physionet_mi.py @@ -16,16 +16,17 @@ class PhysioNet_MI(DataSet): Motor Imagery dataset [1] Reference : - [1] Schalk, G., McFarland, D.J., Hinterberger, T., Birbaumer, N., Wolpaw, J.R. (2004) BCI2000: + [1] Schalk, G., McFarland, D.J., Hinterberger, T., Birbaumer, N., Wolpaw, J.R. (2004) BCI2000: A General-Purpose Brain-Computer Interface (BCI) System. IEEE TBME 51(6):1034-1043. - + [2] Goldberger AL, Amaral LAN, Glass L, Hausdorff JM, Ivanov PCh, Mark RG, Mietus JE, Moody GB, - Peng C-K, Stanley HE. (2000) PhysioBank, PhysioToolkit, and PhysioNet: Components of a New Research + Peng C-K, Stanley HE. (2000) PhysioBank, PhysioToolkit, and PhysioNet: Components of a New Research Resource for Complex Physiologic Signals. Circulation 101(23):e215-e220. """ + def __init__(self): - super().__init__(title='PhysioNet MI', + super().__init__(title='PhysioNet MI', fs=160, doi='https://doi.org/10.1109/TBME.2004.827072' ) @@ -36,10 +37,10 @@ def load_raw(self, path=None, epoch_duration=1, band=[1., 40.], f_order=6): ''' set_log_level(verbose=False) subjects = range(1, 110) - runs = [4, 8 , 12] # imagined hand tasks + runs = [4, 8, 12] # imagined hand tasks #event_id = dict(left=2, right=3) - ev_id = dict(T1=2, T2=3) - subjects = range(1,110) + ev_id = dict(T1=2, T2=3) + subjects = range(1, 110) if epoch_duration.__class__ is list: # t = [[1., 3.],[1.9, 3.9]] t = epoch_duration @@ -49,43 +50,48 @@ def load_raw(self, path=None, epoch_duration=1, band=[1., 40.], f_order=6): low_pass, high_pass = band epochs = [] labels = [] - for subj in subjects: - raw_names = ['{f}/S{s:03d}/S{s:03d}R{r:02d}.edf'.format(f=path, s=subj, r=r) for r in runs] - raw = concatenate_raws([read_raw_edf(f, preload=True) for f in raw_names]) + for subj in subjects: + raw_names = [ + '{f}/S{s:03d}/S{s:03d}R{r:02d}.edf'.format(f=path, s=subj, r=r) for r in runs] + raw = concatenate_raws( + [read_raw_edf(f, preload=True) for f in raw_names]) if raw.info['sfreq'] != 160: - raw.info['sfreq'] = 160 + raw.info['sfreq'] = 160 # strip channel names of "." characters raw.rename_channels(lambda x: x.strip('.')) # Apply band-pass filter - raw.filter(low_pass, high_pass, method='iir', iir_params=dict(order=f_order, ftype='butter')) + raw.filter(low_pass, high_pass, method='iir', + iir_params=dict(order=f_order, ftype='butter')) events, _ = events_from_annotations(raw, event_id=ev_id) picks = pick_types(raw.info, meg=False, eeg=True, stim=False, eog=False, - exclude='bads') + exclude='bads') epos = Epochs(raw, events, ev_id, tmin, tmax, proj=True, picks=picks, - baseline=None, preload=True) - lb = epos.events[:, -1] - 2 + baseline=None, preload=True) + lb = epos.events[:, -1] - 2 # epo = [epos.copy().crop(tmin=f[0], tmax=f[1]) for f in t] epo = np.vstack(epo) epo = np.transpose(epo, (2, 1, 0)) epochs.append(epo) labels.append(np.repeat(lb, 2)) - + self.ch_names = raw.info['ch_names'] - return epochs, labels # list of epochs: len(epochs) = subjects, epochs: samples, channels trials + # list of epochs: len(epochs) = subjects, epochs: samples, channels + # trials + return epochs, labels def generate_set(self, load_path=None, epoch=1, band=[1., 40.], - order=6, save_folder=None): + order=6, save_folder=None): self.epochs, self.y = self.load_raw(load_path, epoch, band, order) self.subjects = self._get_subjects(n_subjects=109) self.paradigm = self._get_paradigm() - + # save dataset if not os.path.isdir(save_folder): os.mkdir(save_folder) fileName = save_folder + '/physionet_mi.pkl' f = gzip.open(fileName, 'wb') - pickle.dump(self, f, protocol=pickle.HIGHEST_PROTOCOL) + pickle.dump(self, f, protocol=pickle.HIGHEST_PROTOCOL) f.close() def load_set(self, fileName=None): @@ -95,15 +101,14 @@ def load_set(self, fileName=None): else: raise FileNotFoundError f.close() - return data + return data def get_path(self): NotImplementedError def _get_subjects(self, n_subjects=0): - return [subject.Subject(id='S'+str(s),gender=None, age=0, handedness='') - for s in range(1, n_subjects+1)] + return [subject.Subject(id='S' + str(s), gender=None, age=0, handedness='') + for s in range(1, n_subjects + 1)] def _get_paradigm(self): return None - diff --git a/aawedha/io/san_diego.py b/aawedha/io/san_diego.py index c14de4f..15bdf46 100644 --- a/aawedha/io/san_diego.py +++ b/aawedha/io/san_diego.py @@ -17,73 +17,84 @@ class SanDiego(DataSet): Detecting Steady-State Visual Evoked Potentials," PLoS One, vol.10, no.10, e140703, 2015. ''' + def __init__(self): super().__init__(title='San_Diego', - ch_names=['PO7','PO3','POz','PO4','PO8','O1','Oz','O2'], - fs=256, - doi = 'http://dx.doi.org/10.1371/journal.pone.0140703' - ) + ch_names=['PO7', 'PO3', 'POz', + 'PO4', 'PO8', 'O1', 'Oz', 'O2'], + fs=256, + doi='http://dx.doi.org/10.1371/journal.pone.0140703' + ) - def load_raw(self, path=None, epoch_duration=1, - band=[5.0, 45.0], order=6, augment=False): - list_of_files = glob.glob(path + 's*.mat') - list_of_files.sort() + def load_raw(self, path=None, epoch_duration=1, + band=[5.0, 45.0], order=6, augment=False): + list_of_files = sorted(glob.glob(path + 's*.mat')) - epoch_duration = np.round(np.array(epoch_duration) * self.fs).astype(int) - onset = 39 # onset in samples - n_subjects = 10 + epoch_duration = np.round( + np.array(epoch_duration) * self.fs).astype(int) + onset = 39 # onset in samples + n_subjects = 10 X = [] Y = [] augmented = 0 for subj in range(n_subjects): data = loadmat(list_of_files[subj]) - eeg = data['eeg'].transpose((2,1,3,0)) #samples, channels, trials, targets + # samples, channels, trials, targets + eeg = data['eeg'].transpose((2, 1, 3, 0)) eeg = bandpass(eeg, band=band, fs=self.fs, order=order) if augment: augmented = 4 - v = [eeg[onset+(stride*self.fs):onset+(stride*self.fs)+epoch_duration, :, :, :] for stride in range(augmented)] + v = [eeg[onset + (stride * self.fs):onset + (stride * self.fs) + + epoch_duration, :, :, :] for stride in range(augmented)] eeg = np.concatenate(v, axis=2) samples, channels, blocks, targets = eeg.shape - y = np.tile(np.arange(1, targets+1), (int(blocks/augmented),1)) - y = np.tile(y, (1,augmented)) - y = y.reshape((1,blocks*targets),order='F') + y = np.tile(np.arange(1, targets + 1), + (int(blocks / augmented), 1)) + y = np.tile(y, (1, augmented)) + y = y.reshape((1, blocks * targets), order='F') del v else: - eeg = eeg[onset:onset+epoch_duration, :, :, :] + eeg = eeg[onset:onset + epoch_duration, :, :, :] samples, channels, blocks, targets = eeg.shape - y = np.tile(np.arange(1, targets+1), (blocks,1)) - y = y.reshape((1,blocks*targets),order='F') - - del data # save some RAM + y = np.tile(np.arange(1, targets + 1), (blocks, 1)) + y = y.reshape((1, blocks * targets), order='F') + + del data # save some RAM - X.append(eeg.reshape((samples, channels, blocks*targets),order='F')) + X.append( + eeg.reshape( + (samples, + channels, + blocks * + targets), + order='F')) Y.append(y) - X = np.array(X) Y = np.array(Y).squeeze() return X, Y - def generate_set(self, load_path=None, epoch=1, - band=[5.0, 45.0], - order=6, save_folder=None, augment=False): - self.epochs, self.y = self.load_raw(load_path, epoch, band, order, augment) + def generate_set(self, load_path=None, epoch=1, + band=[5.0, 45.0], + order=6, save_folder=None, augment=False): + self.epochs, self.y = self.load_raw( + load_path, epoch, band, order, augment) self.subjects = self._get_subjects(n_subjects=10) self.paradigm = self._get_paradigm() self.save_set(save_folder) - + def _get_subjects(self, n_subjects=0): - return [Subject(id='S'+str(s),gender='M',age=0, handedness='') - for s in range(1, n_subjects+1)] + return [Subject(id='S' + str(s), gender='M', age=0, handedness='') + for s in range(1, n_subjects + 1)] def _get_paradigm(self): return SSVEP(title='SSVEP_JFPM', stimulation=4000, break_duration=1000, repetition=15, - stimuli=12, phrase='', - stim_type='ON_OFF', frequencies=[9.25, 11.25, 13.25, 9.75, 11.75, 13.75, 10.25, - 12.25, 14.25, 10.75, 12.75, 14.75], - phase=[0.0, 0.0, 0.0, 0.5*np.pi, 0.5*np.pi, 0.5*np.pi, - np.pi, np.pi, np.pi, 1.5*np.pi, 1.5*np.pi, 1.5*np.pi] ) + stimuli=12, phrase='', + stim_type='ON_OFF', frequencies=[9.25, 11.25, 13.25, 9.75, 11.75, 13.75, 10.25, + 12.25, 14.25, 10.75, 12.75, 14.75], + phase=[0.0, 0.0, 0.0, 0.5 * np.pi, 0.5 * np.pi, 0.5 * np.pi, + np.pi, np.pi, np.pi, 1.5 * np.pi, 1.5 * np.pi, 1.5 * np.pi]) def get_path(self): - NotImplementedError \ No newline at end of file + NotImplementedError diff --git a/aawedha/io/tsinghua.py b/aawedha/io/tsinghua.py index c085266..8711f26 100644 --- a/aawedha/io/tsinghua.py +++ b/aawedha/io/tsinghua.py @@ -21,53 +21,62 @@ class Tsinghua(DataSet): def __init__(self): super().__init__(title='Tsinghua', - ch_names=['FP1','FPZ','FP2','AF3','AF4','F7','F5','F3', - 'F1','FZ','F2','F4','F6','F8','FT7','FC5', - 'FC3','FC1','FCz','FC2','FC4','FC6','FT8','T7', - 'C5','C3','C1','Cz','C2','C4','C6','T8','M1', - 'TP7','CP5','CP3','CP1','CPZ','CP2','CP4','CP6', - 'TP8','M2','P7','P5','P3','P1','PZ','P2', - 'P4','P6','P8','PO7','PO5','PO3','POz','PO4', - 'PO6','PO8','CB1','O1','Oz','O2','CB2'], - fs=250, - doi='http://dx.doi.org/10.1073/pnas.1508080112' - ) + ch_names=['FP1', 'FPZ', 'FP2', 'AF3', 'AF4', 'F7', 'F5', 'F3', + 'F1', 'FZ', 'F2', 'F4', 'F6', 'F8', 'FT7', 'FC5', + 'FC3', 'FC1', 'FCz', 'FC2', 'FC4', 'FC6', 'FT8', 'T7', + 'C5', 'C3', 'C1', 'Cz', 'C2', 'C4', 'C6', 'T8', 'M1', + 'TP7', 'CP5', 'CP3', 'CP1', 'CPZ', 'CP2', 'CP4', 'CP6', + 'TP8', 'M2', 'P7', 'P5', 'P3', 'P1', 'PZ', 'P2', + 'P4', 'P6', 'P8', 'PO7', 'PO5', 'PO3', 'POz', 'PO4', + 'PO6', 'PO8', 'CB1', 'O1', 'Oz', 'O2', 'CB2'], + fs=250, + doi='http://dx.doi.org/10.1073/pnas.1508080112' + ) def load_raw(self, path=None, epoch_duration=1, band=[5.0, 45.0], order=6, augment=False): ''' ''' - list_of_files = np.array(glob.glob(path+'/S*.mat')) - indices = np.array([int(re.findall(r'\d+', n)[0]) for n in list_of_files]) - 1 - onset = int(0.5*self.fs) - epoch_duration = np.round(np.array(epoch_duration+0.5) * self.fs).astype(int) - - n_subjects = 35 - + list_of_files = np.array(glob.glob(path + '/S*.mat')) + indices = np.array([int(re.findall(r'\d+', n)[0]) + for n in list_of_files]) - 1 + onset = int(0.5 * self.fs) + epoch_duration = np.round( + np.array(epoch_duration + 0.5) * self.fs).astype(int) + + n_subjects = 35 + X = [] Y = [] augmented = 0 - + for subj in range(n_subjects): - data = loadmat(list_of_files[indices==subj][0]) - eeg = data['data'].transpose((1,0,2,3)) + data = loadmat(list_of_files[indices == subj][0]) + eeg = data['data'].transpose((1, 0, 2, 3)) del data eeg = bandpass(eeg, band=band, fs=self.fs, order=order) if augment: augmented = 4 tg = eeg.shape[2] - v = [eeg[onset+(stride*self.fs):onset+(stride*self.fs)+epoch_duration, :, :, :] for stride in range(augmented)] + v = [eeg[onset + (stride * self.fs):onset + (stride * self.fs) + + epoch_duration, :, :, :] for stride in range(augmented)] eeg = np.concatenate(v, axis=2) - samples, channels, targets, blocks = eeg.shape - y = np.tile(np.arange(1, tg+1), (1, augmented)) - y = np.tile(y, (1,blocks)) - del v #saving some RAM + samples, channels, targets, blocks = eeg.shape + y = np.tile(np.arange(1, tg + 1), (1, augmented)) + y = np.tile(y, (1, blocks)) + del v # saving some RAM else: - eeg = eeg[onset:epoch_duration, :, :, :] + eeg = eeg[onset:epoch_duration, :, :, :] samples, channels, targets, blocks = eeg.shape - y = np.tile(np.arange(1, targets+1), (1,blocks)) - - X.append(eeg.reshape((samples, channels, blocks*targets),order='F')) + y = np.tile(np.arange(1, targets + 1), (1, blocks)) + + X.append( + eeg.reshape( + (samples, + channels, + blocks * + targets), + order='F')) Y.append(y) del eeg del y @@ -76,16 +85,15 @@ def load_raw(self, path=None, epoch_duration=1, Y = np.array(Y).squeeze() return X, Y - - def generate_set(self, load_path=None, epoch=1, band=[5.0,45.0], - order=6, save_folder=None, augment=False): + def generate_set(self, load_path=None, epoch=1, band=[5.0, 45.0], + order=6, save_folder=None, augment=False): ''' ''' - self.epochs , self.y = self.load_raw(load_path, - epoch, band, order, - augment) + self.epochs, self.y = self.load_raw(load_path, + epoch, band, order, + augment) self.subjects = self._get_subjects(path=load_path) - self.paradigm = self._get_paradigm() + self.paradigm = self._get_paradigm() self.save_set(save_folder) def _get_subjects(self, n_subject=0, path=None): @@ -95,31 +103,28 @@ def _get_subjects(self, n_subject=0, path=None): f = open(sub_file, 'r') info = f.read().split('\n')[2:] f.close() - return [Subject(id=s.split()[0],gender=s.split()[1], - age=s.split()[2], handedness=s.split()[3]) - for s in info if len(s)>0] + return [Subject(id=s.split()[0], gender=s.split()[1], + age=s.split()[2], handedness=s.split()[3]) + for s in info if len(s) > 0] def _get_paradigm(self): ''' ''' return SSVEP(title='SSVEP_JFPM', stimulation=5000, break_duration=500, repetition=6, - stimuli=40, phrase='', - stim_type='Sinusoidal', frequencies=[8. , 9. , 10. , 11. , 12. , 13. , 14. , 15. , 8.2, 9.2, - 10.2, 11.2, 12.2, 13.2, 14.2, 15.2, 8.4, 9.4, 10.4, 11.4, - 12.4, 13.4, 14.4, 15.4, 8.6, 9.6, 10.6, 11.6, 12.6, 13.6, - 14.6, 15.6, 8.8, 9.8, 10.8, 11.8, 12.8, 13.8, 14.8, 15.8], - phase=[0., 1.57079633, 3.14159265, 4.71238898, 0. , - 1.57079633, 3.14159265, 4.71238898, 1.57079633, 3.14159265, - 4.71238898, 0. , 1.57079633, 3.14159265, 4.71238898, - 0. , 3.14159265, 4.71238898, 0. , 1.57079633, - 3.14159265, 4.71238898, 0. , 1.57079633, 4.71238898, - 0. , 1.57079633, 3.14159265, 4.71238898, 0. , - 1.57079633, 3.14159265, 0. , 1.57079633, 3.14159265, - 4.71238898, 0. , 1.57079633, 3.14159265, 4.71238898] - ) - + stimuli=40, phrase='', + stim_type='Sinusoidal', frequencies=[8., 9., 10., 11., 12., 13., 14., 15., 8.2, 9.2, + 10.2, 11.2, 12.2, 13.2, 14.2, 15.2, 8.4, 9.4, 10.4, 11.4, + 12.4, 13.4, 14.4, 15.4, 8.6, 9.6, 10.6, 11.6, 12.6, 13.6, + 14.6, 15.6, 8.8, 9.8, 10.8, 11.8, 12.8, 13.8, 14.8, 15.8], + phase=[0., 1.57079633, 3.14159265, 4.71238898, 0., + 1.57079633, 3.14159265, 4.71238898, 1.57079633, 3.14159265, + 4.71238898, 0., 1.57079633, 3.14159265, 4.71238898, + 0., 3.14159265, 4.71238898, 0., 1.57079633, + 3.14159265, 4.71238898, 0., 1.57079633, 4.71238898, + 0., 1.57079633, 3.14159265, 4.71238898, 0., + 1.57079633, 3.14159265, 0., 1.57079633, 3.14159265, + 4.71238898, 0., 1.57079633, 3.14159265, 4.71238898] + ) def get_path(self): NotImplementedError - - \ No newline at end of file diff --git a/aawedha/models/EEGModels.py b/aawedha/models/EEGModels.py index 7e8cac9..8c0ed4d 100644 --- a/aawedha/models/EEGModels.py +++ b/aawedha/models/EEGModels.py @@ -5,42 +5,42 @@ Requirements: (1) tensorflow-gpu == 1.12.0 (2) 'image_data_format' = 'channels_first' in keras.json config - (3) Data shape = (trials, kernels, channels, samples), which for the + (3) Data shape = (trials, kernels, channels, samples), which for the input layer, will be (trials, 1, channels, samples). - + To run the EEG/MEG ERP classification sample script, you will also need (4) mne >= 0.17.1 (5) PyRiemann >= 0.2.5 (6) scikit-learn >= 0.20.1 (7) matplotlib >= 2.2.3 - + To use: - + (1) Place this file in the PYTHONPATH variable in your IDE (i.e.: Spyder) (2) Import the model as - - from EEGModels import EEGNet - + + from EEGModels import EEGNet + model = EEGNet(nb_classes = ..., Chans = ..., Samples = ...) - + (3) Then compile and fit the model - + model.compile(loss = ..., optimizer = ..., metrics = ...) fitted = model.fit(...) predicted = model.predict(...) Portions of this project are works of the United States Government and are not - subject to domestic copyright protection under 17 USC Sec. 105. Those - portions are released world-wide under the terms of the Creative Commons Zero - 1.0 (CC0) license. - - Other portions of this project are subject to domestic copyright protection - under 17 USC Sec. 105. Those portions are licensed under the Apache 2.0 - license. The complete text of the license governing this material is in - the file labeled LICENSE.TXT that is a part of this project's official - distribution. + subject to domestic copyright protection under 17 USC Sec. 105. Those + portions are released world-wide under the terms of the Creative Commons Zero + 1.0 (CC0) license. + + Other portions of this project are subject to domestic copyright protection + under 17 USC Sec. 105. Those portions are licensed under the Apache 2.0 + license. The complete text of the license governing this material is in + the file labeled LICENSE.TXT that is a part of this project's official + distribution. """ from tensorflow.keras.models import Model @@ -55,9 +55,9 @@ from tensorflow.keras import backend as K -def EEGNet(nb_classes, Chans = 64, Samples = 128, - dropoutRate = 0.5, kernLength = 64, F1 = 8, - D = 2, F2 = 16, norm_rate = 0.25, dropoutType = 'Dropout'): +def EEGNet(nb_classes, Chans=64, Samples=128, + dropoutRate=0.5, kernLength=64, F1=8, + D=2, F2=16, norm_rate=0.25, dropoutType='Dropout'): """ Keras Implementation of EEGNet http://iopscience.iop.org/article/10.1088/1741-2552/aace8c/meta @@ -65,44 +65,44 @@ def EEGNet(nb_classes, Chans = 64, Samples = 128, version (version v1 and v2 on arxiv). We strongly recommend using this architecture as it performs much better and has nicer properties than our earlier version. For example: - - 1. Depthwise Convolutions to learn spatial filters within a - temporal convolution. The use of the depth_multiplier option maps + + 1. Depthwise Convolutions to learn spatial filters within a + temporal convolution. The use of the depth_multiplier option maps exactly to the number of spatial filters learned within a temporal - filter. This matches the setup of algorithms like FBCSP which learn - spatial filters within each filter in a filter-bank. This also limits + filter. This matches the setup of algorithms like FBCSP which learn + spatial filters within each filter in a filter-bank. This also limits the number of free parameters to fit when compared to a fully-connected - convolution. - + convolution. + 2. Separable Convolutions to learn how to optimally combine spatial filters across temporal bands. Separable Convolutions are Depthwise - Convolutions followed by (1x1) Pointwise Convolutions. - - - While the original paper used Dropout, we found that SpatialDropout2D - sometimes produced slightly better results for classification of ERP - signals. However, SpatialDropout2D significantly reduced performance + Convolutions followed by (1x1) Pointwise Convolutions. + + + While the original paper used Dropout, we found that SpatialDropout2D + sometimes produced slightly better results for classification of ERP + signals. However, SpatialDropout2D significantly reduced performance on the Oscillatory dataset (SMR, BCI-IV Dataset 2A). We recommend using the default Dropout in most cases. - + Assumes the input signal is sampled at 128Hz. If you want to use this model for any other sampling rate you will need to modify the lengths of temporal - kernels and average pooling size in blocks 1 and 2 as needed (double the - kernel lengths for double the sampling rate, etc). Note that we haven't - tested the model performance with this rule so this may not work well. - + kernels and average pooling size in blocks 1 and 2 as needed (double the + kernel lengths for double the sampling rate, etc). Note that we haven't + tested the model performance with this rule so this may not work well. + The model with default parameters gives the EEGNet-8,2 model as discussed in the paper. This model should do pretty well in general, although it is - advised to do some model searching to get optimal performance on your - particular dataset. + advised to do some model searching to get optimal performance on your + particular dataset. We set F2 = F1 * D (number of input filters = number of output filters) for the SeparableConv2D layer. We haven't extensively tested other values of this parameter (say, F2 < F1 * D for compressed learning, and F2 > F1 * D for - overcomplete). We believe the main parameters to focus on are F1 and D. + overcomplete). We believe the main parameters to focus on are F1 and D. Inputs: - + nb_classes : int, number of classes to classify Chans, Samples : number of channels and time points in the EEG data dropoutRate : dropout fraction @@ -110,15 +110,15 @@ def EEGNet(nb_classes, Chans = 64, Samples = 128, that setting this to be half the sampling rate worked well in practice. For the SMR dataset in particular since the data was high-passed at 4Hz we used a kernel - length of 32. + length of 32. F1, F2 : number of temporal filters (F1) and number of pointwise - filters (F2) to learn. Default: F1 = 8, F2 = F1 * D. + filters (F2) to learn. Default: F1 = 8, F2 = F1 * D. D : number of spatial filters to learn within each temporal convolution. Default: D = 2 dropoutType : Either SpatialDropout2D or Dropout, passed as a string. """ - + if dropoutType == 'SpatialDropout2D': dropoutType = SpatialDropout2D elif dropoutType == 'Dropout': @@ -126,65 +126,63 @@ def EEGNet(nb_classes, Chans = 64, Samples = 128, else: raise ValueError('dropoutType must be one of SpatialDropout2D ' 'or Dropout, passed as a string.') - - input1 = Input(shape = (1, Chans, Samples)) - ################################################################## - block1 = Conv2D(F1, (1, kernLength), padding = 'same', - input_shape = (1, Chans, Samples), - use_bias = False)(input1) - block1 = BatchNormalization(axis = 1)(block1) - block1 = DepthwiseConv2D((Chans, 1), use_bias = False, - depth_multiplier = D, - depthwise_constraint = max_norm(1.))(block1) - block1 = BatchNormalization(axis = 1)(block1) - block1 = Activation('elu')(block1) - block1 = AveragePooling2D((1, 4))(block1) - block1 = dropoutType(dropoutRate)(block1) - - block2 = SeparableConv2D(F2, (1, 16), - use_bias = False, padding = 'same')(block1) - block2 = BatchNormalization(axis = 1)(block2) - block2 = Activation('elu')(block2) - block2 = AveragePooling2D((1, 8))(block2) - block2 = dropoutType(dropoutRate)(block2) - - flatten = Flatten(name = 'flatten')(block2) - - dense = Dense(nb_classes, name = 'dense', - kernel_constraint = max_norm(norm_rate))(flatten) - softmax = Activation('softmax', name = 'softmax')(dense) - - return Model(inputs=input1, outputs=softmax) + input1 = Input(shape=(1, Chans, Samples)) + ################################################################## + block1 = Conv2D(F1, (1, kernLength), padding='same', + input_shape=(1, Chans, Samples), + use_bias=False)(input1) + block1 = BatchNormalization(axis=1)(block1) + block1 = DepthwiseConv2D((Chans, 1), use_bias=False, + depth_multiplier=D, + depthwise_constraint=max_norm(1.))(block1) + block1 = BatchNormalization(axis=1)(block1) + block1 = Activation('elu')(block1) + block1 = AveragePooling2D((1, 4))(block1) + block1 = dropoutType(dropoutRate)(block1) + + block2 = SeparableConv2D(F2, (1, 16), + use_bias=False, padding='same')(block1) + block2 = BatchNormalization(axis=1)(block2) + block2 = Activation('elu')(block2) + block2 = AveragePooling2D((1, 8))(block2) + block2 = dropoutType(dropoutRate)(block2) + + flatten = Flatten(name='flatten')(block2) + + dense = Dense(nb_classes, name='dense', + kernel_constraint=max_norm(norm_rate))(flatten) + softmax = Activation('softmax', name='softmax')(dense) + return Model(inputs=input1, outputs=softmax) -def EEGNet_SSVEP(nb_classes = 12, Chans = 8, Samples = 256, - dropoutRate = 0.5, kernLength = 256, F1 = 96, - D = 1, F2 = 96, dropoutType = 'Dropout'): - """ SSVEP Variant of EEGNet, as used in [1]. +def EEGNet_SSVEP(nb_classes=12, Chans=8, Samples=256, + dropoutRate=0.5, kernLength=256, F1=96, + D=1, F2=96, dropoutType='Dropout'): + """ SSVEP Variant of EEGNet, as used in [1]. Inputs: - + nb_classes : int, number of classes to classify Chans, Samples : number of channels and time points in the EEG data dropoutRate : dropout fraction kernLength : length of temporal convolution in first layer F1, F2 : number of temporal filters (F1) and number of pointwise - filters (F2) to learn. + filters (F2) to learn. D : number of spatial filters to learn within each temporal convolution. dropoutType : Either SpatialDropout2D or Dropout, passed as a string. - - + + [1]. Waytowich, N. et. al. (2018). Compact Convolutional Neural Networks for Classification of Asynchronous Steady-State Visual Evoked Potentials. - Journal of Neural Engineering vol. 15(6). + Journal of Neural Engineering vol. 15(6). http://iopscience.iop.org/article/10.1088/1741-2552/aae5d8 """ - + if dropoutType == 'SpatialDropout2D': dropoutType = SpatialDropout2D elif dropoutType == 'Dropout': @@ -192,162 +190,160 @@ def EEGNet_SSVEP(nb_classes = 12, Chans = 8, Samples = 256, else: raise ValueError('dropoutType must be one of SpatialDropout2D ' 'or Dropout, passed as a string.') - - input1 = Input(shape = (1, Chans, Samples)) + + input1 = Input(shape=(1, Chans, Samples)) ################################################################## - block1 = Conv2D(F1, (1, kernLength), padding = 'same', - input_shape = (1, Chans, Samples), - use_bias = False)(input1) - block1 = BatchNormalization(axis = 1)(block1) - block1 = DepthwiseConv2D((Chans, 1), use_bias = False, - depth_multiplier = D, - depthwise_constraint = max_norm(1.))(block1) - block1 = BatchNormalization(axis = 1)(block1) - block1 = Activation('elu')(block1) - block1 = AveragePooling2D((1, 4))(block1) - block1 = dropoutType(dropoutRate)(block1) - - block2 = SeparableConv2D(F2, (1, 16), - use_bias = False, padding = 'same')(block1) - block2 = BatchNormalization(axis = 1)(block2) - block2 = Activation('elu')(block2) - block2 = AveragePooling2D((1, 8))(block2) - block2 = dropoutType(dropoutRate)(block2) - - flatten = Flatten(name = 'flatten')(block2) - - dense = Dense(nb_classes, name = 'dense')(flatten) - softmax = Activation('softmax', name = 'softmax')(dense) - - return Model(inputs=input1, outputs=softmax) + block1 = Conv2D(F1, (1, kernLength), padding='same', + input_shape=(1, Chans, Samples), + use_bias=False)(input1) + block1 = BatchNormalization(axis=1)(block1) + block1 = DepthwiseConv2D((Chans, 1), use_bias=False, + depth_multiplier=D, + depthwise_constraint=max_norm(1.))(block1) + block1 = BatchNormalization(axis=1)(block1) + block1 = Activation('elu')(block1) + block1 = AveragePooling2D((1, 4))(block1) + block1 = dropoutType(dropoutRate)(block1) + + block2 = SeparableConv2D(F2, (1, 16), + use_bias=False, padding='same')(block1) + block2 = BatchNormalization(axis=1)(block2) + block2 = Activation('elu')(block2) + block2 = AveragePooling2D((1, 8))(block2) + block2 = dropoutType(dropoutRate)(block2) + + flatten = Flatten(name='flatten')(block2) + + dense = Dense(nb_classes, name='dense')(flatten) + softmax = Activation('softmax', name='softmax')(dense) + return Model(inputs=input1, outputs=softmax) -def EEGNet_old(nb_classes, Chans = 64, Samples = 128, regRate = 0.0001, - dropoutRate = 0.25, kernels = [(2, 32), (8, 4)], strides = (2, 4)): +def EEGNet_old(nb_classes, Chans=64, Samples=128, regRate=0.0001, + dropoutRate=0.25, kernels=[(2, 32), (8, 4)], strides=(2, 4)): """ Keras Implementation of EEGNet_v1 (https://arxiv.org/abs/1611.08024v2) This model is the original EEGNet model proposed on arxiv https://arxiv.org/abs/1611.08024v2 - - with a few modifications: we use striding instead of max-pooling as this - helped slightly in classification performance while also providing a - computational speed-up. - + + with a few modifications: we use striding instead of max-pooling as this + helped slightly in classification performance while also providing a + computational speed-up. + Note that we no longer recommend the use of this architecture, as the new version of EEGNet performs much better overall and has nicer properties. - + Inputs: - + nb_classes : total number of final categories Chans, Samples : number of EEG channels and samples, respectively regRate : regularization rate for L1 and L2 regularizations dropoutRate : dropout fraction - kernels : the 2nd and 3rd layer kernel dimensions (default is + kernels : the 2nd and 3rd layer kernel dimensions (default is the [2, 32] x [8, 4] configuration) strides : the stride size (note that this replaces the max-pool used in the original paper) - + """ # start the model - input_main = Input((1, Chans, Samples)) - layer1 = Conv2D(16, (Chans, 1), input_shape=(1, Chans, Samples), - kernel_regularizer = l1_l2(l1=regRate, l2=regRate))(input_main) - layer1 = BatchNormalization(axis=1)(layer1) - layer1 = Activation('elu')(layer1) - layer1 = Dropout(dropoutRate)(layer1) - + input_main = Input((1, Chans, Samples)) + layer1 = Conv2D(16, (Chans, 1), input_shape=(1, Chans, Samples), + kernel_regularizer=l1_l2(l1=regRate, l2=regRate))(input_main) + layer1 = BatchNormalization(axis=1)(layer1) + layer1 = Activation('elu')(layer1) + layer1 = Dropout(dropoutRate)(layer1) + permute_dims = 2, 1, 3 - permute1 = Permute(permute_dims)(layer1) - - layer2 = Conv2D(4, kernels[0], padding = 'same', - kernel_regularizer=l1_l2(l1=0.0, l2=regRate), - strides = strides)(permute1) - layer2 = BatchNormalization(axis=1)(layer2) - layer2 = Activation('elu')(layer2) - layer2 = Dropout(dropoutRate)(layer2) - - layer3 = Conv2D(4, kernels[1], padding = 'same', - kernel_regularizer=l1_l2(l1=0.0, l2=regRate), - strides = strides)(layer2) - layer3 = BatchNormalization(axis=1)(layer3) - layer3 = Activation('elu')(layer3) - layer3 = Dropout(dropoutRate)(layer3) - - flatten = Flatten(name = 'flatten')(layer3) - - dense = Dense(nb_classes, name = 'dense')(flatten) - softmax = Activation('softmax', name = 'softmax')(dense) - - return Model(inputs=input_main, outputs=softmax) + permute1 = Permute(permute_dims)(layer1) + layer2 = Conv2D(4, kernels[0], padding='same', + kernel_regularizer=l1_l2(l1=0.0, l2=regRate), + strides=strides)(permute1) + layer2 = BatchNormalization(axis=1)(layer2) + layer2 = Activation('elu')(layer2) + layer2 = Dropout(dropoutRate)(layer2) + layer3 = Conv2D(4, kernels[1], padding='same', + kernel_regularizer=l1_l2(l1=0.0, l2=regRate), + strides=strides)(layer2) + layer3 = BatchNormalization(axis=1)(layer3) + layer3 = Activation('elu')(layer3) + layer3 = Dropout(dropoutRate)(layer3) -def DeepConvNet(nb_classes, Chans = 64, Samples = 256, - dropoutRate = 0.5): + flatten = Flatten(name='flatten')(layer3) + + dense = Dense(nb_classes, name='dense')(flatten) + softmax = Activation('softmax', name='softmax')(dense) + + return Model(inputs=input_main, outputs=softmax) + + +def DeepConvNet(nb_classes, Chans=64, Samples=256, + dropoutRate=0.5): """ Keras implementation of the Deep Convolutional Network as described in Schirrmeister et. al. (2017), Human Brain Mapping. - - This implementation assumes the input is a 2-second EEG signal sampled at + + This implementation assumes the input is a 2-second EEG signal sampled at 128Hz, as opposed to signals sampled at 250Hz as described in the original paper. We also perform temporal convolutions of length (1, 5) as opposed - to (1, 10) due to this sampling rate difference. - - Note that we use the max_norm constraint on all convolutional layers, as + to (1, 10) due to this sampling rate difference. + + Note that we use the max_norm constraint on all convolutional layers, as well as the classification layer. We also change the defaults for the - BatchNormalization layer. We used this based on a personal communication + BatchNormalization layer. We used this based on a personal communication with the original authors. - + ours original paper pool_size 1, 2 1, 3 strides 1, 2 1, 3 conv filters 1, 5 1, 10 - - Note that this implementation has not been verified by the original - authors. - + + Note that this implementation has not been verified by the original + authors. + """ # start the model - input_main = Input((1, Chans, Samples)) - block1 = Conv2D(25, (1, 5), - input_shape=(1, Chans, Samples), - kernel_constraint = max_norm(2., axis=(0,1,2)))(input_main) - block1 = Conv2D(25, (Chans, 1), - kernel_constraint = max_norm(2., axis=(0,1,2)))(block1) - block1 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block1) - block1 = Activation('elu')(block1) - block1 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block1) - block1 = Dropout(dropoutRate)(block1) - - block2 = Conv2D(50, (1, 5), - kernel_constraint = max_norm(2., axis=(0,1,2)))(block1) - block2 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block2) - block2 = Activation('elu')(block2) - block2 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block2) - block2 = Dropout(dropoutRate)(block2) - - block3 = Conv2D(100, (1, 5), - kernel_constraint = max_norm(2., axis=(0,1,2)))(block2) - block3 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block3) - block3 = Activation('elu')(block3) - block3 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block3) - block3 = Dropout(dropoutRate)(block3) - - block4 = Conv2D(200, (1, 5), - kernel_constraint = max_norm(2., axis=(0,1,2)))(block3) - block4 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block4) - block4 = Activation('elu')(block4) - block4 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block4) - block4 = Dropout(dropoutRate)(block4) - - flatten = Flatten()(block4) - - dense = Dense(nb_classes, kernel_constraint = max_norm(0.5))(flatten) - softmax = Activation('softmax')(dense) - + input_main = Input((1, Chans, Samples)) + block1 = Conv2D(25, (1, 5), + input_shape=(1, Chans, Samples), + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(input_main) + block1 = Conv2D(25, (Chans, 1), + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(block1) + block1 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block1) + block1 = Activation('elu')(block1) + block1 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block1) + block1 = Dropout(dropoutRate)(block1) + + block2 = Conv2D(50, (1, 5), + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(block1) + block2 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block2) + block2 = Activation('elu')(block2) + block2 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block2) + block2 = Dropout(dropoutRate)(block2) + + block3 = Conv2D(100, (1, 5), + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(block2) + block3 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block3) + block3 = Activation('elu')(block3) + block3 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block3) + block3 = Dropout(dropoutRate)(block3) + + block4 = Conv2D(200, (1, 5), + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(block3) + block4 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block4) + block4 = Activation('elu')(block4) + block4 = MaxPooling2D(pool_size=(1, 2), strides=(1, 2))(block4) + block4 = Dropout(dropoutRate)(block4) + + flatten = Flatten()(block4) + + dense = Dense(nb_classes, kernel_constraint=max_norm(0.5))(flatten) + softmax = Activation('softmax')(dense) + return Model(inputs=input_main, outputs=softmax) @@ -355,51 +351,50 @@ def DeepConvNet(nb_classes, Chans = 64, Samples = 256, def square(x): return K.square(x) + def log(x): - return K.log(K.clip(x, min_value = 1e-7, max_value = 10000)) + return K.log(K.clip(x, min_value=1e-7, max_value=10000)) -def ShallowConvNet(nb_classes, Chans = 64, Samples = 128, dropoutRate = 0.5): +def ShallowConvNet(nb_classes, Chans=64, Samples=128, dropoutRate=0.5): """ Keras implementation of the Shallow Convolutional Network as described in Schirrmeister et. al. (2017), Human Brain Mapping. - - Assumes the input is a 2-second EEG signal sampled at 128Hz. Note that in + + Assumes the input is a 2-second EEG signal sampled at 128Hz. Note that in the original paper, they do temporal convolutions of length 25 for EEG - data sampled at 250Hz. We instead use length 13 since the sampling rate is + data sampled at 250Hz. We instead use length 13 since the sampling rate is roughly half of the 250Hz which the paper used. The pool_size and stride in later layers is also approximately half of what is used in the paper. - - Note that we use the max_norm constraint on all convolutional layers, as + + Note that we use the max_norm constraint on all convolutional layers, as well as the classification layer. We also change the defaults for the - BatchNormalization layer. We used this based on a personal communication + BatchNormalization layer. We used this based on a personal communication with the original authors. - + ours original paper pool_size 1, 35 1, 75 strides 1, 7 1, 15 - conv filters 1, 13 1, 25 - - Note that this implementation has not been verified by the original + conv filters 1, 13 1, 25 + + Note that this implementation has not been verified by the original authors. We do note that this implementation reproduces the results in the - original paper with minor deviations. + original paper with minor deviations. """ # start the model - input_main = Input((1, Chans, Samples)) - block1 = Conv2D(40, (1, 13), - input_shape=(1, Chans, Samples), - kernel_constraint = max_norm(2., axis=(0,1,2)))(input_main) - block1 = Conv2D(40, (Chans, 1), use_bias=False, - kernel_constraint = max_norm(2., axis=(0,1,2)))(block1) - block1 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block1) - block1 = Activation(square)(block1) - block1 = AveragePooling2D(pool_size=(1, 35), strides=(1, 7))(block1) - block1 = Activation(log)(block1) - block1 = Dropout(dropoutRate)(block1) - flatten = Flatten()(block1) - dense = Dense(nb_classes, kernel_constraint = max_norm(0.5))(flatten) - softmax = Activation('softmax')(dense) - - return Model(inputs=input_main, outputs=softmax) - + input_main = Input((1, Chans, Samples)) + block1 = Conv2D(40, (1, 13), + input_shape=(1, Chans, Samples), + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(input_main) + block1 = Conv2D(40, (Chans, 1), use_bias=False, + kernel_constraint=max_norm(2., axis=(0, 1, 2)))(block1) + block1 = BatchNormalization(axis=1, epsilon=1e-05, momentum=0.1)(block1) + block1 = Activation(square)(block1) + block1 = AveragePooling2D(pool_size=(1, 35), strides=(1, 7))(block1) + block1 = Activation(log)(block1) + block1 = Dropout(dropoutRate)(block1) + flatten = Flatten()(block1) + dense = Dense(nb_classes, kernel_constraint=max_norm(0.5))(flatten) + softmax = Activation('softmax')(dense) + return Model(inputs=input_main, outputs=softmax) diff --git a/aawedha/paradigms/base.py b/aawedha/paradigms/base.py index bd16351..dcf4cfa 100644 --- a/aawedha/paradigms/base.py +++ b/aawedha/paradigms/base.py @@ -3,6 +3,7 @@ """ from abc import ABCMeta, abstractmethod + class Paradigm(metaclass=ABCMeta): """ Paradigms: @@ -28,7 +29,7 @@ class Paradigm(metaclass=ABCMeta): number of stimulus presented in the experiment. stim_type: str - stimulus presented to subject (ERP) / type of stimulations used in the experiment (SSVEP) / type of cue presented before task (MI) + stimulus presented to subject (ERP) / type of stimulations used in the experiment (SSVEP) / type of cue presented before task (MI) phrase : str sequence of characters to be spelled by the subject during the experiments @@ -37,8 +38,8 @@ class Paradigm(metaclass=ABCMeta): ------- """ - def __init__(self, title=None, control=None, stimulation=0, - break_duration=0, repetition=0, stimuli=0, stim_type=None, phrase=None): + def __init__(self, title=None, control=None, stimulation=0, + break_duration=0, repetition=0, stimuli=0, stim_type=None, phrase=None): self.title = title self.control = control self.stimulation = stimulation @@ -46,4 +47,4 @@ def __init__(self, title=None, control=None, stimulation=0, self.repetition = repetition self.stimuli = stimuli self.stim_type = stim_type - self.phrase = phrase \ No newline at end of file + self.phrase = phrase diff --git a/aawedha/paradigms/erp.py b/aawedha/paradigms/erp.py index e1edf50..12bc49d 100644 --- a/aawedha/paradigms/erp.py +++ b/aawedha/paradigms/erp.py @@ -1,9 +1,10 @@ """ - ERP paradigm + ERP paradigm """ from aawedha.paradigms.base import Paradigm + class ERP(Paradigm): """ Attributes @@ -38,10 +39,11 @@ class ERP(Paradigm): Methods ------- """ - def __init__(self, title='ERP', control='Sync', stimulation=100, - break_duration=100, repetition=10, stimuli=12, phrase='12345', - stim_type='flash', flashing_mode='SC', speller=[]): + + def __init__(self, title='ERP', control='Sync', stimulation=100, + break_duration=100, repetition=10, stimuli=12, phrase='12345', + stim_type='flash', flashing_mode='SC', speller=[]): super(ERP, self).__init__(title, control, stimulation, break_duration, repetition, - stimuli, stim_type, phrase) + stimuli, stim_type, phrase) self.flashing_mode = flashing_mode - self.speller = speller \ No newline at end of file + self.speller = speller diff --git a/aawedha/paradigms/hybrid.py b/aawedha/paradigms/hybrid.py index 5cf3279..882a4af 100644 --- a/aawedha/paradigms/hybrid.py +++ b/aawedha/paradigms/hybrid.py @@ -1,17 +1,18 @@ """ - Hybrid aradigm + Hybrid aradigm """ # from aawedha.paradigms.base import Paradigm from aawedha.paradigms.erp import ERP from aawedha.paradigms.ssvep import SSVEP + class HybridLARESI(object): - def __init__(self, title='LARESI_HYBRID', erp = None, - ssvep=None, mode='calib' + def __init__(self, title='LARESI_HYBRID', erp=None, + ssvep=None, mode='calib' ): - + self.title = title self.experiment_mode = mode if erp: @@ -24,16 +25,16 @@ def __init__(self, title='LARESI_HYBRID', erp = None, self._set_SSVEP() def _set_ERP(self): - self.erp = ERP(title='ERP_LARESI', stimulation=100, - break_duration=50, repetition=5, - phrase='123456789', - flashing_mode='SC', - speller=['1','2','3','4','5','6','7','8','9']) + self.erp = ERP(title='ERP_LARESI', stimulation=100, + break_duration=50, repetition=5, + phrase='123456789', + flashing_mode='SC', + speller=['1', '2', '3', '4', '5', '6', '7', '8', '9']) def _set_SSVEP(self): self.ssvep = SSVEP(title='SSVEP_LARESI', control='Async', - stimulation=1000, - break_duration=2000, repetition=15, - stimuli=5, phrase='', - stim_type='Sinusoidal', frequencies=['14','12','10','8'], - ) \ No newline at end of file + stimulation=1000, + break_duration=2000, repetition=15, + stimuli=5, phrase='', + stim_type='Sinusoidal', frequencies=['14', '12', '10', '8'], + ) diff --git a/aawedha/paradigms/motor_imagery.py b/aawedha/paradigms/motor_imagery.py index 8c40f5c..203b4c6 100644 --- a/aawedha/paradigms/motor_imagery.py +++ b/aawedha/paradigms/motor_imagery.py @@ -4,6 +4,7 @@ from aawedha.paradigms.base import Paradigm + class MotorImagery(Paradigm): """ Attributes @@ -32,20 +33,20 @@ class MotorImagery(Paradigm): stim_type : str stimulus presented to subject: flash / face / inverted_face ... - + Methods ------- """ - def __init__(self, title='MotorImagery', - control='Sync', - stimulation=3000, - break_duration=2000, - repetition=72, - stimuli=4, - phrase=None, - stim_type='Arrow Cue/Beep' - ): - super().__init__(title, control, stimulation, - break_duration, repetition, - stimuli, stim_type, phrase) - \ No newline at end of file + + def __init__(self, title='MotorImagery', + control='Sync', + stimulation=3000, + break_duration=2000, + repetition=72, + stimuli=4, + phrase=None, + stim_type='Arrow Cue/Beep' + ): + super().__init__(title, control, stimulation, + break_duration, repetition, + stimuli, stim_type, phrase) diff --git a/aawedha/paradigms/ssvep.py b/aawedha/paradigms/ssvep.py index a4eb715..750e713 100644 --- a/aawedha/paradigms/ssvep.py +++ b/aawedha/paradigms/ssvep.py @@ -1,9 +1,10 @@ """ - SSVEP paradigm + SSVEP paradigm """ from aawedha.paradigms.base import Paradigm + class SSVEP(Paradigm): """ Attributes @@ -34,17 +35,18 @@ class SSVEP(Paradigm): frequencies : list frequencies presented at each stimuli in Hz - + phase : list phase of each frequency in rad Methods ------- """ - def __init__(self, title='SSVEP', control='Sync', stimulation=4000, - break_duration=4000, repetition=10, stimuli=4, phrase='1234', - stim_type='ON_OFF', frequencies=['7.5','8.57','10','12'], phase=None): + + def __init__(self, title='SSVEP', control='Sync', stimulation=4000, + break_duration=4000, repetition=10, stimuli=4, phrase='1234', + stim_type='ON_OFF', frequencies=['7.5', '8.57', '10', '12'], phase=None): super().__init__(title, control, stimulation, break_duration, repetition, - stimuli, stim_type, phrase) + stimuli, stim_type, phrase) self.frequencies = frequencies - self.phase = phase \ No newline at end of file + self.phase = phase diff --git a/aawedha/paradigms/subject.py b/aawedha/paradigms/subject.py index f22869c..072f55f 100644 --- a/aawedha/paradigms/subject.py +++ b/aawedha/paradigms/subject.py @@ -2,13 +2,14 @@ Subjects information class """ + class Subject(object): """ Attributes ---------- id : str subject identifier - + gender : str subject gender : M / F @@ -24,8 +25,9 @@ class Subject(object): Methods ------- """ - - def __init__(self, id='S', gender='M', age=18, handedness='R', condition='healthy'): + + def __init__(self, id='S', gender='M', age=18, + handedness='R', condition='healthy'): self.id = id self.gender = gender self.age = age diff --git a/aawedha/process/base.py b/aawedha/process/base.py index cd62599..6324c71 100644 --- a/aawedha/process/base.py +++ b/aawedha/process/base.py @@ -4,12 +4,13 @@ # ?import model # ?import pipeline + class Process(metaclass=ABCMeta): """Process base class Attributes ---------- - data : Dataset object + data : Dataset object model : model object A Keras model to be trained @@ -17,7 +18,7 @@ class Process(metaclass=ABCMeta): pipeline : scikit-learn pipeline object the process steps in a pipeline to be run - params : + params : dict of pipeline parameters Methods @@ -27,12 +28,12 @@ class Process(metaclass=ABCMeta): run() """ + def __init__(self, data=None, model=None, pipeline=None, params=None): self.data = data self.model = model self.pipeline = pipeline - self.params = params - + self.params = params @abstractmethod def split_set(self, n_sets=30, split=None): @@ -51,12 +52,10 @@ def split_set(self, n_sets=30, split=None): """ pass - def run(self): pass - class SingleSubject(Process): """Single subject analysis process """ @@ -64,4 +63,4 @@ class SingleSubject(Process): class MultipleSubjects(Process): """Multiple subject analysis process - """ \ No newline at end of file + """ diff --git a/aawedha/scripts/cross_subj.py b/aawedha/scripts/cross_subj.py index 4a4130e..566e88b 100644 --- a/aawedha/scripts/cross_subj.py +++ b/aawedha/scripts/cross_subj.py @@ -7,19 +7,17 @@ from tf.keras.metrics import AUC import numpy as np - # generate set : convert raw EEG to a multisubject dataset: # RUN ONCE - # aawedha DataSet Obeject: + # aawedha DataSet Obeject: # X : data tensor : (subjects, samples, channels, trials) - # Y : labels : (subjects, trials) - + # Y : labels : (subjects, trials) data_folder = 'drive/My Drive/data/inria_data' save_path = 'drive/inria_data/epoched/' - t= [0. , 1.25] + t = [0., 1.25] ds = Inria_ERN() - ds.generate_set(load_path=data_folder, epoch=t, + ds.generate_set(load_path=data_folder, epoch=t, save_folder=save_path) # load dataset @@ -29,22 +27,22 @@ # select an evaluation type: Multiple Subjects subjects = len(ds1.subjects) - prt = [14,1,1] - evl = CrossSubject(n_subjects=subjects, partition=prt, - dataset=ds1) + prt = [14, 1, 1] + evl = CrossSubject(n_subjects=subjects, partition=prt, + dataset=ds1) evl.generate_split(nfolds=30) subjects, samples, channels, trials = evl.dataset.epochs.shape n_classes = np.unique(evl.dataset.y).size # select a model - evl.model = EEGNet(nb_classes = n_classes, Chans = channels, Samples = samples, - dropoutRate = 0.5, kernLength = 32, F1 = 8, D = 2, F2 = 16, - dropoutType = 'Dropout') - - evl.model.compile(loss='binary_crossentropy', - optimizer='adam', - metrics = ['accuracy', AUC()] - ) + evl.model = EEGNet(nb_classes=n_classes, Chans=channels, Samples=samples, + dropoutRate=0.5, kernLength=32, F1=8, D=2, F2=16, + dropoutType='Dropout') + + evl.model.compile(loss='binary_crossentropy', + optimizer='adam', + metrics=['accuracy', AUC()] + ) # train model evl.run_evaluation() - # save model \ No newline at end of file + # save model diff --git a/aawedha/scripts/single_subj.py b/aawedha/scripts/single_subj.py index 00a7377..967fcc7 100644 --- a/aawedha/scripts/single_subj.py +++ b/aawedha/scripts/single_subj.py @@ -1,7 +1,7 @@ ''' - main script for evaluating a model(s) on a dataset (s) ''' -from aawedha.io.inria_ern import Inria_ERN +from aawedha.io.inria_ern import Inria_ERN from aawedha.evaluation.single_subject import SingleSubject from aawedha.paradigms.erp import ERP from aawedha.paradigms.subject import Subject @@ -14,16 +14,16 @@ K.set_image_data_format('channels_first') # generate set : convert raw EEG to a multisubject dataset: # RUN ONCE -# aawedha DataSet Obeject: +# aawedha DataSet Obeject: # X : data tensor : (subjects, samples, channels, trials) -# Y : labels : (subjects, trials) +# Y : labels : (subjects, trials) ''' data_folder = 'drive/My Drive/data/inria_data' save_path = 'drive/inria_data/epoched/' t= [0. , 1.25] ds = Inria_ERN() -ds.generate_set(load_path=data_folder, epoch=t, +ds.generate_set(load_path=data_folder, epoch=t, save_folder=save_path) ''' @@ -38,19 +38,19 @@ # subject = 3 # subject number # Signle Subject evaluation -prt = [2,1,1] +prt = [2, 1, 1] evl = SingleSubject(n_subjects=subjects, partition=prt, dataset=dt, verbose=2) -evl.generate_split(nfolds=1) #nfolds=10 +evl.generate_split(nfolds=1) # nfolds=10 # select a model -evl.model = EEGNet(nb_classes = n_classes, Chans = channels, Samples = samples, - dropoutRate = 0.5, kernLength = 32, F1 = 8, D = 2, F2 = 16, - dropoutType = 'Dropout' +evl.model = EEGNet(nb_classes=n_classes, Chans=channels, Samples=samples, + dropoutRate=0.5, kernLength=32, F1=8, D=2, F2=16, + dropoutType='Dropout' ) -evl.model.compile(loss='categorical_crossentropy', - optimizer='adam', - metrics = ['accuracy'] +evl.model.compile(loss='categorical_crossentropy', + optimizer='adam', + metrics=['accuracy'] ) # train/test model evl.run_evaluation() @@ -60,4 +60,4 @@ # print(f'Subject N: {subject+1} Acc: {evl.results["acc"]} ') # print(f' Acc: {evl.results["acc"]} ') -# Visualize learning curves \ No newline at end of file +# Visualize learning curves diff --git a/aawedha/visualization/__init__.py b/aawedha/visualization/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/aawedha/visualization/viz.py b/aawedha/visualization/viz.py new file mode 100644 index 0000000..4108ca0 --- /dev/null +++ b/aawedha/visualization/viz.py @@ -0,0 +1,76 @@ +from sklearn.metrics import auc +from scipy import interp +import matplotlib.pyplot as plt +import pandas as pd +import seaborn as sns +import numpy as np + + +def plot_train_val_curve(model_history={}): + ''' + ''' + ky = list(model_history.keys()) + for k in range(len(ky) // 2): + plt.figure(k) + plt.plot(model_history[ky[k]]) + plt.plot(model_history[ky[k + len(ky) // 2]]) + plt.title(ky[k].upper()) + plt.ylabel(ky[k].upper()) + plt.xlabel('Epoch') + plt.legend(['Train', 'val'], loc='upper left') + plt.show() + + +def plot_roc_curve(results={}, nfolds=4): + ''' + ''' + mean_fpr = np.linspace(0, 1, 100) + tprs = [] + # + fg = plt.figure() + for fld in range(nfolds): + fpr = results['fpr'][0][fld] + tpr = results['tpr'][0][fld] + roc_auc = results['auc'][0, fld] # (subject, fold) + tprs.append(interp(mean_fpr, fpr, tpr)) + fg.plot(fpr, tpr, lw=1, alpha=0.3, + label=f'ROC fold {fld} AUC = {roc_auc}') + + fg.plot([0, 1], [0, 1], linestyle='--', lw=2, + color='r', label='Chance', alpha=.8) + mean_tpr = np.mean(tprs, axis=0) + mean_auc = auc(mean_fpr, mean_tpr) + std_auc = np.std(results['auc']) + fg.plot(mean_fpr, mean_tpr, color='b', + label=r'Mean ROC (AUC = %0.2f $\pm$ %0.2f)' % (mean_auc, std_auc), + lw=2, alpha=.8) + + std_tpr = np.std(tprs, axis=0) + tprs_upper = np.minimum(mean_tpr + std_tpr, 1) + tprs_lower = np.maximum(mean_tpr - std_tpr, 0) + fg.fill_between(mean_fpr, tprs_lower, tprs_upper, color='grey', alpha=.2, + label=r'$\pm$ 1 std. dev.') + + fg.xlim([-0.05, 1.05]) + fg.ylim([-0.05, 1.05]) + fg.xlabel('False Positive Rate') + fg.ylabel('True Positive Rate') + fg.title('Receiver operating characteristic example') + fg.legend(loc="lower right") + fg.show() + return fg + + +def plot_confusion_matrix(cm, target_names, title='Confusion matrix', + cmap="Blues"): + ''' + ''' + cm = 100 * cm.astype('float') / cm.sum(axis=1)[:, np.newaxis] + + df = pd.DataFrame(data=cm, columns=target_names, index=target_names) + g = sns.heatmap(df, annot=True, fmt=".1f", linewidths=.5, vmin=0, vmax=100, + cmap=cmap) + g.set_title(title) + g.set_ylabel('True label') + g.set_xlabel('Predicted label') + return g