Kaggle BCI Challenge eegnet not fixed.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Kaggle BCI Challenge @NER 2015\n",
    "\n",
    "**Author**:\n",
    "Yundong Wang"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "This project applied machine learning techniques to the P300 speller classification challenge at Kaggle BCI Challenge @ NER 2015. The P300 speller is a brain-computer interface paradigm that allows one to input text or commands to a computer via brain activity measured by electroencephalography (EEG).\n",
    "\n",
    "xDAWN “Covariance algorithm and Tangent Space algorithm” which uses the spatial correlations of the simultaneously recorded EEG as features for the error detection classification are implemented.\n",
    "\n",
    "We applied multiple machine learning techniques (SVM, Logistic Regression, Random Forest, etc) and deep learning model dedicated for EEG data ([EEGNET](https://arxiv.org/pdf/1611.08024.pdf)). Lastly, we built a StackNet model that achieved the 4th best AUC (0.803) result compared to the candidates pool.\n",
    "\n",
    "The presentation slides can be found [here.](https://docs.google.com/presentation/d/16Orsav1K4bPTFhPx3TsfYNxzgrvTMb2Gcf9RTC6WYmU/edit?ts=5c8aa323#slide=id.g5231037314_0_0)\n",
    "\n",
    "![video](images2/p300.gif)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Dataset\n",
    "Each session include 60 target stimulus, however, the last session of each subject contains 100 target stimulus. Which makes 340 target stimulus for each subject\n",
    "\n",
    "Data is collected at 200 Hz across 26 subjects (16 for training, 10 for testing). Each subject participated in 5 different sessions.\n",
    "\n",
    "0 or 1 for bad or good feedback, respectively. Bad feedback is when the selected item is different from the expected item. Good feedback is when the selected item is similar to the expected item.\n",
    "\n",
    "The Kaggle description and dataset can be downloaded [here.](https://www.kaggle.com/c/inria-bci-challenge)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 1: Setup\n",
    "\n",
    "- To run python code, please download data from https://www.kaggle.com/c/inria-bci-challenge/data and put *train* and *test* directories under *data*. **Python version: 3.6+**\n",
    "\n",
    "- To run *preprocess.py*, install **pyriemann** packages using pip.\n",
    "\n",
    "- To run *EEGNET.py*, install EEGModels at https://github.com/vlawhern/arl-eegmodels and required **TensorFlow** packages.\n",
    "\n",
    "- To run *StackNet.py*, install **pystacknet** at https://github.com/h2oai/pystacknet, install **LightGBM** and **XGBoost** using conda or pip.\n",
    "\n",
    "Before we can begin processing our data, we must first import some useful packages. <br>\n",
    "- **numpy** is used widely for processing numerical data, and supports matrix operations natively. We will be using NumPy arrays to store our data.\n",
    "- **scipy** contains useful functions to calculate filter coefficients and carry out filtering. We will be using SciPy filters to clean our EEG data.\n",
    "- **matplotlib** is used to create plots. We will use it to visually explore our filter and EEG data.\n",
    "- **sklearn** used for machine learning. We will use this to train our classifier and make predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np                                      # for dealing with data\n",
    "from scipy.signal import butter, sosfiltfilt, sosfreqz  # for filtering\n",
    "import matplotlib.pyplot as plt                         # for plotting\n",
    "from scipy import interp\n",
    "from sklearn.discriminant_analysis import LinearDiscriminantAnalysis\n",
    "from sklearn.model_selection import cross_val_score\n",
    "from sklearn.metrics import roc_curve, auc\n",
    "from sklearn.model_selection import StratifiedKFold\n",
    "import os\n",
    "from os import listdir\n",
    "from os.path import isfile, join, isdir\n",
    "import pickle"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 2: Pre-processing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 1: Design a filter**<br>\n",
    "\n",
    "- Filter.......: Butterworth, noncausal filtfilt\n",
    "- Order......: 5nd\n",
    "- High pass: 1 Hz\n",
    "- Low pass.: 40 Hz"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create filtering variables\n",
    "fs = 200.0     # 200 Hz sampling rate\n",
    "lowcut = 1.0   # 0.1 Hz is the lowest frequency we will pass\n",
    "highcut = 40.0 # 30  Hz is the highest frequency we will pass."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def butter_bandpass_filter(raw_data, fs, lowcut = 1.0, highcut =40.0, order = 5):\n",
    "    '''\n",
    "    The filter I want to apply to my raw eeg data.\n",
    "    :raw_data (nparray): data you want to process\n",
    "    :fs (float): sampling rate\n",
    "    :lowcut (float, optional): lowest frequency we will pass\n",
    "    :highcut (float, optional): highest frequency we will pass\n",
    "    :order (int, optional): order of filter\n",
    "    '''\n",
    "    nyq = 0.5 * fs\n",
    "    low = lowcut / nyq\n",
    "    high = highcut / nyq\n",
    "    sos = butter(order, [low, high], analog = False, btype = 'band', output = 'sos')\n",
    "    filted_data = sosfiltfilt(sos, raw_data)\n",
    "    return filted_data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 2: Determine our epoch, baseline, and window of focus**<br>\n",
    "\n",
    "- Epoch...: 0 ~ 700ms\n",
    "- Baseline: 0 ~ 100ms \n",
    "\n",
    "The following is a visual explanation of the baseline, ERP window selection, and sub-sampling routine (on a single channel, using the figure provided by Alessandro D'Amico):\n",
    "<img src=\"images2/img1.png\">\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "epoch_s = 0      # epoch starting time relative to stmulus in miliseconds\n",
    "epoch_e = 700    # epoch ending time relative to stmulus in miliseconds\n",
    "bl_s = 0         # baseline starting time relative to stmulus in miliseconds\n",
    "bl_e = 100       # baseline ending time relative to stmulus in miliseconds\n",
    "\n",
    "\n",
    "# number of mark per epoch\n",
    "epoch_len = int((abs(epoch_s) + abs(epoch_e)) * (fs / 1000)) \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 3: Understand our data**<br>\n",
    "\n",
    "- 16 participants in the training dataset.\n",
    "- 10 participants in the testing dataset.\n",
    "- 5 trials per participant\n",
    "- 60 stimulus per trial in the first 4 trials, and 100 stimulus per trial in the last trial.\n",
    "- Thus, 340 stimulus per trial.\n",
    "- 56 channel EEG data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_subj_num = 16\n",
    "test_subj_num = 10\n",
    "stimulus_per_subj = 340\n",
    "trial_per_subj = 5\n",
    "\n",
    "channels = ['Fp1', 'Fp2', 'AF7', 'AF3', 'AF4', 'AF8', 'F7', 'F5', 'F3', 'F1',\n",
    "    'Fz', 'F2', 'F4', 'F6', 'F8', 'FT7', 'FC5', 'FC3', 'FC1', 'FCz',\n",
    "    'FC2', 'FC4', 'FC6', 'FT8', 'T7', 'C5', 'C3', 'C1', 'Cz', 'C2',\n",
    "    'C4', 'C6', 'T8', 'TP7', 'CP5', 'CP3', 'CP1', 'CPz', 'CP2', 'CP4',\n",
    "    'CP6', 'TP8', 'P7', 'P5', 'P3', 'P1', 'Pz', 'P2', 'P4', 'P6', 'P8',\n",
    "    'PO7', 'POz', 'P08', 'O1', 'O2']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "train_labels = pd.read_csv(\"./data/TrainLabels.csv\")\n",
    "sample_train_data = pd.read_csv(\"./data/train/Data_S02_Sess01.csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>IdFeedBack</th>\n",
       "      <th>Prediction</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>S02_Sess01_FB001</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>S02_Sess01_FB002</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>S02_Sess01_FB003</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>S02_Sess01_FB004</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>S02_Sess01_FB005</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5435</th>\n",
       "      <td>S26_Sess05_FB096</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5436</th>\n",
       "      <td>S26_Sess05_FB097</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5437</th>\n",
       "      <td>S26_Sess05_FB098</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5438</th>\n",
       "      <td>S26_Sess05_FB099</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5439</th>\n",
       "      <td>S26_Sess05_FB100</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5440 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "            IdFeedBack  Prediction\n",
       "0     S02_Sess01_FB001           1\n",
       "1     S02_Sess01_FB002           1\n",
       "2     S02_Sess01_FB003           0\n",
       "3     S02_Sess01_FB004           0\n",
       "4     S02_Sess01_FB005           1\n",
       "...                ...         ...\n",
       "5435  S26_Sess05_FB096           1\n",
       "5436  S26_Sess05_FB097           0\n",
       "5437  S26_Sess05_FB098           0\n",
       "5438  S26_Sess05_FB099           0\n",
       "5439  S26_Sess05_FB100           1\n",
       "\n",
       "[5440 rows x 2 columns]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "train_labels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Time</th>\n",
       "      <th>Fp1</th>\n",
       "      <th>Fp2</th>\n",
       "      <th>AF7</th>\n",
       "      <th>AF3</th>\n",
       "      <th>AF4</th>\n",
       "      <th>AF8</th>\n",
       "      <th>F7</th>\n",
       "      <th>F5</th>\n",
       "      <th>F3</th>\n",
       "      <th>...</th>\n",
       "      <th>P4</th>\n",
       "      <th>P6</th>\n",
       "      <th>P8</th>\n",
       "      <th>PO7</th>\n",
       "      <th>POz</th>\n",
       "      <th>P08</th>\n",
       "      <th>O1</th>\n",
       "      <th>O2</th>\n",
       "      <th>EOG</th>\n",
       "      <th>FeedBackEvent</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.000</td>\n",
       "      <td>509.630351</td>\n",
       "      <td>609.688545</td>\n",
       "      <td>522.032153</td>\n",
       "      <td>454.102593</td>\n",
       "      <td>309.486489</td>\n",
       "      <td>524.804335</td>\n",
       "      <td>484.333383</td>\n",
       "      <td>587.340998</td>\n",
       "      <td>469.651564</td>\n",
       "      <td>...</td>\n",
       "      <td>295.834834</td>\n",
       "      <td>485.530586</td>\n",
       "      <td>411.161892</td>\n",
       "      <td>448.726097</td>\n",
       "      <td>99.151239</td>\n",
       "      <td>613.529369</td>\n",
       "      <td>428.286969</td>\n",
       "      <td>378.933590</td>\n",
       "      <td>-906.668876</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0.005</td>\n",
       "      <td>864.521388</td>\n",
       "      <td>1038.483761</td>\n",
       "      <td>917.023604</td>\n",
       "      <td>775.970207</td>\n",
       "      <td>576.123392</td>\n",
       "      <td>934.962455</td>\n",
       "      <td>833.427809</td>\n",
       "      <td>996.070640</td>\n",
       "      <td>823.715523</td>\n",
       "      <td>...</td>\n",
       "      <td>528.713692</td>\n",
       "      <td>894.541770</td>\n",
       "      <td>735.827300</td>\n",
       "      <td>780.487458</td>\n",
       "      <td>239.407648</td>\n",
       "      <td>1067.349342</td>\n",
       "      <td>782.319229</td>\n",
       "      <td>711.286079</td>\n",
       "      <td>-1484.107119</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.010</td>\n",
       "      <td>754.576625</td>\n",
       "      <td>904.117869</td>\n",
       "      <td>803.235608</td>\n",
       "      <td>678.039830</td>\n",
       "      <td>489.573235</td>\n",
       "      <td>803.064230</td>\n",
       "      <td>706.392172</td>\n",
       "      <td>868.996293</td>\n",
       "      <td>717.565698</td>\n",
       "      <td>...</td>\n",
       "      <td>459.066586</td>\n",
       "      <td>870.217213</td>\n",
       "      <td>691.971456</td>\n",
       "      <td>689.656321</td>\n",
       "      <td>255.876323</td>\n",
       "      <td>933.617965</td>\n",
       "      <td>679.366071</td>\n",
       "      <td>699.014004</td>\n",
       "      <td>-1313.435186</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0.015</td>\n",
       "      <td>801.593761</td>\n",
       "      <td>957.870867</td>\n",
       "      <td>834.026569</td>\n",
       "      <td>716.254450</td>\n",
       "      <td>519.006819</td>\n",
       "      <td>860.058239</td>\n",
       "      <td>769.518284</td>\n",
       "      <td>921.777695</td>\n",
       "      <td>756.664761</td>\n",
       "      <td>...</td>\n",
       "      <td>487.899954</td>\n",
       "      <td>905.564172</td>\n",
       "      <td>723.595463</td>\n",
       "      <td>723.944898</td>\n",
       "      <td>247.990115</td>\n",
       "      <td>992.180045</td>\n",
       "      <td>713.405603</td>\n",
       "      <td>714.734949</td>\n",
       "      <td>-1391.966973</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>0.020</td>\n",
       "      <td>769.054899</td>\n",
       "      <td>920.060832</td>\n",
       "      <td>811.596932</td>\n",
       "      <td>692.750986</td>\n",
       "      <td>494.315305</td>\n",
       "      <td>817.206918</td>\n",
       "      <td>727.500302</td>\n",
       "      <td>892.004312</td>\n",
       "      <td>729.239924</td>\n",
       "      <td>...</td>\n",
       "      <td>464.032652</td>\n",
       "      <td>769.258159</td>\n",
       "      <td>642.962959</td>\n",
       "      <td>687.852792</td>\n",
       "      <td>179.835749</td>\n",
       "      <td>949.832471</td>\n",
       "      <td>680.906506</td>\n",
       "      <td>607.057329</td>\n",
       "      <td>-1347.494166</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>131996</th>\n",
       "      <td>659.980</td>\n",
       "      <td>326.493877</td>\n",
       "      <td>148.619853</td>\n",
       "      <td>367.120885</td>\n",
       "      <td>329.607132</td>\n",
       "      <td>568.760051</td>\n",
       "      <td>228.081535</td>\n",
       "      <td>386.342230</td>\n",
       "      <td>193.539354</td>\n",
       "      <td>263.262863</td>\n",
       "      <td>...</td>\n",
       "      <td>497.935681</td>\n",
       "      <td>136.897186</td>\n",
       "      <td>194.816484</td>\n",
       "      <td>235.516599</td>\n",
       "      <td>698.470304</td>\n",
       "      <td>115.904670</td>\n",
       "      <td>346.770258</td>\n",
       "      <td>-22.426859</td>\n",
       "      <td>-1179.103975</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>131997</th>\n",
       "      <td>659.985</td>\n",
       "      <td>299.398029</td>\n",
       "      <td>122.348935</td>\n",
       "      <td>317.862613</td>\n",
       "      <td>299.081994</td>\n",
       "      <td>540.210149</td>\n",
       "      <td>183.619873</td>\n",
       "      <td>354.533814</td>\n",
       "      <td>148.285289</td>\n",
       "      <td>208.522216</td>\n",
       "      <td>...</td>\n",
       "      <td>483.786225</td>\n",
       "      <td>87.994431</td>\n",
       "      <td>158.168639</td>\n",
       "      <td>220.717339</td>\n",
       "      <td>681.475238</td>\n",
       "      <td>69.722673</td>\n",
       "      <td>304.952629</td>\n",
       "      <td>-71.490262</td>\n",
       "      <td>-1245.844338</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>131998</th>\n",
       "      <td>659.990</td>\n",
       "      <td>313.987327</td>\n",
       "      <td>143.185940</td>\n",
       "      <td>352.162704</td>\n",
       "      <td>321.971403</td>\n",
       "      <td>551.396419</td>\n",
       "      <td>219.448645</td>\n",
       "      <td>345.700353</td>\n",
       "      <td>196.750843</td>\n",
       "      <td>253.721728</td>\n",
       "      <td>...</td>\n",
       "      <td>482.381296</td>\n",
       "      <td>130.181271</td>\n",
       "      <td>228.417928</td>\n",
       "      <td>229.943069</td>\n",
       "      <td>718.838360</td>\n",
       "      <td>110.441383</td>\n",
       "      <td>325.000921</td>\n",
       "      <td>62.368775</td>\n",
       "      <td>-1150.955909</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>131999</th>\n",
       "      <td>659.995</td>\n",
       "      <td>308.194641</td>\n",
       "      <td>121.836127</td>\n",
       "      <td>323.251844</td>\n",
       "      <td>313.754887</td>\n",
       "      <td>575.047001</td>\n",
       "      <td>202.242791</td>\n",
       "      <td>358.987420</td>\n",
       "      <td>165.629556</td>\n",
       "      <td>231.005230</td>\n",
       "      <td>...</td>\n",
       "      <td>511.892613</td>\n",
       "      <td>103.343795</td>\n",
       "      <td>229.027509</td>\n",
       "      <td>233.729526</td>\n",
       "      <td>782.742133</td>\n",
       "      <td>85.836597</td>\n",
       "      <td>327.645835</td>\n",
       "      <td>38.482607</td>\n",
       "      <td>-1304.594442</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>132000</th>\n",
       "      <td>660.000</td>\n",
       "      <td>236.894562</td>\n",
       "      <td>118.063306</td>\n",
       "      <td>273.854465</td>\n",
       "      <td>243.895105</td>\n",
       "      <td>416.279695</td>\n",
       "      <td>188.173085</td>\n",
       "      <td>270.812314</td>\n",
       "      <td>167.005916</td>\n",
       "      <td>215.650455</td>\n",
       "      <td>...</td>\n",
       "      <td>353.931780</td>\n",
       "      <td>122.418483</td>\n",
       "      <td>154.473702</td>\n",
       "      <td>177.953948</td>\n",
       "      <td>500.095611</td>\n",
       "      <td>109.008917</td>\n",
       "      <td>269.794326</td>\n",
       "      <td>18.927431</td>\n",
       "      <td>-805.816250</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>132001 rows × 59 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "           Time         Fp1          Fp2         AF7         AF3         AF4  \\\n",
       "0         0.000  509.630351   609.688545  522.032153  454.102593  309.486489   \n",
       "1         0.005  864.521388  1038.483761  917.023604  775.970207  576.123392   \n",
       "2         0.010  754.576625   904.117869  803.235608  678.039830  489.573235   \n",
       "3         0.015  801.593761   957.870867  834.026569  716.254450  519.006819   \n",
       "4         0.020  769.054899   920.060832  811.596932  692.750986  494.315305   \n",
       "...         ...         ...          ...         ...         ...         ...   \n",
       "131996  659.980  326.493877   148.619853  367.120885  329.607132  568.760051   \n",
       "131997  659.985  299.398029   122.348935  317.862613  299.081994  540.210149   \n",
       "131998  659.990  313.987327   143.185940  352.162704  321.971403  551.396419   \n",
       "131999  659.995  308.194641   121.836127  323.251844  313.754887  575.047001   \n",
       "132000  660.000  236.894562   118.063306  273.854465  243.895105  416.279695   \n",
       "\n",
       "               AF8          F7          F5          F3  ...          P4  \\\n",
       "0       524.804335  484.333383  587.340998  469.651564  ...  295.834834   \n",
       "1       934.962455  833.427809  996.070640  823.715523  ...  528.713692   \n",
       "2       803.064230  706.392172  868.996293  717.565698  ...  459.066586   \n",
       "3       860.058239  769.518284  921.777695  756.664761  ...  487.899954   \n",
       "4       817.206918  727.500302  892.004312  729.239924  ...  464.032652   \n",
       "...            ...         ...         ...         ...  ...         ...   \n",
       "131996  228.081535  386.342230  193.539354  263.262863  ...  497.935681   \n",
       "131997  183.619873  354.533814  148.285289  208.522216  ...  483.786225   \n",
       "131998  219.448645  345.700353  196.750843  253.721728  ...  482.381296   \n",
       "131999  202.242791  358.987420  165.629556  231.005230  ...  511.892613   \n",
       "132000  188.173085  270.812314  167.005916  215.650455  ...  353.931780   \n",
       "\n",
       "                P6          P8         PO7         POz          P08  \\\n",
       "0       485.530586  411.161892  448.726097   99.151239   613.529369   \n",
       "1       894.541770  735.827300  780.487458  239.407648  1067.349342   \n",
       "2       870.217213  691.971456  689.656321  255.876323   933.617965   \n",
       "3       905.564172  723.595463  723.944898  247.990115   992.180045   \n",
       "4       769.258159  642.962959  687.852792  179.835749   949.832471   \n",
       "...            ...         ...         ...         ...          ...   \n",
       "131996  136.897186  194.816484  235.516599  698.470304   115.904670   \n",
       "131997   87.994431  158.168639  220.717339  681.475238    69.722673   \n",
       "131998  130.181271  228.417928  229.943069  718.838360   110.441383   \n",
       "131999  103.343795  229.027509  233.729526  782.742133    85.836597   \n",
       "132000  122.418483  154.473702  177.953948  500.095611   109.008917   \n",
       "\n",
       "                O1          O2          EOG  FeedBackEvent  \n",
       "0       428.286969  378.933590  -906.668876              0  \n",
       "1       782.319229  711.286079 -1484.107119              0  \n",
       "2       679.366071  699.014004 -1313.435186              0  \n",
       "3       713.405603  714.734949 -1391.966973              0  \n",
       "4       680.906506  607.057329 -1347.494166              0  \n",
       "...            ...         ...          ...            ...  \n",
       "131996  346.770258  -22.426859 -1179.103975              0  \n",
       "131997  304.952629  -71.490262 -1245.844338              0  \n",
       "131998  325.000921   62.368775 -1150.955909              0  \n",
       "131999  327.645835   38.482607 -1304.594442              0  \n",
       "132000  269.794326   18.927431  -805.816250              0  \n",
       "\n",
       "[132001 rows x 59 columns]"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sample_train_data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 4: Create placeholder for loading our data**<br>\n",
    "- There are 16 training participants, 5 trials per participant.\n",
    "- There are 10 test participants, 5 trials per participant.\n",
    "- There are 340 trials per participant, 56 channels for each trial, and our epoch length is 260.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(16, 5) (10, 5)\n",
      "(0, 340, 56, 140) (0, 340, 56, 140)\n"
     ]
    }
   ],
   "source": [
    "train_list_arr = np.array(sorted(listdir('./data/train')))\n",
    "train_list_np = np.reshape(\n",
    "    train_list_arr, (train_subj_num, trial_per_subj))\n",
    "test_list_arr = np.array(sorted(listdir('./data/test')))\n",
    "test_list_np = np.reshape(\n",
    "    test_list_arr, (test_subj_num, trial_per_subj))\n",
    "print(train_list_np.shape, test_list_np.shape)\n",
    "\n",
    "train_data_list = np.empty(\n",
    "    (0, stimulus_per_subj, len(channels), epoch_len), float)\n",
    "test_data_list = np.empty(\n",
    "    (0, stimulus_per_subj, len(channels), epoch_len), float)\n",
    "print(train_data_list.shape, test_data_list.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 5: Define our epoch function**<br>\n",
    "We apply this function to each EEG trial."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_epoch(file_path, channels, fs, eeg_filter, stimulus_times=None, baseline=True,  epoch_s=0, epoch_e=700, bl_s=0, bl_e=100):\n",
    "    \"\"\"\n",
    "    :description: Generating epoch given csv file. Make sure the csv file layout meets the requirement.\n",
    "        It should contain 'Time' column that represents timepoints, and the time should start from 0.\n",
    "        If your csv file does not have FeedBackEvent indicating the stimulus, you must pass stumulus_times.\n",
    "        Here we used a butter bandpass filter, but you can change to your favorite one.\n",
    "\n",
    "    :file_path (String): path to your csv file\n",
    "    :channels ([String]): array of channels to epoch\n",
    "    :fs (float): sampling rate\n",
    "    :eeg_filter (function): the filter you want to apply to raw eeg data\n",
    "    :stimulus_times ([float], optional): The time points that stimulus occur\n",
    "    :baseline (boolean, optional): whether you want to apply baseline correction after epoching\n",
    "    :epoch_s (int, optional): epoch starting time relative to stmulus in miliseconds\n",
    "    :epoch_e (int, optional): epoch ending time relative to stmulus in miliseconds\n",
    "    :bl_s (int, optional): baseline starting time relative to stmulus in miliseconds\n",
    "    :bl_e (int, optional): baseline ending time relative to stmulus in miliseconds\n",
    "\n",
    "    :rtype (3d-nparray): epoched data with dimension (stimulus_per_subj, number_of_channels, number_of_time_points)\n",
    "    \"\"\"\n",
    "    # read dataand data selection\n",
    "    train_data = pd.read_csv(file_path)\n",
    "\n",
    "    train_data.loc[:, 'Time'] = train_data.loc[:, 'Time']*1000\n",
    "    raw_eeg = train_data[channels].values.T\n",
    "\n",
    "    time_df = train_data['Time'].values\n",
    "    train_data['index'] = train_data.index.values\n",
    "    if stimulus_times is None:\n",
    "        mark_indices = np.asarray(\n",
    "            train_data[train_data['FeedBackEvent'] == 1].index).flatten()\n",
    "    else:\n",
    "        mark_indices = np.round(np.asarray(\n",
    "            stimulus_times).flatten() * fs).astype(int)\n",
    "\n",
    "    # Define the bounds of our epoch as well as our baseline\n",
    "    # index in epoch_df where our baseline begins\n",
    "    b_s = int((abs(epoch_s) + bl_s) * (fs / 1000))\n",
    "    # index in epoch_df where our baseline ends\n",
    "    b_e = int((abs(epoch_s) + bl_e) * (fs / 1000))\n",
    "    # Let's calculate the length our epoch with our given sampling rate\n",
    "    epoch_len = int((abs(epoch_s) + abs(epoch_e)) * (fs / 1000))\n",
    "\n",
    "    # Let's define some helpful variables to make our extraction easier\n",
    "    # effectively the number of indices before marker we want\n",
    "    e_s = int((epoch_s * (fs / 1000)))\n",
    "    # effectively the number of indices after marker we want\n",
    "    e_e = int((epoch_e * (fs / 1000)))\n",
    "\n",
    "    # Epoch the data\n",
    "    final_epoch = np.empty((mark_indices.shape[0], epoch_len, 0), float)\n",
    "    for channel in channels:\n",
    "        epoch = np.zeros(shape=(int(mark_indices.shape[0]), epoch_len))\n",
    "        raw_eeg = train_data[channel].values\n",
    "\n",
    "        ################# You may want to apply your own filter ################\n",
    "        clean_eeg = eeg_filter(raw_eeg, fs, lowcut, highcut, 5)\n",
    "        ########################################################################\n",
    "\n",
    "        for i, mark_idx in enumerate(mark_indices):\n",
    "            # grab the appropriate samples around the stimulus onset\n",
    "            epoch[i, :] = clean_eeg[mark_idx + e_s: mark_idx + e_e]\n",
    "\n",
    "        # Baseline correction\n",
    "        if baseline:\n",
    "            for i in range(0, int(epoch.shape[0])):\n",
    "                epoch[i, :] = epoch[i, :] - np.mean(epoch[i, b_s:b_e])\n",
    "\n",
    "        # stack epoch of each channel\n",
    "        final_epoch = np.dstack((final_epoch, epoch))\n",
    "    final_epoch = np.swapaxes(final_epoch, 1, 2)\n",
    "    return final_epoch"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 6: Epoch all data using epoch function we wrote**<br>\n",
    "Now we have generated trainable data that can be fed into models. Skip this if we have already saved them in the first run."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "if not isfile(\"./data/train_data.npy\"):\n",
    "    for training_participant_id in range(train_subj_num):\n",
    "        subject_dir_list = train_list_np[training_participant_id]\n",
    "        subject_epoch = np.empty((0, len(channels), epoch_len), float)\n",
    "        for trial_id in range(trial_per_subj):\n",
    "            subject_dir = subject_dir_list[trial_id]\n",
    "            data = generate_epoch('./data/train/'+subject_dir, channels, fs,\n",
    "                butter_bandpass_filter, epoch_s = epoch_s, epoch_e = epoch_e, bl_s = bl_s, bl_e = bl_e)\n",
    "            subject_epoch = np.vstack((subject_epoch, data))\n",
    "        subject_epoch = np.reshape(\n",
    "            subject_epoch, (1, stimulus_per_subj, len(channels), epoch_len))\n",
    "        train_data_list = np.vstack((train_data_list, subject_epoch))\n",
    "\n",
    "    print('Epoched training data shape: ' + str(train_data_list.shape))\n",
    "\n",
    "    for testing_participant_id in range(test_subj_num):\n",
    "        subject_dir_list = test_list_np[testing_participant_id]\n",
    "        subject_epoch = np.empty((0, len(channels), epoch_len), float)\n",
    "        for trial_id in range(trial_per_subj):\n",
    "            subject_dir = subject_dir_list[trial_id]\n",
    "            data = generate_epoch('./data/test/'+subject_dir, channels, fs,\n",
    "                butter_bandpass_filter, epoch_s = epoch_s, epoch_e = epoch_e, bl_s = bl_s, bl_e = bl_e)\n",
    "            subject_epoch = np.vstack((subject_epoch, data))\n",
    "        subject_epoch = np.reshape(\n",
    "            subject_epoch, (1, stimulus_per_subj, len(channels), epoch_len))\n",
    "        test_data_list = np.vstack((test_data_list, subject_epoch))\n",
    "\n",
    "    print('Epoched testing data shape: ' + str(test_data_list.shape))\n",
    "\n",
    "    np.save('./data/train_data.npy', train_data_list)\n",
    "    np.save('./data/test_data.npy', test_data_list)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoched training data shape: (16, 340, 56, 140)\n",
      "Epoched testing data shape: (10, 340, 56, 140)\n"
     ]
    }
   ],
   "source": [
    "train_data_list = np.load('./data/train_data.npy')\n",
    "test_data_list = np.load('./data/test_data.npy')\n",
    "print('Epoched training data shape: ' + str(train_data_list.shape))\n",
    "print('Epoched testing data shape: ' + str(test_data_list.shape))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Samples of Epoched EEG data\n",
    "<table>\n",
    "    <tbody>\n",
    "        <tr>\n",
    "            <td>\n",
    "                <img src=\"./images2/Subject_2_Fz.png\"/>\n",
    "            </td>\n",
    "            <td>\n",
    "                <img src=\"./images2/Subject_2_Pz.png\"/>\n",
    "            </td>\n",
    "        </tr>\n",
    "        <tr>\n",
    "            <td>\n",
    "                <img src=\"./images2/Fz.png\"/>\n",
    "            </td>\n",
    "            <td>\n",
    "                <img src=\"./images2/Pz.png\"/>\n",
    "            </td>\n",
    "        </tr>\n",
    "    </tbody>\n",
    "</table>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 7: Apply xDawn Covariance algorithm and Tangent Space algorithm**<br>\n",
    "After epoching our data, our data is in (subject, trial, channel, stimulus). Last two dimensions correspond to **spatial** and **temporal** dimensions.\n",
    "\n",
    "Therefore, further feature exctraction and dimensionality reduction is necessary for standard machine learning algorithm to work."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pyriemann.estimation import XdawnCovariances\n",
    "from pyriemann.tangentspace import TangentSpace"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "Y_train = pd.read_csv('data/TrainLabels.csv')['Prediction'].values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "if not isfile('./data/X_train.npy'):\n",
    "    XC = XdawnCovariances(nfilter=5) # our transformer\n",
    "    TS = TangentSpace(metric='riemann')\n",
    "    # collapse first two dimension\n",
    "    train_data = np.reshape(train_data_list, \n",
    "    (train_subj_num * stimulus_per_subj, len(channels), epoch_len))\n",
    "    test_data = np.reshape(test_data_list, \n",
    "    (test_subj_num * stimulus_per_subj, len(channels), epoch_len))\n",
    "\n",
    "    # transform our data\n",
    "    X_train = XC.fit_transform(train_data, Y_train)\n",
    "    X_train = TS.fit_transform(X_train)\n",
    "    X_test = XC.transform(test_data)\n",
    "    X_test = TS.transform(X_test)\n",
    "\n",
    "    # save to local\n",
    "    np.save('./data/X_train', X_train)\n",
    "    np.save('./data/X_test', X_test)\n",
    "    np.save('./data/Y_train', Y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 4: Cross Validation and Classification with Machine Learning "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 1: Load our data**<br>\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load our data\n",
    "X_train = np.load('./data/X_train.npy')\n",
    "X_test = np.load('./data/X_test.npy')\n",
    "Y_train = np.load('./data/Y_train.npy')\n",
    "\n",
    "# Only used for scoring\n",
    "Y_test = np.reshape(pd.read_csv('./data/true_labels.csv', header=None).values, 3400)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Transformed training data shape: (5440, 210)\n",
      "Training label shape: (5440,)\n",
      "Transformed testing data shape: (3400, 210)\n",
      "Testing label shape: (3400,)\n"
     ]
    }
   ],
   "source": [
    "print('Transformed training data shape: ' + str(X_train.shape))\n",
    "print('Training label shape: ' + str(Y_train.shape))\n",
    "print('Transformed testing data shape: ' + str(X_test.shape))\n",
    "print('Testing label shape: ' + str(Y_test.shape))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 2a: Use RandomForest for prediction**<br>\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.ensemble import RandomForestClassifier\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "\n",
    "# Create the parameter grid based on the results of random search \n",
    "param_grid = {\n",
    "    'bootstrap': [True],\n",
    "    'max_depth': [80, 90, 100, 110],\n",
    "    'max_features': [2, 3],\n",
    "    'n_estimators': [100, 200, 300, 1000]\n",
    "}\n",
    "\n",
    "# Fit the grid search to the data\n",
    "if not isfile('random_forest_gs'): \n",
    "    rf = RandomForestClassifier()\n",
    "    # Instantiate the grid search model\n",
    "    grid_search = GridSearchCV(estimator = rf, param_grid = param_grid, \n",
    "               cv = 4, n_jobs = -1)\n",
    "    grid_search.fit(X_train, Y_train)\n",
    "    pickle.dump(grid_search, open('random_forest_gs', 'wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'bootstrap': True, 'max_depth': 80, 'max_features': 3, 'n_estimators': 200}\n"
     ]
    }
   ],
   "source": [
    "grid_search = pickle.load(open('random_forest_gs', 'rb'))\n",
    "print(grid_search.best_params_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Step 3a: Evaluate the performance of RandomForest**<br>\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEWCAYAAAB42tAoAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3hU1dbA4d9KSAESWkINJZTQO6EpIipcFAURCyiK9SIWUMGunxUVsV4UxY4iVVRAQQQUQQSk9xogQEJJKGmE9P39sYcYYhIGyGQyyXqfZ57knDlzZs1MctbsLsYYlFJKqfx4uTsApZRSxZsmCqWUUgXSRKGUUqpAmiiUUkoVSBOFUkqpAmmiUEopVSBNFKrQiUh1EVkqIoki8k4xiMeISCN3x+EuIuInIttEpMYFPt7pz1NEeohIVAH3TxSR0RcShzNE5AcRudpV5y+tNFF4OBGJFJHTIpIkIkcc/4gBuY65RER+d/yjx4vITyLSPNcxFUTkfRE54DhXhGM7+ALCGgocAyoYY0blEfNEEUlzPM8JEVkoIk0v4HncKtfrOHMbWMQxRIpIz3McNhRYaow5kuuxviKyo6ALe47H5/t5FiYRqSkiX4jIYcff6w4ReVlEyjt+vyePxzwiImscm2OA11wZY2mkiaJk6GuMCQDaAu2AZ87cISJdgQXAbKAWUB/YCPwlIg0cx/gCvwEtgKuBCsAlwHGg0wXEUw/YZgoezTnWEXMIEA18cQHPUxyMNcYE5LhNP98TiIi3KwLL4X5gUh77nwBinHi8M5/nRRORKsAKoCzQ1RgTCPQCKgENga+BIXk89A7HfRhjVgEVRCTclbGWOsYYvXnwDYgEeubYHgvMzbH9J/BRHo/7BfjG8ft9wFEg4Dye9xJgNRDv+HmJY/9EIB1IA5JyxpbjsROB0Tm2+wCncmw3BH7HJqpjwGSgUq7X/DiwyfH80wH/HPc/ARwGDgH3AAZo5LivIvANEAvsB54HvBz33QX8BbwHxAF7Ha/zLuAg9qJ6Z36vI9drbAb84TjPVqBfrsd9DMwDTgE9AT/gbeCA47OYAJR1HB8M/Ow41wnHZ+qFvfhnAacd7/WTecRR13F/mVz76wPbgWuAqAI+5399no5Y33e8v4ccv/s5ju+R83zYLy7rgETH5zStgPdsNLD5zOeRx/21gQygXq73OQ0IzrHvM+BFd/9vlqSblihKEBGpjf3Hj3Bsl8Ne6L7L4/AZ2G9rYP/55xtjkpx8nirAXGAcEAS8C8wVkSBjzF3YC/uZb9qLznGu8sCtZ2I+sxt4A1sCagbUAV7K9dBbsKWf+kBr7MUcR/30447XFuZ4bTl9gE0WDYDLsd9Q785xf2dsAgoCpmAvbB2BRsDtwIe5q/byeE0+wE/Yklw1YDgwWUSa5DjsNmwVSSCwDHgTaIwtFTbClrRecBw7CogCqgLVgWcBY4y5A5tY+jre67F5hNMK2GuMycjjfXgWm0Tylc/n+RzQxRFrG2yp8/k83gdfYBY2oVXB/h3eWMDT9QR+MMZk5RNLFLAYW4I4YwgwzxhzLMe+7Y64VCHRRFEyzBKRRP751vuiY38V7Gd8OI/HHMZ+UwV7UczrmPxcC+w2xkwyxmQYY6YCO4C+53GOx0UkDvtNsxs5/vmNMRHGmIXGmFRjTCw2EV2e6/HjjDGHjDEnsBflto79twBfGWO2GGNOkSPBOKp4BgLPGGMSjTGRwDucfeHZZ4z5yhiTif0GXAd4xRHLAuy315wN44+LSJzjduZi1QUIAMYYY9KMMb9jSwS35njcbGPMX46LYirwX+AxY8wJY0wi8DowyHFsOlAT+0063RjzpzHG2WqgStj3OJuI3IAtYfzo5DlyG4x9T2Icn8/LnP0entEF8AHed8Q9E1v6zI8zf4dfn3kuEfFyxPJ1rmMSsa9bFRJNFCVDf2Prc3sATfknAZzEVk3UzOMxNbHVOmCrePI6Jj+1sNU2Oe3Hfgt21tvGmEpAKPZbbfa3bRGpJiLTRCRaRBKAb/nnNZ2Rs2E2GXthPhPbwVxxnREM+Obalzvuozl+Pw1gjMm9L2eJ4m1jTCXH7UyMtYCDub4Z536enDFWBcoBa88kHWC+Yz/AW9gS1wIR2SsiT+O8k9hSC5BdghuLLeX8i4hMyNEw/2w+58z9+e937MvruOhcSS33301Ozvwd/gDUFJEu2L/3ctjSbU6B2Go6VUg0UZQgxpgl2Drltx3bp7CNgzfncfgt2AZsgEVAb8dFxBmHsA2cOdXFNkqfF2PMAeAR4H8iUtax+w1su0JrY0wFbJWPOHnKw9hSQM64zjiG/XZeL9f95x33ORwC6ji+8eb3PDkvnsewCahFjqRT0djGfhyln1HGmAbYUttIEbkqj/PkZRPQQETKOLbDsMn5TxE5wj8X3iMiEmqMGWb+aZh/vYDXl/s9PJTHcYeBEBGRXMfmZxFwQ6737SzGmGRgJrbK6Q5gmjEmLddhzbAdNlQh0URR8rwP9BKRM1UxTwN3isgIEQkUkcqOfuxdsVUGYOuQDwLfi0hTEfESkSAReVZE+uTxHPOAxiJym4iUcXQJbY6tXjlvxpiF2AvNUMeuQGzDaZyIhGAbp501A7hLRJo72mjOVMPhqE6aAbzmeC/qASOxJZbC9De2kfpJEfERkR7YC/y0vA52lDw+A94TkWoAIhIiIr0dv18nIo0cF9wEINNxA1sCapBfII56/d3803ttCzaRtnXcznRkaMvZpZyCTAWeF5Gqju7TL5D3e7gC2/g8wvF3MoCCe9G9i+1x97XjsznzPrwrIq1zHPc1tgrxRv5d7QS2mvIXJ1+LcoImihLGUWf8DfB/ju1lQG9gAPYb3n5sT5RuxpjdjmNSsQ2JO4CF2IvRKmxVzd95PMdx4DpsI+tx4EngulwNiufrLeyF1Q+bwNpjezTNxX7rdYox5hdssvwdW13ze65DhmMv4nuxjchTgC8vIu68YkgD+mE7FhwDPgKGGGN2FPCwpxzxrnRUty3in+q4MMd2Evbi+5Ex5g/HfW9gL9pxIvJ4Puf+BEe9vqNN6ciZG7YXVZZjOzOfx+c2GliDLa1sxvZq+tcgOsf7MADb0eAk9uKe72fpaG+6BFvq+9vR7vYb9u8gZ2eHpY590caYs9o8RKQjtgfdKidfi3KCON8mppTyRI7kux64yhhzPp0WPI6IfA98YYyZ5+5YShJNFEoppQqkVU9KKaUKpIlCKaVUgTRRKKWUKlCZcx9SvAQHB5vQ0FB3h6GUUh5l7dq1x4wxVc995L95XKIIDQ1lzZo15z5QKaVUNhEpaFR8gbTqSSmlVIE0USillCqQJgqllFIF0kShlFKqQJoolFJKFUgThVJKqQK5LFGIyJciEiMiW/K5X0RknIhEiMgmEWnvqliUUkpdOFeWKCZi1zTOzzXY6ZPDsOsQfOzCWJRSqtRKz8xzGXKnuWzAnTFmqYiEFnDI9cA3jmUSV4pIJRGpWdKnQVZKqQthjCEiJonN0fEcS0rlaEIqe2OTOByfQmJKBqkZmaSmZ5GVc0ZwY+i2ZRm9d6+4qOd258jsEM5eUSvKse9fiUJEhuJY/axu3YJWUlRKKc+WnpnFjsOJrN1/gi2HEog8dor9J5KJS04jPfOfJODv40VoUHlqVy5HxbI++Pt44VvGC2/HyrMVj0Zz9YTXCFuzhEP1Gl9UTO5MFHmtgZzn4hjGmE+BTwHCw8N1AQ2llEfIyMwiKTWDhNMZJKSkczo9k9NpmZxOzyQlPZP40+nEJKRyNCGFo4mpxCSksP94MqfT7WKDwQF+NKhaniuaVCUowI/QoHJ0qFeZ6hX8CfArw9nLkedgDITfATt3wjvvUGvECPDxueDX4c5EEYVdu/eM2uS9QLtSShUbaRlZnDiVxrGkVGKTUolNTOVIfAqH41M4mmB/njiVSmJKBslp515d1ttLqBrgR/UKftSpUo5LGgbTrm4lOtSrTK1KZc8vuOXLoVUrCAyEzz+H4GCoU+fcjzsHdyaKOcDDIjIN6AzEa/uEUsqdjDGcTE5n//FTHDiRzIHjyew/kczBE8kcS0rlWFIa8afT83xsUHlfalT0p1ZFf1qFVKCCvw+B/j4E+pehQln7s7xvGcr6euHv401ZH28C/MsQVN4Pb698SgbOOn4cnn7aJocXX4SXXoJ27S7unDm4LFGIyFSgBxAsIlHAi4APgDFmAjAP6INdND0ZuNtVsSillDGGpNQMYhNTiXGUAqLjTtvbSfvzUNzpf5UCqlfwo07lcjSpEcilAX4ElfcjONCX4AA/ggP8bGmgoh9+Zbzd8aLgm2/g8cfh5El44gl7K2Su7PV06znuN8BDrnp+pVTpY4whNjGVXUeT2HU0kd0xieyJOcWRhBRiE1Oz6/5zqlLel5BKZWlYtTzdw6oSUrksdauUo15QOepULkdZXzckAGc99RS89RZccglMmGCrnVzA49ajUEqpMwkhIiaJ3TGOpHA0iV0xicQl/1M1VKW8L42qBdCubiWqBfpR9cwtwJ8aFf2oVaks5Xw97DJ4+jScOmXbH+69F8LC7E8v1w2L87B3SClVmmRlGaLjThMRk+RIConZvyekZGQfV8G/DI2rB9KnVU0aVwugcfVAGtcIJDjAz43Ru8D8+fDQQ9C2LXz/PTRpYm8upolCKVVsRMedZvGOGNZEnmB3TBJ7Y0+dVV0UHOBLw6oB9Gtbi0ZVA2hULZBG1QKoXsEv/66iJcGhQ/Doo/DddzYxPPxwkT69JgqllNskpqSzOTqeJbti+WNHLDuPJgJQo4I/jWsE0qVBEI2qBdhb1QAql/d1c8Ru8NtvcMMNkJYGr75qG6v9irakpIlCKeVyaRlZ7Im1bQk7jiSy03GLjjsNQBkvoVP9KjzXoRlXNK1Gw6rlS3YJwRnp6XaQXJs20KcPjB4NjRq5JRRNFEqpQmWMYU9sEr9tj2HLoQR2Hklgb+wpMrLspAplvISGVQPoUK8yt3WuS9MagXSqX4VA/wsfOVyiJCTA//0f/P03/PWXbbSeNs2tIWmiUEpdtKwsw4aoOBZsPcqCrUfYe+wUACGVytK0RiA9m1WnSY1AmtaoQP3g8viW0aVw/sUYmDkTHnkEjhyBBx+E1FQoV87dkWmiUEpdmLSMLFbuPc6vW4+wcNtRYhJTKeMldG0YxN2XhtKreQ1qVPR3d5ieITYW7rwTfvnFjqiePRs6dnR3VNk0USilnJaSnsmSXbHM33KERduPkpiSQVkfb3o0qUrvFjW4okk1KpbTKqTzVqECHDsG779vu7+WKV6X5uIVjVKq2DmVmsHinTH8suUIi3fEkJyWSaVyPvRuUYOrW9SgW1gw/j7FePRycbV0Kbz2mh0PERAAK1e6dNDcxdBEoZT6l7jkNH7fEcP8LUdYsiuW1Iwsgsr70r9dCNe0rEGXBkH4eBfPi1qxd+yY7eI6cSKEhkJkJLRsWWyTBGiiUEo5HDiezIJttr1hzf6TZGYZqgX6MahjHa5pVZOOoVUufpbT0swY+OormyQSEuCZZ+D554tFY/W5aKJQqpTKyjJsjo5n4bajLNx2NHuwW+PqAQy7vAG9mtegdUhFvDQ5FJ5vv4Xmze0Efi1auDsap2miUKoUMcaw/mAcP6yLYuG2oxxNSMVLoGNoFZ6/thm9mlenXlB5d4dZciQnw+uvw7BhULu2bY+oWLFYVzPlRROFUqVATGIKP66L5ru1UUTEJOHv40WPxtXo1bw6VzatVjqnxnC1efNsD6bISAgJgQcegMqV3R3VBdFEoVQJlZaRxe87Ypi59iCLd8aSmWXoUK8yYwa04trWNXUktKtERdkJ/L7/Hpo1gyVLoHt3d0d1UTRRKFXCxCen88WyvUz++wDHT6VRNdCP/17WgJs61KZRtQB3h1fyvfYazJ1rq5xGjQJfzy+taaJQqoSIS07ji2X7mPhXJImpGfRsVp3bOtehe1hVymhXVtdatQrKlrUrzI0ebXs2NWjg7qgKjSYKpTxcXHIan/+5j4nLI0lKzeCaljUYcVUYzWpWcHdoJV98PDz7LHz8MVx3HcyZA0FB9laCaKJQykOdPJXG58v28vXy/SSlZnBtq5oMv6oRTWtognA5Y2D6dHjsMYiJgeHD7VoRJZQmCqU8zLGkVL5cto+vl0eSnJ5Jn1Y1GXFlGE1qBLo7tNLj229hyBAID4eff4YOHdwdkUtpolDKA6SkZ7Jo+1F+XBfNkl2xZBrDta1qMuKqMBpX1wRRJFJTYe9e25PpllsgI8MmC++SP8+VJgqliqmsLMPKfceZtT6aXzYfITE1g+oV/Li3W31uDq9No2qaIIrM4sV2HERyMuzebZcivftud0dVZDRRKFXM7DqayA/ropm9IZrD8SmU9/XmmlY1uaFdCF0aBOl8S0UpJgYefxwmTbK9mD79tMjXqy4ONFEoVQzEJqYye0M0P6yLZtvhBLy9hO5hwTzTpxm9mlWnrG/Jr94odiIioFMnSEqC556zt7Jl3R2VW2iiUMqNElPS+WTJXj5ftpeU9Cza1K7Ii32b07dNLYIDSt8312IhIcEuJNSwIdx7L9xzj22XKMU0USjlBhmZWUxdfZD3F+7i+Kk0rm9bi+FXNtJ2B3c6dQpeeQU++ww2bbKT+L31lrujKhY0UShVhIwxLNoew5hftrMn9hSd6lfhq2ub0bp2JXeHVrr99BM8/DAcOGBLER6wRkRR0kShVBHZFBXHa3O38/e+EzSoWp7PhoTTs1k1RLRx2m0yMmxX1x9/tOtD/PkndOvm7qiKHU0USrmQMYYVe47z2Z97WbwzlqDyvrzavyWDOtbRpUTdyRgQgTJloGZNGDPGjrIuARP4uYImCqVcID0zi583HeKzpfvYdjiBoPK+PNazMfd0C9Xpvd1t5Uq7TsRnn0H79jB+vLsjKvY0UShViOJPpzN11QEm/hXJkYQUGlULYMyAVvRvF4K/j3ZxdauTJ+0Efp98ArVq2W3lFJcmChG5Gvgf4A18bowZk+v+usDXQCXHMU8bY+a5MialXOHgiWS+WLaPGWsOkpyWySUNg3hjQCsub1xV15wuDqZPhxEj4Ngxu6jQyy9DoPYwc5bLEoWIeAPjgV5AFLBaROYYY7blOOx5YIYx5mMRaQ7MA0JdFZNShSkhJZ35m4/w4/poVu47jrcIfdvU4t5u9WkZUtHd4amcduyA0FCYPx/atXN3NB7HlSWKTkCEMWYvgIhMA64HciYKA5yZE7kicMiF8Sh10dIysliyK5ZZ66NZuP0oaRlZhAaVY8SVYQzqVIeaFUvnyN1iJyUF3nzTtkH07WurnJ5/vlRM4OcKrkwUIcDBHNtRQOdcx7wELBCR4UB5oGdeJxKRocBQgLp16xZ6oEoVxBjD2v0n+XF9NHM3HyYuOZ2g8r7c2rEO/duF0LZOJe3iWpwsWgQPPmgn7xs1yiYKH+1AcDFcmSjy+s8xubZvBSYaY94Rka7AJBFpaYzJOutBxnwKfAoQHh6e+xxKuURCSjpfLYtk5rqDHDxxGn8fL/7TvAY3tAuhW1iwdm8tbo4ehZEjYcoUaNQIFiyAXr3cHVWJ4MpEEQXUybFdm39XLd0LXA1gjFkhIv5AMBDjwriUKlBKeiaTVuxn/B8RxCWnc1lYMI9e1ZjeLWsQ4KcdBYuthQth5kx44QV45hnw93d3RCWGK//qVwNhIlIfiAYGAbflOuYAcBUwUUSaAf5ArAtjUipfGZlZzFwbxf9+283h+BS6N67Kk72baMN0cbZxo61iuukmGDwYLr0U6td3d1QljssShTEmQ0QeBn7Fdn390hizVUReAdYYY+YAo4DPROQxbLXUXcYYrVpSRcoYwy9bjvD2gp3sjT1F2zqVeOeWNlzSMNjdoan8JCXBiy/C//5nezP1729HWWuScAmXlqMdYyLm5dr3Qo7ftwGXujIGpQqybPcxxv66g01R8YRVC+CTOzrwn+bVtXG6OJs1C4YPh6goGDoU3njDJgnlMvruqlJp48E4xv66g78ijhNSqSxv39yGG9qF6Opxxd3mzXDDDdCqlR1Ed8kl7o6oVNBEoUqVfcdOMXb+Dn7ZcoQq5X154brmDO5SF78y2r++2EpPt7O6XnmlTRBz59reTNrltchoolClQlxyGv/7bTeTVuzHt4wXj/YM477LGmgvpuJu+XIYNgy2boWdO2231z593B1VqaP/JapES8vI4psVkYz7bTdJqRkM7FiHx3o1plqgdp0s1k6cgKeftjO81qkDP/xgk4RyC00UqkQyxjB/yxHGzN/B/uPJdG9clef6NKNJDZ0IrthLSYG2beHQITuy+qWXICDA3VGVapooVImz8WAco+duY3XkSRpXD2Di3R3p0aSau8NS5xIVZdep9veHV1+1yaJNG3dHpdBEoUqQ6LjTjJ2/g9kbDhEc4MvrN7TilvDalNGpNoq306dtF9c337Qjq/v2hTvvdHdUKgenEoWI+AJ1jTERLo5HqfN28EQyHy/Zw8w1UYjAQ1c0ZNjlDXUlOU+wYIGdwG/PHrj9dujUyd0RqTycM1GIyLXAu4AvUF9E2gIvGmNucHVwShUk8tgpxi+O4Mf10YjAzeF1eOiKRoRU0qm+PcLw4fDhhxAWZmd8veoqd0ek8uFMieIV7PTgiwGMMRtERLsfKLeJiEli/OIIZm+Ixsfbi9u71OP+yxvoWhCeIDPT/vT2hi5dIDgYnnpKJ/Ar5pxJFOnGmLhcUxrofEyqyO04ksAHv0cwb/Nh/Mt4c2+3+vy3ewPt6uop1q2zYyLuuMOWJgYPdndEyknOJIrtInIL4OWYCfYRYKVrw1LqHzuPJPLuwp38uvUoAX5leODyhtzbrT5BAX7uDk05IzHRTv09bhxUrQo1a7o7InWenEkUDwMvAFnAD9jZYJ9xZVBKgW2kfm/hLn7cEE2AbxkeuSqMuy8NpVI5X3eHppy1YAHcc48dEzFsGLz+OlSq5O6o1HlyJlH0NsY8BTx1ZoeIDMAmDaUKXUxCCh8ujmDqqgN4iTC0ewOGdW9I5fKaIDyOry9Uqwbffw+dc6+ErDyFnGv5BxFZZ4xpn2vfWmNMB5dGlo/w8HCzZs0adzy1crH45HQ+WbqHL//aR0amYWDHOoy4KozqFbQNwmOkp8O770JCArz2mt2XlQVeOpbF3RzX7fALeWy+JQoR6Y1dpjRERN7NcVcFbDWUUoXidFomXy3fx4Q/9pCQkkG/NrUY2asxocHl3R2aOh/Llv0zgd/NN/+TIDRJeLyCqp5igC1ACrA1x/5E4GlXBqVKh7SMLKavPsC43yOITUzlyqbVePw/TWheq4K7Q1Pn4/hx28X1iy+gbl346Se47jp3R6UKUb6JwhizHlgvIpONMSlFGJMq4TKzDHM2RvPuwl0cPHGaTqFV+GhwezqGVnF3aOpCHD8O06bBk0/a3k3ltSRY0jjTmB0iIq8BzYHsymJjTGOXRaVKrNjEVO79ejWbouJpXrMCX93dkh6Nq+rSo55m+3aYMcOuW924MRw4AFU00ZdUziSKicBo4G3gGuButI1CXYCok8nc/vnfHE1I5f2BbenXphZeuvSoZ0lOto3Ub71lp/6+914746smiRLNmVamcsaYXwGMMXuMMc8DV7g2LFXSRMQkctPHKzhxKo1v7+tM/3YhmiQ8zfz50LKlHQtx2212xbnatd0dlSoCzpQoUsXWC+wRkWFANKCT+yunbY6K586vVuElwvT7u9KspjZWe5ykJDv1RlAQLF4MPXq4OyJVhJwpUTwGBAAjgEuB/wL3uDIoVXL8vfc4t362krI+3swcpknCo2Rmwrff2p8BAXaG140bNUmUQucsURhj/nb8mgjcASAiWt5U5/T7jqM88O06alcuy7f3ddbZXT3J2rVw//32Z9mycOONutpcKVZgiUJEOopIfxEJdmy3EJFv0EkB1TnMWh/N0G/W0rh6IDPu76pJwlPEx8OIEXYBoeho2+11wAB3R6XcLN9EISJvAJOBwcB8EXkOuybFRkC7xqp8fbZ0L49O30B4aGWm/LezzvLqSW680S4m9OCDsGMHDBwI2nW51Cuo6ul6oI0x5rSIVAEOObZ3Fk1oytNkZRlen7edz5ft49pWNXl3YBv8yni7Oyx1Lnv32um/AwNt11cvL+jY0d1RqWKkoKqnFGPMaQBjzAlghyYJlZ+0jCwem7GBz5ft486u9Rh3aztNEsVdWprt6tqiBYwebfd17qxJQv1LQSWKBiJyZipxAUJzbGOM0YpLBUBSagYPfLuWP3cf44neTXiwR0MdaV3cLV1qJ/Dbvh1uusm2SyiVj4ISxY25tj90ZSDKM8UmpnLPxNVsO5zAWze15ubwOu4OSZ3Le+/ByJEQGgpz50KfPu6OSBVzBU0K+FtRBqI8z6aoOB6esp6YxBQ+G9KBK5tWd3dIKj9ZWXDqlG2HuPZaiI2F55+HcuXcHZnyADpRvDpv6ZlZvL9oFzd8tJy0jCym/LeLJonibOtWuPxyuOsuu924sW2b0CShnOTSRCEiV4vIThGJEJE817AQkVtEZJuIbBWRKa6MR128iJgkbvp4Oe8v2k2/NrX49dHutK9b2d1hqbwkJ8Mzz0DbtrYt4rrr4BwrWiqVF2fmegJARPyMManncbw3MB7oBUQBq0VkjjFmW45jwoBngEuNMSdFROeQKqaysgxfr4hkzC87KOfrzUeD29OnVU13h6Xys369HSgXGQl33w1jx0JwsLujUh7qnIlCRDoBXwAVgboi0ga4zxgz/BwP7QREGGP2Os4zDTs2Y1uOY/4LjDfGnAQwxsSc/0tQrhYdd5onvtvI8j3HubJpNcYMaEU1Xce6eDLGDpCrW9fevv4aund3d1TKwzlTohgHXAfMAjDGbBQRZ6YZDwEO5tiOAjrnOqYxgIj8BXgDLxlj5jtxblUEjDH8uD6aF2dvJcsYxgxoxcCOdbTra3GUkWFHVM+ZAwsX2llelyxxd1SqhHAmUXgZY/bnujhkOvG4vK4muStIywBhQA+gNvCniLQ0xsSddSKRocBQgLp16zrx1OpiHU9K5bkftzB/6xE6hlbmnZvbUjdIGz+LpVWr7JiI9evhmmsgIQEqa7uRKjzOJIkeN+kAACAASURBVIqDjuon42h3GA7scuJxUUDOTvW1sdOA5D5mpTEmHdgnIjuxiWN1zoOMMZ8CnwKEh4dra5wLZWRm8e3K/by7cBcp6Vk8c01T7rusAd66yFDxk5QETz0FH38MNWvCd9/ZuZq0xKcKmTOJ4gFs9VNd4CiwyLHvXFYDYSJSH7vY0SDgtlzHzAJuBSY6ZqhtDOx1LnRV2JbvOcbLc7ax82gi3RoF81K/5jSqFujusFR+fHzgjz9g+HB49VWooGt9KNdwJlFkGGMGne+JjTEZIvIw8Cu2/eFLY8xWEXkFWGOMmeO47z8isg1bnfWEMeb4+T6XujjRcad5fe525m4+TO3KZfnkjg78p3l1bYsojiIi4JVXYPx4O3hu7Vrw144FyrXEnKNftYjsAXYC04EfjDGJRRFYfsLDw82aNWvcGUKJkZKeySdL9vLxkggAHuzRiKHdG+Dvo5P5FTupqbaL62uvga+vnXrjssvcHZXyICKy1hgTfiGPdWaFu4Yicgm26uhlEdkATDPGTLuQJ1TuZ4zh161HePXn7UTHneba1jV5tk8zQirp4kLF0uLF8MADsHOnXR/i3XehVi13R6VKEacG3BljlgPLReQl4H3sgkaaKDxQ1Mlknvp+E39FHKdpjUCm/rcLXRsGuTsslR9jbCkiPR3mz4fevd0dkSqFnBlwF4AdKDcIaAbMBi5xcVzKBfbEJnH753+TlJLBK9e34LZOdSnjrdN9FTtZWfDFF3D11VCnDkyaBJUq2bWrlXIDZ64SW4AuwFhjTCNjzChjzN8ujksVsu2HExj4yQrSM7OYMawrQ7qGapIojjZtgm7dYOhQ+Pxzu69mTU0Syq2cqXpqYIzJcnkkymU2HoxjyJerKOvjzeT/dqZh1QB3h6RyS0qCl1+2a0VUrgwTJ8KQIe6OSimggEQhIu8YY0YB34vIv7pG6Qp3nmHVvhPcM3E1lcv7MOW+LtSpoqOri6WXXoJ33oH77oMxY+wUHEoVEwWVKKY7furKdh5q6a5Yhk5aQ0ilsky+rws1Kmp/+2Ll4EG7mFDTpvD009C/v612UqqYybeS2hizyvFrM2PMbzlv2EZtVYwt2HqE+75eQ2hQeabf31WTRHGSkWG7uDZrBvffb/cFB2uSUMWWM62Z9+Sx797CDkQVntkbonlg8jqa1arAtKFdCA7wc3dI6oyVKyE8HEaNgh497DTgShVzBbVRDMR2ia0vIj/kuCsQiMv7Ucrdpq06wDM/bqZjaBW+vKsjAX5Or02lXG3uXOjb1w6W++EHW9Wk06QoD1DQVWQVcBw76+v4HPsTgfWuDEpdmAlL9jDmlx10b1yVT27vQFlfnYrD7YyBQ4cgJAR69rTzND3yiJ2nSSkPkW+iMMbsA/ZhZ4tVxZgxhjG/7OCTpXu5tnVN3r2lDX5lNEm43a5d8OCD9ue2bRAQAM8/7+6olDpvBVU9LTHGXC4iJzl7wSEBjDGmisujU+eUkZnFsz9uZsaaKG7vUpeX+7XUtSPcLSXFdnF94w07UO7MT6U8VEFVT2eWO9UV2YuplPRMRkxdz4JtRxlxVRiP9QzTqcHd7cgRu0b17t1w6622d1ONGu6OSqmLUlD32DOjsesA3saYTKArcD9QvghiUwVISEnnzi9XsWDbUV7q25yRvRprknCn9HT7s3p1mygWLIApUzRJqBLBme6xs7DLoDYEvsGOoZji0qhUgY4lpXLrpytZu/8k/xvUlrsure/ukEqvrCyYMAEaNoSoKNuL6fPPoVcvd0emVKFxJlFkOda0HgC8b4wZDoS4NiyVn4Mnkrl5wgr2xCbx2Z3hXN9WPwq32bgRLrnErhURFvZPqUKpEsaZRJEhIjcDdwA/O/b5uC4klZ9thxK4acJyjielMvm+zlzRpJq7QyqdjIHHH4cOHWDvXjsN+KJFUF9LdqpkcnZk9hXYacb3ikh9YKprw1K5Ld0Vyy2frMBLhO+GXUKHetrpzG1E4ORJuPdeu+rc7bfrwDlVop1zzWwAESkDNHJsRhhjMlwaVQFK45rZM9Yc5NkfNtOoWgAT7+6k8za5w/79dqDcCy9A+/a2bcJL1/NQnuNi1sw+51+6iFwGRABfAF8Cu0Tk0gt5MnV+jDG8t3AXT87cRNeGQXw3TCf3K3Lp6TB2LDRvDgsX2hIEaJJQpYozEwG9B/QxxmwDEJFmwCTggjKTck56ZhbP/rCZ79ZGcVOH2rwxoBU+uiJd0Vq+3M7uumULXH89jBsHdeu6OyqlipwzicL3TJIAMMZsFxFfF8ZU6iWmpPPg5HX8ufsYj1wVxqM6kM49Fi2C+HiYNcsmCqVKqXO2UYjIRCAVW4oAGAyUM8bc6drQ8lbS2yiOxKdw98TV7DqayBs3tOKWjnXcHVLpYYztwVS1KlxzDaSm2qqnAF06Vnk+l7ZRAMOAPcCTwFPAXuzobFXIlkcc44aP/uLA8VN8eVdHTRJFaccOuPJKuPNO+Ooru8/PT5OEUpyj6klEWgENgR+NMWOLJqTSJzElndfn7WDqqgOEBpVjxrCutKhV0d1hlQ6nT8Prr8Obb0L58vDJJ3bdaqVUtoJmj30Wu5LdOqCjiLxijPmyyCIrJRbvjOHZHzZzNCGFod0b8FjPxrqORFH66ScYPdqOhXj7bTtXk1LqLAWVKAYDrY0xp0SkKjAP2z1WFYL45HRe+Xkb36+LIqxaAB89cAnt6lZ2d1ilw5EjsGEDXH013HwzhIZCp07ujkqpYqugRJFqjDkFYIyJFRHtm1lIFmw9wnOztnDiVBoPX9GI4Vc10oWGikJmpq1aeuYZ8PWFAwfsOhGaJJQqUEGJokGOtbIFaJhz7WxjzACXRlZCnVmutFnNCnx1V0dahmhbRJFYtw6GDYPVq+2SpB99pIsJKeWkghLFjbm2P3RlIKVBUmoGHy2OoEeTqnx6Rzi+ZbSQViT27bOlhuBgu0bEoEE6N5NS56GgNbN/K8pASoPpqw+SkJLBI1eFaZJwNWNg82Zo3drO6vrVV9C3L1Sq5O7IlPI4erUqIumZWXy5bB+d6lfRRmtX27cPrrsO2rWDTZvsvjvu0CSh1AVyaaIQkatFZKeIRIjI0wUcd5OIGBEpsfNHzd10mOi409zfvYG7Qym50tJgzBho0QKWLLHdXZs3d3dUSnk8Z+Z6AkBE/IwxqedxvDcwHugFRAGrRWROznmjHMcFAiOAv509t6cxxvDJ0r00qhagiw25SmamXW1u7VoYMADefx/q6Mh2pQqDM9OMdxKRzcBux3YbEfnAiXN3wq5dsdcYkwZMA/KaWe1VYCyQ4nzYnmVZxDG2H05gaPcGeHlpI2qhSkiwP7294Z577AC677/XJKFUIXKm6mkccB1wHMAYsxG74t25hAAHc2xHkWutbRFpB9QxxvxMAURkqIisEZE1sbGxTjx18WGMYcKSPVQL9OP6trXcHU7JYQxMnAgNGsDs2Xbfgw/atgmlVKFyJlF4GWP259qX6cTj8vrqnD1VrWMA33vAqHOdyBjzqTEm3BgTXrVqVSeeuvj4Ytk+/oo4zv2XN9RBdYVl2zbo0QPuvhuaNoWGDd0dkVIlmjOJ4qCIdAKMiHiLyKPALiceFwXkLP/XBg7l2A4EWgJ/iEgk0AWYU5IatBfvjOH1edu5pmUN7r4k1N3hlAxjx0KbNnYxoc8/h6VLoWVLd0elVInmTKJ4ABgJ1AWOYi/oDzjxuNVAmIjUdyx0NAiYc+ZOY0y8MSbYGBNqjAkFVgL9jDElYrGJiJhERkxZT9MaFXjnljbaNnGxzqybUqMGDB5spwW/915dklSpInDOXk/GmBjsRf68GGMyRORh4FfAG/jSGLNVRF4B1hhj5hR8Bs918lQa9369Bj8fbz67M5xyvk53LlO5HToEjzwCl10GI0bAkCH2ppQqMue8gonIZ+RoWzjDGDP0XI81xszDzjqbc98L+Rzb41zn8wTpmVk8NGUdh+NSmDq0CyGVdD6hC5KZaedjeu45u8rcJZe4OyKlSi1nvuouyvG7P3ADZ/dmUjm88tM2lu85zjs3t6FDPR2BfUE2bLCLB61dC//5j00Y2mCtlNs4U/U0Pee2iEwCFrosIg/29fJIJq3cz/2XN+DGDrXdHY7nio+3VU7Tp9v1InQCP6Xc6kIqz+sD9Qo7EE83aUUkL87ZSs9m1Xmyd1N3h+NZjIHvvoPdu21V0+WXw9694O/v7siUUjg3MvukiJxw3OKwpYlnXR+a5/j8z73832ybJMYPboe39nBy3p490KcPDBxoB86lp9v9miSUKjYKLFGIiABtgGjHrixjzL8atkuzj/6IYOz8nVzTsgb/G9ROpw93VmqqnbRv9Gjw8YH//c+OrC6jPcSUKm4K/K80xhgR+dEY06GoAvIUxhjG/RbBe4t20a9NLd69pQ1lvDVJOO3gQXj1VbtGxPvvQ0jIuR+jlHILZ65sq0Skvcsj8SDGGN5esJP3Fu3ixva1eW9gW00SzoiNhQ8dCyU2amSn4vjuO00SShVz+V7dRORMaaMbNlnsFJF1IrJeRNYVTXjF01u/7mT84j3c2qkOb93UWtskziUrC774ws7LNHIk7Nxp9zfQtTmU8gQFVT2tAtoD/YsoFo+waNtRPvpjD4M61uG1/q10ao5z2bIFHngAli2zo6snTIAmTdwdlVLqPBSUKATAGLOniGIp9mISUnjy+000q1mBl69voUniXNLS7IC5tDT48ku46y4dE6GUByooUVQVkZH53WmMedcF8RRbWVmGUd9tJDktg3GD2uqU4QX5/Xc7FsLXF2bMsFVOwcHujkopdYEKaoH1BgKw04HndStVvvxrH3/uPsb/XdecsOql7uU7JyoKbrwRrroKvvnG7uvWTZOEUh6uoBLFYWPMK0UWSTG2JTqeN+fvoFfz6tzWqa67wyl+MjJsb6b/+z87md8bb9ipwJVSJcI52yhKu+S0DB6Ztp4q5X1588bWiNax/9sdd8C0aXDNNTB+PNSv7+6IlFKFqKBEcVWRRVFMGWN4/sct7D12im/v7UyV8r7uDqn4iIuzo6gDAuChh2yV0403amO1UiVQvm0UxpgTRRlIcfTxkj38sD6aR64K49JGWs8O2An8pk2DZs1sVRPYdoibbtIkoVQJpcOJ8zF/yxHGzt9J3za1eOSqMHeHUzxEREDv3nDrrVC7Ntx+u7sjUkoVAU0UedgSHc9j0zfQtk4l3rpJ2yUAmDIFWraEv/+2DdcrV0IHnQJMqdJAp+rMJSYhhfu+XkPlcj58OqQD/j6lfLxEerqd3TU83FYvjR0LtWq5OyqlVBHSEkUOmVmGoZPWkpCSzud3dqRaYCleEyEmxvZmGjjQbjduDN9+q0lCqVJIE0UO01YfYMPBON4Y0IrmtSq4Oxz3yMqCTz+18zFNnw4tWtixEUqpUkurnhwSUtJ5Z8EuOoVWoV+bUvqtee9e20C9YgX06AEff2yn31BKlWqaKBw++G03J5PTeKFv89LbeF2xoh0f8fXXttqptL4PSqmzaNUTsO/YKSYuj+TmDrVpGVLR3eEUrTlzYMAAW70UFGSnBR8yRJOEUiqbJgrgtbnb8fX24vHepWidhAMHoH9/uP562LULDh+2+730T0IpdbZSf1VYtvsYi7Yf5aErG5WOXk4ZGfD223Zk9YIF8OabsH69HUCnlFJ5KNVtFBmZWbzy81bqVCnLPZeWkonsMjPh88/hyivhgw8gNNTdESmlirlSXaKYuuoAu44m8VyfZiV7YN3Jk/DUU5CYCH5+8Ndftm1Ck4RSygmlNlHEJ6fz7sJddGlQhd4targ7HNcwBiZPtl1c33kHFi+2+4OCtLFaKeW0UpsoPvtzL3Gn03nhuhYlszvsrl3Qq5cdFxEaCmvWQL9+7o5KKeWBSmUbRUJKOl+viOSaljVK7gjsRx+1yeGjj2DoUPAuwVVrSimXKpWJYtKK/SSmZPBgj0buDqVwLVxoq5nq1LGjqv38oEYJrVZTShUZl1Y9icjVIrJTRCJE5Ok87h8pIttEZJOI/CYi9VwZD9ilTb9Yto8eTaqWnMF1R47AbbfBf/5ju7sC1KunSUIpVShclihExBsYD1wDNAduFZHmuQ5bD4QbY1oDM4GxrornjGmrDnLiVBoPX1ECShNZWTBhgi1FfP89vPiiHSOhlFKFyJUlik5AhDFmrzEmDZgGXJ/zAGPMYmNMsmNzJeDSUV+pGZl8unQvnepXITy0iiufqmi88QY88IBdQGjTJnjpJfAvBYMGlVJFypVtFCHAwRzbUUDnAo6/F/glrztEZCgwFKBu3boXHNDMtVEcSUhh7E2tL/gcbpeYCMeOQf36MGyY/XnrrdrdVSnlMq4sUeR15TJ5HihyOxAOvJXX/caYT40x4caY8KpVq15QMCnpmXz4ewTt6lbisrDgCzqHWxkDP/4IzZvbxYSMseMhbrtNk4RSyqVcmSiigDo5tmsDh3IfJCI9geeAfsaYVFcFM/nvAxyOT+GJ3k08b9zE/v12DMSAAVClCowbp8lBKVVkXFn1tBoIE5H6QDQwCLgt5wEi0g74BLjaGBPjqkBOpWbw0eIILm0UxCUNPaw0sWIF9Oxpf3/7bXjkEShTKns1K6XcxGUlCmNMBvAw8CuwHZhhjNkqIq+IyJkhwm8BAcB3IrJBROa4Ipav/trH8VNpPP4fD5pGPCHB/mzfHu65B7Zvh1GjNEkopYqcS686xph5wLxc+17I8XtPVz4/wIHjyYxfvIdezavTrm5lVz/dxTt+HJ5+2k4BvnUrBATYWV6VUspNSvRcT1lZhsdnbqSMl/ByvxbuDqdgxsA339gxEV99ZRustR1CKVUMlOh6jG9WRLJq3wnG3tiaWpXKujuc/MXH29Xm/vgDuna1g+hae3AXXqVUiVJiE0XksVO8OX8nPZpU5ebwYrp6mzG21FChAgQHw6efwr336nKkSqlipcRekV6cs5Uy3sKYAa2LZ3fYX3+1DdVRUTZZfPcd/Pe/miSUUsVOibwqnU7L5K+IYwzuXI8aFYvZlBaHD8OgQXD11ZCcDDEu6xWslFKFokQmig0H48jIMnSqX8x6OY0fbxurZ82Cl1+28zO1b+/uqJRSqkAlso1i7f4TAHSoW8wm/lu7Fjp3tgkjLMzd0SillFNKZIliVeRJGlcPoGI5H/cGkpBgV5pbu9Zuf/SRbZvQJKGU8iAlLlFMWhHJ0l2xXNm0uvuCMAZmzoRmzey8TEuW2P3+/jo2QinlcUpUopix+iD/N3srPZtVZ9R/GrsniH374Lrr4OaboVo1O1fTyJHuiUUppQpBiUkUszdE89QPm7gsLJgPb2uHj7ebXtrkybB0Kbz3HqxebdsklFLKg5WIxuy9sUmMnLGRjqFV+PSOcPx9vIs2gD//hNRUO8vrE0/AXXdB7WI6yE8ppc5TiShRbDmUQGaW4ZXrW1DWtwiTxLFjdmbX7t3hlVfsPj8/TRJKqRKlRJQook7aZbdrVy5XNE9oDEycaEsP8fHw1FPwf/9XNM+tXCY9PZ2oqChSUlLcHYpSF8zf35/atWvj41N4vT5LRKLYG3uKKuV9CfAropczb54tSVx6qZ3Ar2XLonle5VJRUVEEBgYSGhpaPKd9UeocjDEcP36cqKgo6tevX2jn9fiqJ2MMK/Ycp2Ooi0dhJyfDX3/Z3/v0gdmzbaO1JokSIyUlhaCgIE0SymOJCEFBQYVeKvb4RBF5PJnouNN0a+TCJU5/+cUmhGuugbg4OxaiXz+dwK8E0iShPJ0r/oY9/kr34e8R+Hp7cWUzFwywi4624yH69LGN1D/9BJUqFf7zKKVUMebRiWLDwTi+XxfFPd3qE1LYCxPFxEDz5vDzzzB6NGzcCJdfXrjPoVQu3t7etG3blpYtW9K3b1/i4uIK5byRkZG0LKRq0rvuuov69evTtm1b2rZty7hx4wrlvHn5448/WL58eb73z5o1i1fO9DgshowxjBgxgkaNGtG6dWvWrVv3r2MSExOz38u2bdsSHBzMo48+CsCECRNo1aoVbdu2pVu3bmzbtg2AzZs3c9dddxXtC/GkW4cOHYwxxqSmZ5p+H/xpwkcvNIkp6abQREX98/v//mdMREThnVsVa9u2bXN3CKZ8+fLZvw8ZMsSMHj26UM67b98+06JFi0I515133mm+++67C3psRkbGeR3/4osvmrfeeivf+7t27WpiY2OdPl96eiFeK5wwd+5cc/XVV5usrCyzYsUK06lTp3M+pn379mbJkiXGGGPi4+Oz98+ePdv07t07e/uqq64y+/fvz/Mcef0tA2vMBV53PbLXkzGG52dtZmNUPONva184vZ3i4+H55+GTT2DlSjv994gRF39e5ZFe/mkr2w4lFOo5m9eqwIt9nV+7vWvXrmzatAmApKQkrr/+ek6ePEl6ejqjR4/m+uuvJzIykmuuuYZu3bqxfPlyQkJCmD17NmXLlmXt2rXcc889lCtXjm7dumWfNyUlhQceeIA1a9ZQpkwZ3n33Xa644gomTpzIrFmzyMzMZMuWLYwaNYq0tDQmTZqEn58f8+bNo0qV/Gdknjp1Kq+//jrGGK699lrefPNNAAICAhg5ciS//vor77zzDmXLlmXkyJEkJSURHBzMxIkTqVmzJuPGjWPChAmUKVOG5s2bM2bMGCZMmIC3tzfffvstH3zwAZdddln28+3atQs/Pz+Cg2375E8//cTo0aNJS0sjKCiIyZMnU716dV566SUOHTpEZGQkwcHBTJo0iaeffpo//viD1NRUHnroIe6///583+OLMXv2bIYMGYKI0KVLF+Li4jh8+DA1a9bM8/jdu3cTExOT/TorVKiQfd+pU6fOan/o27cv06ZN48knn7yoGJ3hkVVPn/25lxlrohhxZSOubZ33G+40Y2DGDDuB3/jxMGwYNGxYOIEqdYEyMzP57bff6NevH2D7xv/444+sW7eOxYsXM2rUKOyXRHtxeeihh9i6dSuVKlXi+++/B+Duu+9m3LhxrFix4qxzjx8/HrDVF1OnTuXOO+/M7iWzZcsWpkyZwqpVq3juuecoV64c69evp2vXrnzzzTfZ53jiiSeyq0o2b97MoUOHeOqpp/j999/ZsGEDq1evZtasWYC9wLVs2ZK///6bzp07M3z4cGbOnJmdyJ577jkAxowZw/r169m0aRMTJkwgNDSUYcOG8dhjj7Fhw4azkgTAX3/9Rfsc67l069aNlStXsn79egYNGsTYsWOz71u7di2zZ89mypQpfPHFF1SsWJHVq1ezevVqPvvsM/bt21fge5zTwIEDz6oqOnPL+f6cER0dTZ06dbK3a9euTXR0dL6f+9SpUxk4cOBZCWH8+PE0bNiQJ5988qxqvvDwcP788898z1WYPK5EYQy8/esuejarzqM9L3LiP2NgwAC7kFD79jBnDoSHF06gyqOdzzf/wnT69Gnatm1LZGQkHTp0oFevXoAtRT/77LMsXboULy8voqOjOXr0KEB2ewFAhw4diIyMJD4+nri4OC53tKvdcccd/PLLLwAsW7aM4cOHA9C0aVPq1avHrl27ALjiiisIDAwkMDCQihUr0rdvXwBatWqVXboBeOutt7jpppuyt2fPnk2PHj2oWrUqAIMHD2bp0qX0798fb29vbrzxRgB27tzJli1bsl9XZmZm9rfr1q1bM3jwYPr370///v3P+V4dPnw4+/nAjoMZOHAghw8fJi0t7axxBP369aNsWduOuWDBAjZt2sTMmTMBiI+PZ/fu3dSuXTvP97hGjRpnPe/06dPPGdsZeSWagnolTZs2jUmTJp2176GHHuKhhx5iypQpjB49mq+//hqAatWqcejQIadjuRgelyjSMrNIy8yiT6saeHldYDew9HTw8bHdXLt1gyuvhAcfBO8iniNKqVzKli3Lhg0biI+P57rrrmP8+PGMGDGCyZMnExsby9q1a/Hx8SE0NDS7FODn55f9eG9vb06fPo0xJt8LUl4XrzNynsvLyyt728vLi4yMjHwfV9A5/f398Xb8bxljaNGixb9KOQBz585l6dKlzJkzh1dffZWtW7fme06w71V8fHz29vDhwxk5ciT9+vXjjz/+4KWXXsq+r3z58mfF+sEHH9C7d++zzjdx4sR83+OcBg4cyM6dO/+1f+TIkQwZMuSsfbVr1+bgwYPZ21FRUdSqVSvP17Nx40YyMjLo0KFDnvcPGjSIBx54IHs7JSUlO/m5msdVPaVlZAIQGlz+HEfm448/oHVrO2AOYNQoGD5ck4QqVipWrMi4ceN4++23SU9PJz4+nmrVquHj48PixYvZv39/gY+vVKkSFStWZNmyZQBMnjw5+77u3btnb+/atYsDBw7QpEmTi4q3c+fOLFmyhGPHjpGZmcnUqVOzSzM5NWnShNjY2OxEkZ6eztatW8nKyuLgwYNcccUVjB07lri4OJKSkggMDCQxMTHP52zWrBkRERHZ2/Hx8YSEhABkf+vOS+/evfn4449JT08H7Htw6tQpp9/j6dOns2HDhn/dcicJsCWZb775BmMMK1eupGLFivm2T0ydOpVbb731rH27d+/O/n3u3LmE5Vj0bNeuXYXWk+1cPK5EkZqRBUC9Kuc5r1NsLDz+OHzzDdSvD4GBLohOqcLTrl072rRpw7Rp0xg8eDB9+/YlPDyctm3b0rRp03M+/quvvspuzM757fnBBx9k2LBhtGrVijJlyjBx4sSzShIXombNmrzxxhtcccUVGGPo06dPng3Bvr6+zJw5kxEjRhAfH09GRgaPPvoojRs35vbbbyc+Ph5jDI899hiVKlWib9++3HTTTcyePftfjdndu3fPbkcQEV566SVuvvlmQkJC6NKlC/v27csz1vvuu4/IyEjat2+PMYaqVasya9asC3qPz6VPnz7MmzePRo0aUa5cOb766qvs+9q2bcuGDRuyt2fMmMG8efPOevyHH37IKBR+ygAACf5JREFUokWL8PHxoXLlymclwMWLF3PttddedIzOkIKKjMVRrbAWpvKt77Dl5d7Oj0CcOhUeegiSkuxEfs89B+WKaAJB5TG2b99Os2bN3B2GOg+PPPIIffv2pWfPnu4OpUilpqZy+eWXs2zZMsqU+ff3/bz+lkVkrTHmghphPa7qKSPTUL2C//kNU8/IsFNwbNgAr72mSUKpEuLZZ58lOTnZ3WEUuQMHDjBmzJg8k4QreFzVEzix7PSpU/Dqq1C3rm2kvv12e9N5fJQqUapXr57dhbg0CQsLO6u9wtU8rkSRnJZJ9Qr++R/w88/QogW8+SY4uvwhoklCOcXTqmKVys0Vf8MelyjSM7MY0jX033dERdkxEX37Qvnydgrw998v8viU5/L39+f48eOaLJTHMo71KPz9C/gyfQE8suqpYdU8usbu3Qu//gpvvAEjR4Kvb9EHpjxa7dq1iYqKIjY21t2hKHXBzqxwV5g8MlFkW7UKVqyARx6x61YfOABBQe6OSnkoHx+fQl0VTKmSwqVVTyJytYjsFJEIEXk6j/v9RGS64/6/RSTUmfN6xcfZRuouXeDdd23jNWiSUP/f3t3HXFnXcRx/f0RJSMOM2VQUdKCFhERklC0jzCktLMd4GKg0jEkZU6M/Gm3Zwx+muRaJIhlDmxrBtJjpyBmKY9zKvZTHWRIyY3NhRcwZmuKnP36/23M63vc5133HebjP/X1tZzvXda6H7/nunOt3rt91ru8vhFAHdWsoJA0ClgOXAWOBOZLGViy2ADhoezTwE+BHtbZ70uFXGfnpSanK6+LFsGNHuiYRQgihLup5RnEBsMf2Xtv/AX4FVN6qeTnQdavhOmCqatwgMeLQAY4ZeSZs3ZouVpeV4Q0hhHD01fMaxenAX8um9wOf6GkZ229JOgR8APh7+UKSFgIL8+Qbgzo7d9JD4awBZjgVuRrAIhclkYuSyEVJnwt61bOh6O7MoPJ/h0WWwfZKYCWApM6+3obebiIXJZGLkshFSeSiRFJnX9etZ9fTfuCMsukRQGXx9HeWkXQsMAz4Zx1jCiGE0Ev1bCi2AmMknSVpMDAbWF+xzHrg6vx8BvAHx91OIYTQUurW9ZSvOVwHbAAGAats75L0fdIg3+uBXwC/lLSHdCYxu8CmV9Yr5n4oclESuSiJXJRELkr6nIt+V2Y8hBBCY/W7Wk8hhBAaKxqKEEIIVbVsQ1Gv8h/9UYFc3Chpt6Ttkh6XNLIZcTZCrVyULTdDkiW17V8ji+RC0sz82dgl6f5Gx9goBb4jZ0raKOnZ/D2Z1ow4603SKkkHJO3s4XVJWpbztF3SxEIbtt1yD9LF778AZwODgW3A2IplvgasyM9nA2uaHXcTczEFGJqfLxrIucjLnQhsAjqASc2Ou4mfizHAs8D78/QpzY67iblYCSzKz8cC+5odd51y8RlgIrCzh9enAY+S7mGbDDxdZLutekZRl/If/VTNXNjeaLtrPMgO0j0r7ajI5wLgB8AtwOuNDK7BiuTiq8By2wcBbB9ocIyNUiQXBrrq/Qzj3fd0tQXbm6h+L9rlwL1OOoCTJJ1aa7ut2lB0V/7j9J6Wsf0W0FX+o90UyUW5BaRfDO2oZi4kfRQ4w/bDjQysCYp8Ls4BzpG0WVKHpEsbFl1jFcnFTcA8SfuBR4BvNCa0ltPb4wnQuuNRHLXyH22g8PuUNA+YBFxU14iap2ouJB1DqkI8v1EBNVGRz8WxpO6nz5LOMp+SNM72v+ocW6MVycUcYLXt2yR9knT/1jjbb9c/vJbSp+Nmq55RRPmPkiK5QNLFwFJguu03GhRbo9XKxYnAOOAJSftIfbDr2/SCdtHvyG9tv2n7ReBPpIaj3RTJxQLg1wC2twDHkwoGDjSFjieVWrWhiPIfJTVzkbtb7iI1Eu3aDw01cmH7kO3htkfZHkW6XjPddp+LobWwIt+R35D+6ICk4aSuqL0NjbIxiuTiJWAqgKQPkxqKgTjm7Xrgqvzvp8nAIdsv11qpJbueXL/yH/1OwVzcCpwArM3X81+yPb1pQddJwVwMCAVzsQG4RNJu4AjwLdv/aF7U9VEwF98Efi7pBlJXy/x2/GEp6QFSV+PwfD3mu8BxALZXkK7PTAP2AP8GvlJou22YqxBCCEdRq3Y9hRBCaBHRUIQQQqgqGooQQghVRUMRQgihqmgoQgghVBUNRWg5ko5Ieq7sMarKsqN6qpTZy30+kauPbsslL87twzaulXRVfj5f0mllr90taexRjnOrpAkF1rle0tD/d99h4IqGIrSiw7YnlD32NWi/c22fTyo2eWtvV7a9wva9eXI+cFrZa9fY3n1UoizFeQfF4rweiIYi9Fk0FKFfyGcOT0n6Y358qptlzpP0TD4L2S5pTJ4/r2z+XZIG1djdJmB0XndqHsNgR671/548/2aVxgD5cZ53k6QlkmaQam7dl/c5JJ8JTJK0SNItZTHPl/SzPsa5hbKCbpLulNSpNPbE9/K8xaQGa6OkjXneJZK25DyulXRCjf2EAS4aitCKhpR1Oz2U5x0APm97IjALWNbNetcCP7U9gXSg3p/LNcwCLszzjwBza+z/i8AOSccDq4FZtj9CqmSwSNLJwJeB82yPB35YvrLtdUAn6Zf/BNuHy15eB1xRNj0LWNPHOC8llenostT2JGA8cJGk8baXkWr5TLE9JZfy+A5wcc5lJ3Bjjf2EAa4lS3iEAe9wPliWOw64PffJHyHVLaq0BVgqaQTwoO0XJE0FPgZszeVNhpAane7cJ+kwsI9Uhvpc4EXbf86v3wN8HbidNNbF3ZJ+BxQuaW77FUl7c52dF/I+Nuft9ibO95LKVZSPUDZT0kLS9/pU0gA92yvWnZznb877GUzKWwg9ioYi9Bc3AH8DziedCb9rUCLb90t6GvgCsEHSNaSyyvfY/naBfcwtLyAoqdvxTXJtoQtIReZmA9cBn+vFe1kDzASeBx6ybaWjduE4SaO43QwsB66QdBawBPi47YOSVpMK31US8JjtOb2INwxw0fUU+othwMt5/IArSb+m/4eks4G9ubtlPakL5nFghqRT8jInq/iY4s8DoySNztNXAk/mPv1hth8hXSju7p9Hr5LKnnfnQeBLpDES1uR5vYrT9pukLqTJudvqfcBrwCFJHwQu6yGWDuDCrvckaaik7s7OQnhHNBShv7gDuFpSB6nb6bVulpkF7JT0HPAh0pCPu0kH1N9L2g48RuqWqcn266Tqmmsl7QDeBlaQDroP5+09STrbqbQaWNF1MbtiuweB3cBI28/keb2OM1/7uA1YYnsbaXzsXcAqUndWl5XAo5I22n6F9I+sB/J+Oki5CqFHUT02hBBCVXFGEUIIoapoKEIIIVQVDUUIIYSqoqEIIYRQVTQUIYQQqoqGIoQQQlXRUIQQQqjqv70QU1u1Nt+gAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn.metrics import roc_auc_score, roc_curve\n",
    "\n",
    "# Predict the probability of our label score\n",
    "Y_pred_proba = grid_search.predict_proba(X_test)\n",
    "\n",
    "r_roc_auc = roc_auc_score(Y_test, Y_pred_proba[:,1])\n",
    "fpr, tpr, thresholds = roc_curve(Y_test, Y_pred_proba[:,1])\n",
    "plt.figure()\n",
    "plt.plot(fpr, tpr, label='RandomForest (area = %0.2f)' % r_roc_auc)\n",
    "plt.plot([0, 1], [0, 1],'r--')\n",
    "plt.xlim([0.0, 1.0])\n",
    "plt.ylim([0.0, 1.05])\n",
    "plt.xlabel('False Positive Rate')\n",
    "plt.ylabel('True Positive Rate')\n",
    "plt.title('ROC of RandomForest (4-fold CV)')\n",
    "plt.legend(loc=\"lower right\")\n",
    "plt.savefig('images2/Randomforest_ROC')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Extra: StackNet**<br>\n",
    "The following is a visual explanation of the baseline, ERP window selection, and sub-sampling routine (on a single channel, using the figure provided by Alessandro D'Amico):\n",
    "<img src=\"images2/StackNet.png\">"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/anaconda3/lib/python3.7/site-packages/sklearn/externals/joblib/__init__.py:15: FutureWarning: sklearn.externals.joblib is deprecated in 0.21 and will be removed in 0.23. Please import this functionality directly from joblib, which can be installed with: pip install joblib. If this warning is raised when loading pickled models, you may need to re-serialize those models with scikit-learn 0.21+.\n",
      "  warnings.warn(msg, category=FutureWarning)\n",
      "/opt/anaconda3/lib/python3.7/site-packages/lightgbm/__init__.py:48: UserWarning: Starting from version 2.2.1, the library file in distribution wheels for macOS is built by the Apple Clang (Xcode_8.3.3) compiler.\n",
      "This means that in case of installing LightGBM from PyPI via the ``pip install lightgbm`` command, you don't need to install the gcc compiler anymore.\n",
      "Instead of that, you need to install the OpenMP library, which is required for running LightGBM on the system with the Apple Clang compiler.\n",
      "You can install the OpenMP library by the following command: ``brew install libomp``.\n",
      "  \"You can install the OpenMP library by the following command: ``brew install libomp``.\", UserWarning)\n"
     ]
    }
   ],
   "source": [
    "from pystacknet.pystacknet import StackNetClassifier\n",
    "from sklearn.tree import DecisionTreeClassifier\n",
    "from sklearn.ensemble import BaggingClassifier, ExtraTreesClassifier, GradientBoostingClassifier, AdaBoostClassifier\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as LDA\n",
    "from sklearn.model_selection import KFold\n",
    "from lightgbm import LGBMClassifier\n",
    "from xgboost import XGBClassifier\n",
    "from sklearn.gaussian_process import GaussianProcessClassifier\n",
    "from sklearn.neural_network import MLPClassifier, MLPRegressor"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "models = [\n",
    "    #1ST layer #\n",
    "    [GridSearchCV(LDA(solver='lsqr'),  {'shrinkage': (0.01, 0.1, 1.)}, cv=3, n_jobs=-1),\n",
    "     GridSearchCV(SVC(class_weight='balanced', probability=True),  {'kernel': [\n",
    "                  'rbf'], 'gamma': [10**(-2), 10**(-3)], 'C': [100, 10, 1, 0.1]}, cv=3, n_jobs=-1),\n",
    "     GridSearchCV(LogisticRegression(class_weight='balanced',\n",
    "                                     solver='lbfgs', max_iter=100000), {'C': [100, 10, 1, 0.1]}, cv=3, n_jobs=-1),\n",
    "     LGBMClassifier(boosting_type='gbdt', num_leaves=40, max_depth=-1, learning_rate=0.01, n_estimators=1000, class_weight='balanced', subsample_for_bin=1000, objective=\"xentropy\",\n",
    "                    min_split_gain=0.0, min_child_weight=0.01, min_child_samples=10, subsample=0.9, subsample_freq=1, colsample_bytree=0.5, reg_alpha=0.0, reg_lambda=0.0, random_state=1, n_jobs=1),\n",
    "     GaussianProcessClassifier(),\n",
    "     XGBClassifier(max_depth=5, learning_rate=0.3, reg_lambda=0.1, n_estimators=300,\n",
    "                   objective=\"binary:logistic\", n_jobs=1, booster=\"gblinear\", random_state=1, colsample_bytree=0.4),\n",
    "     XGBClassifier(max_depth=20, learning_rate=0.1, n_estimators=300, objective=\"binary:logistic\",\n",
    "                   n_jobs=1, booster=\"gbtree\", random_state=1, colsample_bytree=0.4),\n",
    "     XGBClassifier(max_depth=100, learning_rate=0.1, n_estimators=300, objective=\"rank:pairwise\",\n",
    "                   n_jobs=1, booster=\"gbtree\", random_state=1, colsample_bytree=0.4),\n",
    "     ],\n",
    "\n",
    "    #2ND layer #\n",
    "    [\n",
    "        RandomForestClassifier(max_depth=50, n_estimators=50),\n",
    "        AdaBoostClassifier(base_estimator=DecisionTreeClassifier(\n",
    "            max_depth=5), n_estimators=200),\n",
    "        ExtraTreesClassifier(max_depth=5, n_estimators=50),\n",
    "    ],\n",
    "\n",
    "    #3RD layer #\n",
    "    [\n",
    "        RandomForestClassifier(\n",
    "            max_depth=5, class_weight='balanced', n_estimators=100),\n",
    "        SVC(kernel='linear', class_weight='balanced', probability=True),\n",
    "        MLPClassifier(hidden_layer_sizes=(32, 16), activation=\"relu\", solver=\"sgd\", alpha=0.01,\n",
    "                      batch_size=30, learning_rate=\"adaptive\", learning_rate_init=0.001, power_t=0.5,\n",
    "                      max_iter=100, shuffle=True, random_state=1, tol=0.0001, nesterovs_momentum=True, momentum=0.9, validation_fraction=0.1, early_stopping=True,\n",
    "                      beta_1=0.1, beta_2=0.1, epsilon=0.1)\n",
    "    ]\n",
    "]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "# leave 4 subject out\n",
    "kf = KFold(4)\n",
    "generator = kf.split(X_train, Y_train)\n",
    "\n",
    "# build StackNet\n",
    "stacknet = StackNetClassifier(models, metric=\"auc\", folds=generator, \n",
    "           restacking=False, use_retraining=True, use_proba=True,\n",
    "           n_jobs=-1, verbose=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "# make predictions\n",
    "if not isfile('data/stacknet_prediction.npy'):\n",
    "    stacknet.fit(X_train, Y_train)\n",
    "    y_probs = stacknet.predict_proba(X_test)[:, 1]\n",
    "    np.save('data/stacknet_prediction.npy', y_probs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEWCAYAAAB42tAoAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3hUZfbA8e9JIAmQ0BKQXqRIBxFBwQIiCKxtLWtXxJUVFdgFXVTUVUHF+lN2sWBDFEVFQVRERQEbCIQmvQuhSCe09PP74x1CCCkDZHInM+fzPHmSe+fOvSeXcM+8XVQVY4wxJj8RXgdgjDEmuFmiMMYYUyBLFMYYYwpkicIYY0yBLFEYY4wpkCUKY4wxBbJEYYKSiJwmIj+KyH4RecHDOHqLyM9eXd8XQzMRmXcK7+8kIqtF5ICIXFnIsY+JyPsFvL5BRC4+2VgKuXa0iKwQkaqBOL85eZYozJH//Id9D5JtIjJGRGJzHdNRRH7wPbj3icgXItIs1zHlReQlEdnoO9ca33bCSYTVF9gJlFfVwXnEXEtEPhWRnb54fheR3r7X6omIikipk7iu33xJREXk/lz7k0Sksx/v9zfOYcDzeby/kYikFPRg93kC+J+qxqrqpMLiOhUi0l5EpojIXhHZLSJzROR2EakpIhki0iCP90wUkedVNRV4GxgSyBjNibNEYY64TFVjgTbAmcCDR14QkXOBb4HPgRpAfWAR8IuInO47Jgr4HmgO9ADKAx2BXUD7k4inLrBM8x8R+h6wyXdcPHAr8OdJXOdU7QaGiEj5QJxcRKoDXYC8HvCjgLl+nKYusLQo48qL7+/kB2Am0BD379IP6Kmqm3F/H7fkek9loBfwrm/XB8BtIhId6HiN/yxRmGOo6jbgG1zCOOJZYKyqvqyq+1V1t6o+DMwGHvMdcytQB/irqi5T1SxV3a6qw1R1Sl7X8pVS5vpKBHNFpKNv/xjgNuDfvpJJXlUdZwNjVPWgqmao6gJV/dr32o++73t97z9XRBr4SkS7fKWQcSJSMUcstUXkMxHZ4Tvmf/nE/JyI/CwiFXy7lgOzgH/lc3yEiDwgImt95/3Y93DMM848TtENmK+qKbnOez2wF/fwzZeIrAVOB77wXSNaRGqIyGTfJ/41InJnAe+/RUT+8MU+tKBrAc8B76rqM6q6U51EVf2b7/V3yZUogOuBpar6O4CqJgF7gHMKuZYpRpYozDFEpBbQE1jj2y6LKxl8ksfhH+MeZAAXA1NV9YCf16kMfAWMxH3yfBH4SkTiVbU3MA541lddMi2PU8wGRonI9SJSJ9drF/i+V/S9fxYgwNO4ElFToDa+JCcikcCXwB9APaAmMD5XvBEi8gbQCuiuqvtyvPwI8K8cCSCnAcCVwIW+a+/BlQTyizO3lsDKXLGUx1UnHVcll5uqNgA24isx+qp3PgSSfPFcAzwlIl1zv9dXtfgq7uFeA/fvVCuv6/j+Ts4FJhQQzkQgQUTOy7HvFmBsruOWA60L+91M8bFEYY6YJCL7cdU524H/+PZXxv2dbM3jPVuBI+0P8fkck5+/AKtV9T1fieBDYAVwmZ/vvxb4CfeQXi8iC0Xk7PwOVtU1qvqdqqaq6g5cYrrQ93J73IPwfl8JJUVVczZgl8Y9XCvjHriHcp17Ia5qLq+69X8AQ1U1yfeQfgy45gTaTyoC+3PtGwa8paqb/DxHNhGpDZwHDPH9nguBNzn+kz64JPKlqv7oi/0RICufU1ci/78TAFT1MO4Dx62+WBoBZ+Gqm3Laj/u9TZCwRGGOuFJV44DOQBOOJoA9uIdD9TzeUx3X4AyuLSKvY/JTA/cJPqc/cJ/mC6Wqe1T1AVVtDpwGLMQlO8nreBGpKiLjRWSziCQD73P0d6wN/KGqGflcriFwBfC4qqblc8yjQD8RqZZrf11goq9xdy/u03KmL2Z/7AHicvwebXClt//L62ARWeqrYjogIufncUgNYLeq5kw++d33GrgPDgCo6kHcv3N+ceb3d5LTu8DfRCQGl5ymqur2XMfE4arVTJCwRGGOoaozgTH4etn4Hg6zcJ/gc/sbR+vIpwGXiEg5Py+1BfcQzakOsPkEQ0ZVd/rirYH71J9XA/jTvv2tVLU8cDOuOgrcw7BOAZ/ylwO3A1+LyBn5xLAC+Ax4KNdLm3CNuRVzfMX4Gnf9mbp5MdA4x3ZnXPXYRhHZBtwHXC0i831xNPdVMcWq6k95nG8LUFlE4nLsy+++b8UlUSC7eik+ryB9paxZwNUF/TK+mHbhEu/NHF/tBK5qcFFB5zHFyxKFyctLQDffp1eAB3A9UQaISJyIVBKR4bg66cd9xxzphfSpiDTx1enHi8hDItIrj2tMARqLyI0iUkpErgOa4doKCiUiz4hIC99743C9a9ao6i5gB+7T7ek53hIHHMA1HNcEcnZpnYN7KI4QkXIiEiMinXJez1c19hAwLa8unj6P4xJKzmqT14AnRaSuL+4qInKF77W84sztO6Ct7xM4wGigAa6zQRvf+b8CLingHDl/j03Ar8DTvt+zFXAHrk0otwnApSJynq9X2xMU/Mz4N9BbRO4XkXgAEWktIuNzHTcWeAZ3n77I+YLv36Yyrg3KBAlLFOY4vjr8sbg6aXz19ZcAV+EeqH/gutCep6qrfcek4qpEVuAebsm4B3AC8Fse19gFXIprkN2Fe8hc6isd+KMsrnF0L7AOVzq53HfuQ8CTuO67e0XkHNxDvC2wD/dg/SxHLJm4tpGGuIbfJOC6PGJ+F/ew/EFE6uXx+npcwsxZqnoZmAx862sDmg10KCDO3Of8E9fl9Ioj71HVbUe+cMkvxfdv5q8bcKWSLbh7+B9V/S6Pay8F7sG1IWzFVS8l5XdSVf0VuMj3tU5EduMSW+5eb2NxpZiPfH83Od2I6zmVe7/xkNjCRcYEN1/vo3eB9gWMKynxxI2dWARckEe7hfGQJQpjjDEFsqonY4wxBbJEYYwxpkCWKIwxxhQooLNrBkJCQoLWq1fP6zCMMaZESUxM3KmqVU7mvSUuUdSrV4958056an5jjAlLIpJ7JgS/WdWTMcaYAlmiMMYYUyBLFMYYYwpkicIYY0yBLFEYY4wpkCUKY4wxBQpYohCRt0Vku4gsyed1EZGRvjV7F4tI20DFYowx5uQFskQxBuhRwOs9gUa+r764tXmNMcYUEVUlIzOLlPTMUzpPwAbcqeqPec3Zn8MVwFjftMmzRaSiiFRX1RNZd9kYY0JCcko6B1MzyMhU0jOzOJSWyY4DqagqaRnKhl0HiRQhPSuLjEzlz+QUZq7aQXy5KDKylMws9761Ow5SKkIQgfSMLC5ZNYtLVs86pdi8HJldkxzr8eIWRKlJHouzi0hfXKmDOnXqFEtwxhhzMlLSM9m0+xCpGVls25fC4fRMFmzcS2ZWFgALN+0lIkJYsHEvZaMiyVIlJT3rpK+3ee9hujapSmSEUCoighY1K7DnUDodI/fT/ZVhnD5nJjtPb3JKv5OXiULy2Jfn4hiqOhq3Uhbt2rWzBTSMMZ5SVZL2HOaPXYf4cvEWklPSWbYlmb2H09l7KD3f91UsW5rMTOVweiYXN61KlkLDqrGUjhQOpWVSL74cZaMiKR0ZQalIQRXqVC5L6cgISkcK8bHRlI2KpFSEEBkhiOT1GAVUoV07WLkSXniBhAEDoHTpk/59vUwUSeRYuB2ohVua0RhjAk5V+XnNTnYeSCU9U8nIVHYfTGXFtv3sOpBGTOkIMhWyspSMrCyysmDOht1El4ogNeP4EkCTanGUy1SubluLslGRtKhZgahSEZSPKUW1CmWoGhdN6cgAdzT99Vdo2RLi4uDNNyEhAWrXLvx9hfAyUUwG7vUtvN4B2GftE8aYopaVpUxITOLbZdtYuGkvOw+k+fW+BlXKERtT2n16F/cJvmODeLJUaV6jAmkZWXRqGE9CbDQta1UgulRkgH+TAuzaBQ884JLDf/4Djz0GZ55ZZKcPWKIQkQ+BzkCCiCQB/wFKA6jqa7gF13sBa4BDwO2BisUYEz5S0jN5dcZavvp9K7sPprH74NHEEBUZQb34slzepiYCHEzNoGvT06hVqQylIl0df2x0KcpEefjQPxGqMHYs3Hcf7NkD99/vvopYIHs93VDI6wrcE6jrG2NCl/oagDftOcTyrck8O3UlANuSU8jMOtqMWS4qkts71aNM6Uj+cUEDKpQ9+Xr6oDRkCDz3HHTsCK+95qqdAqDErUdhjAk/63YcYNKCzUxITCIyUti0+3Cex51zemXqJ5QjITaae7o0JKZ0CSkZnIjDh+HgQdf+cMcd0KiR+x4RuPYPSxTGmKCgqsxau4vZ63ZxMC2TQ2mZfDY/CeCYxuMWNcvT4ax4ykVFUiUumhY1K1CrUhkaVo3zKvTiM3Uq3HMPtGkDn34KZ5zhvgLMEoUxptipKn8mp7J8azLLtiazafchxs/ddMwx5WNKEVUqgsNpmdx6bl3Oa5hAx4YJxEaH4WNryxb45z/hk09cYrj33mK9fBjecWNMcdq2L4WFm/ay+2Aa3yzdRnpmFiu27T+mkblKXDT14suSlpHFqJva0ryG61pqgO+/h7/+FdLSYNgw11gdHV2sIViiMMackqwsZeWf+0lJzyQ9U5m/cQ97D6Xz5k/ryMjKe3zsde1q06xGeZpWL0+T6nGUjwmxRuaikJ7uBsm1bg29esHw4dCwoSehWKIwxpywZVuS+WbpNt6b/ccxJYO8DLioIbUrl+WsupUoX6Y08eWi8h9RbCA5GR55BH77DX75xTVajx/vaUiWKIwxBdq+P4WZK3cwdck2lm9NJtPXvnBE61oVaFq9PJe0qEZUZASlIoQm1cqHXlfUQFOFCRNg4EDYtg3uvhtSU6FsWa8js0RhjHHVR/M37uGdXzcQHRnB5r2H2XUwjTXbDxxzXHy5KKrERdOmdkXOa5jADe3rUCrQ01KEgx074Lbb4Ouv3Yjqzz+Hs8/2OqpsliiMCTM79qeyZMs+UtOz2Hsojc/mb2bOht3HHNOgSjkOpWVyYeMqlIoQLm52Gl2bVqVqXIxHUYe48uVh50546SXX/bVUcD2agysaY0xApKRn8vPqnQyd9Psx1UY5XXVmTbo3P41LmlezNoTi8OOP8OSTbjxEbCzMnh3QQXOnwhKFMSFEVdmfmsG6HQd555f1LNi4l9joUizbmpx9TM2KZbjz/Po0rhZHxTJRVCxbmuoVYiw5FJedO10X1zFjoF492LABWrQI2iQBliiMKZEOpWWwfGsya3ccJCU9k237UliyJZkfV+047tiWNStwaavqVImL5pqzatG8RgUPIjaowjvvuCSRnAwPPggPPxwUjdWFsURhTAlwIDWDldv2s3DTXoZ9uSzf4+rGl6VptfK0r1+Z1rUr0LZOJSspBJP334dmzdwEfs2bex2N3yxRGBOEtu9PYcXW/fyx+xBPfrUsz6Uy/3lxIy5oXIWEctFUjo0iplSE9UAKNocOwVNPwV13Qa1arj2iQoWgrmbKiyUKY4JIemYWl/33Z1Zs23/ca/d0acAFjapQo2IZalcO/uqKsDdliuvBtGED1KwJ/fpBpUpeR3VSLFEYEySWbtnHX0b+nL19b5eG9GhRjdPKx1Alrnjn9jGnICnJTeD36afQtCnMnAkXXOB1VKfEEoUxHlJVJi/awrNTV7J5r1tjoUaFGGb+u0vg11c2gfHkk/DVV67KafBgiIryOqJTZonCGI/MWLmd3u/MPWbf45c359Zz61oDdEkzZw6UKeNWmBs+3PVsOv10r6MqMpYojCkmWVnKpIWbGf7V8uMm0pt0Tyfa1K7oUWTmpO3bBw89BK++CpdeCpMnQ3y8+wohliiMCbCfVu9gxsodjJ21gfRMN+32mXUqUj++HL071aNVLUsQJY4qfPQR/OtfsH079O/v1ooIUZYojAmQ9Mwsznvmh2OmzKgbX5a3bmsXHst2hrL334dbb4V27eDLL+Gss7yOKKAsURgTAKpKo6FfZ29/2q8jZ9UtmV0jjU9qKqxb53oy/e1vkJHhkkVkpNeRBZwlCmOK2IadB/lwzsbs7eVP9KBMVOg/TELa9OluHMShQ7B6tVuK9PbbvY6q2FiiMKYIpGZk8uFvG/ll7S6+W/Zn9v7p93W2JFGSbd8O990H773nejGNHl3s61UHA0sUxpyi7ftTaP/k99nbFcuW5rZz69G9+WnUTyjnYWTmlKxZA+3bw4EDMHSo+ypTxuuoPGGJwhg/pWZksmDjXpZuSea7ZdtIz3SrwqkePWbBI92oVK7kD7AKa8nJbiGhBg3gjjugTx/XLhHGLFEY4wdV5YyHpx63/+KmpyECFzauwg3t6xAZYQPlSqyDB+GJJ+CNN2DxYjeJ33PPeR1VULBEYUwBNu46xLbkFP72+qzsfZPv7eTmX4qNJsISQ2j44gu4917YuNGVIkrAGhHFyRKFMXkY99sfDJ245Lj9ix/rTvmY0h5EZAIiI8N1dZ040a0P8dNPcN55XkcVdCxRGOOzaNNeJi/awsQFm7On2Oje7DSuPqsW1crH0KpWBZuDKVSoggiUKgXVq8OIEW6UdQhM4BcIlihMWMrIzGLGyh1M+X0rny3YfNzrMaUjeOm6M+nRopoH0ZmAmj3brRPxxhvQti2MGuV1REHPEoUJeapK4h97WLolmU/nJ7E4ad9xx7SvV5mGp8XSrdlpnNcwwab4DkV79rgJ/F5/HWrUcNvGLwFNFCLSA3gZiATeVNURuV6vA7wLVPQd84CqTglkTCa8HE7LpOmjx/dWOq9hAmfWqcgVbWpQPyHWeiuFuo8+ggEDYOdOt6jQ449DnM235a+AJQoRiQRGAd2AJGCuiExW1Zwrwz8MfKyqr4pIM2AKUC9QMZnw8sWiLfT/cEH29id3ncsZ1eKsMTocrVgB9erB1Klw5pleR1PiBLJE0R5Yo6rrAERkPHAFkDNRKFDe93MFYEsA4zEhTlWZtW4XP67ayYTETew84Bqk7zy/PoO7n0FMaZtKI2ykpMAzz7g2iMsuc1VODz8cFhP4BUIgE0VNYFOO7SSgQ65jHgO+FZH+QDng4rxOJCJ9gb4AderUKfJATWj4bP5mBn+y6Jh9w69swc3n1PUoIuOJadPg7rvd5H2DB7tEUdpKkacikIkir0pfzbV9AzBGVV8QkXOB90SkhapmHfMm1dHAaIB27drlPocxAMxZvxuAd/u057yGCdbuEG7+/BMGDYIPPoCGDeHbb6FbN6+jCgmBTBRJQO0c27U4vmrpDqAHgKrOEpEYIAHYHsC4TAmnqhxIzWDhpr18v3w7K7YlIwhrdhwgITaKCxtX8TpE44XvvoMJE+DRR+HBByEmxuuIQkYgE8VcoJGI1Ac2A9cDN+Y6ZiPQFRgjIk2BGGBHAGMyJdjB1AyenbqCd2f9kefr5zdK4JzTQ2utYlOIRYtcFdM118BNN0GnTlC/vtdRhZyAJQpVzRCRe4FvcF1f31bVpSLyBDBPVScDg4E3RORfuGqp3qpqVUvmGKkZmcxcuYO+7yVm7zuvYQLdmp1Go9Ni6dggwcPojCcOHID//Adeftn1ZrrySjfK2pJEQAR0HIVvTMSUXPsezfHzMqBTIGMwJdvaHQfo+sLMY/YtfLQbFcvaVAtha9Ik6N8fkpKgb194+mmXJEzA2N01Qeun1Tu45a05AFSNi+bzeztRvUJ4LhxjfH7/Hf76V2jZ0g2i69jR64jCgs1TYIJSSnomD038HYAb2tdh1oNdLUmEq/R0+OEH93PLlvDVV5CYaEmiGFmiMEFnzvrdNHlkKpt2H6Z9vco8fVVL6+oarn79Fc46y3VzXbPG7evVy8ZFFDNLFCZoHEzNYNT0NdmLBLWsWYFRN7X1OCrjid27XftDp06wdy989pkbG2E8YW0UxlP7U9IZ99tG9h1O59UZa7P3D+7WmP5dG3kYmfFMSgq0aQNbtriR1Y89BrGxXkcV1ixRGM+s33mQLs/POGZfo6qxvHdHB6pVsMFSYScpya1THRMDw4a5ZNG6tddRGSxRGI+8N/sPHpnklhpNiI1m9oMXUcrWgAhPhw+7Lq7PPONGVl92Gdx2m9dRmRz8ShQiEgXUUdU1AY7HhKiU9Ex+WbOTNdsPMHfDHqYt/xOAvheczoM9m9gSo+Hq22/dBH5r18LNN0P79l5HZPJQaKIQkb8ALwJRQH0RaQP8R1X/GujgTMm3PTmFa1+fxR+7Dh332v2XnME9XayBMmz17w//+x80auRmfO3a1euITD78KVE8gZsefDqAqi4UEfvfbQqlqrR/6vvs7Qd6NqFrk6rUiS9LdClbFyAsZWa675GRcM45kJAAQ4bYBH5Bzp9Eka6qe3NVDdh8TKZQs9btAqBWpTL8MLgzUaWsDSKszZ8Pd90Ft9ziShM33eR1RMZP/vzPXS4ifwMiRKS+iLwEzA5wXKaE274/hRvf+A2Al65rY0kinO3fD//6F5x9NmzcCNWrex2ROUH+/O+9FzgLyAI+A1KAgYEMypRskxZspv2TrsopulQE7epV9jgi45lvv4WmTd0sr//4h1u7+pprvI7KnCB/qp4uUdUhwJAjO0TkKlzSMCbb9uQU7nh3Hr9v3gfA0F5Nua1jPW+DMt6KioKqVeHTT6FD7pWQTUkhhS3/ICLzVbVtrn2JqnpWQCPLR7t27XTevHleXNrkIyMzi/4fLuDrJduy9716U1t6trQqhrCTng4vvgjJyfDkk25fVhZEWNWj13zP7XYn8958SxQicglumdKaIvJijpfK46qhjAHgoYm/ZyeJ+y85gzvOq09MaevVFHZ+/tk1Vi9dCtdeezRBWJIo8QqqetoOLMG1SSzNsX8/8EAggzLBLyU9kzd/Wsfz367K3jdt0AU0rBrnYVTGE7t2uS6ub70FderAF1/ApZd6HZUpQvkmClVdACwQkXGqmlKMMZkgl5KeyTlPf8/eQ+kA9GpZjUcvbW7zM4WrXbtg/Hj497/h0UehXDmvIzJFzJ/G7Joi8iTQDMh+Eqhq44BFZYJWZpby4nerspPEgke6UamcLUsadpYvh48/dutWN27sur1Wtt5tocqfRDEGGA48D/QEbsfaKMKKqvLY5KVMXbqNP5NTs/eve6oXEbagUHg5dMg1Uj/3nJv6+4473IyvliRCmj+JoqyqfiMiz6vqWuBhEfkp0IGZ4LD3UBptnvgue7t3x3qUjhQubnqaJYlwM3Wqm8Bv/Xo3u+tzz0GVKl5HZYqBP4kiVdz8HWtF5C5gM1A1sGGZYHA4LfOYJDHv4YtJiI32MCLjmQMH3NQb8fEwfTp07ux1RKYY+ZMo/gXEAgOAJ4EKQJ9ABmW8t3TLPv4y8ufs7fVP97KpwMNNZiZ8+CHccIOrZpo2DZo0gWj7sBBuCk0Uqvqb78f9wC0AIlIrkEEZb933ySImJCYBEFUqghVP9LAkEW4SE92UG4mJUKYMXH21rTYXxgocCSMiZ4vIlSKS4NtuLiJjsUkBQ9aGnQezk8TYPu1ZOayHtUWEk337YMAAt4DQ5s2u2+tVV3kdlfFYQSOznwauBhbhGrAn4iYDfAa4q3jCM8Vl3+F0/jl+AdNX7gDgrgsbcEFja6gMO1dfDT/8APfcA8OHQ4UKXkdkgkBBVU9XAK1V9bCIVAa2+LZXFk9oprjMWLmd3u/Mzd5+p/fZdGli/RXCxrp1rvdSXJzr+hoR4aYEN8anoESRoqqHAVR1t4issCQRWsbP2cgDn/2evd2kWhxj+7SnankbYR0W0tLg+edh2DBX3fTMMzbDq8lTQYnidBE5MpW4APVybKOqVnFZQiWnpHPHmLnM3bAHgPb1KvP01S1pUCXW48hMsfnxRzeB3/Llbn2IAQO8jsgEsYISxdW5tv8XyEBM8Tn/mensO+ym4Hjlprb0sunAw8v//R8MGgT16sFXX0GvXl5HZIJcQZMCfl+cgZjicevbc7KTxNqnehFpPZrCQ1YWHDzo2iH+8hfYsQMefhjKlvU6MlMC+DPgzoSIZo9O5VBaJgCzH+xqSSJcLF3qqpmOrDTXuDE89ZTXUZkSJKAriohIDxFZKSJrRCTPNSxE5G8iskxElorIB4GMJ1xt25dCvQe+yk4SP/27i00JHg4OHYIHH4Q2bVxbxKWXQiErWhqTF79LFCISraqphR+ZfXwkMAroBiQBc0Vksqouy3FMI+BBoJOq7hER65MZAFeO+iX7598f605cTGkPozHFYsECN1Buwwa4/XZ49llISPA6KlNCFVqiEJH2IvI7sNq33VpE/uvHudsDa1R1naqmAeNxYzNyuhMYpap7AFR1+wlFbwr1Z3IK25JTiIwQNoz4iyWJUHekxFCnjvuaORPeftuShDkl/lQ9jQQuBXYBqOoioIsf76sJbMqxneTbl1NjoLGI/CIis0Wkhx/nNX6at2E3HZ5yfRKuP7u2x9GYgMrIgJdegq5d3WR+8fEuSVxwgdeRmRDgT9VThKr+kWtSuEw/3pdXS2nuCtJSQCOgM1AL+ElEWqjq3mNOJNIX6AtQp04dPy4dvrKylKGTlvDhnI3Z+64/uzZPXNHCw6hMQM2Z4xqrFyyAnj0hORkqVfI6KhNC/EkUm0SkPaC+dof+wCo/3pcE5PwYWws3DUjuY2arajqwXkRW4hLH3JwHqepoYDRAu3btrDUuH7kXGbq4aVUua12DK9rkLsiZkHDgAAwZAq++CtWrwyefuLmabKZfU8T8SRT9cNVPdYA/gWm+fYWZCzQSkfq4xY6uB27Mdcwk4AZgjG+G2sbAOv9CN7mdNXwaAKUjhdkPdiXeFhkKbaVLw4wZ0L+/m4ajfHmvIzIhyp9EkaGq15/oiVU1Q0TuBb4BIoG3VXWpiDwBzFPVyb7XuovIMlx11v2quutEr2Xg1zU7ycxyha0Vw3raGIlQtWYNPPEEjBrlBs8lJkKMdXU2gSVaSL9qEVkLrAQ+Aj5T1f3FEVh+2rVrp/PmzfMyhKBzKC2DZo9+A8Bbt7Wja9PTPI7IFLnUVNfF9cknISrKTb1x/vleR2VKEBFJVNV2J/PeQns9qWoDYDhwFvC7iEwSkRMuYZjAufC5GQDUrlzGkkQomj7drS736KNw5ZWwYoUlCVOs/KeH1RAAACAASURBVBqZraq/quoAoC2QDIwLaFTGb2kZWezY78ZBfvNP6woZclRdKSI9HaZOdSvO1ajhdVQmzBTaRiEisbiBctcDTYHPgY4Bjsv46flv3RIht51bl7JRNnVXSMjKgrfegh49oHZteO89qFjRrV1tjAf8KVEsAc4BnlXVhqo6WFV/C3Bcxg9/Jqcw+kfXSezfPZp4HI0pEosXw3nnQd++8Oabbl/16pYkjKf8+Qh6uqpmBTwSc8KOjLo+v1EC5aKtNFGiHTgAjz/u1oqoVAnGjIFbb/U6KmOAAhKFiLygqoOBT0XkuK5RtsKdd77+fWt2lVONCjG8d4ctX1niPfYYvPAC/P3vMGKEm4LDmCBR0MfQj3zfbWW7ILLvcDr9xs3P3n6r99keRmNOyaZNbjGhJk3ggQdcj6bzzvM6KmOOU9AKd3N8PzZV1WOShW8gna2AV4x27E/lfz+s5t1ZfwDQ/6KGDO5+hsdRmZOSkQEjR7rurmed5SbvS0iwJGGClj8V2304vlRxRx77TABkZGZx3ejZJP6xJ3tfxbKl+efFjT2Mypy02bPdBH6LFrklSf9n/41M8CuojeI6XJfY+iLyWY6X4oC9eb/LFLUP5mzMThJ9OtXn0cuaeRyROWlffQWXXebGQXz2matqsgn8TAlQUIliDm4Nilq4leqO2A8sCGRQ5qhHP18KwLRBF9KwaqzH0ZgTpgpbtkDNmnDxxW6epoED3TxNxpQQBbVRrAfW42aLNcVsQmIS/56wCICqcdGWJEqiVavg7rvd92XLIDYWHn7Y66iMOWH5DrgTkZm+73tEZHeOrz0isrv4Qgw/nyYmcd8ni/BNBsvn93byNiBzYlJSXHfXli1h3jx48EEbMGdKtIKqno4sd2qL7RajRZv2MvgTV5KYMuB8mtWwNQZKlG3b3PKjq1fDDTfAiy9CtWpeR2XMKcm3RJFjNHZtIFJVM4FzgX8A5YohtrCzY38qV4z6BYAz61S0JFGSpKe776ed5hLFt9/CBx9YkjAhwZ+5nibhlkFtAIzFTQz4QUCjCkMHUjN4espyAK5rV5uJd1t1U4mQlQWvvQYNGkBSkuvF9Oab0K2b15EZU2T8GUeRparpInIV8JKqjhQR6/VUxHq89CNJew4TXy6Kp69q6XU4xh+LFsE//gG//QYXXXS0VGFMiPGnRJEhItcCtwBf+vaVDlxI4efzhZtJ2nMYgJ+HXESELWMa3FThvvvcqOp169w04NOmQf36XkdmTED4kyj64Bq2n1XVdSJSH/gwsGGFj9dmrmXg+IUAfPyPcykTFelxRKZQIrBnD9xxB6xcCTffbAPnTEgrtOpJVZeIyACgoYg0Adao6pOBDy20bU9Oof1TR6fL6t2xHu3rV/YwIlOgP/5wA+UefRTatoU33oAIvxaINKbE82eFu/OB94DNgADVROQWVf0l0MGFsqS9rqqpRoUYJvTrSI2K1s8+KKWnuzUiHn/cbV93nUsUliRMGPGnMfv/gF6qugxARJriEke7QAYW6n5evROAp65qaUkiWP36q2usXrIErrjCzfhap47XURlT7PxJFFFHkgSAqi4XkagAxhTyZq/bxYvfrQKgTe2KHkdj8jVtGuzbB5MmuURhTJjyp/w8X0ReF5HzfF+vYpMCnrRZa3dx/ejZAFzUpCoVy1rODRqqMHYsfP212x4yxM3RZEnChDl/EsVdwFrg38AQYB1udLY5QQ9+tpgb3nBJ4tZz6/K2rU4XPFascGMhbrsN3nnH7YuOdhP5GRPmCqx6EpGWQANgoqo+WzwhhZ6U9EwufG46fyanAnBvl4bcd4mtThcUDh+Gp56CZ56BcuXg9dfdutXGmGwFLVz0EG4lu/nA2SLyhKq+XWyRhYjf1u3iOl9VE8BP/+5C7cplPYzIHOOLL2D4cDcW4vnn3VxNxphjFFSiuAlopaoHRaQKMAWwROGng6kZDBy/kGnL/wSgff3KjO3TnpjSNqDOc9u2wcKF0KMHXHst1KsH7dt7HZUxQaugRJGqqgcBVHWHiFjH8RNw1/uJ/OTrAju4W2P6d23kcUSGzExXtfTggxAVBRs3unUiLEkYU6CCEsXpOdbKFqBBzrWzVfWqgEZWgj02eWl2klg1vCdRpSzHem7+fLjrLpg71y1J+sortpiQMX4qKFFcnWv7f4EMJFSs2b6fMb9uAODl69tYkggG69e7UkNCglsj4vrrbW4mY05AQWtmf5/fayZ/r0xfC8CwK5pzRZuaHkcTxlTh99+hVSs3q+s778Bll0FFG+BozImyj7tFaMPOg3y2YDMAV5xpScIz69fDpZfCmWfC4sVu3y23WJIw5iQFNFGISA8RWSkia0TkgQKOu0ZEVERK9PxRA8a7Aev9L2pI+RhbsqPYpaXBiBHQvDnMnOm6uzZr5nVUxpR4/sz1BICIRKtq6gkcHwmMAroBScBcEZmcc94o33FxwADgN3/PHazKx5QmLqYUg7vbYLpil5kJHTtCYiJcdRW89BLUru11VMaEhEJLFCLSXkR+B1b7tluLyH/9OHd73NoV61Q1DRgP5DVpzjDgWSDF/7CDV+PT4rwOIbwkJ7vvkZHQp48bQPfpp5YkjClC/lQ9jQQuBXYBqOoi3Ip3hakJbMqxneTbl01EzgRqq+qXFEBE+orIPBGZt2PHDj8ubUKeKowZA6efDp9/7vbdfbdrmzDGFCl/EkWEqv6Ra1+mH+/Lq/+hZr/oBvD9HzC4sBOp6mhVbaeq7apUqeLHpYtfVpby85qdqGrhB5tTs2wZdO4Mt98OTZpAgwZeR2RMSPMnUWwSkfaAikikiPwTWOXH+5KAnOX/WsCWHNtxQAtghohsAM4BJpfUBu3HvlgKQOlI60gWUM8+C61bu8WE3nwTfvwRWrTwOipjQpo/T7V+wCCgDvAn7oHez4/3zQUaiUh930JH1wOTj7yoqvtUNUFV66lqPWA2cLmqzjvB38Fzh9IyGDvLFbpe+Ftrj6MJUUdKatWqwU03uWnB77jDliQ1phgU2utJVbfjHvInRFUzRORe4BsgEnhbVZeKyBPAPFWdXPAZSo6hE5cAcOf59alVyWaGLVJbtsDAgXD++TBgANx6q/syxhSbQhOFiLxBjraFI1S1b2HvVdUpuFlnc+57NJ9jOxd2vmCjqpz/7HSS9hwG4O7ODT2OKIRkZrr5mIYOhfR01/XVGOMJf8ZRTMvxcwzwV47tzRS2tu9PzU4SUwacT6VytqxpkVi40C0elJgI3bu7hGEN1sZ4xp+qp49ybovIe8B3AYuoBHniCzd28LlrWtGsRnmPowkh+/a5KqePPnLrRdgEfsZ4yu+R2TnUB+oWdSAlTXpmFl/9vhWA7s2qeRxNCacKn3wCq1e7qqYLL4R16yAmxuvIjDH4NzJ7j4js9n3txZUmHgp8aMErJT2T/h+4eZ16NK9GhbI2r9NJW7sWevWC665zA+fS091+SxLGBI0CSxQiIkBrYLNvV5aG8YiyP5NTeGrKcj5feHQ4yD1drAH7pKSmukn7hg+H0qXh5ZfdyOpSJ1PINcYEUoH/K1VVRWSiqp5VXAEFq+SUdDo8dXSJjta1KjDm9vbWgH2yNm2CYcPcGhEvvQQ1bVp2Y4KVPx/f5ohIW1WdH/BogtRHczcy5NPfAWhduyIT+3UkIsIaWE/Yjh2ugfree6FhQzcVx+mnex2VMaYQ+bZRiMiRJHIeLlmsFJH5IrJARMImaajqMUnik3+ca0niRGVlwVtvuXmZBg2ClSvdfksSxpQIBZUo5gBtgSuLKZagk5aRxZBP3Qpp5zdK4L07OngcUQm0ZAn06wc//+xGV7/2Gpxh63UYU5IUlCgEQFXXFlMsQWXT7kOc/+z07O0RV7fyMJoSKi3NDZhLS4O334bevW1MhDElUEGJooqIDMrvRVV9MQDxBIUjU3MAxJSOYOb9XTitvHXX9NsPP7ixEFFR8PHHrsopIcHrqIwxJ6mgcRSRQCxuOvC8vkLWyO/XAFCrUhlWDOtpScJfSUlw9dXQtSuMHev2nXeeJQljSriCShRbVfWJYoskSNz+zhymr3Sr6I37u7VJ+CUjA/73P3jkETeZ39NPu6nAjTEhodA2inAxdck27no/MXt74t0dqRtfzsOISpBbboHx46FnTxg1CurX9zoiY0wRKihRdC22KILAw5PcmhLt6lbi2WtacXqVWI8jCnJ797pR1LGxcM89rsrp6qutsdqYEJRvG4Wq7i7OQLz0Z3IKOw+kEhddign9OlqSKIiqKz00beqqmsC1Q1xzjSUJY0KUrSMJXPf6LADu72H9+wu0Zg1ccgnccAPUqgU33+x1RMaYYhD2iWLGyu1s2HUIgGvPqu1xNEHsgw+gRQv47TfXcD17NpwV9lOAGRMWwn6qznG/bQTgwzvPoUxUpMfRBKH0dDe7a7t2rnrp2WehRg2vozLGFKOwL1EsTtpLk2pxnNsg3utQgsv27a4303XXue3GjeH99y1JGBOGwjpRvPnTOv5MTmXngVSvQwkeWVkwerSbj+mjj6B5czc2whgTtsK26ikjM4vhXy0HYOQNZ3ocTZBYt841UM+aBZ07w6uvuuk3jDFhLWwTxUMT3dThV7apQccGNsUEABUquPER777rqp2su6sxhjBMFDv2p9JpxA+kZWYBcFfnBh5H5LHJk2HMGPjkE4iPd9OCR4R1jaQxJpeweyJMXbI1O0l8P/hCmlQr73FEHtm4Ea68Eq64Alatgq1b3X5LEsaYXMLuqfDI50sBmDv0YhqE4wjsjAx4/nk3svrbb+GZZ2DBAjeAzhhj8hBWVU+z1u7K/jkhNsrDSDyUmQlvvgkXXQT//S/Uq+d1RMaYIBdWJYrRP7rF+j7+x7lIODXU7tkDQ4bA/v0QHQ2//OLaJixJGGP8EDaJYs763dnrTLStU9HjaIqJKowb57q4vvACTPct7Rofbz2ajDF+C/mqJ1XlylG/sChpHwAP9WpCqcgwyI+rVsHdd8P330P79vDNN9CmjddRGWNKoJBPFC98uyo7Sbz4t9Zc1TZMGm3/+U+YNw9eeQX69oVIm8fKGHNyQjpRrN95kP9Nd+tfL338EspFh/SvC99956qZatd2o6qjo6FaNa+jMsaUcAGtgxGRHiKyUkTWiMgDebw+SESWichiEfleROoW1bX3HU6ny/MzALj+7NqhnSS2bYMbb4Tu3V13V4C6dS1JGGOKRMAShYhEAqOAnkAz4AYRaZbrsAVAO1VtBUwAni2q629PTgHgrLqVGHF1q6I6bXDJyoLXXnOliE8/hf/8x42RMMaYIhTIEkV7YI2qrlPVNGA8cEXOA1R1uqoe8m3OBoqsAWGHb0bY2zvVK6pTBp+nn4Z+/dwCQosXw2OPQUyM11EZY0JMIOtjagKbcmwnAR0KOP4O4Ou8XhCRvkBfgDp16vh18Z9X7wSgWvkQe3Du3w87d0L9+nDXXe77DTdYd1djTMAEskSR15NL8zxQ5GagHfBcXq+r6mhVbaeq7apUqeLXxY90gW1Xr7Jfxwc9VZg4EZo1c4sJqbrxEDfeaEnCGBNQgUwUSUDORahrAVtyHyQiFwNDgctVtchWEPpmybaiOpX3/vgDLr8crroKKleGkSMtORhjik0gq57mAo1EpD6wGbgeuDHnASJyJvA60ENVtxflxY/MEFvizZoFF1/sfn7+eRg4EEqFcA8uY0zQCViJQlUzgHuBb4DlwMequlREnhCRy32HPQfEAp+IyEIRmVxU1xeBS1tVL6rTFb/kZPe9bVvo0weWL4fBgy1JGGOKXUCfOqo6BZiSa9+jOX6+OJDXL5F27YIHHnBTgC9dCrGxbpZXY4zxSBhMelRCqMLYsW5MxDvvuAZra4cwxgSBkK3HSE0vQW0U+/a51eZmzIBzz3WD6FqF6CBBY0yJE3Ilig07D1Lvga/YvPcwUcE+S6z6eguXLw8JCTB6NPz8syUJY0xQCfIn6Ym76tVfs3/+V7fGHkZSiG++cQ3VSUmuiumTT+DOO23NamNM0Ampp9KXi7ew+2Aa9eLLsu6pXtSuXNbrkI63dStcfz306AGHDsH2Iu0VbIwxRS6kEsWCjXsBeOqqlkREBGFD8KhRrrF60iR4/HE3P1Pbtl5HZYwxBQqpxmwBykVF0rFBgteh5C0xETp0cAmjUSOvozHGGL+ETIlCVRk/dxPpmXlOJ+WN5GS30lxiott+5RXXNmFJwhhTgoRMiWJR0j4OpGZ4HYaj6taHGDjQtUnUqeOmArcpwI0xJVDIlCiuHPULAOP+XtBM5sVg/Xq49FK49lqoWtXN1TRokLcxGWPMKQiJRHFkNTuAjg3iPYwEGDcOfvwR/u//YO5c1yZhjDElWImvelJVnpqyHID7LzkD8WLai59+gtRUN8vr/fdD795Qq8gW6zPGGE+V+BLF4qR9TFrolrm4um0xP5x37nQzu15wATzxhNsXHW1JwhgTUkp8ieKLRS5JvHpTW6pVKKbGYlUYM8aVHvbtgyFD4JFHiufaJmSkp6eTlJRESkpK4Qcb46eYmBhq1apF6dKli+ycJT5RHEzLBOD8xv4tkVokpkxxJYlOndwEfi1aFN+1TchISkoiLi6OevXqeVNlakKOqrJr1y6SkpKoX79+kZ23xFc9/bR6B/HlooiNDnDOO3QIfnE9q+jVCz7/3DVaW5IwJyklJYX4+HhLEqbIiAjx8fFFXkot0YniQGoGSXsOs/tQWmAv9PXXLiH07Al797pJ/C6/3CbwM6fMkoQpaoH4myrRT7qJ85MA6Hdhg8BcYPNmNx6iVy/XSP3FF1CxYmCuZYwxQarEJoqsLOWRz5cCcEP7OkV/ge3boVkz+PJLGD4cFi2CCy8s+usY46Enn3yS5s2b06pVK9q0acNvv/0GwEsvvcShQ4dO6pxjxozh3nvvPW7/Y489RtmyZdmeY8bk2NjYQs/31FNP5fuaqnLRRReRfGSN+SCUmJhIy5YtadiwIQMGDED1+GmGnnvuOdq0aUObNm1o0aIFkZGR7N69G4C9e/dyzTXX0KRJE5o2bcqsWbMAuO+++/jhhx+K5XcosYmiz7tzAahYtnTRTie+ebP7XrUqDBsGS5bA0KEQFVV01zAmCMyaNYsvv/yS+fPns3jxYqZNm0bt2rWBU0sUBUlISOCFF144ofcUlCimTJlC69atKV++vN/ny8zMPKHrn6p+/foxevRoVq9ezerVq5k6depxx9x///0sXLiQhQsX8vTTT3PhhRdSuXJlAAYOHEiPHj1YsWIFixYtomnTpgD079+fESNGFMvvUGJ7Pc1YuQOAmfd3KZoT7tsHDz8Mr78Os2e76b8HDCiacxtTiMe/WMqyLUX7qbhZjfL857Lm+b6+detWEhISiI6OBtxDHGDkyJFs2bKFLl26kJCQwPTp0+nXrx9z587l8OHDXHPNNTz++OMAzJ07l4EDB3Lw4EGio6P5/vvvj7nGV199xfDhw/niiy8A6NOnD2PGjGHIkCHZD8Ij3n//fUaOHElaWhodOnTglVdeYejQoRw+fJg2bdrQvHlzxo0bd8x7xo0bR9++fbO3r7zySjZt2kRKSgoDBw7Mfi02NpZBgwbxzTff8MILL1CmTBkGDRrEgQMHSEhIYMyYMVSvXp033niD0aNHk5aWRsOGDXnvvfcoW/bkP4hu3bqV5ORkzj33XABuvfVWJk2aRM+ePfN9z4cffsgNN9wAQHJyMj/++CNjxowBICoqiijfh9a6deuya9cutm3bRrVq1U46Rn+UyBLF3A2uSHZjhzpUKHOKfYVV4eOPoWlTN/33XXdBgwC1eRgTRLp3786mTZto3Lgxd999NzNnzgRgwIAB1KhRg+nTpzN9+nTAVVHNmzePxYsXM3PmTBYvXkxaWhrXXXcdL7/8MosWLWLatGmUKVMm+/wTJ05kxIgRTJkyJTsJxcbG0qdPH15++eVjYlm+fDkfffQRv/zyCwsXLiQyMpJx48YxYsQIypQpw8KFC49LEgC//PILZ511Vvb222+/TWJiIvPmzWPkyJHs2rULgIMHD9KiRQt+++03OnToQP/+/ZkwYQKJiYn06dOHoUOHAnDVVVcxd+7c7E/ub7311nHXnD59enY1Uc6vjh07Hnfs5s2bqZVjAG6tWrXYfKTWIg+HDh1i6tSpXH311QCsW7eOKlWqcPvtt3PmmWfy97//nYMHD2Yf37ZtW3450hszgEpkieK5b1YCcEGjU1x3QhWuusotJNS2LUyeDO3aFUGExpyYgj75B0psbCyJiYn89NNPTJ8+neuuu44RI0bQu3fv4479+OOPGT16NBkZGWzdupVly5YhIlSvXp2zzz4b4Jjqn+nTpzNv3jy+/fbb46qFBgwYQJs2bRg8eHD2vu+//57ExMTscx0+fJiqVasW+jvs3r2buLi47O2RI0cyceJEADZt2sTq1auJj48nMjIy++G7cuVKlixZQrdu3QBXFVW9enUAlixZwsMPP8zevXs5cOAAl1xyyXHX7NKlCwsXLiw0NiDP9oiCeiV98cUXdOrUKbu0lZGRwfz58/nvf/9Lhw4dGDhwICNGjGDYsGEAVK1alS1btvgVy6kokYlizvrdVCxbmh4tqp/cCdLToXRp1831vPPgoovg7rshMrJoAzUmyEVGRtK5c2c6d+5My5Yteffdd49LFOvXr+f5559n7ty5VKpUid69e5OSkoKq5vvQO/3001m3bh2rVq2iXa4PXxUrVuTGG2/klVdeyd6nqtx22208/fTTJxR/qVKlyMrKIiIighkzZjBt2jRmzZpF2bJl6dy5c/Z4gpiYGCJ9/79VlebNm2c3CufUu3dvJk2aROvWrRkzZgwzZsw47pjp06fzr3/967j9ZcuW5ddffz1mX61atUhKSsreTkpKokaNGvn+PuPHj8+udjry/lq1atHBN7noNddcc0y7REpKyjGluEApcVVPaZlZAFx0RuGfNvI0Ywa0auUGzAEMHgz9+1uSMGFn5cqVrF69Ont74cKF1K1bF4C4uDj2798PuHrycuXKUaFCBf7880++/vprAJo0acKWLVuYO9d1LNm/fz8ZGW5NmLp16/LZZ59x6623snTp0uOuPWjQIF5//fXs47t27cqECROye0Tt3r2bP/74A4DSpUuTnp6e5+9wxhlnsG7dOgD27dtHpUqVKFu2LCtWrGD27Nn5vmfHjh3ZiSI9PT07xv3791O9enXS09PzrOqCoyWK3F+5kwRA9erViYuLY/bs2agqY8eO5YorrsjzvPv27WPmzJnHvF6tWjVq167NypWuFuX777+nWbNm2a+vWrWKFsUw6LfEJYoU35Qd7etXLuTIXHbsgNtugy5d3EyvOYqrxoSjAwcOcNttt9GsWTNatWrFsmXLeOyxxwDo27cvPXv2pEuXLrRu3ZozzzyT5s2b06dPHzp16gS4htWPPvqI/v3707p1a7p163bMiOAzzjiDcePGce2117J27dpjrp2QkMBf//pXUlNTAWjWrBnDhw+ne/futGrVim7durF169bsWFq1asVNN9103O/wl7/8JftTf48ePcjIyKBVq1Y88sgjnHPOOXn+3lFRUUyYMIEhQ4bQunVr2rRpk/2QHzZsGB06dKBbt240adLk5G9uDq+++ip///vfadiwIQ0aNMhuyH7ttdd47bXXso+bOHEi3bt3p1y5cse8/7///S833XQTrVq1YuHChTz00EOAS3Br1qw5rsQWCJJXHVowi67eSKvf9hLf/PMCzqjm58P+ww/hnnvgwAE3kd/QoXAKPRmMKQrLly/P7upoTs7WrVu59dZb+e6777wOpdhNnDiR+fPnZ7dX5JTX35aIJKrqSWWVEtlGAdCwauEDdbJlZLgpOF57zQ2iM8aEhOrVq3PnnXeSnJx8QmMpQkFGRsYxHQICqUSWKKb9+CvnNypgttiDB91guTp1XCP1kd/R5tUxQcRKFCZQirpEUeLaKACqlS9g3Ykvv4TmzeGZZ2DVKrdPxJKECUol7YOaCX6B+JsqkYkiT0lJbkzEZZdBuXJuCvCXXvI6KmPyFRMTw65duyxZmCJzZD2KmJiiXcStxLZRHGfdOvjmG3j6aRg0yOZmMkHvSB/7HTt2eB2KCSFHVrgrSiU7UcyZA7NmwcCBbt3qjRshPt7rqIzxS+nSpYt0FTJjAiWgVU8i0kNEVorIGhF5II/Xo0XkI9/rv4lIPX/OWzn9kGukPuccePFF13gNliSMMSYAApYoRCQSGAX0BJoBN4hI7r6pdwB7VLUh8H/AM4WdNyHtIPHtWrtZXgcMgN9/d20SxhhjAiKQJYr2wBpVXaeqacB4IPfY9SuAd30/TwC6SiHr+FXbvQ1q14a5c11jdZj1nTbGmOIWyDaKmsCmHNtJQIf8jlHVDBHZB8QDO3MeJCJ9gSOTzqfKvHlLyDG1cBhLINe9CmN2L46ye3GU3YujzjjZNwYyUeRVMsjdD9CfY1DV0cBoABGZd7KDRkKN3Yuj7F4cZffiKLsXR4nIvJN9byCrnpKA2jm2awG5J07PPkZESgEVgN0BjMkYY8wJCmSimAs0EpH6IhIFXA9MznXMZOA238/XAD+ojT4yxpigErCqJ1+bw73AN0Ak8LaqLhWRJ4B5qjoZeAt4T0TW4EoS1/tx6tGBirkEsntxlN2Lo+xeHGX34qiTvhclblJAY4wxxSt05noyxhgTEJYojDHGFChoE0Wgpv8oify4F4NEZJmILBaR70WkrhdxFofC7kWO464RERWRkO0a6c+9EJG/+f42lorIB8UdY3Hx4/9IHRGZLiILfP9PenkRZ6CJyNsisl1EluTzuojISN99Wiwibf06saoG3Reu8XstcDoQBSwCmuU65m7gNd/P1wMfeR23h/eiC1DW93O/cL4XvuPigB+B2UA7r+P28O+iEbAAqOTbrup13B7ei9FAP9/PzYANXscdOkhd2wAABjpJREFUoHtxAdAWWJLP672Ar3Fj2M4BfvPnvMFaogjI9B8lVKH3QlWnq+oh3+Zs3JiVUOTP3wXAMOBZIKU4gytm/tyLO4FRqroHQFW3F3OMxcWfe6HAkfl+KnD8mK6QoKo/UvBYtCuAserMBiqKSPXCzhusiSKv6T9q5neMqmYAR6b/CDX+3Iuc7sB9YghFhd4LETkTqK2qXxZnYB7w5++iMdBYRH4Rkdki0qPYoite/tyLx4CbRSQJmAL0L57Qgs6JPk+A4F2Posim/wgBfv+eInIz0A64MKAReafAeyEiEbhZiHsXV0Ae8ufvohSu+qkzrpT5k4i0UNW9AY6tuPlzL24AxqjqCyJyLm78VgtVzQp8eEHlpJ6bwVqisOk/jvLnXiAiFwNDgctVNbWYYituhd2LOKAFMENENuDqYCeHaIO2v/9HPlfVdFVdD6zEJY5Q48+9uAP4GEBVZwExuAkDw41fz5PcgjVR2PQfRxV6L3zVLa/jkkSo1kNDIfdCVfepaoKq1lPVerj2mstV9aQnQwti/vwfmYTr6ICIJOCqotYVa5TFw597sRH4//buN0SqKozj+PdHWK1JgkSRBFkYVpIuZSH5IsyS/lCUiFuYtoGEUoSWvQiDCnohWS8ys60kNDBZVhSkP5iIWYimEv6pRTJMIpCSEAnbILanF+ds3rbZuzPbmrPO7wMDO+f+Oc8c2PvMPXd4zjQASdeREkUjrkG7CZibf/00GTgZEcf6O6gup57izJX/GHKqHItlwAigIz/P/yEi7j9rQZ8hVY5FQ6hyLDYD0yV1At3AsxHxy9mL+syociyeAd6VtIg01dJ6Ln6xlLSONNV4SX4e8wIwDCAi2kjPZ+4BvgN+Ax6r6rzn4FiZmdkgqtepJzMzqxNOFGZmVsqJwszMSjlRmJlZKScKMzMr5URhdUdSt6R9hdeYkn3H9FUps8Y+P8vVR/fnkhfjBnCO+ZLm5r9bJY0ubFsl6fpBjnOPpOYqjlkoafh/7dsalxOF1aOuiGguvI7+T/3OjoiJpGKTy2o9OCLaIuL9/LYVGF3YNi8iOgclytNxrqS6OBcCThQ2YE4UNiTkO4cvJH2VX7dW2Ge8pN35LuSApGty+yOF9rclnddPd58DY/Ox0/IaBgdzrf8LcvtSnV4D5NXc9qKkxZJmkmpurc19NuU7gUmSFkh6pRBzq6Q3BhjnTgoF3SS9JWmv0toTL+W2p0gJa5ukbbltuqSdeRw7JI3opx9rcE4UVo+aCtNOG3Pbz8CdEXEj0AIsr3DcfOD1iGgmXah/zOUaWoApub0bmN1P//cBByVdCKwGWiLiBlIlgwWSRgEPAuMjYgLwcvHgiFgP7CV982+OiK7C5vXAjML7FqB9gHHeRSrT0WNJREwCJgC3SZoQEctJtXymRsTUXMrjeeCOPJZ7gaf76ccaXF2W8LCG15UvlkXDgBV5Tr6bVLeot53AEklXABsi4rCkacBNwJ5c3qSJlHQqWSupCzhKKkM9Dvg+Ir7N29cATwArSGtdrJL0EVB1SfOIOC7pSK6zczj3sSOft5Y4LyKVqyiuUDZL0uOk/+vLSQv0HOh17OTcviP3cz5p3Mz65ERhQ8Ui4CdgIulO+F+LEkXEB5K+BO4FNkuaRyqrvCYinquij9nFAoKSKq5vkmsL3UIqMvcQ8CRwew2fpR2YBRwCNkZEKF21q46TtIrbUuBNYIakq4DFwM0RcULSalLhu94EbImIh2uI1xqcp55sqBgJHMvrB8whfZv+B0lXA0fydMsm0hTMVmCmpEvzPqNU/Zrih4Axksbm93OA7XlOf2REfEx6UFzpl0e/ksqeV7IBeIC0RkJ7bqspzoj4gzSFNDlPW10MnAJOSroMuLuPWHYBU3o+k6ThkirdnZn9zYnChoqVwKOSdpGmnU5V2KcF+FrSPuBa0pKPnaQL6qeSDgBbSNMy/YqI30nVNTskHQT+BNpIF90P8/m2k+52elsNtPU8zO513hNAJ3BlROzObTXHmZ99vAYsjoj9pPWxvwHeI01n9XgH+ETStog4TvpF1rrczy7SWJn1ydVjzcyslO8ozMyslBOFmZmVcqIwM7NSThRmZlbKicLMzEo5UZiZWSknCjMzK/UXZYybmOQQS/IAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "y_probs = np.load('data/stacknet_prediction.npy')\n",
    "r_roc_auc = roc_auc_score(Y_test, y_probs)\n",
    "fpr, tpr, thresholds = roc_curve(Y_test, y_probs)\n",
    "plt.figure()\n",
    "plt.plot(fpr, tpr, label='StackNet (area = %0.2f)' % r_roc_auc)\n",
    "plt.plot([0, 1], [0, 1],'r--')\n",
    "plt.xlim([0.0, 1.0])\n",
    "plt.ylim([0.0, 1.05])\n",
    "plt.xlabel('False Positive Rate')\n",
    "plt.ylabel('True Positive Rate')\n",
    "plt.title('ROC of StackNet (4-fold CV)')\n",
    "plt.legend(loc=\"lower right\")\n",
    "plt.savefig('images2/StackNet_ROC')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 5: Deep Learning in BCI: EEGNET"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "import tensorflow as tf\n",
    "from tensorflow.keras.models import Model\n",
    "from tensorflow.keras.layers import Dense, Activation, Permute, Dropout\n",
    "from tensorflow.keras.layers import Conv2D, MaxPooling2D, AveragePooling2D\n",
    "from tensorflow.keras.layers import SeparableConv2D, DepthwiseConv2D\n",
    "from tensorflow.keras.layers import BatchNormalization\n",
    "from tensorflow.keras.layers import SpatialDropout2D\n",
    "from tensorflow.keras.regularizers import l1_l2\n",
    "from tensorflow.keras.layers import Input, Flatten\n",
    "from tensorflow.keras.constraints import max_norm\n",
    "from tensorflow.keras import backend as K\n",
    "from tensorflow.keras import optimizers\n",
    "from tensorflow.keras import utils as np_utils\n",
    "from tensorflow.keras.callbacks import ModelCheckpoint"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "def EEGNet(nb_classes, Chans = 64, Samples = 128, \n",
    "             dropoutRate = 0.5, kernLength = 64, F1 = 8, \n",
    "             D = 2, F2 = 16, norm_rate = 0.25, dropoutType = 'Dropout'):\n",
    "    \"\"\" Keras Implementation of EEGNet\n",
    "    http://iopscience.iop.org/article/10.1088/1741-2552/aace8c/meta\n",
    "\n",
    "    Note that this implements the newest version of EEGNet and NOT the earlier\n",
    "    version (version v1 and v2 on arxiv). We strongly recommend using this\n",
    "    architecture as it performs much better and has nicer properties than\n",
    "    our earlier version. For example:\n",
    "        \n",
    "        1. Depthwise Convolutions to learn spatial filters within a \n",
    "        temporal convolution. The use of the depth_multiplier option maps \n",
    "        exactly to the number of spatial filters learned within a temporal\n",
    "        filter. This matches the setup of algorithms like FBCSP which learn \n",
    "        spatial filters within each filter in a filter-bank. This also limits \n",
    "        the number of free parameters to fit when compared to a fully-connected\n",
    "        convolution. \n",
    "        \n",
    "        2. Separable Convolutions to learn how to optimally combine spatial\n",
    "        filters across temporal bands. Separable Convolutions are Depthwise\n",
    "        Convolutions followed by (1x1) Pointwise Convolutions. \n",
    "        \n",
    "    \n",
    "    While the original paper used Dropout, we found that SpatialDropout2D \n",
    "    sometimes produced slightly better results for classification of ERP \n",
    "    signals. However, SpatialDropout2D significantly reduced performance \n",
    "    on the Oscillatory dataset (SMR, BCI-IV Dataset 2A). We recommend using\n",
    "    the default Dropout in most cases.\n",
    "        \n",
    "    Assumes the input signal is sampled at 128Hz. If you want to use this model\n",
    "    for any other sampling rate you will need to modify the lengths of temporal\n",
    "    kernels and average pooling size in blocks 1 and 2 as needed (double the \n",
    "    kernel lengths for double the sampling rate, etc). Note that we haven't \n",
    "    tested the model performance with this rule so this may not work well. \n",
    "    \n",
    "    The model with default parameters gives the EEGNet-8,2 model as discussed\n",
    "    in the paper. This model should do pretty well in general, although it is\n",
    "\tadvised to do some model searching to get optimal performance on your\n",
    "\tparticular dataset.\n",
    "\n",
    "    We set F2 = F1 * D (number of input filters = number of output filters) for\n",
    "    the SeparableConv2D layer. We haven't extensively tested other values of this\n",
    "    parameter (say, F2 < F1 * D for compressed learning, and F2 > F1 * D for\n",
    "    overcomplete). We believe the main parameters to focus on are F1 and D. \n",
    "\n",
    "    Inputs:\n",
    "        \n",
    "      nb_classes      : int, number of classes to classify\n",
    "      Chans, Samples  : number of channels and time points in the EEG data\n",
    "      dropoutRate     : dropout fraction\n",
    "      kernLength      : length of temporal convolution in first layer. We found\n",
    "                        that setting this to be half the sampling rate worked\n",
    "                        well in practice. For the SMR dataset in particular\n",
    "                        since the data was high-passed at 4Hz we used a kernel\n",
    "                        length of 32.     \n",
    "      F1, F2          : number of temporal filters (F1) and number of pointwise\n",
    "                        filters (F2) to learn. Default: F1 = 8, F2 = F1 * D. \n",
    "      D               : number of spatial filters to learn within each temporal\n",
    "                        convolution. Default: D = 2\n",
    "      dropoutType     : Either SpatialDropout2D or Dropout, passed as a string.\n",
    "\n",
    "    \"\"\"\n",
    "    \n",
    "    if dropoutType == 'SpatialDropout2D':\n",
    "        dropoutType = SpatialDropout2D\n",
    "    elif dropoutType == 'Dropout':\n",
    "        dropoutType = Dropout\n",
    "    else:\n",
    "        raise ValueError('dropoutType must be one of SpatialDropout2D '\n",
    "                         'or Dropout, passed as a string.')\n",
    "    \n",
    "    input1   = Input(shape = (1, Chans, Samples))\n",
    "\n",
    "    ##################################################################\n",
    "    block1       = Conv2D(F1, (1, kernLength), padding = 'same',\n",
    "                                   input_shape = (1, Chans, Samples),\n",
    "                                   use_bias = False)(input1)\n",
    "    block1       = BatchNormalization(axis = 1)(block1)\n",
    "    block1       = DepthwiseConv2D((Chans, 1), use_bias = False, \n",
    "                                   depth_multiplier = D,\n",
    "                                   depthwise_constraint = max_norm(1.),\n",
    "                                  data_format='channels_first')(block1)\n",
    "    block1       = BatchNormalization(axis = 1)(block1)\n",
    "    block1       = Activation('elu')(block1)\n",
    "    block1       = AveragePooling2D((1, 4), data_format='channels_first')(block1)\n",
    "    block1       = dropoutType(dropoutRate)(block1)\n",
    "    \n",
    "    block2       = SeparableConv2D(F2, (1, 16),\n",
    "                                   use_bias = False, padding = 'same')(block1)\n",
    "    block2       = BatchNormalization(axis = 1)(block2)\n",
    "    block2       = Activation('elu')(block2)\n",
    "    block2       = AveragePooling2D((1, 8), data_format='channels_first')(block2)\n",
    "    block2       = dropoutType(dropoutRate)(block2)\n",
    "        \n",
    "    flatten      = Flatten(name = 'flatten')(block2)\n",
    "    \n",
    "    dense        = Dense(nb_classes, name = 'dense', \n",
    "                         kernel_constraint = max_norm(norm_rate))(flatten)\n",
    "    softmax      = Activation('softmax', name = 'softmax')(dense)\n",
    "    \n",
    "    return Model(inputs=input1, outputs=softmax)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "training data shape: (4080, 1, 56, 140)\n",
      "validation data shape: (1360, 1, 56, 140)\n",
      "testing data shape: (3400, 1, 56, 140)\n",
      "(4080,) (1360,)\n"
     ]
    }
   ],
   "source": [
    "subj, trial, channel, sample = train_data_list.shape\n",
    "X_train_valid = np.reshape(train_data_list, (-1, 1, channel, sample))\n",
    "X_test = np.reshape(test_data_list, (-1, 1, channel, sample))\n",
    "\n",
    "Y_train_valid = pd.read_csv('data/TrainLabels.csv')['Prediction'].values\n",
    "split_thres = int(subj * trial * 0.25)\n",
    "X_train, X_valid = X_train_valid[split_thres:, :], X_train_valid[:split_thres, :]\n",
    "Y_train, Y_valid = Y_train_valid[split_thres:], Y_train_valid[:split_thres]\n",
    "\n",
    "print('training data shape: ' + str(X_train.shape))\n",
    "print('validation data shape: ' + str(X_valid.shape))\n",
    "print('testing data shape: ' + str(X_test.shape))\n",
    "print(Y_train.shape, Y_valid.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "# configure EEGNET model\n",
    "model = EEGNet(nb_classes=2, Chans=channel, Samples=sample)\n",
    "\n",
    "# compile the model and set the optimizers\n",
    "#sgd = optimizers.SGD(lr=0.001, decay=1e-6, momentum=0.9, nesterov=True)\n",
    "\n",
    "model.compile(loss='sparse_categorical_crossentropy', optimizer='adam',\n",
    "              metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train on 4080 samples, validate on 1360 samples\n",
      "Epoch 1/100\n",
      "WARNING:tensorflow:Can save best model only with val_loss available, skipping.\n"
     ]
    },
    {
     "ename": "AbortedError",
     "evalue": " Operation received an exception:Status: 3, message: could not create a reorder primitive descriptor, in file tensorflow/core/kernels/mkl_conv_ops.cc:868\n\t [[node model/separable_conv2d/separable_conv2d/depthwise (defined at /opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/framework/ops.py:1751) ]] [Op:__inference_distributed_function_2021]\n\nFunction call stack:\ndistributed_function\n",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mAbortedError\u001b[0m                              Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-30-78d5523b4714>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m     11\u001b[0m fittedModel = model.fit(X_train, Y_train, batch_size=32, epochs=100,\n\u001b[1;32m     12\u001b[0m         \u001b[0mverbose\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvalidation_data\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mX_valid\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mY_valid\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 13\u001b[0;31m         callbacks=[checkpointer], class_weight=class_weights)\n\u001b[0m",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/keras/engine/training.py\u001b[0m in \u001b[0;36mfit\u001b[0;34m(self, x, y, batch_size, epochs, verbose, callbacks, validation_split, validation_data, shuffle, class_weight, sample_weight, initial_epoch, steps_per_epoch, validation_steps, validation_freq, max_queue_size, workers, use_multiprocessing, **kwargs)\u001b[0m\n\u001b[1;32m    726\u001b[0m         \u001b[0mmax_queue_size\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mmax_queue_size\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    727\u001b[0m         \u001b[0mworkers\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mworkers\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 728\u001b[0;31m         use_multiprocessing=use_multiprocessing)\n\u001b[0m\u001b[1;32m    729\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    730\u001b[0m   def evaluate(self,\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/keras/engine/training_v2.py\u001b[0m in \u001b[0;36mfit\u001b[0;34m(self, model, x, y, batch_size, epochs, verbose, callbacks, validation_split, validation_data, shuffle, class_weight, sample_weight, initial_epoch, steps_per_epoch, validation_steps, validation_freq, **kwargs)\u001b[0m\n\u001b[1;32m    322\u001b[0m                 \u001b[0mmode\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mModeKeys\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mTRAIN\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    323\u001b[0m                 \u001b[0mtraining_context\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mtraining_context\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 324\u001b[0;31m                 total_epochs=epochs)\n\u001b[0m\u001b[1;32m    325\u001b[0m             \u001b[0mcbks\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmake_logs\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mepoch_logs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtraining_result\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mModeKeys\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mTRAIN\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    326\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/keras/engine/training_v2.py\u001b[0m in \u001b[0;36mrun_one_epoch\u001b[0;34m(model, iterator, execution_function, dataset_size, batch_size, strategy, steps_per_epoch, num_samples, mode, training_context, total_epochs)\u001b[0m\n\u001b[1;32m    121\u001b[0m         step=step, mode=mode, size=current_batch_size) as batch_logs:\n\u001b[1;32m    122\u001b[0m       \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 123\u001b[0;31m         \u001b[0mbatch_outs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mexecution_function\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0miterator\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    124\u001b[0m       \u001b[0;32mexcept\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mStopIteration\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merrors\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mOutOfRangeError\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    125\u001b[0m         \u001b[0;31m# TODO(kaftan): File bug about tf function and errors.OutOfRangeError?\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/keras/engine/training_v2_utils.py\u001b[0m in \u001b[0;36mexecution_function\u001b[0;34m(input_fn)\u001b[0m\n\u001b[1;32m     84\u001b[0m     \u001b[0;31m# `numpy` translates Tensors to values in Eager mode.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     85\u001b[0m     return nest.map_structure(_non_none_constant_value,\n\u001b[0;32m---> 86\u001b[0;31m                               distributed_function(input_fn))\n\u001b[0m\u001b[1;32m     87\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     88\u001b[0m   \u001b[0;32mreturn\u001b[0m \u001b[0mexecution_function\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/def_function.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwds)\u001b[0m\n\u001b[1;32m    455\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    456\u001b[0m     \u001b[0mtracing_count\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_get_tracing_count\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 457\u001b[0;31m     \u001b[0mresult\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_call\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwds\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    458\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mtracing_count\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_get_tracing_count\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    459\u001b[0m       \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_call_counter\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcalled_without_tracing\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/def_function.py\u001b[0m in \u001b[0;36m_call\u001b[0;34m(self, *args, **kwds)\u001b[0m\n\u001b[1;32m    518\u001b[0m         \u001b[0;31m# Lifting succeeded, so variables are initialized and we can run the\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    519\u001b[0m         \u001b[0;31m# stateless function.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 520\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_stateless_fn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwds\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    521\u001b[0m     \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    522\u001b[0m       \u001b[0mcanon_args\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcanon_kwds\u001b[0m \u001b[0;34m=\u001b[0m\u001b[0;31m \u001b[0m\u001b[0;31m\\\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/function.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m   1821\u001b[0m     \u001b[0;34m\"\"\"Calls a graph function specialized to the inputs.\"\"\"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1822\u001b[0m     \u001b[0mgraph_function\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_maybe_define_function\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1823\u001b[0;31m     \u001b[0;32mreturn\u001b[0m \u001b[0mgraph_function\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_filtered_call\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m  \u001b[0;31m# pylint: disable=protected-access\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   1824\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1825\u001b[0m   \u001b[0;34m@\u001b[0m\u001b[0mproperty\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/function.py\u001b[0m in \u001b[0;36m_filtered_call\u001b[0;34m(self, args, kwargs)\u001b[0m\n\u001b[1;32m   1139\u001b[0m          if isinstance(t, (ops.Tensor,\n\u001b[1;32m   1140\u001b[0m                            resource_variable_ops.BaseResourceVariable))),\n\u001b[0;32m-> 1141\u001b[0;31m         self.captured_inputs)\n\u001b[0m\u001b[1;32m   1142\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1143\u001b[0m   \u001b[0;32mdef\u001b[0m \u001b[0m_call_flat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcaptured_inputs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcancellation_manager\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/function.py\u001b[0m in \u001b[0;36m_call_flat\u001b[0;34m(self, args, captured_inputs, cancellation_manager)\u001b[0m\n\u001b[1;32m   1222\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mexecuting_eagerly\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1223\u001b[0m       flat_outputs = forward_function.call(\n\u001b[0;32m-> 1224\u001b[0;31m           ctx, args, cancellation_manager=cancellation_manager)\n\u001b[0m\u001b[1;32m   1225\u001b[0m     \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1226\u001b[0m       \u001b[0mgradient_name\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_delayed_rewrite_functions\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mregister\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/function.py\u001b[0m in \u001b[0;36mcall\u001b[0;34m(self, ctx, args, cancellation_manager)\u001b[0m\n\u001b[1;32m    509\u001b[0m               \u001b[0minputs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    510\u001b[0m               \u001b[0mattrs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"executor_type\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mexecutor_type\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"config_proto\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mconfig\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 511\u001b[0;31m               ctx=ctx)\n\u001b[0m\u001b[1;32m    512\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    513\u001b[0m           outputs = execute.execute_with_cancellation(\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/eager/execute.py\u001b[0m in \u001b[0;36mquick_execute\u001b[0;34m(op_name, num_outputs, inputs, attrs, ctx, name)\u001b[0m\n\u001b[1;32m     65\u001b[0m     \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     66\u001b[0m       \u001b[0mmessage\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0me\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmessage\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 67\u001b[0;31m     \u001b[0msix\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mraise_from\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcore\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_status_to_exception\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0me\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcode\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmessage\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     68\u001b[0m   \u001b[0;32mexcept\u001b[0m \u001b[0mTypeError\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0me\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     69\u001b[0m     keras_symbolic_tensors = [\n",
      "\u001b[0;32m/opt/anaconda3/lib/python3.7/site-packages/six.py\u001b[0m in \u001b[0;36mraise_from\u001b[0;34m(value, from_value)\u001b[0m\n",
      "\u001b[0;31mAbortedError\u001b[0m:  Operation received an exception:Status: 3, message: could not create a reorder primitive descriptor, in file tensorflow/core/kernels/mkl_conv_ops.cc:868\n\t [[node model/separable_conv2d/separable_conv2d/depthwise (defined at /opt/anaconda3/lib/python3.7/site-packages/tensorflow_core/python/framework/ops.py:1751) ]] [Op:__inference_distributed_function_2021]\n\nFunction call stack:\ndistributed_function\n"
     ]
    }
   ],
   "source": [
    "# set a valid path for your system to record model checkpoints\n",
    "checkpointer = ModelCheckpoint(filepath='/tmp/checkpoint.h5', verbose=1,\n",
    "                               save_best_only=True)\n",
    "# Weighted loss\n",
    "weight_0 = 1/(len([y for y in Y_train_valid if y == 0]))\n",
    "weight_1 = 1/(len([y for y in Y_train_valid if y == 1]))\n",
    "class_weights = {0: weight_0, 1: weight_1}\n",
    "\n",
    "#\n",
    "# # fit the model\n",
    "fittedModel = model.fit(X_train, Y_train, batch_size=32, epochs=100,\n",
    "        verbose=2, validation_data=(X_valid, Y_valid),\n",
    "        callbacks=[checkpointer], class_weight=class_weights)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sess = tf.Session(config=tf.ConfigProto(log_device_placement=True))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.python.client import device_lib\n",
    "\n",
    "def get_available_gpus():\n",
    "    local_device_protos = device_lib.list_local_devices()\n",
    "    return [x.name for x in local_device_protos if x.device_type == 'GPU']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "get_available_gpus()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(tf.__version__)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}