{"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"}],"dockerImageVersionId":30698,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Библиотеки","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom functools import cache\nfrom typing import Literal\nimport matplotlib.pyplot as plt\nimport mne\n\nimport keras_cv\nimport keras\nfrom keras import ops\nimport tensorflow as tf\n\nimport os\nimport glob\nfrom pathlib import Path\n\nimport pywt\nimport torch","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-05-16T05:05:02.872850Z","iopub.execute_input":"2024-05-16T05:05:02.873312Z","iopub.status.idle":"2024-05-16T05:05:02.879612Z","shell.execute_reply.started":"2024-05-16T05:05:02.873260Z","shell.execute_reply":"2024-05-16T05:05:02.878237Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Импорт данных","metadata":{}},{"cell_type":"code","source":"path = '/kaggle/input/hms-harmful-brain-activity-classification'\n\ndef load_eeg(eeg_id, data_type:Literal['train','test']='train'):\n    eeg_df = pd.read_parquet(f'{path}/{data_type}_eegs/{eeg_id}.parquet')\n    return eeg_df","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:05:02.881529Z","iopub.execute_input":"2024-05-16T05:05:02.882216Z","iopub.status.idle":"2024-05-16T05:05:02.896797Z","shell.execute_reply.started":"2024-05-16T05:05:02.882160Z","shell.execute_reply":"2024-05-16T05:05:02.893881Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df = pd.read_csv(f'{path}/train.csv')\ntrain_df","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:05:02.899886Z","iopub.execute_input":"2024-05-16T05:05:02.900293Z","iopub.status.idle":"2024-05-16T05:05:03.032296Z","shell.execute_reply.started":"2024-05-16T05:05:02.900265Z","shell.execute_reply":"2024-05-16T05:05:03.031075Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"eeg_df = load_eeg(train_df.iloc[0,0])\neeg_df\n","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:05:03.033621Z","iopub.execute_input":"2024-05-16T05:05:03.033974Z","iopub.status.idle":"2024-05-16T05:05:03.066357Z","shell.execute_reply.started":"2024-05-16T05:05:03.033951Z","shell.execute_reply":"2024-05-16T05:05:03.065127Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"info = mne.create_info(\n    eeg_df.columns.to_list(),\n    ch_types=([\"eeg\"]*(len(eeg_df.columns)-1))+['ecg'],\n    sfreq=200\n)\n\ninfo.set_montage(\"standard_1020\")\n\ndata_values = eeg_df.values.T\nraw = mne.io.RawArray(data_values*1e-6, info)\n\nraw.plot_sensors(show_names=True)\nplt.tight_layout()\n\n# график для простмотра и отчистки данных\nraw.plot(show_scrollbars=False, show_scalebars=False)\nplt.tight_layout()","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:09.008621Z","iopub.execute_input":"2024-05-16T05:07:09.008968Z","iopub.status.idle":"2024-05-16T05:07:11.053694Z","shell.execute_reply.started":"2024-05-16T05:07:09.008946Z","shell.execute_reply":"2024-05-16T05:07:11.052772Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"mne.set_config('MNE_BROWSE_RAW_SIZE','16,8')\nraw.plot() #start=, duration=  если поскролить хотим то вписываем это\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:19.497282Z","iopub.execute_input":"2024-05-16T05:07:19.497632Z","iopub.status.idle":"2024-05-16T05:07:20.317058Z","shell.execute_reply.started":"2024-05-16T05:07:19.497609Z","shell.execute_reply":"2024-05-16T05:07:20.315459Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"raw_filtered = raw.notch_filter(50).filter(0.1, 45)\nraw_filtered.plot(start=20, duration=10)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:25.444689Z","iopub.execute_input":"2024-05-16T05:07:25.445130Z","iopub.status.idle":"2024-05-16T05:07:26.369293Z","shell.execute_reply.started":"2024-05-16T05:07:25.445103Z","shell.execute_reply":"2024-05-16T05:07:26.368009Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"raw_filtered = raw.copy().filter(l_freq=1, h_freq=70,).notch_filter(60, picks='eeg')\nraw_filtered.plot(start=20, duration=10)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:31.262841Z","iopub.execute_input":"2024-05-16T05:07:31.263225Z","iopub.status.idle":"2024-05-16T05:07:32.140307Z","shell.execute_reply.started":"2024-05-16T05:07:31.263195Z","shell.execute_reply":"2024-05-16T05:07:32.138672Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ica = mne.preprocessing.ICA(n_components=0.95)\nica.fit(raw)","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:37.296178Z","iopub.execute_input":"2024-05-16T05:07:37.296528Z","iopub.status.idle":"2024-05-16T05:07:37.853260Z","shell.execute_reply.started":"2024-05-16T05:07:37.296498Z","shell.execute_reply":"2024-05-16T05:07:37.852440Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ica.plot_components() # зачем то плотит несколько раз","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:40.428473Z","iopub.execute_input":"2024-05-16T05:07:40.429168Z","iopub.status.idle":"2024-05-16T05:07:44.462540Z","shell.execute_reply.started":"2024-05-16T05:07:40.429136Z","shell.execute_reply":"2024-05-16T05:07:44.460918Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ica.plot_sources(raw) ","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:47.048543Z","iopub.execute_input":"2024-05-16T05:07:47.048971Z","iopub.status.idle":"2024-05-16T05:07:48.077401Z","shell.execute_reply.started":"2024-05-16T05:07:47.048934Z","shell.execute_reply":"2024-05-16T05:07:48.076259Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"raw_reconstructed = raw.copy()\nica.exclude = [4,'EKG']\nica.apply(raw_reconstructed)\nraw_reconstructed.plot()\n","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:52.129002Z","iopub.execute_input":"2024-05-16T05:07:52.129413Z","iopub.status.idle":"2024-05-16T05:07:53.356668Z","shell.execute_reply.started":"2024-05-16T05:07:52.129384Z","shell.execute_reply":"2024-05-16T05:07:53.355457Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ica.plot_components()","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:07:58.129142Z","iopub.execute_input":"2024-05-16T05:07:58.129514Z","iopub.status.idle":"2024-05-16T05:08:02.126878Z","shell.execute_reply.started":"2024-05-16T05:07:58.129486Z","shell.execute_reply":"2024-05-16T05:08:02.125466Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Хоть данный метод и увеличивает качество данных, проблема состоит в том, что мы должны руками перебиреть все файлы. Для нас не особо подходит. :(","metadata":{}},{"cell_type":"code","source":"bipolar = [\n    ['Fp1', 'F7'], ['F7', 'T3'], ['T3', 'T5'], ['T5', 'O1'],    # Left Temporal\n    ['Fp2', 'F8'], ['F8', 'T4'], ['T4', 'T6'], ['T6', 'O2'],    # Right Temporal\n    ['Fp1', 'F3'], ['F3', 'C3'], ['C3', 'P3'], ['P3', 'O1'],    # Left Parasagittal\n    ['Fp2', 'F4'], ['F4', 'C4'], ['C4', 'P4'], ['P4', 'O2'],    # Right Parasagittal\n    ['Fz', 'Cz'], ['Cz', 'Pz'],   # Central\n]\n\nanode, cathode = list(map(list,zip(*bipolar)))\n\nraw_bip_ref = mne.set_bipolar_reference(raw_filtered, anode=anode, cathode=cathode)\nraw_bip_ref.plot(start=20, duration=10)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:08:02.129395Z","iopub.execute_input":"2024-05-16T05:08:02.130561Z","iopub.status.idle":"2024-05-16T05:08:02.920227Z","shell.execute_reply.started":"2024-05-16T05:08:02.130523Z","shell.execute_reply":"2024-05-16T05:08:02.919169Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"raw_bip_ref.time_as_index(1)","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:08:06.477543Z","iopub.execute_input":"2024-05-16T05:08:06.478000Z","iopub.status.idle":"2024-05-16T05:08:06.486193Z","shell.execute_reply.started":"2024-05-16T05:08:06.477968Z","shell.execute_reply":"2024-05-16T05:08:06.484652Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# взято отсюда https://www.kaggle.com/code/kimbyungchun/preprocess-with-mne-for-human/notebook\n\nmne.set_config('MNE_BROWSE_RAW_SIZE','16,8')\n\ndef load_eeg(eeg_id, data_type:Literal['train','test']='train'):\n    BASE_PATH = '/kaggle/input/hms-harmful-brain-activity-classification'\n    eeg_df = pd.read_parquet(f'{BASE_PATH}/{data_type}_eegs/{eeg_id}.parquet')\n    return eeg_df\n\ndef get_filtered_bipolar(eeg_id, data_type:Literal['train','test']='train'):\n    bipolar = [\n        ['Fp1', 'F7'], ['F7', 'T3'], ['T3', 'T5'], ['T5', 'O1'],    # Left Temporal\n        ['Fp2', 'F8'], ['F8', 'T4'], ['T4', 'T6'], ['T6', 'O2'],    # Right Temporal\n        ['Fp1', 'F3'], ['F3', 'C3'], ['C3', 'P3'], ['P3', 'O1'],    # Left Parasagittal\n        ['Fp2', 'F4'], ['F4', 'C4'], ['C4', 'P4'], ['P4', 'O2'],    # Right Parasagittal\n        ['Fz', 'Cz'], ['Cz', 'Pz'],   # Central\n    ]\n    anode, cathode = list(map(list,zip(*bipolar)))\n\n    eeg_df = load_eeg(eeg_id, data_type)\n    \n    info = mne.create_info(\n        eeg_df.columns.to_list(),\n        ch_types=([\"eeg\"]*(len(eeg_df.columns)-1))+['ecg'],\n        sfreq=200\n    )\n    \n    info.set_montage(\"standard_1020\")\n    \n    raw = mne.io.RawArray(\n        eeg_df.to_numpy().T*1e-6,    # µV to V\n        info\n    ).filter(l_freq=1, h_freq=70,).notch_filter(60, picks='eeg')\n    \n    return mne.set_bipolar_reference(raw, anode=anode, cathode=cathode)\n\ndef bipolar(eeg_id, offset=0., duration=10.0, data_type:Literal['train','test']='train', plot=True):\n    bip_ref = get_filtered_bipolar(eeg_id, data_type)\n    \n    start = 25.0 + offset - duration/2\n    if plot:\n        bip_ref.plot(start=start, duration=duration)\n    \n    start_idx = bip_ref.time_as_index(start).item()\n    stop_idx = bip_ref.time_as_index(start + duration).item()\n    return bip_ref.get_data(start=start_idx, stop=stop_idx)","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:08:09.364581Z","iopub.execute_input":"2024-05-16T05:08:09.364953Z","iopub.status.idle":"2024-05-16T05:08:09.377749Z","shell.execute_reply.started":"2024-05-16T05:08:09.364930Z","shell.execute_reply":"2024-05-16T05:08:09.376506Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data = bipolar(1628180742, 0, 10)\ndata.shape","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:10:06.314503Z","iopub.execute_input":"2024-05-16T05:10:06.314906Z","iopub.status.idle":"2024-05-16T05:10:07.228623Z","shell.execute_reply.started":"2024-05-16T05:10:06.314882Z","shell.execute_reply":"2024-05-16T05:10:07.227439Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Хоть данные и становятся чище, мы все нарушаем изначальную структуру данных из-за чего быстрое применение и тестироване моделей может оказатся затруднительно. тем не менее она может быть использована для преобразования данных непосредственно перед передачей в модель","metadata":{}},{"cell_type":"markdown","source":"# linear denoise","metadata":{}},{"cell_type":"code","source":"EEG_SAMPLING_TIME = 50  #second\nEEG_SAMPLING_RATE = 200 #Hz\nEEG_DURATION = EEG_SAMPLING_RATE * EEG_SAMPLING_TIME\nN_INTERACTIVE = 100","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:13:21.267528Z","iopub.execute_input":"2024-05-16T05:13:21.267942Z","iopub.status.idle":"2024-05-16T05:13:21.274098Z","shell.execute_reply.started":"2024-05-16T05:13:21.267916Z","shell.execute_reply":"2024-05-16T05:13:21.272276Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_eeg(df, moving_avg=1):\n    fig, axs = plt.subplots(20, 1, figsize=(30, 15), sharex=True)\n    for i, ax in enumerate(axs):\n        ax.plot(df.iloc[:,i], color=\"black\")\n        for vline in df[df.iloc[:,i].isna()].index:\n            line_min = df.iloc[:,i].min()\n            line_max = df.iloc[:,i].max()\n            ax.vlines(vline, line_min, line_max, color='red')\n        ax.set_ylabel(df.columns[i], rotation=0)\n        ax.set_yticklabels([])\n        ax.set_yticks([])\n        ax.set_xticks([])\n        ax.spines[[\"top\", \"bottom\", \"left\", \"right\"]].set_visible(False)\n\ndef get_train_eeg(q:dict) -> pd.DataFrame:\n    parquet_df = pd.read_parquet(path/\"train_eegs\"/f\"train_eegs/{q['eeg_id']}.parquet\")\n    eeg_start_index = int(EEG_SAMPLING_RATE * q[\"eeg_label_offset_seconds\"])\n    return parquet_df.iloc[eeg_start_index:eeg_start_index+EEG_DURATION]\n\ndef denoise(x, wavelet='db8', level=1): # dmey for seizure patient denoise, db8 for healthy patient denoise\n    # paper: https://jart.icat.unam.mx/index.php/jart/article/view/339/336\n    def _maddest(d, axis=None):\n        return np.mean(np.absolute(d - np.mean(d, axis)), axis)\n    ret = {key:[] for key in x.columns}\n    for pos in x.columns:\n        coeff = pywt.wavedec(x[pos], wavelet, mode=\"per\")\n        sigma = (1/0.6745) * _maddest(coeff[-level])\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        ret[pos]=pywt.waverec(coeff, wavelet, mode='per')\n    return pd.DataFrame(ret)\n\ndef interpolate(raw_df):\n    df = raw_df.copy()\n    df = df.interpolate(\n        method='linear',\n        axis=0,\n        limit=1, # ref to 1 value\n        limit_direction=\"both\", # interpolate from pre and post values\n        limit_area='inside',\n    )\n    return df\n\ndef replace_outlier(series, bias=1.5, upper=0.95, lower=0.05):\n    lower_clip = series.quantile(lower)\n    upper_clip = series.quantile(upper)\n    iqr = upper_clip - lower_clip\n    \n    outlier_min = lower_clip - (iqr) * bias\n    outlier_max = upper_clip + (iqr) * bias\n\n    series = series.clip(outlier_min, outlier_max)\n    return series","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:13:25.274390Z","iopub.execute_input":"2024-05-16T05:13:25.274919Z","iopub.status.idle":"2024-05-16T05:13:25.290189Z","shell.execute_reply.started":"2024-05-16T05:13:25.274894Z","shell.execute_reply":"2024-05-16T05:13:25.288870Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"BASE_DIR = Path(\"/kaggle/input/hms-harmful-brain-activity-classification\")\nTRAIN_EEG_DIR = BASE_DIR/\"train_eegs\"\nWAVELET_DECOMPOSITION_TREE_LEVEL = 1\nWAVELET_TYPE = \"dmey\"\n\ntrain_df = pd.read_csv(BASE_DIR/\"train.csv\")\neeg_list = glob.glob(str(TRAIN_EEG_DIR/\"*.parquet\"))\nif True or os.environ.get('KAGGLE_KERNEL_RUN_TYPE','') == 'Interactive':\n    print(\"Running on development environment\")\n#     train_df = train_df.sample(N_INTERACTIVE)\n    eeg_list = eeg_list[:N_INTERACTIVE]\nelif os.environ.get('KAGGLE_KERNEL_RUN_TYPE','') == 'Batch':\n    print(\"Running on production environment\")","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:14:29.389488Z","iopub.execute_input":"2024-05-16T05:14:29.389907Z","iopub.status.idle":"2024-05-16T05:14:29.537803Z","shell.execute_reply.started":"2024-05-16T05:14:29.389880Z","shell.execute_reply":"2024-05-16T05:14:29.536590Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"","metadata":{}},{"cell_type":"code","source":"preprocessed_eeg_dir = \"train_eegs\"\nos.makedirs(preprocessed_eeg_dir, exist_ok=True)\n\nignore_parquet_num = 0\nignore_parquet_eeg_ids = list()\nfor p_eeg in eeg_list:\n    eeg_id = int(p_eeg.split(\"/\")[-1].split(\".\")[0])\n    eeg = pd.read_parquet(p_eeg)\n    eeg = interpolate(eeg)\n    if eeg.isna().sum().sum(): # discard EEG(still contail missing value)\n        ignore_parquet_num+=1\n        ignore_parquet_eeg_ids.append(eeg_id)\n        continue\n    eeg = denoise(eeg, wavelet=WAVELET_TYPE)\n    eeg.to_parquet(os.path.join(preprocessed_eeg_dir, f\"{eeg_id}.parquet\"))\nprint(f\"{ignore_parquet_num=}/{len(eeg_list)}\")","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:14:32.743306Z","iopub.execute_input":"2024-05-16T05:14:32.743774Z","iopub.status.idle":"2024-05-16T05:14:43.426811Z","shell.execute_reply.started":"2024-05-16T05:14:32.743742Z","shell.execute_reply":"2024-05-16T05:14:43.425388Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"preprocessed_traindf = train_df.loc[~train_df[\"eeg_id\"].isin(ignore_parquet_eeg_ids)].copy()\npreprocessed_traindf.to_csv(\"train.csv\")","metadata":{"execution":{"iopub.status.busy":"2024-05-16T05:15:08.651148Z","iopub.execute_input":"2024-05-16T05:15:08.651517Z","iopub.status.idle":"2024-05-16T05:15:09.188941Z","shell.execute_reply.started":"2024-05-16T05:15:08.651493Z","shell.execute_reply":"2024-05-16T05:15:09.187879Z"},"trusted":true},"execution_count":null,"outputs":[]}]}