Rohit Garg
468 linhas
16 KiB
Python
468 linhas
16 KiB
Python
#!/usr/bin/env python
|
||
# coding: utf-8
|
||
|
||
# In[ ]:
|
||
|
||
|
||
# -*- coding: utf-8 -*-
|
||
"""feature_extraction_25GB_RAM.ipynb
|
||
|
||
Automatically generated by Colaboratory.
|
||
|
||
Original file is located at
|
||
https://colab.research.google.com/drive/1QnVj7GyyJhLPrYF4vBTppMwynXqmOTEJ
|
||
"""
|
||
|
||
# Commented out IPython magic to ensure Python compatibility.
|
||
|
||
|
||
import EEGExtract as eeg
|
||
from scipy import io,signal
|
||
import numpy as np
|
||
import pandas as pd
|
||
from sklearn import preprocessing
|
||
import pandas as pd
|
||
import pickle
|
||
|
||
|
||
class load_data:
|
||
'''
|
||
Load the preprocessed data here, store the paramters
|
||
'''
|
||
def __init__(self,name):
|
||
self.name = name #name of dataset
|
||
self.X = None
|
||
self.Y = None
|
||
self.Z = None
|
||
self.freq = None #(in Hz) is same for all datasets
|
||
self.channels = None
|
||
self.ch_type = 'eeg'
|
||
self.eegData = None
|
||
self.use_autoreject = 'n'
|
||
def load_arrays(self):
|
||
if self.name == 'DREAMER':
|
||
array = np.load('original_data/DREAMER.npz')
|
||
self.freq = 128
|
||
self.channels = ['AF3','F7','F3','FC5','T7','P7','O1','O2','P8','T8','FC6','F4','F8','AF4']
|
||
if self.name == 'DEAP':
|
||
array = np.load('original_data/DEAP.npz')
|
||
self.freq = 128
|
||
# 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39
|
||
self.channels = ['F1', 'AF3', 'F3', 'F7', 'FC5', 'FC1', 'C3', 'T7', 'CP5', 'CP1', 'P3', 'P7', 'PO3', 'O1', 'Oz', 'Pz', 'Fp2', 'AF4', 'Fz', 'F4', 'F8', 'FC6', 'FC2', 'Cz', 'C4', 'T8', 'CP6', 'CP2', 'P4', 'P8', 'PO4', 'O2', 'hEOG','vEOG', 'zEMG','tEMG','GSR','Respiration belt','Plethysmograph','Temperature']
|
||
if self.name == 'OASIS':
|
||
#array = np.load('original_data/OASIS.npz')
|
||
if self.use_autoreject == 'y':
|
||
with open('preprocessed_data/oasis/with_autoreject.p','rb') as file:
|
||
self.X = pickle.load(file)
|
||
self.channels = ['AF3', 'F7', 'F3', 'FC5', 'T7', 'P7', 'O1', 'O2', 'P8', 'T8', 'FC6', 'F4', 'F8', 'AF4']
|
||
self.freq = 128
|
||
self.X ,self.Y= merge_dictionary(self.X)
|
||
(a,b,c) = self.X.shape
|
||
self.X = np.reshape(self.X,(a,c,b))
|
||
else:
|
||
array = np.load('preprocessed_data/oasis/without_autoreject.npz')
|
||
self.freq = 128
|
||
self.channels = ['AF3','F7','F3','FC5','T7','P7','O1','O2','P8','T8','FC6','F4','F8','AF4']
|
||
self.X = array['X']
|
||
self.Y = array['Y']
|
||
(a,b,c) = self.X.shape
|
||
self.X = np.reshape(self.X,(a,c,b))
|
||
|
||
else:
|
||
self.X = array['X']
|
||
if self.name == 'DEAP':
|
||
self.X = self.X[:,:,[1,3,2,4,7,11,13,31,29,25,21,19,20,17]] # To maintain uniformity across all datasets, only 14 electrodes selected
|
||
self.channels = ['AF3', 'F7', 'F3', 'FC5', 'T7', 'P7', 'O1', 'O2', 'P8', 'T8', 'FC6', 'F4', 'F8', 'AF4']
|
||
if self.name != 'OASIS':
|
||
self.Y = array['Y']
|
||
#self.Z = array['Z']
|
||
self.reshape_data()
|
||
def reshape_data(self):
|
||
'''
|
||
reshapes data in the format EEGExtract module expects i.e channels x timepoints x epochs
|
||
'''
|
||
|
||
(epochs,timepoints,channels) = self.X.shape
|
||
self.eegData = np.reshape(self.X,(channels,timepoints,epochs))
|
||
|
||
class features:
|
||
############################ Complexity Features #############################
|
||
#1>
|
||
@staticmethod
|
||
def ShannonRes(eegData,**args):
|
||
#Shannon Entropy
|
||
ShannonRes = eeg.shannonEntropy(eegData, bin_min=-200, bin_max=200, binWidth=2)
|
||
return ShannonRes
|
||
#2>
|
||
@staticmethod
|
||
def ShannonRes_sub_band_delta(eegData,fs):
|
||
# Subband Information Quantity
|
||
# delta (0.5–4 Hz)
|
||
eegData_delta = eeg.filt_data(eegData, 0.5, 4, fs)
|
||
ShannonRes_delta = eeg.shannonEntropy(eegData_delta, bin_min=-200, bin_max=200, binWidth=2)
|
||
return ShannonRes_delta
|
||
#3>
|
||
@staticmethod
|
||
def ShannonRes_sub_band_theta(eegData,fs):
|
||
# theta (4–8 Hz)
|
||
eegData_theta = eeg.filt_data(eegData, 4, 8, fs)
|
||
ShannonRes_theta = eeg.shannonEntropy(eegData_theta, bin_min=-200, bin_max=200, binWidth=2)
|
||
return ShannonRes_theta
|
||
|
||
#4>
|
||
@staticmethod
|
||
def ShannonRes_sub_band_alpha(eegData,fs):
|
||
# alpha (8–12 Hz)
|
||
eegData_alpha = eeg.filt_data(eegData, 8, 12, fs)
|
||
ShannonRes_alpha = eeg.shannonEntropy(eegData_alpha, bin_min=-200, bin_max=200, binWidth=2)
|
||
return ShannonRes_alpha
|
||
|
||
#5>
|
||
@staticmethod
|
||
def ShannonRes_sub_band_beta(eegData,fs):
|
||
# beta (12–30 Hz)
|
||
eegData_beta = eeg.filt_data(eegData, 12, 30, fs)
|
||
ShannonRes_beta = eeg.shannonEntropy(eegData_beta, bin_min=-200, bin_max=200, binWidth=2)
|
||
return ShannonRes_beta
|
||
|
||
#6>
|
||
@staticmethod
|
||
def ShannonRes_sub_band_gamma(eegData,fs):
|
||
# gamma (30–100 Hz)
|
||
eegData_gamma = eeg.filt_data(eegData, 30, 63, fs)
|
||
ShannonRes_gamma = eeg.shannonEntropy(eegData_gamma, bin_min=-200, bin_max=200, binWidth=2)
|
||
return ShannonRes_gamma
|
||
|
||
|
||
#7>
|
||
@staticmethod
|
||
def Hojorth_Mobility(eegData,**args):
|
||
# Hjorth Mobility
|
||
# Hjorth Complexity
|
||
HjorthMob, HjorthComp = eeg.hjorthParameters(eegData)
|
||
return HjorthMob
|
||
#8>
|
||
@staticmethod
|
||
def Hojorth_Complexity(eegData,**args):
|
||
# Hjorth Mobility
|
||
# Hjorth Complexity
|
||
HjorthMob, HjorthComp = eeg.hjorthParameters(eegData)
|
||
return HjorthComp
|
||
#9>
|
||
@staticmethod
|
||
def False_Nearest_Neighbour(eegData,**args):
|
||
# False Nearest Neighbor
|
||
FalseNnRes = eeg.falseNearestNeighbor(eegData)
|
||
return FalseNnRes
|
||
##############################################################################
|
||
|
||
########################Category Features#####################################
|
||
#10>
|
||
@staticmethod
|
||
def median_frequency(eegData,fs):
|
||
#fs-sampling frequency
|
||
# Median Frequency
|
||
medianFreqRes = eeg.medianFreq(eegData,fs)
|
||
return medianFreqRes
|
||
|
||
#11>
|
||
@staticmethod
|
||
def band_power_delta(eegData,fs):
|
||
#fs - sampling frequency
|
||
# δ band Power
|
||
bandPwr_delta = eeg.bandPower(eegData, 0.5, 4, fs)
|
||
return bandPwr_delta
|
||
#12>
|
||
@staticmethod
|
||
def band_power_theta(eegData,fs):
|
||
#fs - sampling frequency
|
||
# θ band Power
|
||
bandPwr_theta = eeg.bandPower(eegData, 4, 8, fs)
|
||
return bandPwr_theta
|
||
|
||
#13>
|
||
@staticmethod
|
||
def band_power_alpha(eegData,fs):
|
||
#fs - sampling frequency
|
||
# α band Power
|
||
bandPwr_alpha = eeg.bandPower(eegData, 8, 12, fs)
|
||
return bandPwr_alpha
|
||
|
||
#14>
|
||
@staticmethod
|
||
def band_power_beta(eegData,fs):
|
||
#fs - sampling frequency
|
||
# β band Power
|
||
bandPwr_beta = eeg.bandPower(eegData, 12, 30, fs)
|
||
return bandPwr_beta
|
||
|
||
#15>
|
||
@staticmethod
|
||
def band_power_gamma(eegData,fs):
|
||
#fs - sampling frequency
|
||
# γ band Power
|
||
bandPwr_gamma = eeg.bandPower(eegData, 30, 63, fs)
|
||
return bandPwr_gamma
|
||
|
||
#16>
|
||
@staticmethod
|
||
def standard_deviation(eegData,**args):
|
||
# Standard Deviation
|
||
std_res = eeg.eegStd(eegData)
|
||
return std_res
|
||
|
||
#17>
|
||
@staticmethod
|
||
def regularity(eegData,fs):
|
||
# Regularity (burst-suppression)
|
||
regularity_res = eeg.eegRegularity(eegData,fs)
|
||
return regularity_res
|
||
|
||
|
||
#18>
|
||
@staticmethod
|
||
def Diffuse_slowing(eegData,**args):
|
||
# Diffuse Slowing
|
||
df_res = eeg.diffuseSlowing(eegData)
|
||
return df_res
|
||
|
||
#19>
|
||
@staticmethod
|
||
def Spikes(eegData,fs,**args):
|
||
# Spikes
|
||
minNumSamples = int(70*fs/1000)
|
||
spikeNum_res = eeg.spikeNum(eegData,minNumSamples)
|
||
return spikeNum_res
|
||
|
||
#20>
|
||
@staticmethod
|
||
def delta_burst_after_spike(eegData,fs):
|
||
# Delta burst after Spike
|
||
eegData_delta = eeg.filt_data(eegData, 0.5, 4, fs)
|
||
deltaBurst_res = eeg.burstAfterSpike(eegData,eegData_delta,minNumSamples=7,stdAway = 3)
|
||
return deltaBurst_res
|
||
|
||
#21>
|
||
@staticmethod
|
||
def Sharp_spike(eegData,fs):
|
||
minNumSamples = int(70*fs/1000)
|
||
# Sharp spike
|
||
sharpSpike_res = eeg.shortSpikeNum(eegData,minNumSamples)
|
||
return sharpSpike_res
|
||
|
||
#22>
|
||
@staticmethod
|
||
def Number_of_Burst(eegData,fs):
|
||
# Number of Bursts
|
||
numBursts_res = eeg.numBursts(eegData,fs)
|
||
return numBursts_res
|
||
|
||
#23>
|
||
@staticmethod
|
||
def Burst_length_u_and_sigma_mean(eegData,fs):
|
||
# Burst length μ and σ
|
||
burstLenMean_res,burstLenStd_res = eeg.burstLengthStats(eegData,fs)
|
||
return burstLenMean_res
|
||
|
||
#24>
|
||
@staticmethod
|
||
def Burst_length_u_and_sigma_std(eegData,fs):
|
||
burstLenMean_res,burstLenStd_res = eeg.burstLengthStats(eegData,fs)
|
||
return burstLenStd_res
|
||
|
||
|
||
#25>
|
||
@staticmethod
|
||
def no_of_suprression(eegData,fs):
|
||
# Number of Suppressions
|
||
numSupps_res = eeg.numSuppressions(eegData,fs)
|
||
return numSupps_res
|
||
|
||
#26>
|
||
@staticmethod
|
||
def Suppression_length_u_and_sigma_mean(eegData,fs):
|
||
# Suppression length μ and σ
|
||
suppLenMean_res,suppLenStd_res = eeg.suppressionLengthStats(eegData,fs)
|
||
return suppLenMean_res
|
||
|
||
#27>
|
||
@staticmethod
|
||
def Suppression_length_u_and_sigma_std(eegData,fs):
|
||
# Suppression length μ and σ
|
||
suppLenMean_res,suppLenStd_res = eeg.suppressionLengthStats(eegData,fs)
|
||
return suppLenStd_res
|
||
|
||
##############################################################################
|
||
|
||
def merge_dictionary(dictionary):
|
||
'''
|
||
merge all trial data to form one array
|
||
'''
|
||
no_of_trials = len(list(dictionary.keys()))
|
||
no_of_channels = dictionary[1][0].shape[1]
|
||
length_of_segment = dictionary[1][0].shape[2]
|
||
no_of_epochs_per_trial = dictionary[1][0].shape[0]
|
||
X = np.empty((0,no_of_channels,length_of_segment))
|
||
Y = np.empty((0,2))
|
||
for trial,lst in dictionary.items():
|
||
array = dictionary[trial][0]
|
||
score = dictionary[trial][3]
|
||
X = np.append(X,array,axis = 0)
|
||
for epoch in range(no_of_epochs_per_trial):
|
||
Y = np.append(Y,np.expand_dims(score,axis =0),axis = 0)
|
||
|
||
return X,Y
|
||
|
||
def epoch_data(X,Y, window, stride, sfreq):
|
||
|
||
(channels,timepoints,trials )= X.shape
|
||
X = np.reshape(X,(trials,channels,timepoints))
|
||
segment = int(window*sfreq)
|
||
step = int(stride*sfreq)
|
||
epochPerTrial = int((timepoints-segment)/step + 1)
|
||
count = 0
|
||
X_new = np.empty((trials*epochPerTrial,channels,segment))
|
||
Y_new = np.empty((trials*epochPerTrial,2))
|
||
for trial in range(trials):
|
||
for epoch in range(epochPerTrial):
|
||
X_new[count,:,:] = X[trial,:,epoch*step:(epoch*step)+segment]
|
||
Y_new[count,:] = Y[trial,:2]
|
||
count+=1
|
||
(trials,channels,timepoints) = X_new.shape
|
||
X_new = np.reshape(X_new,(channels,timepoints,trials))
|
||
|
||
return X_new,Y_new
|
||
|
||
def driver_code():
|
||
dataset_dictionary = {0:'DEAP',1:'OASIS',2:'DREAMER'}
|
||
print(dataset_dictionary)
|
||
print('enter number for loading dataset')
|
||
mapping = int(input())
|
||
print('plz wait loading dataset preprocessed arrays')
|
||
dataset = load_data(dataset_dictionary[mapping])
|
||
if mapping == 1:
|
||
print('do you want to use with autoreject data? if yes press y')
|
||
boolean = input()
|
||
if boolean == 'y':
|
||
dataset.use_autoreject = 'y'
|
||
|
||
dataset.load_arrays()
|
||
print('loading complete')
|
||
print('shape of data we will use to make features:',dataset.eegData.shape)
|
||
print('do you want to segment the data before calculating feature values? y/n')
|
||
boolean = input()
|
||
if boolean == 'y':
|
||
window = float(input('enter window size'))
|
||
stride = float(input('enter stride size'))
|
||
dataset.eegData,dataset.Y = epoch_data(dataset.eegData,dataset.Y,window,stride,dataset.freq)
|
||
print('new shapes of X and Y:',dataset.eegData.shape,' ',dataset.Y.shape)
|
||
else:
|
||
window = 0
|
||
stride = 0
|
||
print('features available')
|
||
featuresDict = {0:'shannonEntropy',
|
||
1:'ShannonRes_sub_bands_alpha',
|
||
2:'ShannonRes_sub_bands_beta',
|
||
3:'ShannonRes_sub_bands_delta',
|
||
4:'ShannonRes_sub_bands_theta',
|
||
5:'ShannonRes_sub_bands_gamma',
|
||
6:'Hjorth_mobilty',
|
||
7:'Hjorth_complexity',
|
||
8:'falseNearestNeighbor',
|
||
9:'medianFreq',
|
||
10:'bandPwr_alpha',
|
||
11:'bandPwr_beta',
|
||
12:'bandPwr_gamma',
|
||
13:'bandPwr_theta',
|
||
14:'bandPwr_delta',
|
||
15:'stdDev',
|
||
16:'diffuseSlowing',
|
||
17:'spikeNum',
|
||
18:'deltaBurstAfterSpike',
|
||
19:'shortSpikeNum',
|
||
20:'Sharp spike',
|
||
21:'numBursts',
|
||
22:'burstLen_u_and_sigma_mean',
|
||
23:'burstLen_u_and_sigma_std',
|
||
24:'numSuppressions',
|
||
25:'suppressionLen_u_and_sigma_mean',
|
||
26:'suppressionLen_u_and_sigma_std',
|
||
}
|
||
featureMethod={0:features.ShannonRes,
|
||
1:features.ShannonRes_sub_band_alpha,
|
||
2:features.ShannonRes_sub_band_beta,
|
||
3:features.ShannonRes_sub_band_delta,
|
||
4:features.ShannonRes_sub_band_theta,
|
||
5:features.ShannonRes_sub_band_gamma,
|
||
6:features.Hojorth_Mobility,
|
||
7:features.Hojorth_Complexity,
|
||
8:features.False_Nearest_Neighbour,
|
||
9:features.median_frequency,
|
||
10:features.band_power_alpha,
|
||
11:features.band_power_beta,
|
||
12:features.band_power_gamma,
|
||
13:features.band_power_theta,
|
||
14:features.band_power_delta,
|
||
15:features.standard_deviation,
|
||
16:features.regularity,
|
||
17:features.Diffuse_slowing,
|
||
18:features.Spikes,
|
||
19:features.delta_burst_after_spike,
|
||
20:features.Sharp_spike,
|
||
21:features.Number_of_Burst,
|
||
22:features.Burst_length_u_and_sigma_mean,
|
||
23:features.Burst_length_u_and_sigma_std,
|
||
24:features.no_of_suprression,
|
||
25:features.Suppression_length_u_and_sigma_mean,
|
||
26:features.Suppression_length_u_and_sigma_std,
|
||
}
|
||
|
||
|
||
print(featuresDict)
|
||
|
||
#define path for saving before hand in np.savez line below
|
||
path = 'features/'
|
||
#os.mkdir('features/'+window+'_'+stride)
|
||
if dataset.name == 'DEAP':
|
||
path = path +'deap/'
|
||
elif dataset.name == 'DREAMER':
|
||
path = path + 'dreamer/'
|
||
else:
|
||
if dataset.use_autoreject == 'y':
|
||
path = path +'oasis/with_autoreject/'
|
||
else:
|
||
path = path +'oasis/without_autoreject/'
|
||
boolean = input('do you want to individually make features? y/n')
|
||
if boolean =='n':
|
||
for key in featureMethod.keys():
|
||
feature_matrix = featureMethod[key](eegData = dataset.eegData,fs=dataset.freq)
|
||
filename = featuresDict[key]
|
||
print('saving ---',filename)
|
||
np.savez(path+filename+'_'+str(int(window))+'_'+str(int(stride)),features = feature_matrix , Y = dataset.Y)
|
||
else:
|
||
found_features = False
|
||
while not found_features:
|
||
print('enter feature no')
|
||
key = int(input())
|
||
feature_matrix = featureMethod[key](eegData = dataset.eegData,fs=dataset.freq)
|
||
filename = featuresDict[key]
|
||
print('saving ---',filename)
|
||
np.savez(path+filename+'_'+str(int(window))+'_'+str(int(stride)),features = feature_matrix , Y = dataset.Y)
|
||
boolean = input('do you want to find more features? y/n ')
|
||
if boolean =='n':
|
||
found_features = True
|
||
|
||
|
||
print('feature extraction done!!!!')
|
||
|
||
def __main__():
|
||
driver_code()
|
||
|
||
__main__()
|
||
|
||
if __name__ == 'main':
|
||
driver_code()
|
||
|
||
|
||
|