{"cells":[{"metadata":{},"cell_type":"markdown","source":"## Some words on this notebook\n\nUnder the pretext of 'analysing' a bit the EEG files, I created this notebook which in reality is more a playground for me to get familiar with the data. For this prupose I wanted to investigate how to extract the information from the .mat files, how to plot the data of the EEG files, etc.\nI also produced some artificial data and tried to find patterns that could be used to make a prediction.\nI also wanted to see if it was possible to use the image of an espectogram for the same purposes.\nAll this without understanding really anything about EEG.\n\nAfter a bit of research I realized it could have been interesting to know what kind of montage was used to produce these EEG values, so that we could train an algorithm to identify which electrodes (and the value differences between them) are more likely to yield a good prediction. However, after looking at the EEG plot (below) it seems that the EEGs were produced using one common reference for all electrode measurements.\n\n\nSo, I don't want to spend too much time on getting educated on the vast world of neuro science interpretations. Since I wasn't able to -very- easily find something visually useful, I will stop this notebook and instead I will go directly into trying to train an RNN algorithm for predictions."},{"metadata":{"trusted":true},"cell_type":"code","source":"!rm -r '/kaggle/working/train/'","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#Visualizing an example:\nfrom scipy.io import loadmat\nimport os\nimport numpy as np\nimport matplotlib.pyplot as plt\ninterictal_tst = '/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_interictal_segment_0001.mat'\ninterictal_data = loadmat(interictal_tst)          # loads the file as a dictionary\nfor item in interictal_data.items():\n    print(item)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"segment_tst = '/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_test_segment_0001.mat'\nsegment_data = loadmat(segment_tst)          # loads the file as a dictionary\nfor item in segment_data.items():\n    print(item)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"files= [x for x in os.listdir('/kaggle/working/Patient1/train/preictal/')]\nfiles_i=[x for x in os.listdir('/kaggle/working/Patient1/train/interictal/')]\nlen(files),len(files_i)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"interictal_data.get('interictal_segment_1')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Check for data inconsistencies, like for example, that there are as much data samples needed as data seconds and that these two correspond to the data samples per seconds."},{"metadata":{"trusted":true},"cell_type":"code","source":"# search_key = \"_segment_\"       # search key string\n\n# for dirname, _, filenames in os.walk('/kaggle/input/seizure-prediction/'):\n#     for filename in filenames:\n#         if filename != \"sampleSubmission.csv\":\n#             try:\n#                 segment_dict = loadmat(os.path.join(dirname, filename))        # loads the .mat file as a dictionary\n#                 segment = dict(filter(lambda item: search_key in item[0], segment_dict.items()))  # retrieves only the item which contain the data of interes\n#                 segment_list_of_values = list(segment.values())\n#                 total_num_dataSamples = np.array(segment_list_of_values[0][0][0][0].shape[1], dtype=np.float)\n#                 datarows_duration_secs = np.array(segment_list_of_values[0][0][0][1], dtype=np.float)\n#                 datasamples_perSecond = np.array(segment_list_of_values[0][0][0][2], dtype=np.float)\n\n#                 if total_num_dataSamples != np.around((datarows_duration_secs * datasamples_perSecond)):\n#                     print(\"Data inconsistency in file: %s\" % os.path.join(dirname, filename))\n#             except:\n#                 e = sys.exc_info()[0]\n#                 print(\"Exception %s caught on file %s \" % (e, os.path.join(dirname, filename)))\n#                 break\n# print(\"Check done.\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Does all the files of the same patient or dog (i.e. Patient_1 or Patient_2 or Dog_1, etc.) have the same number of electrodes?\nAnd do they all have the same names?"},{"metadata":{"trusted":true},"cell_type":"code","source":"# search_key = \"_segment_\"       # search key string\n# subdirs = [\n#     '/kaggle/input/seizure-prediction/Patient_1/Patient_1',\n#     '/kaggle/input/seizure-prediction/Patient_2/Patient_2',\n#     '/kaggle/input/seizure-prediction/Dog_1/Dog_1',\n#     '/kaggle/input/seizure-prediction/Dog_2/Dog_2',\n#     '/kaggle/input/seizure-prediction/Dog_3/Dog_3',\n#     '/kaggle/input/seizure-prediction/Dog_4/Dog_4',\n#     '/kaggle/input/seizure-prediction/Dog_5/Dog_5',\n# ]\n\n# for dirname in subdirs:\n#     for dirname, _, filenames in os.walk(dirname):\n#         # auxiliar variables needed to monitor both the\n#         # number of electrodes and electrode names between\n#         # the last file analyzed and the current file analyzed\n#         numElectrodes_lastFile = 0\n#         firstComparisonDone = False\n#         lastFileName = \"\"\n#         electrodeNames_lastFile = []\n#         for filename in filenames:\n#             matFile = scipy.io.loadmat(os.path.join(dirname, filename))\n#             segment = dict(filter(lambda item: search_key in item[0], matFile.items()))  # retrieves only the item which contain the data of interes\n#             segment_list_values = list(segment.values())\n#             # number of electrodes on the currently analyzed file\n#             numElectrodes_currentFile = segment_list_values[0][0][0][3][0].shape[0]\n#             # next two lines of code are relevant only if during the\n#             # first cycle of this loop\n#             numElectrodes_lastFile = numElectrodes_currentFile if True is not firstComparisonDone else numElectrodes_lastFile\n#             firstComparisonDone = True if True != firstComparisonDone else firstComparisonDone\n#             # check if number of electrodes between the current file and the last\n#             # file analyzed are equal\n#             if numElectrodes_currentFile != numElectrodes_lastFile:\n#                 print(\"Discrepancy in number of electrodes between files:\")\n#                 print(\"     Last file    %s with %d electrodes\" % (lastFileName, numElectrodes_lastFile))\n#                 print(\"     Current file %s with %d electrodes\" % (filename, numElectrodes_currentFile))\n#             numElectrodes_lastFile = numElectrodes_currentFile\n#             # collect the name of the electrodes on the current file being analyzed\n#             electrode_names = [name[0] for name in segment_list_values[0][0][0][3][0]]\n#             if len(electrodeNames_lastFile) > 0:\n#                 # check if the electrode names of the last file analyzed and \n#                 # the current file analyzed are equal\n#                 commonNames = set(electrodeNames_lastFile) & set(electrode_names)\n#                 if len(commonNames) != len(electrode_names):\n#                     print(\"Different electrode names on files: \")\n#                     print(\"     %s (last faile name)\" % lastFileName)\n#                     print(\"         Electrode names: \", electrodeNames_lastFile)\n#                     print(\"     %s \" % filename)\n#                     print(\"         Electrode names: \", electrode_names)\n#             electrodeNames_lastFile = electrode_names\n#             lastFileName = filename\n# print(\"Check done\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Here's an example of how an eeg looks like."},{"metadata":{"trusted":true},"cell_type":"code","source":"# Code in this cell was partially taken from:\n# https://matplotlib.org/3.3.1/gallery/specialty_plots/mri_with_eeg.html\n\nimport matplotlib.pyplot as plt\nfrom matplotlib.collections import LineCollection\nimport numpy as np\n\nfile_to_inspect = loadmat('/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_interictal_segment_0001.mat')\n# retrieves only the item which contain the data of interese\nsearch_key = '_segment_'\nsegment = dict(filter(lambda item: search_key in item[0], file_to_inspect.items()))\nsegment = list(segment.values())\ndata = segment[0][0][0][0]\nnum_electrodes = data.shape[0]             # rows (i.e. electrodes) in the data matrix\nn_samples = data.shape[1]                  # number of samples on each row (i.e. electrode's samples)\nelectrode_names = segment[0][0][0][3][0]   # name or labels of the electrodes\nt = 10 * np.arange(n_samples) / n_samples\n\n# create the 'figure'\nfig = plt.figure(\"EEG samples\",figsize=(30,15))\n\nticklocs = []\nax = fig.add_subplot(1, 1, 1)\nax.set_xlim(0, 10)\nax.set_xticks(np.arange(10))\ndata_min = data.min()\ndata_max = data.max()\ndr = (data_max - data_min) * 0.7   #crowd it a bit\ny0 = data_min\ny1 = (num_electrodes -1) * dr + data_max\nax.set_ylim(y0, y1)\n\nsegs = []\nfor electrode in range(num_electrodes):\n    segs.append(np.column_stack((t, data[electrode, :])))\n    ticklocs.append(electrode * dr)\n    \noffsets = np.zeros((num_electrodes, 2), dtype=float)\noffsets[:, 1] = ticklocs\n\nlines = LineCollection(segs, offsets=offsets, transOffset=None)\nax.add_collection(lines)\n\n# Set the yticks to use axes coordinates on the y axis\nax.set_yticks(ticklocs)\nax.set_yticklabels(electrode_names)\n\nax.set_xlabel('Time (minutes)')\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport copy\nimport scipy.signal as signal\nimport scipy.stats as stats\nimport scipy.io as sio\nimport tqdm\n\nclass Dataset:\n    def __init__(self,path):\n        self.path = path\n        if self.path[-1] != '/':\n            self.path += '/'\n        self.files = np.concatenate([[self.path+'/interictal/'+x for x in os.listdir(self.path+'/interictal/')], [self.path+'/preictal/'+x for x in os.listdir(self.path+'/preictal/')]])\n        self.N_interictal = 15000\n        self.N_preictal = 5400\n        self.targets=np.array([0 if x < self.N_interictal else 1 for x in range(len(self.files))])\n\n    def __len__(self):\n        return len(self.files)\n\n    def __getitem__(self, item):\n        data = np.load(self.files[item])\n        target = self.targets[item]\n        data = np.expand_dims(data,axis=0)\n        return data,target\n\n    def split_random(self,N_valid):\n        sample = np.random.choice(len(self), len(self), replace=False)\n        self.files, self.targets = self.files[sample], self.targets[sample]\n        train = copy.deepcopy(self)\n        valid = copy.deepcopy(self)\n        train.files, train.targets = train.files[N_valid:], train.targets[N_valid:]\n        valid.files, valid.targets = valid.files[:N_valid], valid.targets[:N_valid]\n        return train,valid","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def epoch_raw(raw_path, secs=30):\n    raw_files_interictal = sorted([f for f in os.listdir(raw_path) if 'interictal' in f])\n    raw_files_preictal = sorted([f for f in os.listdir(raw_path) if 'preictal' in f])\n    raw_files = {'interictal': raw_files_interictal, 'preictal': raw_files_preictal}\n    pat_dir = '/kaggle/working/Patient1/'\n    train_dir = pat_dir+'/train/'\n    train_dir_interictal = train_dir+'/interictal/'\n    train_dir_preictal = train_dir+'/preictal/'\n    if not os.path.isdir(pat_dir):\n        os.mkdir(pat_dir)\n    if not os.path.isdir(train_dir):\n        os.mkdir(train_dir)\n    if not os.path.isdir(train_dir_interictal):\n        os.mkdir(train_dir_interictal)\n    if not os.path.isdir(train_dir_preictal):\n        os.mkdir(train_dir_preictal)\n    train_dirs = {'interictal': train_dir_interictal, 'preictal': train_dir_preictal}\n    \n    for category in ['preictal']:\n        for _,f in enumerate(raw_files[category]):\n            search_key = \"_segment_\"       # search key string\n            try:\n                matFile = loadmat(os.path.join(raw_path,f))\n            except OSError as error:\n                print(error)\n            segment = dict(filter(lambda item: search_key in item[0], matFile.items()))\n            segment_list_values = list(segment.values())\n            x = np.array(segment_list_values[0][0][0][0], dtype=np.float)\n            Fs = np.array(segment_list_values[0][0][0][2], dtype=np.float)[0][0]\n            L = int(Fs*secs)\n            for electrode in range(x.shape[0]):\n                for i in range(L,x.shape[1]+1,L):\n                    digit = f'_{electrode}_{int(i//L)}'\n                    save_filepath = train_dirs[category]+f[:-4]+digit+'.npy'\n                    spect,_,_,_ = plt.specgram(x[electrode, i-L:i], NFFT=1024, Fs=5000, noverlap=128)\n                    spect = stats.zscore(spect[:200,:],axis=1)\n                    try:\n                        np.save(save_filepath, spect)\n                    except OSError as e:\n                        print(e)\n            print(f\"Finished file {f}\")\n        print(f\"Finished category {category}\")   \n    return\n\npatient_dir = '/kaggle/input/seizure-prediction/Patient_1/Patient_1'\nepoch_raw(patient_dir, 30)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import torch.nn as nn\nimport torch.nn.functional as F\nimport torch.optim as optim\n\nclass NN(nn.Module):\n    def __init__(self,NFILT=256,NOUT=2):\n        super(NN,self).__init__()\n        self.conv0 = nn.Conv2d(1,NFILT,kernel_size=(200,3),padding=(0,1),bias=False)\n        self.bn0 = nn.BatchNorm2d(NFILT)\n        self.gru = nn.GRU(input_size=NFILT,hidden_size=128,num_layers=1,batch_first=True,bidirectional=False)\n        self.fc1 = nn.Linear(128,NOUT)\n\n\n\n    def forward(self, x):\n        x = F.relu(self.bn0(self.conv0(x)))\n        x = x.squeeze().permute(0,2,1)\n        x,_ = self.gru(x)\n        x = F.dropout(x,p=0.5,training=self.training)\n        x = self.fc1(x)\n        return x","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import torch\nimport numpy as np\nfrom scipy.special import softmax,expit\nfrom sklearn.metrics import f1_score,confusion_matrix,cohen_kappa_score,roc_curve,roc_auc_score,average_precision_score\n\n\nclass Statistics(object):\n    def __init__(self):\n        self.target = []\n        self.logits = []\n\n    def reset(self):\n        self.target = []\n        self.logits = []\n\n    @staticmethod\n    def idx2onehot(idx_array):\n        y = np.zeros((idx_array.shape[0], idx_array.max() + 1))\n        y[np.arange(y.shape[0]), idx_array] = 1\n        return y\n\n    @staticmethod\n    def F1(conf):\n        x0 = np.sum(conf, 0)\n        x1 = np.sum(conf, 1)\n        dg = np.diag(conf)\n        f1 = 2 * dg / (x0 + x1)\n        return f1\n\n    @staticmethod\n    def Kappa(conf):\n        x0 = np.sum(conf, 0)\n        x1 = np.sum(conf, 1)\n        N = np.sum(np.sum(conf))\n        ef = np.sum(x0 * x1 / N)\n        dg = np.sum(np.diag(conf))\n        K = (dg - ef) / (N - ef)\n        return K\n\n    def append(self,target,logits):\n        self.logits.append(logits.data.cpu().numpy())\n        self.target.append(target.data.cpu().numpy())\n\n    @staticmethod\n    def random_auprc(target):\n        y_chance = np.zeros((target.max()+1,))\n        for i in range(target.max()+1):\n            y_chance[i] = len(target[target==i]) / len(target)\n\n        return y_chance\n\n\n    def evaluate(self):\n        self.logits = np.concatenate(self.logits)\n        self.target = np.concatenate(self.target).astype('int32')\n\n        self.probs = softmax(self.logits,axis=1)\n        self.argmax = np.argmax(self.probs,axis=1)\n\n        CONF = np.array(confusion_matrix(y_true=self.target,y_pred=self.argmax))\n        F1 = Statistics.F1(CONF)\n        KPS = Statistics.Kappa(CONF)\n        AUROC = roc_auc_score(y_true=Statistics.idx2onehot(self.target),y_score=self.probs,average=None)\n        AUPRC = average_precision_score(y_true=Statistics.idx2onehot(self.target),y_score=self.probs,average=None)\n        AUPRC_chance = self.random_auprc(self.target)\n\n        print(CONF)\n        print(F1)\n        print(KPS)\n        print(AUROC,np.mean(AUROC))\n        print(AUPRC,np.mean(AUPRC))\n        print(AUPRC_chance)\n\n        self.reset()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import torch\nfrom torch.utils.data import DataLoader\n\nN = 50*300+18*300\ntrain_dataset, valid_dataset = Dataset('/kaggle/working/Patient1/train/').split_random(int(0.3*N))\nprint(\"N Train, valid: \", len(train_dataset),len(valid_dataset))\nNWORKERS = 24\nDEVICE = 'cuda' if torch.cuda.is_available() else 'cpu'\n\nTRAIN = DataLoader(dataset=train_dataset,\n                   batch_size=32,\n                   shuffle=True,\n                   drop_last=False,\n                   num_workers=NWORKERS)\n\nVALID = DataLoader(dataset=valid_dataset,\n                   batch_size=32,\n                   shuffle=True,\n                   drop_last=False,\n                   num_workers=NWORKERS)\n\n\nmodel = NN(2).to(DEVICE)\noptimizer = optim.Adam(model.parameters(),lr=1e-3,weight_decay=1e-4)\nloss = nn.CrossEntropyLoss()\nstatistics = Statistics()\nN_epochs = 30\nepoch_loss={'train':np.zeros((N_epochs)),'val':np.zeros((N_epochs))}\nfor epoch in range(N_epochs):\n    model.train()\n    for i,(x,t) in enumerate(TRAIN):\n        optimizer.zero_grad()\n        x = x.to(DEVICE).float()\n        t = t.to(DEVICE).long()\n        y = model(x)\n        J = loss(input=y[:,-1,:],target=t)\n        epoch_loss['train'][epoch] += J.data.cpu().numpy()\n        J.backward()\n        optimizer.step()\n\n        if i%50==0:\n            print('EPOCH:{}\\tITER:{}\\tLOSS:{}'.format(str(epoch).zfill(2),\n                                                      str(i).zfill(5),\n                                                      J.data.cpu().numpy()))\n    epoch_loss['train'][epoch]/=len(TRAIN)\n        \n    # evaluate results for validation set\n    model.eval()\n    for i,(x,t) in enumerate(VALID):\n        x = x.to(DEVICE).float()\n        t = t.to(DEVICE).long()\n        y = model(x)\n        val_loss = loss(input=y[:,-1,:],target=t)\n        epoch_loss['val'][epoch] += val_loss.data.cpu().numpy()\n        statistics.append(target=t,logits=y[:,-1,:])\n    epoch_loss['val'][epoch]/=len(VALID)\n    statistics.evaluate()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"fig=plt.figure()\nplt.title(\"Epoch Loss\")\nplt.plot(epoch_loss['train'])\nplt.plot(epoch_loss['val'])\nplt.legend(['train', 'val'])\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"torch.save(model, \"/kaggle/working/Patient1/cnngru.pth\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def saveSpectrograms(file, electrode):\n    search_key = \"_segment_\"       # search key string\n    matFile = loadmat(os.path.join(file))\n    segment = dict(filter(lambda item: search_key in item[0], matFile.items()))\n    segment_list_values = list(segment.values())\n    x = np.array(segment_list_values[0][0][0][0][electrode], dtype=np.float)\n    Fs = np.array(segment_list_values[0][0][0][2], dtype=np.float)[0][0]\n    L = int(Fs*3)\n    NFFT = 1024  # the length of the windowing segments\n    for i in np.arange(L,x.shape[0]+1,L):\n        _,_, data = signal.spectrogram(x[i-L:i], fs=Fs,nperseg=256,noverlap=128,nfft=NFFT)\n        np.save()\n\n\ndef plotSpectograms(file, electrode=0):\n    search_key = \"_segment_\"       # search key string\n    matFile = loadmat(os.path.join(file))\n    segment = dict(filter(lambda item: search_key in item[0], matFile.items()))\n    segment_list_values = list(segment.values())\n\n    x = np.array(segment_list_values[0][0][0][0][electrode], dtype=np.float)\n    x=x[:15000]\n    #secs = np.array(segment_list_values[0][0][0][1], dtype=np.float)\n    Fs = np.array(segment_list_values[0][0][0][2], dtype=np.float)[0][0]\n    print(\"Fs,\",Fs)\n    NFFT = 1024  # the length of the windowing segments\n    dt = 1/Fs\n    t = np.arange(0.0, 3.0, dt)\n\n    fig, (ax1, ax2) = plt.subplots(nrows=2, figsize=(15,10))\n    ax1.plot(t, x)\n    Pxx, freqs, bins, im = ax2.specgram(x, NFFT=NFFT, Fs=Fs, noverlap=128)\n    print(Pxx.shape)\n    print(freqs.shape)\n    # The `specgram` method returns 4 objects. They are:\n    # - Pxx: the periodogram\n    # - freqs: the frequency vector\n    # - bins: the centers of the time bins\n    # - im: the .image.AxesImage instance representing the data in the plot\n    #plt.show()\n    \nfile = [\n    '/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_interictal_segment_0001.mat',\n    '/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_preictal_segment_0001.mat',\n        '/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_preictal_segment_0011.mat',\n    '/kaggle/input/seizure-prediction/Patient_1/Patient_1/Patient_1_interictal_segment_0023.mat'\n\n    #'/kaggle/input/seizure-prediction/Dog_3/Dog_3/Dog_3_interictal_segment_0005.mat',\n   # '/kaggle/input/seizure-prediction/Dog_3/Dog_3/Dog_3_preictal_segment_0005.mat'\n]\nfor f in file:\n    plotSpectograms(f, electrode=5)\n# 513, 3348","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"search_key = \"_segment_\"       # search key string\nsubdirs = [\n    '/kaggle/input/seizure-prediction/Patient_1/Patient_1'\n]\n\nfor dirname in subdirs:\n    for dirname, _, filenames in os.walk(dirname):\n        print(len(filenames))\n        print(\"N test:\", len([f for f in filenames if \"test\" in f]))\n        print(\"N interictal:\", len([f for f in filenames if \"interictal\" in f]))\n        print(\"N preictal:\", len([f for f in filenames if \"preictal\" in f]))\n        numElectrodes_lastFile = 0\n        firstComparisonDone = False\n        lastFileName = \"\"\n        electrodeNames_lastFile = []\n        for filename in filenames:\n            print(filename)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pwd","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}