{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":59093,"databundleVersionId":7469972,"sourceType":"competition"},{"sourceId":7392733,"sourceType":"datasetVersion","datasetId":4297749}],"dockerImageVersionId":30646,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import pandas as pd, numpy as np, os\nimport matplotlib.pyplot as plt, gc\n\ntrain = pd.read_csv('/kaggle/input/hms-harmful-brain-activity-classification/train.csv')\nprint('Train shape', train.shape )\ndisplay( train.head() )","metadata":{"execution":{"iopub.status.busy":"2024-02-12T05:23:41.643503Z","iopub.execute_input":"2024-02-12T05:23:41.644278Z","iopub.status.idle":"2024-02-12T05:23:43.488094Z","shell.execute_reply.started":"2024-02-12T05:23:41.644240Z","shell.execute_reply":"2024-02-12T05:23:43.486637Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"NAMES = ['FL','FR','PL','PR','C']\n\nFEATS = [['Fp1','F7','F3','C3','T3'],\n             ['Fp2','F8','F4','C4','T4'],\n             ['O1','T5','P3','C3','T3'],\n             ['O2','T6','P4','C4','T4'],\n             ['Fz','Cz','Pz'],]\n\ndirectory_path_add = '/second/path/EEG_ADD_Spectrograms/'\nif not os.path.exists(directory_path_add):\n    os.makedirs(directory_path_add)   ","metadata":{"execution":{"iopub.status.busy":"2024-02-12T05:28:11.283148Z","iopub.execute_input":"2024-02-12T05:28:11.283551Z","iopub.status.idle":"2024-02-12T05:28:11.291312Z","shell.execute_reply.started":"2024-02-12T05:28:11.283524Z","shell.execute_reply":"2024-02-12T05:28:11.289941Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pywt\nprint(\"The wavelet functions we can use:\")\nprint(pywt.wavelist())\n\nUSE_WAVELET = None","metadata":{"execution":{"iopub.status.busy":"2024-02-12T05:28:14.047696Z","iopub.execute_input":"2024-02-12T05:28:14.048160Z","iopub.status.idle":"2024-02-12T05:28:14.054946Z","shell.execute_reply.started":"2024-02-12T05:28:14.048124Z","shell.execute_reply":"2024-02-12T05:28:14.053701Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# DENOISE FUNCTION\ndef maddest(d, axis=None):\n    return np.mean(np.absolute(d - np.mean(d, axis)), axis)\n\ndef denoise(x, wavelet='haar', level=1):    \n    coeff = pywt.wavedec(x, wavelet, mode=\"per\")\n    sigma = (1/0.6745) * maddest(coeff[-level])\n\n    uthresh = sigma * np.sqrt(2*np.log(len(x)))\n    coeff[1:] = (pywt.threshold(i, value=uthresh, mode='hard') for i in coeff[1:])\n\n    ret=pywt.waverec(coeff, wavelet, mode='per')\n    \n    return ret","metadata":{"execution":{"iopub.status.busy":"2024-02-12T05:28:14.771121Z","iopub.execute_input":"2024-02-12T05:28:14.771550Z","iopub.status.idle":"2024-02-12T05:28:14.779978Z","shell.execute_reply.started":"2024-02-12T05:28:14.771517Z","shell.execute_reply":"2024-02-12T05:28:14.778957Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import librosa\n\ndef spectrogram_from_eeg_add(parquet_path, display=False):\n    \n    # LOAD MIDDLE 50 SECONDS OF EEG SERIES\n    eeg = pd.read_parquet(parquet_path)\n    middle = (len(eeg)-10_000)//2\n    eeg = eeg.iloc[middle:middle+10_000]\n    \n    # VARIABLE TO HOLD SPECTROGRAM\n    img = np.zeros((128,256,5),dtype='float32')\n    \n    if display: plt.figure(figsize=(10,7))\n    signals = []\n    for k in range(5):\n        COLS = FEATS[k]\n        if k < 4:\n            for kk in range(4):\n        \n                # COMPUTE PAIR DIFFERENCES\n                x = eeg[COLS[kk]].values - eeg[COLS[kk+1]].values\n\n                # FILL NANS\n                m = np.nanmean(x)\n                if np.isnan(x).mean()<1: x = np.nan_to_num(x,nan=m)\n                else: x[:] = 0\n\n                # DENOISE\n                if USE_WAVELET:\n                    x = denoise(x, wavelet=USE_WAVELET)\n                signals.append(x)\n\n                # RAW SPECTROGRAM\n                mel_spec = librosa.feature.melspectrogram(y=x, sr=200, hop_length=len(x)//256, \n                     n_fft=1024, n_mels=128, fmin=0, fmax=20, win_length=128)\n\n                # LOG TRANSFORM\n                width = (mel_spec.shape[1]//32)*32\n                mel_spec_db = librosa.power_to_db(mel_spec, ref=np.max).astype(np.float32)[:,:width]\n\n                # STANDARDIZE TO -1 TO 1\n                mel_spec_db = (mel_spec_db+40)/40 \n                img[:,:,k] += mel_spec_db\n                \n        # AVERAGE THE 4 MONTAGE DIFFERENCES\n            img[:,:,k] /= 4.0\n            \n        if k == 4:\n            for kk in range(2):\n        \n                # COMPUTE PAIR DIFFERENCES\n                x = eeg[COLS[kk]].values - eeg[COLS[kk+1]].values\n\n                # FILL NANS\n                m = np.nanmean(x)\n                if np.isnan(x).mean()<1: x = np.nan_to_num(x,nan=m)\n                else: x[:] = 0\n\n                # DENOISE\n                if USE_WAVELET:\n                    x = denoise(x, wavelet=USE_WAVELET)\n                signals.append(x)\n\n                # RAW SPECTROGRAM\n                mel_spec = librosa.feature.melspectrogram(y=x, sr=200, hop_length=len(x)//256, \n                     n_fft=1024, n_mels=128, fmin=0, fmax=20, win_length=128)\n\n                # LOG TRANSFORM\n                width = (mel_spec.shape[1]//32)*32\n                mel_spec_db = librosa.power_to_db(mel_spec, ref=np.max).astype(np.float32)[:,:width]\n\n                # STANDARDIZE TO -1 TO 1\n                mel_spec_db = (mel_spec_db+40)/40 \n                img[:,:,k] += mel_spec_db\n                \n        # AVERAGE THE 4 MONTAGE DIFFERENCES\n            img[:,:,k] /= 2.0    \n            \n            \n        \n        if display:\n            plt.subplot(3,2,k+1)\n            plt.imshow(img[:,:,k],aspect='auto',origin='lower')\n            plt.title(f'EEG {eeg_id} - Spectrogram {NAMES[k]}')\n            \n    if display: \n        plt.show()\n        plt.figure(figsize=(10,5))\n        offset = 0\n        for k in range(5):\n            if k>0: offset -= signals[4-k].min()\n            plt.plot(range(10_000),signals[k]+offset,label=NAMES[4-k])\n            offset += signals[4-k].max()\n        plt.legend()\n        plt.title(f'EEG {eeg_id} Signals')\n        plt.show()\n        print(); print('#'*25); print()\n        \n    return img","metadata":{"execution":{"iopub.status.busy":"2024-02-12T05:28:15.955507Z","iopub.execute_input":"2024-02-12T05:28:15.955996Z","iopub.status.idle":"2024-02-12T05:28:15.980174Z","shell.execute_reply.started":"2024-02-12T05:28:15.955958Z","shell.execute_reply":"2024-02-12T05:28:15.978740Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nPATH = '/kaggle/input/hms-harmful-brain-activity-classification/train_eegs/'\nDISPLAY = 4\nEEG_IDS = train.eeg_id.unique()\nall_eegs_add = {}\n\nfor i,eeg_id in enumerate(EEG_IDS):\n    if (i%100==0)&(i!=0): print(i,', ',end='')\n        \n    # CREATE SPECTROGRAM FROM EEG PARQUET\n    img = spectrogram_from_eeg_add(f'{PATH}{eeg_id}.parquet', i<DISPLAY)\n    \n    # SAVE TO DISK\n    if i==DISPLAY:\n        print(f'Creating and writing {len(EEG_IDS)} spectrograms to disk... ',end='')\n    np.save(f'{directory_path_add}{eeg_id}',img)\n    all_eegs_add[eeg_id] = img\n   \nnp.save('eeg_specs_add',all_eegs_add) ","metadata":{"execution":{"iopub.status.busy":"2024-02-12T05:28:16.886883Z","iopub.execute_input":"2024-02-12T05:28:16.887277Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}