LSTM-AE 2nd SWAT.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "cbe2e728",
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "from typing import Dict\n",
    "\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "d6b60065",
   "metadata": {},
   "outputs": [],
   "source": [
    "seed_value = 42\n",
    "\n",
    "import os\n",
    "os.environ['PYTHONHASHSEED'] = str(seed_value)\n",
    "\n",
    "import random\n",
    "random.seed(seed_value)\n",
    "\n",
    "import numpy as np\n",
    "np.random.seed(seed_value)\n",
    "\n",
    "import tensorflow as tf\n",
    "tf.random.set_seed(seed_value)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "96ffde0c",
   "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>Timestamp</th>\n",
       "      <th>FIT101</th>\n",
       "      <th>LIT101</th>\n",
       "      <th>MV101</th>\n",
       "      <th>P101</th>\n",
       "      <th>P102</th>\n",
       "      <th>AIT201</th>\n",
       "      <th>AIT202</th>\n",
       "      <th>AIT203</th>\n",
       "      <th>FIT201</th>\n",
       "      <th>...</th>\n",
       "      <th>P501</th>\n",
       "      <th>P502</th>\n",
       "      <th>PIT501</th>\n",
       "      <th>PIT502</th>\n",
       "      <th>PIT503</th>\n",
       "      <th>FIT601</th>\n",
       "      <th>P601</th>\n",
       "      <th>P602</th>\n",
       "      <th>P603</th>\n",
       "      <th>Normal/Attack</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>28/12/2015 7:50:00 AM</td>\n",
       "      <td>2.462287</td>\n",
       "      <td>541.3740</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>260.7344</td>\n",
       "      <td>8.409895</td>\n",
       "      <td>324.8398</td>\n",
       "      <td>0.000513</td>\n",
       "      <td>...</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>251.4580</td>\n",
       "      <td>1.297536</td>\n",
       "      <td>190.1436</td>\n",
       "      <td>0.000064</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>Normal</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>28/12/2015 7:50:01 AM</td>\n",
       "      <td>2.441469</td>\n",
       "      <td>541.4918</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>260.7344</td>\n",
       "      <td>8.412138</td>\n",
       "      <td>324.3015</td>\n",
       "      <td>0.000256</td>\n",
       "      <td>...</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>251.4580</td>\n",
       "      <td>1.297536</td>\n",
       "      <td>190.1436</td>\n",
       "      <td>0.000064</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>Normal</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>28/12/2015 7:50:02 AM</td>\n",
       "      <td>2.428338</td>\n",
       "      <td>541.2955</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>260.7344</td>\n",
       "      <td>8.414701</td>\n",
       "      <td>323.7888</td>\n",
       "      <td>0.000256</td>\n",
       "      <td>...</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>251.4580</td>\n",
       "      <td>1.297536</td>\n",
       "      <td>190.1436</td>\n",
       "      <td>0.000064</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>Normal</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>28/12/2015 7:50:03 AM</td>\n",
       "      <td>2.425776</td>\n",
       "      <td>542.5124</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>260.7344</td>\n",
       "      <td>8.420149</td>\n",
       "      <td>323.3530</td>\n",
       "      <td>0.000256</td>\n",
       "      <td>...</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>251.4580</td>\n",
       "      <td>1.441707</td>\n",
       "      <td>190.1436</td>\n",
       "      <td>0.000064</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>Normal</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>28/12/2015 7:50:04 AM</td>\n",
       "      <td>2.424174</td>\n",
       "      <td>542.7871</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>260.7344</td>\n",
       "      <td>8.423674</td>\n",
       "      <td>322.8916</td>\n",
       "      <td>0.000256</td>\n",
       "      <td>...</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>251.2817</td>\n",
       "      <td>1.441707</td>\n",
       "      <td>190.1436</td>\n",
       "      <td>0.000064</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>Normal</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 53 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                Timestamp    FIT101    LIT101  MV101  P101  P102    AIT201  \\\n",
       "0   28/12/2015 7:50:00 AM  2.462287  541.3740      2     1     1  260.7344   \n",
       "1   28/12/2015 7:50:01 AM  2.441469  541.4918      2     1     1  260.7344   \n",
       "2   28/12/2015 7:50:02 AM  2.428338  541.2955      2     1     1  260.7344   \n",
       "3   28/12/2015 7:50:03 AM  2.425776  542.5124      2     1     1  260.7344   \n",
       "4   28/12/2015 7:50:04 AM  2.424174  542.7871      2     1     1  260.7344   \n",
       "\n",
       "     AIT202    AIT203    FIT201  ...  P501  P502    PIT501    PIT502  \\\n",
       "0  8.409895  324.8398  0.000513  ...     2     1  251.4580  1.297536   \n",
       "1  8.412138  324.3015  0.000256  ...     2     1  251.4580  1.297536   \n",
       "2  8.414701  323.7888  0.000256  ...     2     1  251.4580  1.297536   \n",
       "3  8.420149  323.3530  0.000256  ...     2     1  251.4580  1.441707   \n",
       "4  8.423674  322.8916  0.000256  ...     2     1  251.2817  1.441707   \n",
       "\n",
       "     PIT503    FIT601  P601  P602  P603  Normal/Attack  \n",
       "0  190.1436  0.000064     1     1     1         Normal  \n",
       "1  190.1436  0.000064     1     1     1         Normal  \n",
       "2  190.1436  0.000064     1     1     1         Normal  \n",
       "3  190.1436  0.000064     1     1     1         Normal  \n",
       "4  190.1436  0.000064     1     1     1         Normal  \n",
       "\n",
       "[5 rows x 53 columns]"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pd.read_excel(\"./data/SWAT.xlsx\")\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "87317792",
   "metadata": {},
   "outputs": [],
   "source": [
    "df['Normal/Attack'].replace({'Normal':0, 'Attack':1}, inplace=True)\n",
    "df[\"anomaly\"] = df[\"Normal/Attack\"]\n",
    "\n",
    "anomaly = df.anomaly\n",
    "df.drop([\"Normal/Attack\"], axis=1, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "62f0185b",
   "metadata": {},
   "outputs": [],
   "source": [
    "important_features = [\n",
    "    'LIT101', 'MV301', 'MV302', ' MV303', 'MV304', 'P302', 'AIT401', 'AIT504', 'FIT504', 'PIT502', 'P602', 'anomaly']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "54358dd3",
   "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>LIT101</th>\n",
       "      <th>MV301</th>\n",
       "      <th>MV302</th>\n",
       "      <th>MV303</th>\n",
       "      <th>MV304</th>\n",
       "      <th>P302</th>\n",
       "      <th>AIT401</th>\n",
       "      <th>AIT504</th>\n",
       "      <th>FIT504</th>\n",
       "      <th>PIT502</th>\n",
       "      <th>P602</th>\n",
       "      <th>anomaly</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>541.3740</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>148.808</td>\n",
       "      <td>12.03538</td>\n",
       "      <td>0.306761</td>\n",
       "      <td>1.297536</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>541.4918</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>148.808</td>\n",
       "      <td>12.03538</td>\n",
       "      <td>0.306761</td>\n",
       "      <td>1.297536</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>541.2955</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>148.808</td>\n",
       "      <td>12.03538</td>\n",
       "      <td>0.306761</td>\n",
       "      <td>1.297536</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>542.5124</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>148.808</td>\n",
       "      <td>12.03538</td>\n",
       "      <td>0.306761</td>\n",
       "      <td>1.441707</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>542.7871</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>148.808</td>\n",
       "      <td>12.03538</td>\n",
       "      <td>0.306761</td>\n",
       "      <td>1.441707</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     LIT101  MV301  MV302   MV303  MV304  P302   AIT401    AIT504    FIT504  \\\n",
       "0  541.3740      1      1       1      1     1  148.808  12.03538  0.306761   \n",
       "1  541.4918      1      1       1      1     1  148.808  12.03538  0.306761   \n",
       "2  541.2955      1      1       1      1     1  148.808  12.03538  0.306761   \n",
       "3  542.5124      1      1       1      1     1  148.808  12.03538  0.306761   \n",
       "4  542.7871      1      1       1      1     1  148.808  12.03538  0.306761   \n",
       "\n",
       "     PIT502  P602  anomaly  \n",
       "0  1.297536     1        0  \n",
       "1  1.297536     1        0  \n",
       "2  1.297536     1        0  \n",
       "3  1.441707     1        0  \n",
       "4  1.441707     1        0  "
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = df[important_features]\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "61ecdceb",
   "metadata": {},
   "outputs": [],
   "source": [
    "first_anomaly_index = anomaly[anomaly == 1].index[0]\n",
    "train_size = int(first_anomaly_index * 0.8)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "77086959",
   "metadata": {},
   "source": [
    "### Autoencoder prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "40afe34c",
   "metadata": {},
   "outputs": [],
   "source": [
    "from utils.data import create_sequences\n",
    "from sklearn.preprocessing import StandardScaler"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "5134b6d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "from utils.pipeline import Pipeline, Predictor\n",
    "from models.lstm_ae import LSTM_AE\n",
    "from tensorflow import keras\n",
    "\n",
    "optimizer = keras.optimizers.Adam(learning_rate=0.0001)\n",
    "class Config:\n",
    "    TIME_STEPS = 50\n",
    "    epochs = 3\n",
    "    batch_size = 10 \n",
    "    Q = 0.99\n",
    "    \n",
    "conf = Config\n",
    "conf.train_size = train_size\n",
    "\n",
    "class LSTM_AE_2ndPredictor:\n",
    "    def fit(self, df_train):\n",
    "        ss = StandardScaler()\n",
    "        self.ss = ss\n",
    "        x_train = df_train[:train_size].drop([\"anomaly\"], axis=1)\n",
    "        x_train = np.array(ss.fit_transform(x_train))\n",
    "        x_train = create_sequences(x_train, conf.TIME_STEPS)\n",
    "\n",
    "        self.model = LSTM_AE(optimizer=optimizer, loss='mae')\n",
    "        self.model.fit(x_train, conf.epochs, conf.batch_size, conf.TIME_STEPS)\n",
    "        \n",
    "    def predict(self, dfs_test: Dict[str, pd.DataFrame]) -> Dict[str, pd.Series]:\n",
    "        predictions = {}\n",
    "        for file_name, df in dfs_test.items():\n",
    "            y_test = df.anomaly\n",
    "            df = df[important_features]\n",
    "            df.drop([\"anomaly\"], axis=1, inplace=True)\n",
    "            x_test = np.array(self.ss.transform(df))\n",
    "            x_test = create_sequences(x_test, conf.TIME_STEPS)\n",
    "            r = pd.Series(np.sum(np.mean(np.abs(self.model.predict(x_test) - x_test), axis=1), axis=1)) \n",
    "            \n",
    "            smoothed = r.rolling(window = 10).mean().fillna(r.iloc[0])\n",
    "            a = smoothed.diff(10).abs()\n",
    "            \n",
    "            IQR = a.quantile(0.75) - a.quantile(0.75)\n",
    "            UCL = a.quantile(0.75) + 1.5 * IQR\n",
    "            \n",
    "            s = a > UCL\n",
    "            w = 10\n",
    "            is_anomaly = s.rolling(window = w).mean().fillna(s.iloc[0])\n",
    "            anomalous_data_indices = pd.Series(is_anomaly.index[is_anomaly > 0]) #+ conf.TIME_STEPS - 1\n",
    "\n",
    "            prediction = pd.Series(data=0, index=df.index)\n",
    "            prediction.iloc[anomalous_data_indices] = 1\n",
    "            predictions[file_name] = prediction\n",
    "            \n",
    "        return predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "id": "355fafb4",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoch 1/3\n",
      "1005/1005 [==============================] - 239s 228ms/step - loss: 0.3111\n",
      "Epoch 2/3\n",
      "1005/1005 [==============================] - 252s 251ms/step - loss: 0.2050\n",
      "Epoch 3/3\n",
      "1005/1005 [==============================] - 202s 201ms/step - loss: 0.1802\n",
      "CPU times: user 22min 19s, sys: 1min 28s, total: 23min 48s\n",
      "Wall time: 12min 33s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "predictor = LSTM_AE_2ndPredictor()\n",
    "predictor.fit(df[:train_size])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "id": "9b67284f",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 1min 1s, sys: 4.77 s, total: 1min 6s\n",
      "Wall time: 38.6 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "df_by_names = {\"swat\": df[train_size:]}\n",
    "y_test_by_names = {name: df.anomaly for name, df in df_by_names.items()}\n",
    "predictions_by_names = predictor.predict(df_by_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "id": "77d08925",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEICAYAAABPgw/pAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAqUElEQVR4nO3de7wVdb3/8deHzU0Bb4CmgkFHTEkBt3j74a0sBTXN0mOmJZGX7Fh0ORYezsPMc/In6clLkkm/zKyU1KOGhkWldjFR8A4Igoa6vSAgIDdF2J/fHzMbZu+91t6z9ppZa82s9/Px2DBrZtbM9/td3/VZM9/vfGfM3RERkezrVu0EiIhIMhTQRURyQgFdRCQnFNBFRHJCAV1EJCcU0EVEckIBPacs8HMzW2Vmj1c7PWmrt/yKFKKAnkFmNsTM3My6d7DaEcAngEHufoiZ9TSzu8xsafjeY1JK2y5mdo+ZrTezl83scx2s+w0ze8nM3jGz183smrZ5MrOJZvbPcHvPm9k+RTbXKr8JZqlk4efzkJltMLOFZvbxDta92swWm9nacN0vRJYNMLNHzGylma02s0fNbExlcpGe8Md3SpivleG0FVn3GDNrNrN1kb9zwmW9zOxnYT1ba2ZPm9m4yuamtiig59cHgaXuvj4y7+/A2cCbKe53KrAJ2A04C7jRzD5SZN0ZQKO77wDsD4wEvtay0MzOBb4EnAj0BU4CVhTZVqH8xhIGmCS/C7cDTwH9gcnAXWY2sMi664FPAjsC5wDXmdn/CZetAyYAA4GdgSnAfZ38kKci4X2eD3yK4PMeQZD/CzpY/3V37xv5+0U4vzvwKnA0Qfn9J3CHmQ1JMK3Z4u76S+EP+A7wGrAWWAQcC/QGNgIDwnUmA5uBHcLX/wVcG06fSBAU3iGotJdFtv0K4ARf+HXA4W32/SXgXWBLuPx7bZY3AcekkOc+BMF8n8i8XwJXxnhvf+BPwI/D193CfB8b473t8ksQAO8HlgOrwulBkfc8DHwfeCT8TPZOqAz2Ad4D+kXm/Q34csz3zwC+VWB+N4LA58CuRd57C/Bj4IGwHB4BPgBcG5bBQuDAyPqTgBfDOroAODWybHz4/muAlcB/J1hP/gGc3+bzm11k3WOAphK2/SzwmaTrdlb+qp6APP4BHw6D0R7h6yHAv4TTf22pcMCs8As1LrLs1HD6GOCA8Is8AlgGfCqyPQe6d5CG8cDfiyzrNKCHgWF1kb9ni7znQGBDm3n/DtzXwX4+R/Cj5WHwHRnO3yucNzEsy38SBOpucfJL8APxGWB7oB9wJ3BvZPnDBD+MHyE40utRYJv3d1AG9xdJx6nA823m3QD8KEa92Q54AxjbZv6zBD+UDvy0g/ffQnAGcxDBwcODYbl9AWgA/ht4KLL+6cAeYR07g+BsYfdIeW4GvhqWz3ZFPrti5bMa2KtIOtcAh0ZejwbWFln3mDDvy8K8XAP0KbLubgQ/7PtW4ntei39qcknHFqAXMNzMerj7Und/MVz2F+Do8BR2BHB9+Lo3cDBBUMfdH3b359y92d2fJTiNP7pSGXD3r7j7TkX+RhR5W1+C4By1hiCgFtvPbR40uewD/ITgiwswKPz/OIIfto8CZxIczcVJ/0p3/1933+DuawmOxtuW3y3uPt/dN7v7+wW2cVIHZXBSkV33DfMc1WEZRPwEeAb4Q5t0jAB2IAigf+9kG/e4+xPu/i5wD/Cuu9/q7luA3xD86LZs9053fz2sY78BFgPR/ofX3f1HYflsbLuj8LMrVj47ufsrRdLYtozWAH2LtKMvBEYBuwMfI/ix+mHblcysB/Br4BfuvrBY4eSdAnoK3H0J8HXgMuAtM5tuZnuEi/9CcNTRCDwH/JEg0BwGLHH3lQBmdmjYsbbczNYAXwYGVDIfXbCOIPBE7UBwSt8hd18MzCc4M4CgGQTgB+6+2t2XAjcBJ8RJiJltb2Y3hR1m7xD8UO5kZg2R1V6Ns60SdakMzOwqgn6Ef/XwcDPK3d9199uBSWY2soNNLYtMbyzwum9kn18IOxJXm9nqcP/ROpZG+UD7MtoBWFck32+6+4LwR+efwLcJzry2Cvs/fklwJH9RSmnOBAX0lIRHL0cQdNY5QYcWBO2HHyY4Nf+Luy8gaF44gSDYt7iNoD11sLvvSHD01nIEk/otMs3sJ22uLIj+zS/ytheA7mY2LDJvJEGgjqM78C/h9CK2NTO0KCXf3yIo50PDM4CjwvnRo8AOt2dmD3RQBg8Uedt84ENmFj0i77AMzOx7wDjgOHdve4bTVg/gQ52s0ykz+yDwU4IA2N/ddwLmUVr5nNVB+awzs72KvHU+QZm0KKWOOJG4FR7V/4ygueUzhc606okCegrM7MNm9jEz60XQprcRaAZw9w3AE8C/sS2A/4PgCDwa0PsBb7v7u2Z2CMHpdovl4fZK+mKHl3n1Dl/2NLPexS4Xc/cve+srC6J/Ba9a8eAKk7uBy82sT3iJ3SkER0+F0nOume0aTg8HLgH+HG5rA0ETwbfNrJ+ZDSK4OuL+mNntR1Duq81sF+C7Md8Xzc+4Dsqg4OVx7v4C8DTw3bB8TyVoWvvfQuub2SUEn+3HW87OIssOM7MjLLjkdDsz+w5B4Hqs1LwU0Idt/RaY2RcJjtBjc/dfd1A+fTtocrkV+KaZ7RmeuX6LoP2/HTP7qJl9MLwSaTBwJfDbyCo3AvsBnyzULFRvFNDT0Yug4q0guERwV4Jg1eIvBEdaj0de9yNsPw99hSAwrgUuBe5oWRAGu+8Dj4Sny4fFTNcigiC3J0E77UaCM4gkfYWgc+8tgnb/C919PoCZHWlm6yLrjgGeM7P1wMzw7z8iyy8iOD1/HXiU4Kzl5pjpuDZMxwpgNvD7LuanKz5L0NG3iqAenObuLYHzrDZnOFcQnKEtiRzZtpRBL4LLQFcSXDF1AnCiu79ebgLDM8P/ISjXZQT9FI+Uu92YbgLuI2hynAf8LpwHQFgGR4YvDyQ44Fkf/v8c4aWt4VnGBQRt7G9Gyu+sCuWj5liBZisREckgHaGLiOSEArqISE4ooIuI5IQCuohITlT8Jj8tBgwY4EOGDKnW7kVEMumJJ55Y4e4Fb/ZWtYA+ZMgQ5s6dW63di4hkkpm9XGyZmlxERHJCAV1EJCcU0EVEckIBXUQkJxTQRURyotOAbmY3m9lbZjavyHIzs+vNbImZPWtmjcknU0REOhPnCP0WYGwHy8cBw8K/8wluZykiIhXW6XXo7v7XTp6ifQpwa/i0kdlmtpOZ7e7ubySVyKIWzoQVL8ChF0CP7ZLf/oLfwpvzYPeRsF+xJ47lzDO/gZVLeG31Rl59e0Onqw/aZXveXvceGzZtKXlXDd17cMApE+m98x6dr1zAc01r+OOCN7v03lJs36s7X9y3mV4L7gR3WPQA7DEK+u4Gr8yGHXaH/sPgqH+Hbg2dbq9kT94KDb1g5BnJbnf1K/D4NFj5Euz2Eei/d/L7SNszv4GX/w59P1DtlMSz3ydh92JPcCxfEgOL9qT1o6qawnntArqZnU9wFM9eexV7mEkJpp8Z/L/HKPjQMeVvr637JsLGVcEXtx4CenMz3HMB4OyOsXucOys3wR5duANzNwve9M9HBjD0pItL3wBww0OL+cP8ZRR+REcyWu4u/ck3/8yez/9s24Jlz7Vfeb9Pwm7Dk03Ae2thxleD6QNOS/YH46lfwT9+FEwv+h1gcMDp0C0jXWvNzXDP+ZEZKVaERDisXAKn/zy1PVR0pKi7TwOmAYwePTq5G7E3l350GG+7zeluvyY5HPMfnP78kWzXo4FfnXto0TXH//xx5r32DivWvccPThvBv44eHHsvT77wMo23jcDLKNstzTB89x2YOfHIzlfuoidefpvP3PhokM5eO8K3FsIVuxde2VOoJ2nWvei29z0JFsZ9GFStiISQD46BL86sXlLimHpoOnUkIomf4teA6Dd5UDhPREQqKImAPgP4Qni1y2HAmoq0n4uISCudNrmY2e3AMcAAM2sieNhuDwB3/wnBcyBPAJYAG4AvppVYEREpLs5VLmd2stwJnmBfRWk9F9Xb/J932/IZ/1mz5ZVNec+0reDnEiedqTyfN7LNxLdfaHsZquueZtmkJOV0ZqQ7W2pVVr5HIvVAAV3aC68D7OxyQCsyHWsX3ZK5xCzNSxbDPRSczJ30CzJdmUh/+mlUQBcRyQkFdBGRnFBAj6MOG4or0yWaNyqNylJ5t5WPgJ7aRS51VmG6kN/yrlIpT2V33VzJnW3TKpMJZ7jVtq3APEmernKRiovXeWORjigruVMqmaqXdl9Y6+1noeOtizLRqZhxFShjBXQRkZxQQBcRyQkFdBGRnMhJQNfQ/2REh/6X+o6u7rLrnY0V/VTi7CztDkUN/W8tix24Gvov1RKnozOL3ykpRJ2ieaCALu3F/G6XNfQ/oR5/SzkQtb7IJcdBL895qxm6ykVERGJSQBcRyQkF9DjqsKE49tD/eiubDpsm6qwsqk7l3VY+AnpaQaXeglVXhv5XYZ/b3lrJz6dKdUFD/6UE+QjokjCL/NvBWlZ4Ot4usnE/9NadtznuOFSnaPo09F9EROJSQBcRyQkF9FjqsF2xYkNFc0Ttz5Wl8m4nJwE97aH/9aIKnaJlbKGyQ/+rdD/0aC419D/7NPRfqkX9ZPVEH3YeKKBLe7EjefQBF6XuIqEHXCSylZjbz/MvXJ7zVjN0lYuIiMSkgC4ikhMK6HHUYT+Rhv53hcqislTebeUjoGvofzLa5DdOi191h/6Xu/NSdlbBfbXar4b+54uucpGKi9d502rof4kdPol1D6Xcmdd68znuOFSnaPoqUMQK6CIiOaGALiKSE7ECupmNNbNFZrbEzCYVWL6XmT1kZk+Z2bNmdkLySa2m+mtXjD3yv/6KpjiVRWWp8rXTaUA3swZgKjAOGA6caWbD26z2n8Ad7n4g8Fngx0kntGMa+p+Mrgz9L6+MyrlKprKfTtV6RSOTGvqfeTUw9P8QYIm7v+Tum4DpwClt1nFgh3B6R+D15JIo1ZLUg5wlC/RZ50GcgL4n8GrkdVM4L+oy4GwzawJmAl8ttCEzO9/M5prZ3OXLl3chuVIRMQN5q0c/lDr0v1tCD7hIZCsdbb+cp3hkSJ7zVjOyM/T/TOAWdx8EnAD80grcrMPdp7n7aHcfPXDgwIR2LVJJCnxSu+IE9NeAwZHXg8J5UV8C7gBw90eB3sCAJBJYE+qw8yVu23gdFk0HVBiVpfJuK05AnwMMM7OhZtaToNNzRpt1XgGOBTCz/QgCutpUREQqqNOA7u6bgYuAPwDPE1zNMt/MLjezk8PVvgWcZ2bPALcD472SN/nQ0P9kZG7ofwU/n2rVBQ39z5l0y7d7rCS4zyTo7IzOuzQyvQAYk2zSpHpKH/pf+h4Suh96yk3aGvovialAGWukqIhITiigi4jkhAJ6LPXXrhi7KbX+iqY4tT9Xlsq7nZwEdA39T0abTtEYTX7lDv3PiurlM82h/53sT5JXA0P/RST31CmaBwro0l7sof/b1iv5vi8ZGfrfemc5Dno5zlrt0FUuIiISkwJ6HHXY+aL7oXeFCqOyVN5t5SOga6RoMrqQ33JLyL25jPeWufOa3Vmx/aY5UrSDeZIgdYpK1ahhtX7os84DBXQpoPSh/6WHg4Q6RVPuqNTQf0mMhv6LiEhcCuix1F+7YuyBompz3UZlUVkq73YU0EVEciInAV1D/5PRlaH/ye6ztHdW9DKXCu6ryH4TPyIttL16q/MVpqH/IjUkt52Hec1XfVFAl/biDv23wtPx3puNof/R2xvkN5iT77zVDF3lIiIiMSmgx1GHvelxr16pw6LpgAqjslTebeUjoGvofzK68JDopPdZobfW+M6K7VdD/4vLSpOROkVFJHVZCYjFZPWHKFkK6FJAF+6HXmJAMEum6qXdl6eh/5IYDf0XEZG4FNBj0emcxJDZ9ueMUnm3o4AuIpITOQnoGvqfjNKH/ie9z5LeWdGPR0P/JQEa+i8i6VOnaB4ooEt7cQ/Na2LofwUfcJHnK0HynLeaoatcakMddr7UYZYToEKrLJV3WwroIiI5ESugm9lYM1tkZkvMbFKRdf7VzBaY2Xwzuy3ZZHZCQ/+T0W7ofwVOw8sZ+l/JIzQN/a9xWWkySrd8u3e2gpk1AFOBTwBNwBwzm+HuCyLrDAMuAca4+yoz2zWtBItIGrISEIvJ6g9RsuIcoR8CLHH3l9x9EzAdOKXNOucBU919FYC7v5VsMqWy4g79L/UdkfW7JRRA0h76X1YuM0SdoumrkaH/ewKvRl43hfOi9gH2MbNHzGy2mY0ttCEzO9/M5prZ3OXLl3ctxSIiUlBSnaLdgWHAMcCZwE/NbKe2K7n7NHcf7e6jBw4cmNCuK6H+Tucq+7zOLOngKCuz7c8ZpfJuJ05Afw0YHHk9KJwX1QTMcPf33f2fwAsEAV4yTGfhkh2qrBCjUxSYAwwzs6EEgfyzwOfarHMvwZH5z81sAEETzEsJprMTGvqfjGrkNyND/6t2NJi9of/vv/8+TU1NvPvuu2Vvq0PNzXD8HcF0997w/PPp7q9cI/4zOEqKmc7evXszaNAgevToEXsXnQZ0d99sZhcBfwAagJvdfb6ZXQ7MdfcZ4bLjzGwBsAW42N1Xxk6FiFRZcke4TU1N9OvXjyFDhiQ2IrigLZth2aZgumdfGFDjjQLLu4E1wIC9O13V3Vm5ciVNTU0MHTo09i7iHKHj7jOBmW3mXRqZduCb4Z9kXcwvYfTLWvLQ/4QCSNon2hr6X7p33303/WCeSfHLw8zo378/pV48opGicdRh50sdZjkBKrQWCubl60oZKqCLiOREPgK6DieT0XbofyUOssr46Cr7qdfC0P8KbLuOvktLly7lttu6fpeSK664IsHUJCMfAV1EylR/TSQK6FInSh/6X2pAMINmLz+IpH0WUU4eMyVHbd6XXnop11577dbXkydP5rrrrmu33qRJk/jb3/7GqFGjuOaaa9iyZQsXX3wxBx98MCNGjOCmm24C4I033uCoo45i1KhR7L///vztb39j0qRJbNy4kVGjRnHWWWdVKmudinWVi9TPaWiL+stxAuqouSKu7903nwWvv5PoNofvsQPf/eRHii6fMGECn/70p/n6179Oc3Mz06dP5/HHH2+33pVXXsnVV1/N/fffD8C0adPYcccdmTNnDu+99x5jxozhuOOO4+677+b4449n8uTJbNmyhQ0bNnDkkUdyww038PTTTyeat3IpoItIrgwZMoT+/fvz1FNPsWzZMg488ED69+/f6ftmzZrFs88+y1133QXAmjVrWLx4MQcffDATJkzg/fff51Of+hSjRo1KOQddp4AuRVXkfuiSax0dSafp3HPP5ZZbbuHNN99kwoQJsd7j7vzoRz/i+OOPb7fsr3/9K7/73e8YP3483/zmN/nCF76QdJITkZM29BROdXX6XBHuzWW8Obl0dL6vvA/9z9eP96mnnsrvf/975syZUzBAA/Tr14+1a9dufX388cdz44038v777wPwwgsvsH79el5++WV22203zjvvPM4991yefPJJAHr06LF13fiq/IALqUOxR4qW/JZt65NM1a7sQ6JT3VV15ahTFKBnz5589KMfZaeddqKhoaHgOiNGjKChoYGRI0cyfvx4Jk6cyNKlS2lsbMTdGThwIPfeey8PP/wwV111FT169KBv377ceuutAJx//vmMGDGCxsZGfv3rX1cye0UpoItI7jQ3NzN79mzuvPPOouv06NGDBx98sNW8K664ot3liOeccw7nnHNOu/dPmTKFKVOmJJPghOSkySVlddj84nWY5/KpzCqrcHkvWLCAvffem2OPPZZhw2r8hl0J0xG6FJevs3CpE8OHD+ell7bdvfu5557j85//fKt1evXqxWOPPVbppKUuHwE9jaPJejxCrUqey7gfOk7lfnVqYeh/wmmok6H/BxxwQM1dL54WNbmIiOSEAroUUPrQ/1KPk83AEzi6Tv/ijDa5zNnVIFvlNV+1pAJFrIAeS/5OQztTfzlOQA6bKyRbFNBFRHJCAV2K0km41LuHH36Yk046CYAZM2Zw5ZVXFl139erV/PjHPy55H5dddhlXX311l9MYpYBeVD2ePlchz+U84EJD/5Pb9taf7/qo91u2bCn5PSeffDKTJk0quryrAT1JCujSXpceEl1yt2gmOkX1kOjsWbp0Kfvuuy9nnXUW++23H6eddhobNmxgyJAhfOc736GxsZE777yTWbNmcfjhh9PY2Mjpp5/OunXrAPj973/PvvvuS2NjI3fffffW7d5yyy1cdNFFACxbtoxTTz2VkSNHMnLkSP7xj38wadIkXnzxRUaNGsXFF18MwFVXXbX1/urfnXL91m19//vfZ5999uGII45g0aJFieU9H9ehp60eO7vqMMvlU6G188AkePO5ZLf5gQNgXPGmD4BFixbxs5/9jDFjxjBhwoStR879+/fnySefZMWKFXz605/mT3/6E3369GHKlCn88Ic/5Nvf/jbnnXceDz74IHvvvTdnnHFGwe1/7Wtf4+ijj+aee+5hy5YtrFu3jiuvvJJ58+ZtveZ91qxZLF68mMcffxx35+Sxx/LXf8yhz57vMH36dJ5++mk2b95MY2MjBx10UCJFo4AuIrkzePBgxowZA8DZZ5/N9dcHR8ctAXr27NksWLBg6zqbNm3i8MMPZ+HChQwdOnTrLQPOPvtspk2b1m77Dz744NabdDU0NLDjjjuyatWqVuvMmjWLWbNmceCBBwKwbs0qFr+0lLUvvcWpp57K9ttvDwRNOUlRQJeiSm9GEWmjkyPptLStuy2v+/TpAwT3KvrEJz7B7bff3mq9JEeUujuXXHIJF1xwQTBjxQuAce2vfpfYPtrKRxu6hv4noyp57vr90CubWg39z5JXXnmFRx99FIDbbruNI444otXyww47jEceeYQlS5YAsH79el544QX23Xdfli5dyosvvgjQLuC3OPbYY7nxxhuBoIN1zZo1Be+vfvPNN29tm3/tjWW8tXwlRx11FPfeey8bN25k7dq13HfffYnlOx8BPVU6SpV6kK96/uEPf5ipU6ey3377sWrVKi688MJWywcOHMgtt9zCmWeeyYgRI7Y2t/Tu3Ztp06Zx4okn0tjYyK677lpw+9dddx0PPfQQBxxwAAcddBALFiygf//+jBkzhv3335+LL76Y4447js997nMcfvjhHHDAAZw24WusXbeexsZGzjjjDEaOHMm4ceM4+OCDE8u3mlw6Y5abo5aklT/0P4k0pPyAiw5e5UrOmte6d+/Or371q1bzli5d2ur1xz72MebMmdPuvWPHjmXhwoXt5o8fP57x48cDsNtuu/Hb3/623Tq33XZbq9cTJ05k4sSJwYsVi7fOnzx5MpMnT46TlZLoCD2W+gvo9ZfjBOiHX6pMAV2Kytcxm9SLIUOGMG/evGonoyoU0EUkcXriVfm6UoY5CehpVJ6WbdbTcWo1hv6X8YCLyo79r+C+iiUhG0P/e/fuzcqVKxXUy+DurFy5kt69e5f0vlidomY2FrgOaAD+n7sXvLjUzD4D3AUc7O5zS0pJrUqq9y5L4naQRVYrtU/NyMb90Ftdz5yzjsNWEszboEGDaGpqYvny5Ylts6DmLfDOW8F093dg+eZ091eudWFaY6azd+/eDBo0qKRddBrQzawBmAp8AmgC5pjZDHdf0Ga9fsBEIH8P6qtDOrrqCpUZQI8ePRg6dGj6O1r7JvxPMNKTD46BL85Mf5/l+Pm/B/9/sboDiw4Blrj7S+6+CZgOnFJgvf8CpgDvJpg+ERGJKU5A3xN4NfK6KZy3lZk1AoPdvcOfHjM738zmmtnc1E/HpGx5bmGQvFFlhQQ6Rc2sG/BD4Fudrevu09x9tLuPHjhwYLm7jm44uW2luc1aV408l3M/9ORSEWNnVeok19D/nEm3fOME9NeAwZHXg8J5LfoB+wMPm9lS4DBghpmNTiqR1aVffqkHWa/n+iGCeAF9DjDMzIaaWU/gs8CMloXuvsbdB7j7EHcfAswGTs7VVS5QZ0cuFuvrER12X/JVLmZkIYiUNPQ/y3Uki+1rWSvvCpRxpwHd3TcDFwF/AJ4H7nD3+WZ2uZkldyNfEREpS6zr0N19JjCzzbxLi6x7TPnJklqQwWM2qVuqrZCbkaIiIpKTgK6h/8nI2AMuKnuZSyV3Vni/GRn6Xx0ZSXfKlTYnAT1F9dgpahYru61GxZf4wxcM/S9f2o/Ja7X5TveV4TqSxU7RLJd3ShTQRURyQgFditJDoiU7VFdBAV1EJDfyEdA19D8ZVRn6X8b90BNMRm3tLLpfDf3PF3WKVplO5aQeZL2e64cIFNA7Z1m/nKsrDI+RXyv6IsYeLKEHXJS9hc623+oyl45XznIVyWJ/SdbOJmph6L+IiGSDAroUlcFjNqlbqq2Qm4Ce5kjRepKxPNfFQ6IrNVK0o3mSGI0UrTb98ks9yHo91w8RZDmgV+roTEP/O1pt23Spu8CS6RRNOQ5p6H8ty1p5q1NUpLZkMvBJvVBAl+IUuyQzVFlBAV1EJDfyEdA19D8ZGvrfwc6qVB9SHfrf2f4kebrKpbCKVbx6HCmq73WXZLrQMthkkeaPXUZlN6BXSl12gsXLc3RYfKm32g2G/pevsp9OjutCjrNWMzT0X6qp1KcQiVSP6ioooIuI5EZOArqG/iej8nmOc1fHou+tjxuiRyY19D/zNPS/2nQqJ/Ug6/VcP0SQ6YBeqaH/LburowoTs/OmnKH/kND90FPuaNLQ/1qWtfJWp6iIiMSkgC5FZfKgTeqUKisooIuI5EY+ArqG/icjc0P/K5jeXA79L7C9zNb7rKRbV7kUpqH/qfLMfrGrKNNllsEmi0yXdzpiBXQzG2tmi8xsiZlNKrD8m2a2wMyeNbM/m9kHk09qldRlQ3IXrnIpsZiCof8JXOVS9hY62X7ry1xS3lsV1WU9r7BaGPpvZg3AVGAcMBw408yGt1ntKWC0u48A7gJ+kHRCpfL0FZfsUG2FeEfohwBL3P0ld98ETAdOia7g7g+5+4bw5WxgULLJFBGRzsQJ6HsCr0ZeN4XzivkS8EChBWZ2vpnNNbO5y5cvj5/KTmnofzIqn2fT0P/4+9XQ/+zL0tB/MzsbGA1cVWi5u09z99HuPnrgwIFl7q3CnaJ11gFTX7ktRUen9lkutSw2Weh+6G11j7HOa8DgyOtB4bxWzOzjwGTgaHd/L5nk1YB67CyKnefI/dBLDAhm5RybR7eTwEY62n4ld1ZNec5bzaiBTlFgDjDMzIaaWU/gs8CM6ApmdiBwE3Cyu7+VfDKlGvQdl+xQZYUYAd3dNwMXAX8AngfucPf5Zna5mZ0crnYV0Be408yeNrMZRTYnIiIpidPkgrvPBGa2mXdpZPrjCadLRERKlN2RolEa+p+MrA39r2RyNfS/xmUl3Rm6yqWiNPQ/VZn9XldTpgstg23QmS7vdGQ3oFdKXfYMWvhvx3kva+h/5N/yVPABF1kMenFlvp5nIP21MPRfRESyQQFdRCQnchLQNfQ/GZXPczn3NK9sajX0XxKQpaH/laWh/2mq6MMjciPLZZaBNuh2NPS/rQwH9ArJfGdRF8TMczl3Cg/uh16+9If+l9HzmyV5zlvNUKeoVJG+45IdqqyggC4ikhsK6CIiOZHdgO5p9v4X2mYddLpE8lyxPuCyhv5X8DOJs6/U66GG/reSZtmkRle5VJna5qQeqJ7ngQJ6Z+qyZzAcmN9J1q2My1wMwxMIIml/Ohr6nxUZSL+G/ouISFwK6CKlyPyRrORZTgK6hv4noxp5zko5a+i/JEBD/4vR0P801Vduk5LlUtOZRx5kOKBXSD2eYm/Ncyf3Q48s7+ze6YV2kUinaOpD/4u9yJl6rOcVp05RERGJSQFdRCQnFNBFRHIiuwG9YkP/6+gh0Zkb+p9gOjphNTH0P81tZ7DzX0P/28luQK+Uuuwsink/9DIfEp3I/dDT7miql17RuqznFaaRolJN+o5LdqiyggK6iEhuKKCLiOREhgN6hTtEstRZlIiK9YqW8c5KfiZx9pXyLSgqMvQ/SzKYfg39rza1zYlINiigd6Ye43nYG9pZ1su4HTpkZuh/GZfyZEme81YzdJWLiIjEFCugm9lYM1tkZkvMbFKB5b3M7Dfh8sfMbEjiKRURkQ51GtDNrAGYCowDhgNnmtnwNqt9CVjl7nsD1wBTkk6oiIh0rHuMdQ4Blrj7SwBmNh04BVgQWecU4LJw+i7gBjMzT+Gx7HPuvo6B835KN5y9wnkrZ/43ax+4NtH9NNDMYOCNNe+xO/DKDw6jOectVD38ffYE/u8Dz/P2+tEV2efur9zP0sv379J7p25x+mxsgKnbJZyqbXZudmb1XMcuy1fyQrc9+bcf/oU/Fll3+Z3fYL1NTnT/vXwTu4fTr009kfetR2LbHti8gj7h9HV/XsJE4LUbTkh0H2lqqa8ATzStY9IP/1LV9HTm8g0rOWjz87x++f6sPOjrHHTiuYnvI05A3xN4NfK6CTi02DruvtnM1gD9gRXRlczsfOB8gL322ouu6N63P29vPxSAtxjGds3rea9b7y5tqzNv8mHm7PBxDn7nTzSwOZV91Jom2593BhzKuEG7csqoPTtc9+RRe7Bi3Sa279nA8D12KGk/A/r04k+Dz2HA6mfLSS477rw99OtV1jY60gBs3vIOCzdv4bk+hzNsh378b68v00w39tswl56+ibUNO7HPxqd4dbuPpJKG1VsG0WwNidfztxlKMw283WM33tzpGB5bsZgevinRfaTtnS1BHX24/+cZtl3fKqemY0+vO5lea4PPsGffXVLZh3V2EG1mpwFj3f3c8PXngUPd/aLIOvPCdZrC1y+G66wotE2A0aNH+9y5cxPIgohI/TCzJ9y94Cl0nDaE14DBkdeDwnkF1zGz7sCOwMrSkyoiIl0VJ6DPAYaZ2VAz6wl8FpjRZp0ZwDnh9GnAg2m0n4uISHGdtqGHbeIXAX8gaFK82d3nm9nlwFx3nwH8DPilmS0B3iYI+iIiUkFxOkVx95nAzDbzLo1MvwucnmzSRESkFPm+Dk9EpI4ooIuI5IQCuohITiigi4jkRKcDi1Lbsdly4OUS3jKANiNPZSuVTXEqm+JUNsXVctl80N0HFlpQtYBeKjObW2x0VL1T2RSnsilOZVNcVstGTS4iIjmhgC4ikhNZCujTqp2AGqayKU5lU5zKprhMlk1m2tBFRKRjWTpCFxGRDiigi4jkRMUDupndbGZvhQ/FaJm3i5n90cwWh//vHM43M7s+fPj0s2bWGHnPOeH6i83snMj8g8zsufA915uZVTaHXVOkXK4ys4Vh3u8xs50iyy4J87jIzI6PzC/4QO/w9sePhfN/E94KOTMKlU9k2bfMzM1sQPi6buoNFC8bM/tqWH/mm9kPIvPrpu4U+V6NMrPZZva0mc01s0PC+dmvN+5e0T/gKKARmBeZ9wNgUjg9CZgSTp8APAAYcBjwWDh/F+Cl8P+dw+mdw2WPh+ta+N5xlc5jguVyHNA9nJ4SKZfhwDNAL2Ao8CLBrY0bwukPAT3DdYaH77kD+Gw4/RPgwmrnudzyCecPJri188vAgHqrNx3UnY8CfwJ6ha93rce6U6RsZrV8vmFdeTgv9abiR+ju/leCe6ZHnQL8Ipz+BfCpyPxbPTAb2MnMdgeOB/7o7m+7+yrgj8DYcNkO7j7bg9K+NbKtmlaoXNx9lru3PMx0NsHToiAol+nu/p67/xNYQvAw760P9Hb3TcB04JTwqOFjBA/whtZlnAlF6g3ANcC3gWjvft3UGyhaNhcCV7r7e+E6b4Xz66ruFCkbB1oegrsj8Ho4nfl6Uytt6Lu5+xvh9JvAbuF0oQdU79nJ/KYC8/NgAsERAJReLv2B1ZEfh1yUi5mdArzm7s+0WaR6A/sAR4ZNJX8xs4PD+ao78HXgKjN7FbgauCScn/l6UysBfavwl07XUkaY2WRgM/DraqelVpjZ9sB/AJd2tm6d6k7QRHAYcDFwR9Xbd2vHhcA33H0w8A2CJ67lQq0E9GXh6Qvh/y2nh8UeUN3R/EEF5meWmY0HTgLOCn/soPRyWUlw+ti9zfws+xeCNuBnzGwpQZ6eNLMPoHoDwdHi3WHzweNAM8ENp1R3gucf3x1O30nQ3AR5qDdV6qgYQutOiqto3Sn6g3D6RFp3Ujwe6aT4J0EHxc7h9C5FOilOqGYnRZnlMhZYAAxss95HaN2x9RJBp1b3cHoo2zq2PhK+505ad2x9pdr5Lbd82ixbyrZO0bqqN0XqzpeBy8PpfQiaDKwe606BsnkeOCacPhZ4Ii/1phqFezvwBvA+wVHElwja6f4MLCbomW8pLAOmEvS+PweMjmxnAkGHzhLgi5H5o4F54XtuIBwNW+t/RcplSfhFfDr8+0lk/clhHhcR6Vkn6Kl/IVw2OTL/Q2HlWxJ+QXtVO8/llk+b5UvZFtDrpt50UHd6Ar8K8/Qk8LF6rDtFyuYI4AmCH63HgIPyUm809F9EJCdqpQ1dRETKpIAuIpITCugiIjmhgC4ikhMK6CIiOaGALiKSEwroIiI58f8BM+W2HmxfVbsAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from utils.metrics import f1_score, far, mar\n",
    "metrics= {}\n",
    "predictions = []\n",
    "\n",
    "draw = 1\n",
    "\n",
    "all_prediction = {\"swat\": []}\n",
    "all_test = {\"swat\": []}\n",
    "\n",
    "for file_name, prediction in predictions_by_names.items():\n",
    "    y_test = y_test_by_names[file_name]\n",
    "    f1 = f1_score(y_test, prediction)\n",
    "    far_score = far(y_test, prediction)\n",
    "    mar_score = mar(y_test, prediction)\n",
    "    all_prediction[\"swat\"].append(prediction)\n",
    "    all_test[\"swat\"].append(prediction)\n",
    "            \n",
    "    metrics[file_name] = (f1, far_score, mar_score)\n",
    "    \n",
    "    if draw:\n",
    "        title = file_name + \" f1 = {:3.2f} far = {:3.2f} mar = {:3.2f}\".format(f1, far_score, mar_score)\n",
    "        plt.title(title)\n",
    "        plt.plot(y_test)    \n",
    "        plt.plot(prediction)\n",
    "        plt.legend([\"y_test\", \"predicted\"])\n",
    "        plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "6b6af5af",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "metric | swat\n",
       "--- | ---\n",
       "F1 |0.36\n",
       "FAR |0.23\n",
       "MAR |0.52"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from utils.data import show_score_table\n",
    "\n",
    "show_score_table(metrics)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "32f7a2a3",
   "metadata": {},
   "source": [
    "### SKAB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "id": "5e96cbe6",
   "metadata": {},
   "outputs": [],
   "source": [
    "# libraries importing\n",
    "import numpy as np\n",
    "from tensorflow import keras\n",
    "from tensorflow.keras import layers\n",
    "from tensorflow.keras.callbacks import EarlyStopping, ModelCheckpoint, ReduceLROnPlateau\n",
    "from sklearn.metrics import mean_absolute_error, mean_squared_error\n",
    "from sklearn.preprocessing import StandardScaler"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "id": "e81c5125",
   "metadata": {},
   "outputs": [],
   "source": [
    "def arch(data):\n",
    "    EPOCHS = 3\n",
    "    BATCH_SIZE = 16\n",
    "    VAL_SPLIT = 0.1\n",
    "    \n",
    "    # model defining\n",
    "    # define encoder\n",
    "    inputs = keras.Input(shape=(data.shape[1], data.shape[2]))\n",
    "    encoded = layers.LSTM(100, activation='relu')(inputs)\n",
    "\n",
    "    # define reconstruct decoder\n",
    "    decoded = layers.RepeatVector(data.shape[1])(encoded)\n",
    "    decoded = layers.LSTM(100, activation='relu', return_sequences=True)(decoded)\n",
    "    decoded = layers.TimeDistributed(layers.Dense(data.shape[2]))(decoded)\n",
    "\n",
    "    # tie it together\n",
    "    model = keras.Model(inputs, decoded)\n",
    "    encoder = keras.Model(inputs, encoded)\n",
    "\n",
    "    model.compile(optimizer='adam', loss='mae', metrics=[\"mse\"])\n",
    "    \n",
    "    # callbacks defining\n",
    "    early_stopping = EarlyStopping(patience=5, verbose=0)\n",
    "#     reduce_lr = ReduceLROnPlateau(factor=0.1, patience=5, min_lr=0.0001, verbose=0)\n",
    "    \n",
    "    # model fitting\n",
    "    history = model.fit(data, data,\n",
    "                        validation_split=VAL_SPLIT,\n",
    "                        epochs=EPOCHS,\n",
    "                        batch_size=BATCH_SIZE,\n",
    "                        verbose=0,\n",
    "                        shuffle=False,\n",
    "                        callbacks=[early_stopping]#, reduce_lr]\n",
    "                        )\n",
    "    return history, model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "id": "f4ea2f1d",
   "metadata": {},
   "outputs": [],
   "source": [
    "from utils.pipeline import Pipeline, Predictor\n",
    "from tensorflow import keras\n",
    "\n",
    "optimizer = keras.optimizers.Adam(learning_rate=0.0001)\n",
    "class Config:\n",
    "    TIME_STEPS = 50\n",
    "    epochs = 3\n",
    "    batch_size = 10 \n",
    "    Q = 0.99\n",
    "    \n",
    "conf = Config\n",
    "\n",
    "class LSTM_AE_SKAB_Predictor:\n",
    "    def fit(self, df_train):\n",
    "        ss = StandardScaler()\n",
    "        self.ss = ss\n",
    "        x_train = df_train[:train_size].drop([\"anomaly\"], axis=1)\n",
    "        x_train = np.array(ss.fit_transform(x_train))\n",
    "        x_train = create_sequences(x_train, conf.TIME_STEPS)\n",
    "        \n",
    "        _, model = arch(x_train)# ~ fit\n",
    "        self.model = model\n",
    "        \n",
    "        healthy_r = pd.Series(np.sum(np.mean(np.abs(self.model.predict(x_train) - x_train), axis=1), axis=1)) \n",
    "        self.UCL = healthy_r.quantile(conf.Q)\n",
    "            \n",
    "    def predict(self, dfs_test: Dict[str, pd.DataFrame]) -> Dict[str, pd.Series]:\n",
    "        predictions = {}\n",
    "        for file_name, df in dfs_test.items():\n",
    "            y_test = df.anomaly\n",
    "            df = df[important_features]\n",
    "            df.drop([\"anomaly\"], axis=1, inplace=True)\n",
    "            x_test = np.array(self.ss.transform(df))\n",
    "            x_test = create_sequences(x_test, conf.TIME_STEPS)\n",
    "            r = pd.Series(np.sum(np.mean(np.abs(self.model.predict(x_test) - x_test), axis=1), axis=1)) \n",
    "            \n",
    "            smoothed = r.rolling(window = 10).mean().fillna(r.iloc[0])\n",
    "            a = smoothed.diff(10).abs()\n",
    "            \n",
    "            IQR = a.quantile(0.75) - a.quantile(0.75)\n",
    "            UCL = a.quantile(0.75) + 1.5 * IQR\n",
    "            \n",
    "            s = a > UCL\n",
    "            w = 10\n",
    "            is_anomaly = s.rolling(window = w).mean().fillna(s.iloc[0])\n",
    "            anomalous_data_indices = pd.Series(is_anomaly.index[is_anomaly > 0]) #+ conf.TIME_STEPS - 1\n",
    "\n",
    "            prediction = pd.Series(data=0, index=df.index)\n",
    "            prediction.iloc[anomalous_data_indices] = 1\n",
    "            predictions[file_name] = prediction\n",
    "            \n",
    "        return predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "id": "4869f5f6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 7min 53s, sys: 34.5 s, total: 8min 27s\n",
      "Wall time: 4min 4s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "skab_predictor = LSTM_AE_SKAB_Predictor()\n",
    "skab_predictor.fit(df[:train_size])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "id": "0ae186f3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 29.4 s, sys: 2.83 s, total: 32.3 s\n",
      "Wall time: 16.5 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "predictions_by_names = skab_predictor.predict(df_by_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "id": "67123a9f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEICAYAAABPgw/pAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAwe0lEQVR4nO2deZwV1Zm/n5fekH1rFwQFI2hcEBtQDC6JJgLRiBoTlyyiURIzTnSSMYPBUeOMRtRfEk3cmIkxJipRx4UYE4lRo1ER0CAiKqBCaBVpdpqt6e7390dVQ/Vd695b996quu/z+TRU1dnec+rUe6vOqW8dUVUMwzCM6NOl3AYYhmEYwWAO3TAMIyaYQzcMw4gJ5tANwzBigjl0wzCMmGAO3TAMIyaYQ48p4vBrEVkvInPLbU+xEZG9ROQFEdksIv+v3PYYRjkwhx5BRGSIiKiIVGeIdizwBWCQqh4lIrUi8oiILHfTfrZItvUTkcdEZIuIrBCR83ykqRWRt0WkMeH4SBF5TUS2uv+PzJDNFGAN0EtVf1BYLQpDRE4SkXdcu58Tkf3TxNtTRB4UkY9EZKOIvCQiRyfEqReRB9zw9SJyf2lqUTyC6iMiMlxEnhCRJhFZJyJPi8hBxbU+3JhDjy/7A8tVdYvn2N+BrwOrilju7UALsBfwNeBOETk0S5orgCbvARGpBZ4Afgf0BX4DPOEeT8X+wGLNQymX5Ycx17wGAI8C/wn0A+YDv08TvQcwDxjlxv0N8EcR6eGJ8yjO+doP2BO4JShbc0FEqgLMLpA+AvQBZgEHuXnNxekzlYuq2l8R/oD/AD4ENgPvAicBXYFtwAA3zjSgFeeuEuC/gJ+726cA/wA2ASuBaz15/xNQoNn9Oyah7G8B24E2N/zHCeGNwGeLUOfuOBfqcM+x3wI3ZkgzFHgbmAg0eo6f7LafJNR7Qoo87gV2umU3A58HjgJeATYAHwO/BGo9aRT4F2Ap8EGAbTAFeDmhTbYBB/tMvwkY5WmD5UCVz7TLcRzfQmAL8CscR/cntx8+A/T1xH8Y58diI/ACcGhCm94JPOXm9fmw9ZEU8fq557V/0H07Kn92h14E3Me+S4ExqtoTGI9zt7wd547sBDfqCcAKYJxn/2/u9hbgmzh3IacAl4jI6W7Y8e7/fVS1h6q+4i1fVX8FfAd4xQ2/Jo863CEiG9L8LUyTbDjQqqpLPMfeADLdff0C+BGO0/NyKLBQ3SvVZWGqvFR1MnA/cJNb32dwfsz+DRgAHIPzg/rdhKSnA0cDh6QyLEP9N4jI1DT1ORSnzh22bQHeS2V3ivJGArXAMvfQWJybgd+IyFoRmSciJ6RL7/JlnKG24cCXcJz5j4B6nCfy73ni/gkYhnPn/zpOG3o5D7ge6InzdJdob7n7SCLHA6tUdW2WeLElsEdNoxNtQB1wiIg0qepyT9jfgBNE5AlgBPATd/85YAzOnRKq+rwnzUIReRDH4T9edOud8r9LsgPMRg+cO0wvG3EcQhIicgbO3edjKcb0e7hpfeWViKq+5tldLiJ347Tfzz3Hf6Kq6zLk0cdPWQn0IHloIKvdItIL5071x6raUe9BOHfpFwEX4DjrJ0TkQFVdkyarX6jqJ26eLwKrVfUf7v5jOD9sAKjqPZ7yrwXWi0hvT/lPqOpL7vb2xIJC0Ee88QbhDOV8P0d7YoXdoRcBVV0GXA5cC6wWkZkiMtAN/hvwWaABeBP4C46jGQss67i7EJGj3Qm1JhHZiHPHPaCU9ciDZqBXwrFeOI/7nRCR7sBNdL5jzCuvVLgTZk+KyCoR2QTcQHL7rfSTV47kbLeI7AH8AZijqj/xBG3DebL7laruVNWZODaPS5WPyycJ6RP3e7hlVonIjSLynts+y9043jYqa/v46CMd8eqB2cAdqvpgQHZGEnPoRUJVH1DVY3Em6xSY7ga9jDOJcwbwN1VdjDPh9UV2D7cAPIAz4TNYVXsDdwHSkX2x7ReRu0SkOc3fW2mSLQGqRWSY59gRQKr4w4AhwIsisgpn8m8f1wEPcdOMEBHxpBmRJq9U3Am8AwxT1V44j+ySECdjO2aof7OI/ChNsrdw6tyRR3fgU+nsFpE6nKeuRuDbCcELU9gY1Lk/D5iEM9/QG+dcQOc2ytY+5e4jiEhfHGc+S1Wvz2RvJWAOvQiIyEEicqJ7sW7HuTNqB1DVrcBrOBNyHQ78ZZw7cK9D7wmsU9XtInIUzgXYQZOb3wE52lUnIl3d3VoR6ZrgMHehqt9xx6NT/aUc73THix8FrhOR7iIyDsdp/DZF9EXAYGCk+3cRzt3kSJw7w+dxhq6+59p9qZvuWZ/V7YnzaN8sIgcDl/hM561Puvr3UNUb0iR7DDhMRL7stvXVOHMB7yRGFJEa4BGc/nG+qranyKuviJzv3lGfhTMM8xKF0xPYAawFuuE8weREufuIO0z1NPCSqqab06gozKEXhzrgRpz3olfhTDpd6Qn/G1CD85pVx35P3PFzl+/idPrNOE7hoY4A90fheuAldwJqrE+73sVxHvviXAjbcJ4gguS7wB7AauBB4BJVfQtARI4TkWYAVW1V1VUdf8A6oN3db1PVFpxJy2/ivKlyIXC6e9wP/47zI7gZ+B/SvzoYKKrahDPWfT2wHmfS9ZyOcPeu9i539zPAqTjj5Bs8d7fHuXmtA05z67IRmApMyjB+ngv34UzIfwgsBuYEkKdfAukjOE+5Y4ALEp4O9ithXUKFdH6JwDAMw4gqdoduGIYRE8yhG4ZhxARz6IZhGDHBHLphGEZMKJtSdMCAATpkyJByFW8YhhFJXnvttTWqWp8qrGwOfciQIcyfP79cxRuGYUQSEVmRLsyGXAzDMGKCOXTDMIyYYA7dMAwjJphDNwzDiAnm0A3DMGJCVocuIveIyGoRWZQmXETkNhFZJiILRaQheDMNwzCMbPi5Q78XmJAhfCLOd4uH4ayneGfhZhmGYRi5kvU9dFV9oeNj8mmYBNznrv04R0T6iMg+qvpxUEampL0NXr0LavaA0RcGm/fSZ2De/8KWJhg0GnoPhrGXQJeAFj5//T7Y+3AYeGR+6d98BHoNhP0/E4w9HXS06bYNfLhhGyvXbc2aZFC/bqxr3sHWlraci6uqruHwSZfRte/A7JE72dkOr97Fym21PNx6HAD7bnydftv/yZt7nb4r2qfWPk8XbWPpgJOSsui7dTkHr5nN2m5DGbb2Ofbd/AYtVd1YsPdZLNjnbAC6tLdy5Me/Z0O/w5h4UB9qP5oHmz+GjxbAQRM7Z7jyVfjUiXDs5c7+/Htgn5Gwb5YH1vY2ePVuOPDzUD88OXxHM8ydAfUHwcGnpM9n+yZ47V4Y8VXouXfmMjtY/ASs8jx41/WAo78D1XXO/puPwB594cDk9isZqk799zoMhrgLNf3957D873DiVbBjE3zwotM+a5Y69u7YDG0t0LwKeiS0xY5NUJe4YFJgxsKCB+HwL0NVXeoon/4S7DOiSOX7/Hyu69CfVNXDUoQ9ibNi99/d/b8C/6GqSaohEZmCcxfPfvvtN2rFirTvx2fnk7fgTteh/fsy6JFSOJUf04fCtoSlJr/9YjAnYvsmuHEwdK+HK5Zlj5+Ka3u7/ycuuVkgnyyGO48BoB3xtzaOkNcaOl3ESfTB6KsYeuoVuSVe/Tbc4XwC/qDt99IitbxX+zW6iDJyxww2Oqus8UGdsybI0B0PJGVxXfU9fKPqmZTZd8QfKct4rPZqNmk3qvsMpNvGxPOVYgGpaaugZQvc/Cl/53jVIrhrHAw+Gr41Ozn83T/Bg+7n1DOd74UPwaMXw9h/gQk+16q4cX/YvoFOJ3HyU7sdZ7H6WS6sex9uOxK61MDVa3ZfPwBDjnPa+qPXs2SSaqGvlOu6FEi2/BUOPRO+8uuCShGR11R1dKqwkipFVXUGMANg9OjRhX2Ivd1zR6i53x1mZGeKxcWDKqMjny2J6wiHgA7bzr6frzzfjz1qqvjdRUenjT7513NZ9OEm1jTv4KazRvDV0YN9F/X6khU0PDACbc+jXT1pPr13Tx6//PNwrdOdFlx1EnR3l8W81vnvg5+kuLP9w2xn3agU7Ir/z1fhHuglW9na3gaHfRkW/Z8TNu0TqHEXf9q2HqYPcbZVd9vn5xx3tPn65anD/bZPm7vux/YN/uIDaDuM/S5M+Ilzx3vvKcFfS4XS7i7i1L7T+d9rX8uW9PaedA389cdQ3RWucpdVnf2f8PJtMPQEOH9W8LbO+xX88fsw6gL40s+Tw28/uujtG8RbLh/iLBPVwSD3mGEYhlFCgnDos4Bvum+7jAU2Fn383DAMw0gi65CLiDwIfBYYICKNwDU462GiqncBT+GsWL8M2ApcUCxjDcMwjPT4ecvl3CzhirOCfYnxDMGXYl3UoMoI8xquHtv8rzVbWH3yW9N2dxpJLN9vfr7iaZrtbGlyqFNWO4rYX1KVHbr+meX8ZrM3bPUpsj2mFDUKImzXi2EAIMV4iyX8mEM3knEvhmzXhKTZ9lVEl2AuuNJet+kKi7LzkIT/44JbH28HkRTHSkrxyzWHbhiGERPMoRuGYcQEc+iGB02x5TdF1AiL5WGxo4OQ2VPoJE2FTfJE16FrPm8gFFRgCcqIHvm9pVJwoZkCi1NOLm/PBNkmRW3fOPTpNHUoyzh5rm9OBU90HbpRRPxdDOK5aCTnCyiYrlfS6zZdYVF+o6LsE4XFIlW9SjQBXMZ+Yg7dMAwjJphDNwzDiAnm0A3DMGJChB16qaX/QeUT4omoTtJ/n0kKLrM9n0S7tkz6XwCxkP5nSx6y+pj03ygXfiY6w3a9GPkSs0nR2E3y+sMcupGMz4uhIOl/QBeclNQRxVH6H1dM+m8YhmFEGHPohmEYMcEcuuEhD+l/ZAfRw2J3WOzoIGT2mPQ/J6Lr0CMr/Y9XByv8LZc8ctAMb7kUSfovlSb9D60j9Hm+bVLUMDoQz78ZYqWYb/JfRFCToiUkXYWj7DzKPlFYLMoo/U+HSf8NwzAMv5hDNwzDiAnm0I3d5PW52OKYUnTCYndY7OggdGPnNimaCxF26KWW/gdURsw6WOG1ySeHUk3mlUD6n5MNAZOprcLaT5Ok/yG1Mx0m/TfKRezmyYwMxOxkV2jnNYduJOP7YvAucJFrEQEtcBFILoWWFmXnEWXbM2HSf8MwDCPCmEM3DMOICebQDQ8m/XeCSlmnsLVfyOwx6X9ORNeha9qdUhRYvnxK2EH9jPhVivQ/pzRRlP4njSuH1RGGWPqf86IpwRNdh24UEZ/fQ/fON+U44RPY5VbKC7eMq7kXjSjbnpFySv/T9ZMiF4s5dMMwjNhgDt0wDCMm+HLoIjJBRN4VkWUiMjVF+H4i8pyI/ENEForIF4M31Sg6+Sj/wzrUmo2wqCTD1oAhMydKc05hIKtDF5Eq4HZgInAIcK6IHJIQ7SrgIVU9EjgHuCNoQ5OpUOl/yDqoFnjB5feWTIZJ0WJJ/33na9L/omLS/4z4uUM/Climqu+ragswE5iUEEeBXu52b+Cj4Ew0ykVQCzkbUSBm57pC+64fh74vsNKz3+ge83It8HURaQSeAv41VUYiMkVE5ovI/KampjzMNUqCz4uh0/sDuUr/u0RwgQuT/kcPk/7nxbnAvao6CPgi8FtJ8bEOVZ2hqqNVdXR9fX1ARRuGYRjgz6F/CAz27A9yj3n5FvAQgKq+AnQFBgRhoFFKvEpRf2N9URvC3E0mw00pGhpiNudUbPw49HnAMBEZKiK1OJOesxLi/BM4CUBEPo3j0G1MxTCM4lGh4+SZyOrQVbUVuBR4Gngb522Wt0TkOhE5zY32A+BiEXkDeBCYrMX+yEen7CtI+l/CO6jwSv8zWVAs2X0O73Ga9L+IJL55FCLpvy+K267VvkxQfQpnstN77GrP9mJgXLCmGeUjd+l/7iUE9D30Ul63KWXkpTYiYKJse0bKKf1PQwna2pSihmEYMcEcumEYRkwwh27sxjNem5MwMoqERSUZtrcwwmaPSf9zIsIO3aT/xcbPkF+h0v/88Er/E4PKLf1PSBdmOtVJMoSFCJP+ZyTCDt0wDCMNsZ3szYw5dCMZn9eCd1GLnL/7YtL/kBBl2zOR4i0Xk/4bhmEYUcEcuuEh90nRqA1h7sak/6kJmT0RmnMKA9F16EVViqbKr/KUon4ouDbank+iDBaYUjSHzPMMKydFPN8lwSZFjbIR1/FVI4m4nWqbFDV2U5mdYTe5S/9zb7GAJkVLeeGmnRONcH+Jsu0ZKaP0P12bmvTfMAzD8Is5dGM3XqWo7yRRG8N0CZtSNKjJv9hNIobNnnBjDt0wjPgSuh+o4hJhh27S/2LjT/pfKPnk4JX+F1MKnqGcjGnysCFtY/vMKx/RTMVK/8tYN5P+G4Zh5EhsJ3szYw7dSMbnxZBKVe2/CJP+h4Mo256JDNL/stXZ3nIxDMMwfGIO3fDglf77G+sL61BrdsKikgxbA4bMngjNOYWB6Dr0kkuiK0/6X5IH07wWic40WVkk2X0uH7eJovQ/iotEp9wPOzYpahiGkRs2KWrspjI7w258Top6v4eeY5uJBNP1SnrdppsFjrLziLLtGckk/S8TJv03DMMw/GIO3dhNJU0ghU36H7d8AsMmRXPBHLphGNEktsNF+RNhh27S/2JTmuulQOl20ksPxZH+5/C5sgDLL0Z+3qwjLv3P+40ik/4bhmFEhwq9ezeHbiTj92IIhfS/pK+55Hg8CkTZ9kyY9N+oeHJfJDq6hEwpWnCDB5xPWIh/RwwUc+iGYcSXCvtB8OXQRWSCiLwrIstEZGqaOF8VkcUi8paIPBCsmSnIZ0V2/5n7PBZU3qVM75+SDGeEWvqfT5o8J+rSDUHlmldOQ1km/S89xbW3OlsEEakCbge+ADQC80Rklqou9sQZBlwJjFPV9SKyZ7EMNgzDyEo5JkVD8DTg5w79KGCZqr6vqi3ATGBSQpyLgdtVdT2Aqq4O1sxSE9eJIr/4lf7nmsITv0tAbVwu6b+f41EgyrZnpJzS//L1Ez8OfV9gpWe/0T3mZTgwXEReEpE5IjIhVUYiMkVE5ovI/KampvwsNgzDMFIS1KRoNTAM+CxwLvA/ItInMZKqzlDV0ao6ur6+PqCijcDwPDJq5MYmc8Sk/8XNJzCiI8QLA34c+ofAYM/+IPeYl0ZglqruVNUPgCU4Dt6IMLF9GjdSEMGTbR00iayTosA8YJiIDMVx5OcA5yXEeRznzvzXIjIAZwjm/QDtTIFJ/+NBodJtH6vAq+Z58efTxzTZpoIo0/kuoJ/t3LmTxsZGtm/fXrgdrdUw/iFn++23ob11935VHWg7tO9MkbC/G0+cdAC9jnOO1eyx+1iQ1I108q/rmTr/EVc5/dBn2V27dmXQoEHU1NT4NiGrQ1fVVhG5FHgaqALuUdW3ROQ6YL6qznLDThaRxUAbcIWqrvVthWEYsaGxsZGePXsyZMiQwhXBLVtgTbuzPfDT0LoDVrc6+zXdHAff1pKcrs/+sKEKECcdwOZVsPljqOsF/T9VmF2paG6CTTXQbQD0GZwc3tQFpAoGHJg1K1Vl7dq1NDY2MnToUN8m+LlDR1WfAp5KOHa1Z1uB77t/RtQRwc+dofdizVn6H9Ajfmkfuk3674ft27cH48xjh//2EBH69+9Pri+PmFLU8GDS/+xhQRMWwVqw+ZgzL5x82tAcumEYRkyIrkM36X/RKclNVr7yepdk6X/m+EUtpyNNkI83xXpUSsw3MtL/4Fi+fDkPPJD/V0puuOGGAK0Jhug6dMMwjAIwh14xVPr4X+7S/1zbTATatfB2LulQbToZeZTHi6NsexquvvlOfv7Lu3btT5s2jVtvvTUp3tSpU3nxxRcZOXIkP/vZz2hra+OKK65gzJgxjBgxgrvvvhuAjz/+mOOPP56RI0dy2GGH8eKLLzJ16lS2bdvGyJEj+drXvlayumXD11suRoXQSSkac8KmFA1KnxAyncOP//AWiz/alH8G2gY7tznbta+AtnNIn51cc3zvtEkuPGcSZ377R1z+9Ym0t7czc+ZM5s6dC3R+X/3GG2/klltu4cknnwRgxowZ9O7dm3nz5rFjxw7GjRvHySefzKOPPsr48eOZNm0abW1tbN26leOOO45f/vKXLFiwIP+6FQFz6IZhxIohgwfSv19f/rHoHT5pWseRRx5J//79nffQMzB79mwWLlzII488AsDGjRtZunQpY8aM4cILL2Tnzp2cfvrpjBw5sgS1yA9z6EZaSru8W4UT1NBH3vkU51xf86VDC8ugZQusWeJsDzzSFRYtzpwGuGjyN7j3oT+wavVaLvzOZb6KUlV+8YtfMH78+KSwF154gT/+8Y9MnjyZ73//+3zzm9/MqRqlIsJj6Cb9jwOq7fmk2rWV5IZStk++bRZj6X+2+kS8n50x6VT+/NzLzHvjrZQOGqBnz55s3rx51/748eO588472bnTGZpZsmQJW7ZsYcWKFey1115cfPHFXHTRRbz++usA1NTU7IrrnzIvcGFUID7v8lKuv+u3CILp2qFYJDrSE4tRtj09tbW1fO4zo+nTuydVVVUp44wYMYKqqiqOOOIIJk+ezGWXXcby5ctpaGhAVamvr+fxxx/n+eef5+abb6ampoYePXpw3333ATBlyhRGjBhBQ0MD999/fymrlxZz6IZhxI729nbmvP4mD999U9o4NTU1PPvss52O3XDDDUmvI55//vmcf/75SemnT5/O9OnTgzE4ICI85GIEj1f6H+1H7uyY9L+4+ZSPxUve58DDRnHSsUcx7ID9y21OSbE7dCM98XwaN1IR6WGjzhwy/ADeX/wPWL8cgDfffJNvfOMb7pcZd0KXKuq69eTVV18tr6FFILoO3aT/MSGP+mQ892m+h54P+fSxqEj/s9YnPv3s8MMPd94XL/bnc0OADbkYhmHEBHPoKYnP42d+5C79z7XFREADaOfySf9jQhzrVG7K+Nl8c+jGbkz6nz2slHZEOR+jLJhDNwzDiAnm0I202MN4JWFnOxXPP/88p556KgCzZs3ixhtvTBt3w4aN3HHvQzmXce2113LLLbfkbaOXCDt0k/7Hgryqk2HhiWJJ/33nYdL/4uLXrsw/UG1tbTmXfNpppzF16tS04Rs2buSO+x7OOd8gibBDN4qGb+m/pNz2WUgwk6IF5xDe0kpD/Oq0fOVHHDzyKL526TQ+fcKZnHXWWWzdupUhh43hP66/lYbPncbDDz/M7NmzOeaYY2hoaOArX/kKzc3NAPz5z3/m4IMPpqGhgUcffXRXvvfeey+XXnopAJ988glnnHEGRxxxBEcccQQvv/wyU6/5b95b0cjIYydwxRVXAHDzzTfv+r76NdNv25XX9ddfz/Dhwzn22GN59913A6t7dN9DN4qAptyMJyFTihZ8RxxwPkHxp6mw6s3802sb7NzqbNf2BG2HXvvAZ/41Y7J3lyzlVzddybgxR3LhVbdxxx13ANC/b29ef24Wa7Q3Z555Js888wzdu3dn+vTp/PSnP+WHP/whF198Mc8++ywHHnggZ599dsr8v/e973HCCSfw2GOP0dbWRnNzMzf++CoWLXqTBX//M/QezOzZs1m6dClz585FVTltwkm88PI8uu+7iZkzZ7JgwQJaW1tpaGhg1KhR+beRB3PohmHEjsGD9mXcmJEAfP3rX+e225y747NPOxmAOXPmsHjxYsaNGwdAS0sLxxxzDO+88w5Dhw5l2LBhu9LOmDEjKf9nn31210e6qqqq6N27N+sbO8eZPXs2s2fP5sgjjwSgeeN6lr6/nM3vr+aMM86gW7dugDOUExTm0I205D6MYuRNub+HXqxzPTH9JKIvdjTD2qXO9sAjoXU7rH47a7LEvtux373bHoDzraIvfOELPPjgg53iBbkCkapy5ZVX8u1vf9s5sGYJIPz8d38MrIxEojuGHlnpf6GExY6gyON76GWQ/ovfPEz6Hwr+ubKRV+a/AcADDzzAscce2yl87NixvPTSSyxbtgyALVu2sGTJEg4++GCWL1/Oe++9B5Dk8Ds46aSTuPPOOwFngnXjxo307NGDzc1bdsUZP34899xzz66x+Q8//oTVTWs5/vjjefzxx9m2bRubN2/mD3/4Q2D1jq5DNwzDSMNBw4dx+28e4tMnnMn69eu55JJLOoXX19dz7733cu655zJixIhdwy1du3ZlxowZnHLKKTQ0NLDnnnumzP/WW2/lueee4/DDD2fUqFEsXryY/v37MW7MSA475vNcccUVnHzyyZx33nkcc8wxHH744Zx14ffY3LyFhoYGzj77bI444ggmTpzImDFjAqu3DbmkpNKHGkol/S+cki5wUciKHmElJtVIpLq6it/94npAYOBIAJYvmg+bP9oV58QTT2TevHlJaSdMmMA777yTdHzy5MlMnjwZgL322osnnniic4TmJh64/QboPgB6Dwbgsssu47LL3CXw1izdFXXatGlMmzYt/wqmwe7Qjd2Y9D97WCntiHI+Rlkwh26kJaY3b0ZK4nO2hwweyKL5r5TbjLJgDt0wjMCJ/4pXxSefNoywQy+19D+ofEz634m86uOV/vvJz6T/ydkWT/rftWtX1q5dWxynnvRSU8yuBxdVZe3atXTt2jWndL4mRUVkAnArUAX8r6qmfLlURL4MPAKMUdX5OVlihAe/E34FzBEKEfweeoyGJXYTfJ0GDRpEY2MjTU1NhWfWugOaVzvbG992lpDb7O5X1TlK0vbW5HRrFLY0OR1ko/ve+o5NsG0D1GyG1S2F25bIjs2wbT3UbYM9mpPDO+rRlMLeFHTt2pVBgwblZEJWhy4iVcDtwBeARmCeiMxS1cUJ8XoClwHxW6ivYqikRaIzEcVFooOicHtqamoYOnRoALYAK16B//uqs33tRlizDB5x9wc2OE5yU2NyunMehEfPdZz+f7qO9KVb4S9Xw/CJcN7MYOzz8urd8PQP4agp8MWbk8N//e/O/xeUV1h0FLBMVd9X1RZgJjApRbz/AqYD2wO0zzAMIzVxeXU0QPw49H2BlZ79RvfYLkSkARisqhl/ekRkiojMF5H5gTyOGUXFrpcKwk52LCh4UlREugA/BX6QLa6qzlDV0ao6ur6+vrCCIyv9D9nX8MpNXnOipZf+55QmFtL/sJL4MkTU6lFce/049A+BwZ79Qe6xDnoChwHPi8hyYCwwS0RGB2WkYRhGTpTjiSME805+HPo8YJiIDBWRWuAcYFZHoKpuVNUBqjpEVYcAc4DTov2WS4U+fu7qkEnrAKXEK7vP+S0XEcrazvkoRYvhJDqKCup11lzzSaxTCJxSZ/K1x61Xp/qVqr+lKacEPzJZHbqqtgKXAk8DbwMPqepbInKdiAT3IV/DMIygCd0PVHHx9R66qj4FPJVw7Oo0cT9buFlGGKjQ55TyUO7voUfxbNtEbhIRVooahmEYXiLs0Est/Q/J1+xi9wiZxwIXnaT/Ce2Rciw5AOm/73Y36X9R0YRzUsr5jiAocrtG2KEbweN2NhFf/a7T58FzfGQXgnFV+S+Tl0/pxXAS5V4kOqSOr4O82yXFpGhYnXyAmEM3DCO+hPVJo0iYQzfSYotEl5ByT4pG8VxH0eYiYw7dMAwjJkTXoWvanaAzD7iMynoEzEq+8vrdO4mByXECkP77k1oRI+l/WPupSf8zEV2HbhiGkY4KHY4xh56SyuwMu28eBPVxJyFpd7IjEtACF/kmLFj6H1AfKffrsHGX/mc9VkLCIP03DMMIJRV6F54Jc+hGWuxyqSTsbMeBCDv0YipFA/ymdtD5hO6RuBz4WSQ6YKWoLRIdDvwqRcOKKUUNwzBypEIfOMyhp6RCe0NB0v/cEKS8k6IZ73zThZVZ+p8xTlDS/5Dd8cZK+m+TooZhGPkTtSGZAjGHbqSn3Dc0lYRJ/3MnijYXGXPohmEYMSG6Dj2Q712nzdznsaDyLmX6kJGX9D/tDinHkgOQ/vvOw6T/Rcak/5mIrkM3gsfjPGI/9JiPUrSYdhQ6KZrvItG5lFEWAlLAloTyt5059JRU+ticv/p7F7XI9VO7jvS/cEp6puI4ZhvHOgFllf6na1OT/hvlJNdViIwCKPekaBSppLr6xBy6YRhGTIiwQzfpfxzw81XHVKk6SLtIdNmk/4WUV0JM+l8eTPpvGIZh+MEcekoqdWyuQ/rvL3Yh0n+ntACk/3mPo+Yj/S8GQd35B/Xd8JDd8Qb1nfd0x0qKTYoahmGkpuwOOnyYQzfSYtdLBWEnOxaYQzcMw4gJ0XXoJv2PB3lJ/zOd++JI/33fv5r0v8gEIf0v8qIhfj5zXCSi69CN4Okk/Q/rBR0QJv3PvYyyEDZ7wo0vhy4iE0TkXRFZJiJTU4R/X0QWi8hCEfmriOwfvKlG6fAp/S9g7QBH+l/OBS7yKSxMiyUERFzqkUQm6X+R6tzRlmGW/otIFXA7MBE4BDhXRA5JiPYPYLSqjgAeAW4K2lCj9MT1Ug8lZZf+29mOA37u0I8Clqnq+6raAswEJnkjqOpzqrrV3Z0DDArWTMMwDCMbfhz6vsBKz36jeywd3wL+lCpARKaIyHwRmd/U1OTfypRoys1ASDWOaNL/opAk3fdFhsnKYkn/fbd70N/oLvIEXr7h5SJxQjysdqYjStJ/Efk6MBq4OVW4qs5Q1dGqOrq+vj7Ioo1A0BRbcSVkSlFfF7oPmwt2GCE781Fz2GWm2kecD4HBnv1B7rFOiMjngWnACaq6IxjzjLLgexzW8z30HMdgRfK7N0/OJ4BM/JdWysJKRBzrhEn/MzAPGCYiQ0WkFjgHmOWNICJHAncDp6nq6uDNNMpB2ft/JVHuSVE72bEgq0NX1VbgUuBp4G3gIVV9S0SuE5HT3Gg3Az2Ah0VkgYjMSpOdYRiGUST8DLmgqk8BTyUcu9qz/fmA7TIMwzByJLpKUZP+x4NCpf9J6Ysj/ffd7qoBv+Ri0v/ORED6X46yXaLr0I3g6ST9L6MdpSBs0v9C44bltdrACZs94cYcupECcf/NPFFWkPTf829hlHAyL44Th1JkOXzZKIP0P2PZhEP6bxiGYUQDc+iGYRgxIcIOPYBJr7RZm/S/VGjB0v+E9IFK/3fjXwZl0v+iYtL/jETYoRvB45X+R+xCyRmT/udeRhmImsMuM+bQjWR8Tt5Imm2/RZj0PwzEsU6Y9N8wEil7/68kTPpvBIA5dMMwjJhgDt0wDCMmRNehm/Q/eMqhFC1U+p/UHj6l/7kuuOzXTs3xzYtsizv7tSHnRaKzlReWycg05xfCJ/3HT9vZWy6GYRiGD8yhGylwhflZ5smkgNdcBEEDmPUv7UsuBXzrIJe8S5lPRUn//YQFUaZJ/w3DMIwCMYduGIYREyLs0E36Hw+Clv4Xkm/6cvznZ9L/omLS/4xE2KEbwaMptuJKyKT/BccNy1tYARMyc8KOOXQjGZ8TZd7vpWf7dnqqIgKZFDXpf4G4dYqbUtSk/4ZhGEaUMYduGIYRE8yhG4ZhxIToOnST/gdPJKX/SYHJcQKQ/vse/YyM9D9LPqGZjcwk/U8VXk4K/aZ94UTXoRtFxOf30AtcJDqQ76EHkIf/wso9qVYE4lgnIPMi0cUu2pSiRgiJ7bUeRsr9PfRYvsFTeZhDNwzDiAnm0A3DMGJChB26Sf+Dp5gTzX7KzD2NJM2ZpZrUy3WS228e6dLmof5Me1792pDjpGi28sIiqU+0I0n6X1JrCsek/4ZhGIYfzKEbybgTa9mmyQr4HDqY9D8h63JNipr0v3TYWy6GYRiGT3w5dBGZICLvisgyEZmaIrxORH7vhr8qIkMCt9QwDMPISFaHLiJVwO3AROAQ4FwROSQh2reA9ap6IPAzYHrQhhqGYRiZEc0y6yoixwDXqup4d/9KAFX9iSfO026cV0SkGlgF1GuGzEePHq3z58/P2eB5j95K/aL/oZtuY09dA8Aq2ZPtUpdzXukY0r4y6ViT9GeLdCs47zptYR/9BIDlXQbnnL6GVvZt/zjv9JnoaNOLu9/GX9cN4Nyj9uP6Mw5PG/+aJxbxm1dWAPD7KWM5+oD+vstavXk7VbccSI20s65Lv7zsBGiq3of6vr2h6R0nsM9+UNMN2ttg7VLnWP8DoUt150w2/BN2bk1dQL8DoKoWtm+EzR/vOnx/7Vf5csvjdKWFL/Sa1SnJXzadBsBHsjdCO/voaiD7OequW6nXtWnj9tZN9NWNAPyzyyDa04zDdsTbTi2ruuyVsUyAKtoZ3P4hv+w6hSdqT6V/+1pmNl/AWunLZulR1H6WC95zvbzL4E7XTwvVVNNGlxSvulze7UZ+vnUq66U3X+35WwBOafkzl2+/g79Xj+XH3X4UuK2TWp7k0u0z+H3tGfxv1wuSwq/b+t+Mav0HH3XZm7WjLmfUKRflVY6IvKaqo1OFVac6mMC+gNfDNQJHp4ujqq0ishHoD6xJMGQKMAVgv/3282V8ksE9+rOu21DWASvpQptUU60teeWVjjUcSBfa6N26lif7X0DD5r/RvX1jYPl339bCqtr9aa7qlVf6rtvb2VTdj/XV9YHZBLAOeLdLL/aoH87EgXVMGrlvxvinjRzImuYWutVWccjA3OoyoHsdzww+nwEbFuZl58fayk6pY7++ddCzDnruDRtWwt6H7Y7YtRe0t0LfIcmZ1B8EUgXaBlvXwZqlsOfB0LIVeu2zK5quXcY7O+rZ0iqs6HsSN8kEDt06l2F9e3TK7snq8xnV/DyravcHoNu2Vj6qHcq2qu5Z69Jj62ss2eMI2iT5clwHrNSd7JCuVNGWMZ/uW+awuPvYjOV5WcVBNPU/jmG1PRDtxnNrzqR36+5LtktLLdu6dGdD9QDfeQbNOqCtZQVravZhp9QCULNdaJNqPqw7AEVorDuQwTuW8nHtEAbs/IhPavdD+hzOS01fZGH3YxjWwzlXa3eOY+7aJcztfQrD9uiRodT8WNJ+Jn9Zt5E5/b7JsC57JIUvaD6Nus1dAajtkdtNjF/83KGfBUxQ1Yvc/W8AR6vqpZ44i9w4je7+e26cNanyhPzv0A3DMCqZTHfofiZFPwS8z1yD3GMp47hDLr2BtbmbahiGYeSLH4c+DxgmIkNFpBY4B5iVEGcWcL67fRbwbKbxc8MwDCN4so6hu2PilwJPA1XAPar6lohcB8xX1VnAr4DfisgynGGvc4pptGEYhpGMn0lRVPUp4KmEY1d7trcDXwnWNMMwDCMXTClqGIYRE8yhG4ZhxARz6IZhGDHBHLphGEZMyCosKlrBIk3AihySDCBBeWrswtomPdY26bG2SU+Y22Z/VU0pEy+bQ88VEZmfTh1V6VjbpMfaJj3WNumJatvYkIthGEZMMIduGIYRE6Lk0GeU24AQY22THmub9FjbpCeSbROZMXTDMAwjM1G6QzcMwzAyYA7dMAwjJpTcoYvIPSKy2l0Uo+NYPxH5i4gsdf/v6x4XEbnNXXx6oYg0eNKc78ZfKiLne46PEpE33TS3iUjqNbtCRpp2uVlE3nHr/piI9PGEXenW8V0RGe85nnJBb/fzx6+6x3/vfgo5MqRqH0/YD0RERWSAu18x/QbSt42I/Kvbf94SkZs8xyum76S5rkaKyBwRWSAi80XkKPd49PuNqpb0DzgeaAAWeY7dBEx1t6cC093tLwJ/AgQYC7zqHu8HvO/+39fd7uuGzXXjipt2YqnrGGC7nAxUu9vTPe1yCPAGUAcMBd7D+bRxlbt9AFDrxjnETfMQcI67fRdwSbnrXGj7uMcH43zaeQUwoNL6TYa+8zngGaDO3d+zEvtOmraZ3XF+3b7yfFz6Tcnv0FX1BZxvpnuZBPzG3f4NcLrn+H3qMAfoIyL7AOOBv6jqOlVdD/wFmOCG9VLVOeq09n2evEJNqnZR1dmq2uruzsFZLQqcdpmpqjtU9QNgGXCU+7dMVd9X1RZgJjDJvWs4EXjETe9t40iQpt8A/Az4IXRaKbhi+g2kbZtLgBtVdYcbZ7V7vKL6Tpq2UaBjEdzewEfuduT7TVjG0PdS1Y7l1VcBHcuWp1qget8sxxtTHI8DF+LcAUDu7dIf2OD5cYhFu4jIJOBDVX0jIcj6DQwHjnOHSv4mImPc49Z34HLgZhFZCdwCXOkej3y/CYtD34X7S2fvUnoQkWlAK3B/uW0JCyLSDfgRcHW2uBVKNc4QwVjgCuChso/vhodLgH9T1cHAv+GsuBYLwuLQP3EfX3D/73g8TLdAdabjg1IcjywiMhk4Ffia+2MHubfLWpzHx+qE41HmUzhjwG+IyHKcOr0uIntj/Qacu8VH3eGDuUA7zgenrO846x8/6m4/jDPcBHHoN2WaqBhC50mKm+k8KXqTu30KnScp5nomKT7AmaDo6273SzNJ8cVyTlIU2C4TgMVAfUK8Q+k8sfU+zqRWtbs9lN0TW4e6aR6m88TWd8td30LbJyFsObsnRSuq36TpO98BrnO3h+MMGUgl9p0UbfM28Fl3+yTgtbj0m3I07oPAx8BOnLuIb+GM0/0VWIozM9/RWALcjjP7/iYw2pPPhTgTOsuACzzHRwOL3DS/xFXDhv0vTbsscy/EBe7fXZ7409w6votnZh1npn6JGzbNc/wAt/Mtcy/QunLXudD2SQhfzm6HXjH9JkPfqQV+59bpdeDESuw7adrmWOA1nB+tV4FRcek3Jv03DMOICWEZQzcMwzAKxBy6YRhGTDCHbhiGERPMoRuGYcQEc+iGYRgxwRy6YRhGTDCHbhiGERP+P39mK6PSg/CYAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from utils.metrics import f1_score, far, mar\n",
    "metrics= {}\n",
    "predictions = []\n",
    "\n",
    "draw = 1\n",
    "\n",
    "all_prediction = {\"swat\": []}\n",
    "all_test = {\"swat\": []}\n",
    "\n",
    "for file_name, prediction in predictions_by_names.items():\n",
    "    y_test = y_test_by_names[file_name]\n",
    "    f1 = f1_score(y_test, prediction)\n",
    "    far_score = far(y_test, prediction)\n",
    "    mar_score = mar(y_test, prediction)\n",
    "    all_prediction[\"swat\"].append(prediction)\n",
    "    all_test[\"swat\"].append(prediction)\n",
    "            \n",
    "    metrics[file_name] = (f1, far_score, mar_score)\n",
    "    \n",
    "    if draw:\n",
    "        title = file_name + \" f1 = {:3.2f} far = {:3.2f} mar = {:3.2f}\".format(f1, far_score, mar_score)\n",
    "        plt.title(title)\n",
    "        plt.plot(y_test)    \n",
    "        plt.plot(prediction)\n",
    "        plt.legend([\"y_test\", \"predicted\"])\n",
    "        plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "id": "12c562fd",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "metric | swat\n",
       "--- | ---\n",
       "F1 |0.4\n",
       "FAR |0.26\n",
       "MAR |0.42"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from utils.data import show_score_table\n",
    "\n",
    "show_score_table(metrics)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}