confident_predictions.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "92d70b84-bf4f-4d83-a88d-0e3c4a7fa959",
   "metadata": {},
   "source": [
    "# Assessment of model performance for high confidence predictions\n",
    "Required inputs: `activity_summary_stats_and_metadata.txt`\n",
    "\n",
    "Output: `predictedObservedClasses_{confident,random}.svg`, violin plots of the observed activity for high confidence predictions, stratified by their predicted activity. Bonus plot is generated for the randomly sampled sequences, but this figure is not used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "428bc5f9-cde3-4896-85ac-cab1ce8f7940",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/ref/bclab/software/dlyon/spack/opt/spack/linux-rocky8-x86_64/gcc-8.5.0/miniconda3-4.10.3-2vnv4wet2ju7hlsej5cgneo6trlw7lly/envs/active-learning/lib/python3.9/site-packages/tqdm/auto.py:22: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html\n",
      "  from .autonotebook import tqdm as notebook_tqdm\n"
     ]
    }
   ],
   "source": [
    "import os\n",
    "import sys\n",
    "\n",
    "import numpy as np\n",
    "from scipy import stats\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.metrics import f1_score, confusion_matrix\n",
    "\n",
    "from src.mpra_tools import loaders, plot_utils, modeling\n",
    "from src.selene_files.enhancer_model import EnhancerModel\n",
    "\n",
    "data_dir = \"Data\"\n",
    "figures_dir = \"Figures\"\n",
    "plot_utils.set_manuscript_params()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6741e47f-e825-4cf9-9181-25a22823dfec",
   "metadata": {},
   "source": [
    "Load in the data and extract the high confidence and randomly sampled sequences. Also load in the CNN trained through Round 3b to obtain predicted activities. This represents the learning round at which the high confidence predictions were made. Both the high confidence and randomly sampled sequences are scored against this CNN to enable a direct comparison between the two groups."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "63fabf17-d61b-403a-9e21-c2b719bd635b",
   "metadata": {},
   "outputs": [],
   "source": [
    "activity_df = loaders.load_data(\n",
    "    os.path.join(data_dir, \"activity_summary_stats_and_metadata.txt\")\n",
    ")\n",
    "activity_df = activity_df[\n",
    "    activity_df[\"data_batch_name\"].isin([\"HighConfidence\", \"Round3c\"])\n",
    "]\n",
    "model = loaders.load_cnn(\n",
    "    os.path.join(\"ModelFitting\", \"CNN_Clf\", \"Round3bFull\", \"Round3b\", \"1\", \"best_model.pth.tar\"),\n",
    "    EnhancerModel(164, 4)\n",
    ")\n",
    "# Define the sort order for the four activity bins\n",
    "class_sort_order = [\"Silencer\", \"Inactive\", \"WeakEnhancer\", \"StrongEnhancer\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c9c8fcd4-563e-451b-b6e3-9cc118ede8ba",
   "metadata": {},
   "source": [
    "## Make predictions\n",
    "For each group, make predictions and return the predicted class, the probability, and the entropy (uncertainty) of the predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "722e1758-f332-40b5-ba3f-fcbc620d2e0c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def predict_group(df):\n",
    "    seqs = modeling.one_hot_encode(df[\"sequence\"])\n",
    "    seqs_batched = modeling.batch(seqs)\n",
    "    preds = modeling.cnn_predict_batches(model, seqs_batched)\n",
    "    preds_df = pd.DataFrame(preds, index=df.index, columns=class_sort_order)\n",
    "    # Go from log prob to probability\n",
    "    preds_df = preds_df.apply(np.exp)\n",
    "    results = pd.concat({\n",
    "        \"pred\": preds_df.idxmax(axis=1),\n",
    "        \"prob\": preds_df.max(axis=1),\n",
    "        \"entropy\": preds_df.apply(stats.entropy, base=2, axis=1)\n",
    "    }, axis=1)\n",
    "    return results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "de0d3d75-ed60-488f-ad19-3d428e0e891a",
   "metadata": {},
   "outputs": [],
   "source": [
    "grouper = activity_df.groupby(\"data_batch_name\", group_keys=True)\n",
    "predictions = grouper.apply(predict_group).rename(\n",
    "    index={\"HighConfidence\": \"Confident\", \"Round3c\": \"Random\"},\n",
    "    level=\"data_batch_name\"\n",
    ")\n",
    "loaders.write_data(predictions, os.path.join(data_dir, \"high_confidence_predictions.txt\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6dcc714a-206e-40bc-8a50-16b92d48bae9",
   "metadata": {},
   "source": [
    "Show differences in the entropy of the two groups."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "04bddd53-557e-4faf-9a6d-990601885b64",
   "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>data_batch_name</th>\n",
       "      <th>Confident</th>\n",
       "      <th>Random</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>2055.000000</td>\n",
       "      <td>6335.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>1.286283</td>\n",
       "      <td>1.488641</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>0.711146</td>\n",
       "      <td>0.570795</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>0.002692</td>\n",
       "      <td>0.019447</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>1.435779</td>\n",
       "      <td>1.120372</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>1.638233</td>\n",
       "      <td>1.775729</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>1.725487</td>\n",
       "      <td>1.932399</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>1.984030</td>\n",
       "      <td>1.999961</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "data_batch_name    Confident       Random\n",
       "count            2055.000000  6335.000000\n",
       "mean                1.286283     1.488641\n",
       "std                 0.711146     0.570795\n",
       "min                 0.002692     0.019447\n",
       "25%                 1.435779     1.120372\n",
       "50%                 1.638233     1.775729\n",
       "75%                 1.725487     1.932399\n",
       "max                 1.984030     1.999961"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, -4.588190444444445, 'Data batch')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAKAAAADRCAYAAABLsblJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAgUUlEQVR4nO2deVQUV9rGn4YWAdlXWVxAjTpuI6KCIEKMglFRQcVdHIMyUSexo0Y0cScYjcyo42SimYAn6DE6KtHRMcEF0IhLUFxA44ZABGQREJSloe/3B1/XQNMFTdNdVdD3d06fpG9X33qrfbh1l+e+JSKEEFAoPKHHdwAU3YYKkMIrVIAUXhHzHYAmkMlkyM3NhampKUQiEd/hUJRACEF5eTkcHR2hp/e/dq9DCDA3NxfdunXjOwyKCuTk5MDZ2Zl53yEEaGpqCqD+4szMzHiOhqKM169fo1u3bsy/lZwOIUD5bdfMzIwKUOAodpHoIITCK1SAFF7RigBzcnLg5+eHUaNGYfv27Ux5dXU1pkyZAi8vLxw4cIC1jKI7aEWAu3fvxqZNm3D16lUkJCSgvLwcAHD06FEEBATg8uXLOHbsGKqqqpSWtUR1dTVev37d6EVRj7dv3+LWrVt4+/YtL+fXyiBk3bp1MDc3BwDU1dVBLK4/TWpqKhYvXgw9PT0MGjQIDx48UFo2dOjQZuuPiorC5s2btRF6h6G4uBjx8fEtHpeVlYWtW7fi888/R48ePVo8ftasWejSpYsGIqxHKwK0srICAOzfvx9ubm4wMjICAJSXl8PExAQA0KVLF1RUVCgta4mIiAhIJBLmvXyIr23evn2Ly5cvY8iQIejatavWz6cuhBBERkaqNDH/8uVL5r8NJ4iVUVdXh+3bt2Pr1q0ai1Vrg5C4uDjEx8c36gOamJjgzZs3AIA3b97AzMxMaVlLdO7cmZly4XLq5eHDhwgICMC2bds4OZ+6xMbGQiwWQ09PDyKRqNlXQ1o6ViwW482bNzh16pTGYtWKAK9fv46DBw/i+PHjMDAwYMrd3NyQnJwMQgjS0tLQt29fpWVCRf6Hwld/SRUqKipw48YNjd4mG2Jubo6ff/4ZdXV1GqlPKwLctm0b8vPzMWHCBPj6+iIpKQnR0dEICQlBQkICRowYgeDgYBgaGiotEyrXr18HAJSVlUEqlfIcjXKOHj3KdIHYKC4uRm5uLnJzc1FUVAQAKCoqYsqKi4ub/b6hoSF++uknjcSrlT7g6dOnm5SNGTMGAHDy5MlG5YaGhk3KhMrTp08BAEZGRjh//jwmTJjAc0RNyczMbHTXUaS4uBh79+5tUn7ixIlG71esWAFra2uldZiamiItLQ3vv/9+24JFB1mK4wJCCNMydOnSBWlpaYIUYEVFBSwsLFg/r66uBgAEBQXBxsYGtbW1KC0thYWFBcRiMYqKinDixAnmODbkU2tthQpQRZ4+fYqGuxfkt672io2NDRwdHQEA3bt35y0OuhSnIr/88kuj0XZ5eXmLrQQfdO7cuV2dhwpQRfLy8pgJdQAwNjbGjRs3eIxIOS0NQDSBTCaDjY2NRuqiAlQR+YStHAsLC1y7do2naNjx8/NrcRTbVl6+fImpU6dqpC4qQBVISUmBTCZrVCYSiZCVlaWx+TBNMXz4cOjr62v1HJaWlo1czW2BCrAFCCE4duyY0ltbly5d8K9//YuHqJpn9OjRGhulKlJcXIxJkyZprD4qwBaIjo5u1PdriLGxMe7du4eUlBSOo2qeadOmaa1lNjY2hpeXl8bqowJkgRCCr776Ci9evGh2WcvGxgaHDh3CxYsXOYyueUQiEYKDg/Hq1SuN1ltQUIDQ0FCN1kkFqISKigp8+umnKCgoYGxlzWFra4uTJ0/iu+++4yA61fD29oaJiQk0lXmltrYWzs7O6Nevn0bqk6NVAQYHB+P58+fM++joaPj6+sLX1xfm5ua4desWXr9+DScnJ/j6+mLs2LHaDEclHjx4gDVr1sDQ0LBVC/o2NjbIyMjAZ599Jpj5wWXLljUZvdfW1gJovPbb8CWfYJcfJ6e4uBjLly/XeIxaWQmpqalBSEgIbt682ahcIpFAIpEgPT0dkZGRcHNzw7Vr1xAeHo7PP/9cG6G0inPnzuH06dOwt7dX6/umpqaQSqWQSCRYv349s9LAF46Ojujfvz/y8/OZiePS0lIATdd+FSktLWVWSCoqKuDp6dlkS6Um0IoAq6ursXLlStZbUmRkJONoTk9Px7lz55CQkIBFixZh0aJFKtXfsJXRhCU/Li4ON27cgJ2dXZvq6dSpE2xtbbF161YsX74cAwYMaHNsbSE8PBwrV65krku+TixfC1ZEvhbccD25pqYG8+bN00p8WrkFm5qawsfHR+lnFRUVKC4uRp8+fQAAPXv2xI4dO3D+/HkcPHgQhYWFLdYfFRUFc3Nz5tVWN3RcXBx+/fVXVvdHaxGJROjatSv27duH3377TSN1qouBgQFGjx7NOM3lI3r5WrDiSy5K+XElJSWYMmWK1lKecD4IOXv2LIKCgpj3Hh4e8PT0hIGBATw8PJCZmdliHRERESgrK2NeOTk5asdz5coVXL9+HZaWlko/b4t3zt7eHrt37+bdwBoSEoKamhq1vmtgYIB3331XwxH9D87dMElJSVi2bBnzfsOGDRg/fjzGjRuH1NRUrF69usU6OnfurJHFcEIIjh49yrquqQnvnJWVFb7++mt88sknbY5XXUQiEby9vVu9dl1SUtKosdAGnAjw3r17SEhIgEQiQWZmZqNlHIlEgnnz5mHLli1YsGCBxm6DqnDx4sVml6004Z3r1KkTnjx5AplM1uKmH20SHByMpKSkVn1HLBazdqU0hVYFGBsby/z/oEGDANTfghvi5OSES5cuaTMMVu7fv9+seVOOJrxz2dnZ6Nmzp1rf1QR6enro168f0tLSVDq+qqoKvr6+Wo0J0PGJ6E6dOnFiJpDJZMzWUz6ZO3cuSkpKVDr2zZs3mD59upYj0nEBzpgxAwUFBVo9ByEE5ubmGvPPtQVLS0uVWnwAsLOzQ6dOnbQbEHTckm9ra4uJEyfi/PnzSt0uDVcNlMG2atCQly9fYtOmTW0PVkO4uLiodJy8y6RtdFqAADB58mRUVFQonYpRZ9WgIS9fvoREIlF7ZUUbjB49WqXj3N3dtRxJPTovQACYPXs2SktLkZmZ2Wj9V51VAzmFhYUIDw9nJtyFgqrLado2tcqhAvx/wsPDsWrVqkZliqsGbCj6BaurqzFgwAAMHjxY84FqGa6fW6TTg5CGiEQijBw5UiOrFq9evdK4b07TsM1dcp3qjgqwAYGBgSgrK2tzPfb29kxGMKHCZtlXd8lOXagAG2BoaIi+ffu2yc9XVFSEwMBADUalHdj6eMbGxpzGQQWowPLly1FeXt5kF5wqlJeXY8CAARg2bJgWItMsbGZbbWXVYoMKUIFOnTph27ZtKCoqalWHvLy8HHZ2dggPD9didJqDrYvAddeBU0u+Mvu9EJOUm5ubY9u2bSpv6qmqqoKVlZVKTh6hoExohBDhCLAtndGamhpMmzaNyacnJyMjA+Hh4UhMTMSFCxcAKE9c3hJcJCm3tLREWFiYSseKRCJERES0q+fUde7cuUkLX1VVxfmSIasAR44ciRUrVuDWrVutrlRuyX/vvfcalcvt9z4+PoiJiQFQn7jc29u7UZLyltC0I5oNW1tbAGhxZBwWFtauxAfU+xQV/9grKio438fCKsDbt28jKCgIe/fuxZgxY/C3v/1N5ZwjbJZ8ZfZ7dZOUa8oRrQoWFhZKByXyllfoUy7KsLOzQ2VlZaOy2tpaTv2YQAt9QE9PT4wbNw4mJib4z3/+g7lz52LXrl1qn0yZ/b49JCkPDg5W6prhwi2iLSwtLZt0s+rq6jh/1h6rAENDQ+Hu7o5Hjx5h3759OH/+PM6dO4dDhw6pfbINGzYgISEBMpkMqamp6NWrV7tIUm5tbd1kt9zbt2/xhz/8gaeI2o6JiYnSVr259L7agHUteM6cOY0czXIuX77c6pM0tOQr2u9DQkIwe/ZsxMTEYNGiRYJKUi5flrt16xaMjIxw//59VFZWwsLCAqWlpZg8eTLPEaqPsbFxEwGy5cDRJqxnNDQ0hIeHB4qKimBtbY19+/bB3d29VROVyiz5ivZ7IScpf/jwIQCwjoblNi1tbNjWNsq6D4ISoEQiwdGjR+Hq6oonT55g/vz5gssCpW3kSRj79esHY2NjbNiwAWfOnMG0adPwzjvvYObMmTA1NRWc5UpVFDdJ8bFpilWAFhYWcHV1BQD07t2b8zVCIWBjY4MPPviAed+3b1+cOXMG+vr6CAsLQ69evXiMru0oCo4rD2BDmhXg1KlT4ePjg19//RVVVVXYsmULgPrBhC4i9/fp6em1e/EBTQXIx1wmqwDlHWyRSAR/f3/OAhIy/fv3BwCVUra1B5Q9K45rmh0F79+/H/fv30efPn2wbNkyzh4BIFTknXQhbLHUBEIQIGuvc+nSpcjNzcWUKVNQWFioUtYqXaE9jnqVIYTlQ9YW8NmzZ0x6tYCAAHh7e3MWlNBxcnLiOwSNoChAPkbBrGesrq5m7EglJSVqGTQ7KkLaZqlJBNUH3Lp1Kzw8PGBjY4Pi4mL8/e9/5zIuQdNR+oBCgFWAmZmZePToEQoLC2FjYyOI/oJQ4Hq9VFsIYRqG9RYcFxcHQghsbW2p+BTgM82aJhH0KPjVq1fo3r07fH194efnp1aWTEVLfk5ODvz8/DBq1Chs374dgHKbvtBpzzas5hBUHzAmJoZxBAP1+e1UhS1L/u7du7Fp0yaMGTMGY8eOxbJlyxibfmuy5GsjSbkq9OvXD6mpqRp/VgZfCPIW/OjRIyQlJWHJkiXIyspCVlYWMjMz8dFHH6lcKZslf926dcx0Tl1dHcRisVKbfktwZclXxNjYGG5ubh1mXVyQAiwsLERsbCzy8vIQExODmJgYfP/991i6dKnKlbJZ8q2srKCvr4/9+/fDzc0NRkZGamXJ59qS31ERQh+wyS3Yy8sLXl5euHPnDoYMGaLxE8bFxSE+Ph7x8fEA6m36RkZG0NPTY2z6DW/9ytBUknJdR9AT0T/++CNsbW3h6uoKFxcXxprVFq5fv46DBw/i+PHjzFSGMps+hRuEIEDWQcjJkyfx+++/a6SlkVvyL126hPz8fEyYMAEAcOjQIV6z5Os6QpheYxVg79692xygoiVfIpEoPY6vLPm6jqAd0bm5uXBxccE777wDoP6vRUjPxKW0HUEL8PDhw1zGQdFRWAUoFouxZs0aFBYWIigoCEOGDEGPHj24jI2iZYQwDcPa5i5ZsgTLly+HVCrFuHHjsHLlSi7jomiJhnud8/PzkZ2djbt37yI7Oxu5ubkq5ebRJKwtYFVVFTw9PSESidCrVy9BbRinqE9Le53ltjuuXN+sAjQ3N0dsbCwqKyvx448/qvyEHYqwabjXOS4uDoWFhThx4gSCgoLg7OyMhQsXcrrXmfUW/N133+HBgwewsrJCcnIyvv32W04ComgX+V5nb2/vRg+otrGxgZOTE9zc3DjdaN/svuAvv/ySs0AouknHcFZS1ELQo2BKx0fQa8FZWVn44YcfGqVxbW1KjuDgYOzatYt5UHN1dTVmzpyJoqIihIaGIiwsTGkZhRuEsBbMKvkZM2ZAKpWiR48ezEtV2JKUK0tILtQk5bqAYpJyrp8TB7QwCFm/fr1alcod0fKN7XJSU1OxePHiRgnJlZUNHTq02fqjoqKwefNmtWKj/A9B9wFdXFywZ88eJCUlITk5GcnJySpXyuaIVpaQvD0kKe+oCNqMUFNTg9u3byMtLQ2EEIhEIqWiag3KEpKrm6ScOqI7BqySj4mJwSeffAJ/f3+sWLGiye1UHZQlJG8PScop2oNVgHv27MHy5ctx+/ZtfPTRR4iKilL7JPfu3UN0dDRCQkKQkJCAESNGIDg4GIaGhkrLKNwg6FvwkSNHcOXKFejp6aGurg6enp6IiIhoVeXKkpQrJiQXcpLyjg4fo15FWCVfV1fH/EXo6+vzkj+Yol2EsC+YtQWcMWMG/Pz84OnpiZSUFAQFBXEZF4UHBHULXrVqFQICAvDw4UPMmTMHAwcO5DIuCgcIYSWkiQB37NiBNWvWYNGiRUyAZ8+eBQCNjIQpwkGQApRnxA8NDW1U3jAZEKVjIAQzQpMzOjs7Izs7G5999hlcXFzg4uKC7t276+yzQToygmwBf/nlF/z1r3/F/fv3ERoaCkII9PX1m2S6orR/hNACNhFgYGAgAgMD8d///hcDBw5Et27dcOPGDYwYMYLz4CjaRQgtIKvkT548iSNHjgAAjh07hj//+c+cBUXhB0G5Ye7evYvVq1cDAHbu3ImMjAzOgqJwgxCW4ljPKBKJmPzOOTk5gli2oWgWxRZPUIbUXbt2Ydq0aSgtLYWNjU2rnhPCZrOPjo7GqVOnAAC3b9/GpUuX0Lt3b/Tv3x99+vSBvr4+Lly40MZLoqiKoB3Ro0aNwu3bt9WqVG6zX7p0KQICAjB//nwYGhpCIpFAIpEgPT0dkZGRcHNzw7Vr19pNkvKOjiBuwXPmzAEAJiuqOhlSU1NT4e3t3chm35DIyEjGUt+ekpR3NITQB2zSAsrTsmVmZqpdaXM2+4qKChQXFzO77+VJyocPH47x48dj0qRJLeaIjoiIaJTs8vXr11SE7ZQmAvTz81M6HBeJRCr3z5qz2Z89e7aRs4YmKRcOgpiGiY2NRUxMDBwcHPDxxx/j2LFj+PTTT1t1C27OZp+UlITRo0cz72mSct2miQDle4Bzc3MRGBgIa2tr+Pv74/HjxypXqmizf/z4MaKjowHU39qdnZ2ZYyUSCbZv347Ro0dj+vTpNEk5hwhhao11FGxiYoKoqCi4u7vj6tWrcHBwULlSZTZ7uSVfbu2S4+TkRJOU84QQ5gFZhz1HjhyBmZkZ4uPjYW5uju+//57LuCg8wMdDyVkFaGhoCAMDA1hbW2Pw4MEoKSnhMi4KBygKTlACXLJkCYqKipCQkICamhrMnz+fy7goPCCIUbCc58+fIyIiAoaGhggICFApZQalfaHY5xOUAGUyGZ49ewaRSISCgoIO+5BmXUbQa8G7d+/G3Llz8fDhQ0ydOhV79uzhMi4KD/DRB2QVYHx8PFJSUriMhcIxfAhOEdZb8OXLl5nlNErHRPGWK6gW8PHjx7CysoKDgwNEIhFEIhGePXvGZWwULSPoPqDcDU3RHQSxEnLs2DG4urpi8ODBSE1NVavS6upqTJkyBV5eXjhw4ABT/vr1azg5OcHX1xdjx45t9liK9lG85QpCgNHR0UhNTUV8fLzaeZjZEo9nZGQgPDwciYmJjLVLnSTlFM0ghD5gEwEaGRnB0tISrq6uag9C2BzRytzPLbmnlUGz5HccmvVgq9skszmi5e7n8+fP4+DBgygsLFQrSTm15GsHQQxCHjx4gD/96U8ghDD/L0fV7Fhsjmhl7md1kpRTS75mEEIfsIkA5dkQgKYZslRF7ogeMGAA0tLSEBkZCaDe/Tx+/HiMGzcOqampWL16NeuxzUEt+dpBEBlSx4wZ0+ZKQ0JCMHv2bMTExGDRokV4/PgxEhISIJFIMG/ePGzZsgULFiyAtbV1k2NpknLuEMIghHUesC0054hWdD/TJOXCQVBuGIruIYhpGAqFS6gAKbxCBUhhoH1ACqcIOkMqpeMjhI3pVIA6jBCSlFMB6jBCSM9GBajD0D4ghVeE8LRMKkAdRvERvHw8klcrAmSz2efk5MDPzw+jRo3C9u3bASi36VO4QbHF6zAtIJvNfvfu3di0aROuXr2KhIQElJeXK7XptwR1RGsGxVuwWKwVb0rzMWijUjab/bp16+Dt7Q2g/onsYrGYJinnkYaCq62thbGxMecxaEWAbDZ7Kysr6OvrY//+/XBzc4ORkZFSm35LREREoKysjHnl5ORo4zI6PA1bwNraWubfjEu00uY2Z7OPi4tDfHw84uPjAdAk5XxiaGjIWLBqa2tV2g6habTSArIlKb9+/ToOHjyI48ePw8DAAABNUs4nZmZmqK2tBVAvQCsrK85j0IoA2ZKUb9u2Dfn5+ZgwYQJ8fX3x4sULmqScRywtLVFTUwOgvk/ORwvIqSW/4U62htAk5fxgbW0NqVQKoN6Y0GFuwZT2gYWFBdMHlMlkVIAUbmkoQD09PV6y4FIB6jCmpqaMJ5CPSWiAClCnEYlETKvHVw5wKkAdR97y0RaQwgvylk8+L8s1VIA6Dr0FU3hFLjy+bsH8nLUdU1paitLSUohEIhBCYGFhAQsLC77DUhv5rZevWzAVYAs0FFxeXh6ysrJACEFlZSW6du2K7t27t2sBylu+DiXA6upqzJw5E0VFRQgNDUVYWBhrOduxfHP9+nVkZGSgsLAQVVVVSvfQZmdnIyMjA/fv38eAAQMwcOBAHiJtG3LhmZqa8nJ+Th3RysqFmKQ8KysLiYmJyM7ORmVlZbMbuCsqKvDw4UNcvHgRWVlZeP78OUpLS7kLto107doVANC9e3dezs+pI1pZuRCTlOfl5bU6a4BUKkVeXh4AtCsBhoWFITU1FVOnTuXl/Jw6opWVCzFJuYODA5ycnFplUXd1dYWDgwMAtKs+obGxMdzc3Hix4wMcO6KVlQsxSXmPHj3oA7o5glNHtLJytmObo3PnzjAzM2v0orRPOHVEK5YbGhoqLaPoDiIihBxdbeT169cwNzdHWVkZbQ0FCtu/EV2Ko/BKh1gJkTfiNEOCcJH/2yjecDuEAMvLywGAZkhoB5SXl8Pc3Jx53yH6gDKZDLm5uTA1NdVqgh35dE9OTk6H6WtydU2EEJSXl8PR0bFRRoYO0QLq6enB2dmZs/N1xKkfLq6pYcsnhw5CKLxCBUjhFSrAVtC5c2ds3LixQyVG4vuaOsQghNJ+oS0ghVeoACm8QgVI4RUqQAqvUAFSeIUKkMIrOifAmJgY+Pj4wNPTE+vWrVPpO69evYK7uzuioqKYB+zI2bRpExITE1Wu5/z5860NWSViY2PRr18/+Pr6wtvbG/7+/qiurm5VHa25Fk2hUwK8desWfv75ZyQmJiIlJQX5+fk4c+ZMi99LT09HQEAAIiIisHbtWrXPf/fuXVy5ckXt77fE2rVrkZiYiCtXrsDT0xNJSUlaO5em6BBmBFU5ceIEPvjgA8aNsXfvXjx58gReXl6oq6vDxx9/DA8PDyxZsgR6enooKyvDqVOnsGXLFjx+/Biurq5ITk7GN998g+nTp+Pt27eoqamBr68v0tPTsWzZMtTU1GDVqlUICgqCp6cn7O3tkZWVhR9++AH79+/H1atXERwcjEGDBmn1Wl+8eIFnz57B398flZWV+OMf/4g9e/bA398fVlZWSE9Pxz/+8Q8MHz68ybXcuXMHH374YaPfZPHixRCLxTAzM0OvXr1w5swZ7Ny5EwEBAW0LlOgQYWFhJCMjo1HZlClTSFZWFpFKpcTDw4M8efKEjBw5khBCyBdffEEOHTpELl26RDZu3EgyMzPJwoULyb///W+yc+dOQgghM2fOJJcuXSLTpk0jOTk5RCqVknfffZfIZDLi7OxMpFIpOXz4MImMjGTq0QYxMTGkb9++xMPDgzg6OpJt27aR06dPk6dPnxKZTEaGDx9OCCGkb9++pKSkhFy9epWEhYUpvZbmfpNJkyaRixcvkjt37pAlS5a0OW6dagEdHByQk5OD/v37A6jPgFBWVsZkBejZsycKCgqYnXkODg7MYwwa8uzZMwwZMgRA/U4/AHj69CnmzZsHACgsLERRURFcXV0hFovh4OCA3377TevXt3btWoSGhuIvf/kLTExM4OjoiPXr18PY2Jh5ApW9vT0sLCyYa1N2Lc39JtbW1nBxcQGAVvcxlaFTfcCpU6fi22+/hUwmAyEEGzduRFJSErKzs1FbW4snT57A2tq6RVNrnz59cOvWLQDAvXv3mLITJ07g4sWLmDVrFqysrJQ+jZJwsPT+5Zdf4ptvvkFERASioqKwZ88eSKVSEEKaxKTsWkxNTVv9m6iLTglw6NChGDNmDLy8vODh4QFnZ2ekpqZi9uzZ8PT0xKpVq1TKEhUYGIjbt2/Dz8+PScexYcMGBAcHw93dHWKxWOmzd11dXXHixAncvHlT49fWECMjI0RERCA5ORkTJ07EhAkT0LVrV7x8+VKla9m8eXOrfxN1oW4YCq/oVAtIER5UgBReoQKk8AoVIIVXqAApvEIFSOEVKsA2kJiYiK5du8LX1xc+Pj7w9vbGTz/91Ox3Dhw4oHLdoaGhKh2bk5PT7Hl79uypUj18QAXYRgICApCYmIjk5GTEx8dDIpEgPz+f9fjIyEiNx3DhwgWkpKRovF4uoALUIDY2NggJCcGZM2dQWlqK6dOnY/z48XBzc8PJkycRGxuL/Px8zJ07V+nnimRkZGD8+PEYPnw44uPjAQBHjhyBn58ffHx8MHHiRFRXV2P79u2IjY1FQkIC4uPjMXz4cAwbNgyff/45gPrUx7NmzcKoUaOwYsUKLn+SlmmznUGHuXTpElm4cGGjsn/+858kKiqK3Lx5k5w6dYoQQkhqaip57733CCGE9OjRgxBCWD9vWLebmxuRSqWkpKSE9OnTh9TV1ZHIyEhSXV1NCCEkMDCQXLlyhcTExJCNGzcSqVRKXF1dyatXrwghhKxevZqUl5cTsVhMcnNzCSGE9O/fnxQUFGjl91AHnXLDcEF2djb69+8PR0dH7N27F8ePH4dIJIJUKm10XEufA8DIkSMhFothYWEBKysrFBcXw97eHgsWLECXLl3w/PnzRt8rKiqCnZ0dLC0tAQA7duwAUO/qkWfwt7Ozw5s3b2Bra6utn6BV0FuwBiksLMTx48fx/vvv46uvvoKnpydiY2MxYcIExgWjp6cHmUzG+nlD7t69C5lMhuLiYpSWlsLAwAAbNmzA4cOHsW/fPojFYhBCmDrt7OxQVFTEJIOcPXs2fv/990bp0IQGbQHbyLlz5+Dr6wt9fX3IZDJ8/fXXsLKywsSJE/Hhhx8iLi4Ozs7OePXqFQDA29sbkydPhkQiUfp5Q8RiMfz9/VFaWoo9e/bAzMwMHh4eGDZsGIyMjGBubo68vDwMGjQIX3zxBdzd3REdHQ1/f38A9QMkLtPWqQN1w1B4RbhtM0UnoAKk8AoVIIVXqAApvEIFSOEVKkAKr/wfc3eqdaaSOIcAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 171.26x205.511 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "grouper = predictions.groupby(\"data_batch_name\")[\"entropy\"]\n",
    "display(grouper.describe().T)\n",
    "fig = plot_utils.violin_plot_groupby(\n",
    "    grouper,\n",
    "    \"Prediction entropy\",\n",
    "    figax=plt.subplots(figsize=plot_utils.get_figsize(0.25, 1.2)),\n",
    ")\n",
    "ax = fig.axes[0]\n",
    "ax.set_xlabel(\"Data batch\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1730da42-61a5-4c09-a0b1-25a30f7ef3ac",
   "metadata": {},
   "source": [
    "## Now report the accuracy of those predictions"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "10dac93b-3c52-4c7c-835b-2caf3ae8b694",
   "metadata": {},
   "source": [
    "Helper function that takes the predicted labels and observed data, makes a violin plot and reports various measures of accuracy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "26ec703c-2685-4602-9734-dbdb54fe368f",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_report_accuracy(df, preds, title):\n",
    "    print(title)\n",
    "    preds = pd.Categorical(preds, categories=class_sort_order, ordered=True)\n",
    "    grouper = df.groupby(preds)\n",
    "    \n",
    "    # fig, ax = plt.subplots(figsize=plot_utils.get_figsize(0.33, 1))\n",
    "    # fig, ax = plt.subplots(figsize=plot_utils.get_figsize(0.33, 0.8))\n",
    "    fig, ax = plt.subplots(figsize=plot_utils.get_figsize(0.27, 1))\n",
    "    fig = plot_utils.violin_plot_groupby(\n",
    "        grouper[\"expression_log2\"],\n",
    "        \"Observed activity\",\n",
    "        class_names=[\"Silencer\", \"Inactive\", \"Weak\\nenh.\", \"Strong\\nenh.\"],\n",
    "        class_colors=plot_utils.get_class_colors(),\n",
    "        figax=(fig, ax),\n",
    "    )\n",
    "    \n",
    "    # Calculate main accuracy metric\n",
    "    observed_bin = df[\"activity_bin\"]\n",
    "    score = f1_score(\n",
    "        observed_bin,\n",
    "        preds,\n",
    "        average=\"weighted\",\n",
    "    )\n",
    "    \n",
    "    ax.set_xlabel(\"Predicted activity\")\n",
    "    for y in [-1, 1, loaders.get_strong_cutoff()]:\n",
    "        ax.axhline(y, color=\"k\", lw=1, linestyle=\"--\")\n",
    "\n",
    "    # Check a few other metrics of accuracy\n",
    "    print(\"Micro F1 score:\")\n",
    "    print(f1_score(\n",
    "        observed_bin,\n",
    "        preds,\n",
    "        average=\"micro\"\n",
    "    ))\n",
    "    print(\"Weighted F1 score:\")\n",
    "    print(score)\n",
    "    print(\"Macro F1 score:\")\n",
    "    print(f1_score(\n",
    "        observed_bin,\n",
    "        preds,\n",
    "        average=\"macro\"\n",
    "    ))\n",
    "    print(\"Per-class F1 score:\")\n",
    "    print(f1_score(\n",
    "        observed_bin,\n",
    "        preds,\n",
    "        average=None,\n",
    "        labels=class_sort_order,\n",
    "    ))\n",
    "    \n",
    "    confusion = confusion_matrix(observed_bin, preds, normalize=\"all\")\n",
    "    print(\"Confusion matrix:\")\n",
    "    print(confusion)\n",
    "    print(f\"{np.diagonal(confusion).sum() * 100}% are correct (should be same as micro F1).\")\n",
    "    off_by_one = np.diagonal(confusion, offset=1).sum() + np.diagonal(confusion, offset=-1).sum()\n",
    "    print(f\"{off_by_one * 100}% are off by one class.\")\n",
    "\n",
    "    ax.text(0.02, 0.98, f\"wF1 ={score:.2f}\", ha=\"left\", va=\"top\", transform=ax.transAxes)\n",
    "    \n",
    "    # Make sure all figures have the same axis\n",
    "    ax.set_ylim((-7.9, 7.3))\n",
    "\n",
    "    if title == \"Confident\":\n",
    "        ax.set_title(\"High Confidence Predictions\")\n",
    "    else:\n",
    "        ax.set_title(\"Random Predictions\")\n",
    "    \n",
    "    return fig, ax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "c9d592ec-9583-49f0-8ba0-453b2b1977de",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Confident\n",
      "Micro F1 score:\n",
      "0.7187347931873479\n",
      "Weighted F1 score:\n",
      "0.7177627465325237\n",
      "Macro F1 score:\n",
      "0.5392817661658966\n",
      "Per-class F1 score:\n",
      "[0.92584655 0.2735426  0.29370629 0.66403162]\n",
      "Confusion matrix:\n",
      "[[2.96836983e-02 7.78588808e-03 2.77372263e-02 7.00729927e-02]\n",
      " [4.08759124e-02 5.25547445e-01 6.81265207e-03 2.57907543e-02]\n",
      " [4.86618005e-04 9.73236010e-04 1.22627737e-01 5.83941606e-03]\n",
      " [1.07055961e-02 1.94647202e-03 8.22384428e-02 4.08759124e-02]]\n",
      "71.8734793187348% are correct (should be same as micro F1).\n",
      "14.452554744525548% are off by one class.\n",
      "Random\n",
      "Micro F1 score:\n",
      "0.3414364640883978\n",
      "Weighted F1 score:\n",
      "0.28767693712262504\n",
      "Macro F1 score:\n",
      "0.2863857464827759\n",
      "Per-class F1 score:\n",
      "[0.24216112 0.29789903 0.11062591 0.49485693]\n",
      "Confusion matrix:\n",
      "[[0.07498027 0.04009471 0.16464088 0.04388319]\n",
      " [0.074191   0.03962115 0.08176796 0.03172849]\n",
      " [0.00520916 0.00331492 0.20883978 0.00789266]\n",
      " [0.02541436 0.01689029 0.16353591 0.01799526]]\n",
      "34.143646408839786% are correct (should be same as micro F1).\n",
      "37.07971586424625% are off by one class.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAALQAAADLCAYAAADHubkmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA8kklEQVR4nO2dd1xTV//HPwmQMEWGIEOWoqCiIIIsJUgRte5RZxWr1da21Gp9XK1gXTz9Vau1VmtrxQdtrfZxPioWK2gdWMWC4EC0oMgG2YHM8/sjJhKSAJlovO/XKy/IPfeM3PvNN+ee8x00QggBBYWeQO/sAVBQaBJKoCn0CkqgKfQKSqAp9Aq9E+i0tDTExMRIHWOxWCgoKEBiYiKSk5Pl1isoKACLxWqzbTabjdjYWAwbNgyhoaGYMWMG6uvrlR7jihUr4Ovriw0bNiAxMVGm3M3NTek2VcHNzQ0sFgsREREICQnB/PnzwefzlW6nI9e3sLAQ586dAwBMnz5drXG3haHWWn4JaS3oyhIbG4vg4GB88803AIAtW7bg448/xk8//aRUO0eOHEFubi6MjIzUGo8mSEtLk/w/bdo0nD9/HiNHjlSprbau7x9//IGCggJER0fj0KFDKrXfEfROQ7dFfHw8EhMTUVZWhuHDhyMiIgJz586V3IiSkhKMHz8egYGBSEhIkKrL4XCQkpKC+fPnS459/PHH2LZtGwDg4MGDCAwMRGhoKJYvXy7p7+2338bIkSMxYMAA3Lx5EwkJCSguLsbIkSORmJiI+Ph4CAQCzJ07FyEhIViwYIGk/aysLIkGnTdvHng8HhITEzF16lSMHj0affv2xZkzZwAA33//PQYPHgw/Pz/s2rULAPD1118jJCQEoaGhOHz4cJvXhsfjoaGhARYWFoiJicGbb76J8PBwcDgcxMTEIDw8HBEREcjOzgYAfPPNNxg0aBDGjRuH0tJSqevb2NiIKVOmIDg4GMHBwcjNzUVCQgISExORkpIi+QX6+++/MXToUISHh2PixImorq5GWloaRo8ejXHjxqFfv3748ccfAQCfffYZQkNDERAQgIMHDyr+IETPSE1NJfb29iQ8PFzysrS0JPn5+SQuLo7s27ePLF++nCQmJhJCCNm7dy+ZO3cuyc/PJ/b29qS+vp5wOBzi4OAg1W5xcTEJCAiQ22dlZSXp27cvYbPZhBBCpk+fTk6dOkXi4uLIxx9/TAghJCkpiSxatIgQQoirqyshhJB9+/aRuLg4curUKbJgwQJCCCGZmZnEzc2NEELIkCFDSF5eHiGEkI0bN5Lvv/+e7Nu3j4wfP54QQsiff/5JoqOjSWlpKfH29ibNzc2Ew+GQJUuWkJycHBIREUEEAgHhcrkkICCAVFdXS43b1dWVhIeHExaLRSIiIshXX31FCCFk7ty55LvvviOEELJz507y+eefE0IIefLkCRk6dCipqakh/fv3J1wulzQ0NBAnJyep6/t///d/JC4ujhBCyOXLl0lSUpLks7b8/IGBgeTevXuEEEJ2795Nli1bRlJTU4mfnx8RCASksLCQ9OnThxBCiLu7OyksLCQ1NTUkKSlJ4f3XyymHWPuJaT03vnPnjkTThoSE4NKlSwCA3r17w9zcHADAYDCk6tjY2KCsrEzqGIfDwalTp+Di4gIfHx+YmJgAAIYNG4acnBwAwMCBAwEAzs7OaG5uljveO3fuYPDgwZLzzczMAAD379+XaOympiZERUWhV69eMm3m5+dj4MCBYDKZAESa+fDhw3j06BGGDx8OAGhsbERBQQF8fX2l+m455WhJnz59AAA5OTm4cuWK5BpVVVXh/v376N+/P4yMjGBkZCTT5sOHDzFx4kQAQGhoKEJDQ+U+K1RUVMDLy0tyzY4dO4YxY8bAx8cHdDpd6polJiZi6dKlKC8vb3MO/lpNOcT06dMH169fBwCkp6dLjtPpii8Hg8FAZGSk1I3Zvn07Tp8+DQ8PD+Tk5Egu/qVLl9CrVy8AAI1Ga3c8np6euHbtGgAgNzcXjY2NAAAvLy8cOnQIaWlpiI+Pl3wxW7fp4eGBu3fvgsfjQSAQYOTIkXBxcYGfnx/S0tKQmpqKGTNmwN3dvd2xtL4Wffr0QUxMDNLS0nD06FHMnDkTbm5uuH37NrhcLjgcDu7cuSNV18vLCzdu3AAAXLt2DcuXLwedTodQKJQ6z8bGBg8ePAAAXLx4UeE143A4OHLkCH799Vf8/vvvWLduHQQCgdxx66WGbo+VK1fi7bffxv79+2FpaQlra+sO1du+fTtiY2Oxd+9e8Hg89OrVC7t374a5uTlWrFiB8PBwEEIQGBiISZMmSbR0e0yaNAmXLl1CcHAw+vTpI/mV2LFjB6ZPnw4ejwczMzMkJSXh6dOnMvXt7Ozw0UcfYejQoSCEYMGCBQgKCsKQIUMQFhYGNpuNsWPHwtLSsuMX6TmLFi3C/PnzwWKxUFdXh7Vr18Le3h6rVq1CUFAQHBwcYGNjI1Vn4cKFmDdvHlgsFggh+PHHH9HQ0IBNmzZJfokAYPfu3Zg/fz4IIejSpQsSExNx9+5dmTEwmUxYW1sjODgYDAYDsbGxMDAwkDteGiGvny3H6dOn4eTkBF9fXyQmJqKgoADx8fGdPSwKDfBaamhXV1e88847MDY2hpmZGfbv39/ZQ6LQEK+lhqbQX17Lh0IK/YUSaAq9ghJoCr1CLx4KhUIhiouLYWFh0aF1Xwr9gBCC+vp6ODo6StbN9UKgi4uL0aNHj84eBkUnUVhYCGdnZwB6ItAWFhYARB+sS5cunTwaCl1RV1eHHj16SO4/oCcCLZ5mdOnShRLo15CW00zqoZBCr3jtBDotLQ3du3cHi8WSvE6ePCkp37Ztm1Lb4Lm5uRKLsnXr1smUz5kzR9JPjx49sGTJEhBC8MEHH0jse1v2T6EmCg1LXyFqa2sJAFJbW9vuuampqWTu3Lkyx7lcLomJiSG9evWS2O12hDFjxpDr168ToVBIRo0aRW7fvi33vLKyMuLv70+qq6tJcnIymTVrFiFEZEvds2fPDvdH8QJ5910vNfTgwYNRW1uL6upqGBoaoqKiArW1tQgNDVVYh8vlYurUqVizZo3c8l9//VVKq7NYLBQVFeHOnTsIDAwEjUZDdHQ0Lly4ILf+559/jtWrV6Nr164YNmwYdu7cCUA0/zM01ItHmZcCvbySb775Js6fPw+hUIi+ffsiNTVVchwAkpOTJbbF3bp1w5EjR2BmZobRo0fLNUQHRP5206ZNa7NfCwsLFBcXyxwvLS3FzZs3sXv3bgCAiYkJTExMwGazMX36dMTFxan4SSlao5cCPX78eOzatQs0Gg3r16/H2bNn0dTUhJUrV6KsrEzGo6Uj/PrrrxJfPTEHDx6UesKur6+Xa3N8+PBhzJo1S+rc8vJyjB8/HjExMZgxY4ZyH5BCIXop0IMGDcLdu3dhamqKcePGYfPmzTAyMoK3t7eMG1VHUaShvby8cPPmTfj7++PcuXNYv369zDkXLlzApk2bJO8bGxsRFRWF9evXY9y4cSqNh0I+ejmHBgA/Pz94eHiARqPBw8OjzfmzOmzZsgWffPIJAgMD4efnB39/f1RUVGDWrFmSc/Ly8uDh4SF5v3PnThQVFWHr1q2S+bgilyIK5dALe+i6ujpYWlqitraW2lh5hWCz2bh//z68vLxgamqqdH15911vNTTFy8/9+/fh7++PtEt/aqxNSqApOg0ejwcAOJ16RWNtUgJN0WnczBRFYSp81qixNimBpug0rt4WxeQoqGSrFCRSHpRAU3QKAoEA90pFkVvLmE749cQZjbRLCTRFp3Dw6P9QzHACANBNLHDimmyAGVWgBJpC5xBC8OulbNCMXxjm32q0xMUr19VuWycCXVdXBycnJ7BYLERGRkqOczgcjB8/HqGhofjhhx863B6Hw0FdXZ3Ui+LVIfHwceTQXKSO8bu6YMdv59VuWycCfffuXbz33ntIS0vDH3/8ITl++PBhjBw5En/++SeOHDmiMDpnazZv3gxLS0vJi/InfHXgcDjYn5oNmpmNTFkGzxGHjqs3l9aJQN+5cwfJyckYNmwY9u3bJzmekZGBsLAw0Ol0+Pj44N69ex1qb9WqVaitrZW8CgsLtTV0Cg2zedd/8I9ZP7llxMIOe5JvgMvlqty+TgTazc0NX375Jc6fP4/9+/ejoqICgMg6TRxp08zMDA0NDR1qj8lkSvwHKT/CV4cnT4tw4k4V6Axjhec8NOuHhF2qxxrUiUAHBQVJQqEGBQUhPz8fAGBubi6JhdzY2EgJpp4T/90BPLPp3+Y5dIYJjmVXorikVKU+dCLQa9euRUpKCoRCITIyMtCzZ08AIjPPS5cugRCCzMxMSdR4Cv3j2s1buPLMtEOBgJ7Z9sf63W3kUWkDnQj00qVLkZCQgKFDh2LKlCkoLi7G1q1bMW3aNKSkpCAwMBCTJ0+GsbHinyKKV5sdh1PAsfJo/0QANBodl4qEePhPgdL96MTA38nJSeIGJcbHxwcAcOzYMV0MgaITycq+ixvVxoBtx+s02Hpj+4ET2LH2Y6X6ojZWKLTOvhN/gGvTU6k6NBoNfxXWKb3iQQn0aw6bzcatW7fAZrO11sft4jqVgmgWm7jjZPIf7Z/YAr30KaToOL+npGDihAnIyMjAoEGDNN5+aWkpHpbVgtf0UKaMV1Uo9bc1dIYJbtyvxxQl3C4pgX7Nycy+0/5JanAq+XcUHv2/Ns+pOrVFYdnDD+THSVGEQoH+7bffMH78+JciHzWF9uDwRc654v0ATfP4qShOic3YZTCykTZRIHwu+LVlMLS0B81QOtEpr6oQVae2oL5RuamQQoHOysrC5s2bERoainnz5sHPz0+phileDXgCkY90qYrhHTqKkU0PMLv3ki1w7qvRfhQ+FK5fvx4ZGRmYPHkyEhIS4Ovri2+//bbDBkQUrwa859ldS0q1I9CWz00bVIVhoNzDpEINzWazcfToUSQlJUEoFOLTTz8Fn8/H2LFjkZKSotYgX0aar6WDm5kFhu9AGAcHdfZwdAb/uYZ+VlOrlfYdunUsS68iTIyUe8xTePagQYMwdepU7Nq1SypIirzUtfoANzMLhN0IbmbW6yXQzzU0XyBs50zV6GYrayaqDKYM+SmQFaFwyvHOO+9g/fr1EmHesGEDAODLL79UY3gvLwzfgaCZmoHhO7Czh6JTBEKRhuYLtRNvSJX84i0xUVKgZTT0gQMHsHfvXmRlZSE5ORmAyKGxvr4en332mVqDe5kxDg7qFM2sbvQgdRAIBBAQ0RyVpyWBZjAY7Z/UBkxDNQX6rbfeQlhYGP79739j1apVIISATqfDzs5OrYHpgszMTIl9NQBYWVnB3d0dzc3NMlMlzu1s9GtqBsN3IB5bW8ksW7m5ucHa2hoVFRUyDgQWFhbw9PSEQCBAVlaWzDh8fHxgZGSER48eobZWem7q5OQEe3t7VFdXIz8/HxkZGVi4cCESEhIwcuRIDBwo+oX4+++/0TpKm7e3N0xMTPD48WNUVVVJldnb28PJyQn19fXIy8uTKjMyMpLYzmRnZ0sCvFRWVqKuXuS+xhcKUVRUJBPM0sbGBq6urmhqapJxwKDRaJLVr3v37qGpqUmq3N3dXXKs6dFNqQ0UmpEJDLvYAkIBeNWyIYjJ8y9YXVUFbt26JVXWo0cPdOvWDc+ePZOpJxPBf//+/YQQQuLj48m6deukXi8r4kjurV/iKPl5eXlyy2u+2kJqv9tNgoKCZMqSkpIIIYR8++23MmUjRoxos9/y8nJCCCFjx46VKduyZQshhJDDhw/LlLm6u0g+E4PBkCnPyckhhBAyf/58mbKVK1cSQkQZClqXOTk5Sdp1cnKSO+b3/7WWrFy5Uub4/PnzCSGE5OTkyJQxGAxJu35+fjLlhw8fJgcOHJDbX0dfTs49ZI59++23hBBC9uzZIxPBX0ZD29vbAxBpKE1RWFiIOXPmgMPhYNy4cVi5ciUAkfOst7c3PD09YWBgIOVvqAoXL16U0dAA4OzsjIyMDKlzObezQXuuoROHJ8rV0IDoFys4OFiqTJxGzMzMTKZdAOjatSsA4Ouvv5bJ1+LkJHLdf+ONN5CRkYHFG9/D9aM34LbQGZHWb0jOS09Pl9HQ4ueZzz//HIsXL5YqE983f39/mTG13Bw7e/asREMf/99ZVNLNsevzJeDwBXj/4/cxdepUqbo2NjaSvlu329I+4+DBg3I1tPiX0XLobBh2tX9RtwMa+tnprRg/eSrmz5klVSb2IY2KipKpJyPQ0dHRAICTJ0/inXfewahRoyRZOlVl+/btiI+PR3h4OCIjI/HBBx/AwsJC4jz7+eefK9Ueh8MBh8ORvBd7ffv6+sr1ejE2Npa1U2jxvi23gm7duqFbt25yywwMDNq0fxA7MsjDysoKxsbGaOoi2gkzcWSiHi+819vayHJ1dYWrq6vcMgsLizbHJJ56AEDqX5kw7uIAALB06IHyykoMVlDXxMSkzXa9vb0V1gMAk56D5W+sAGA6yt4BTqnI9qNHDxeF/Vpbyy4JKpTUuLg4XLhwAYMHD8bKlSvx4MEDRae2y+rVqxEWFgZA9CAizimiyHm2PfTF6/v0H/8Du9uLX4bHzfk6C8lACEFt84vwWzYOPXD972yd9K0MxkzlHioVCvSAAQOwZcsWJCcn459//oG/v7/Kg7K2toaBgQH27NmDQYMGSb61ipxn20OTXt91e39CxZwY1O39SeU2VOVa3lUYmb14iue4svHb2SM66fvin5dh2+uF9zWNRkN5XZPMNKezsbK0aP+kFigU6N9++w0TJ07ElClTEBkZiaKiIrUGduDAARw/fhwJCQmSY4qcZ9tDk17fnIuXQBobwbl4SeU2VIEQgrsV0pZuBiYG+KtA/ehBHeFa1l10c5KetnT16ItLV67qpP+O0t1O/nRPEQoF+urVq9i4cSMuXbqERYsWqSU0169fx/79+/Hf//5Xal1SkfOsLmGGDwPNzAzM8GE67TftahpKrJ7KHL9Xc0fr9jIPHj6C0EJ2Gba7iwcu3pBdhuxMxA+lHUXmofA///kP5syZA0tLS/z2229SZWvXrlVpUBs2bEBpaSlGjRoFAFi3bh0yMjKwdOlSzJ49G1988QXmzJmj9OA1QZf57wDz39F5vyevHYehnQG4BdLHa3tU4+CJA5g/bYHW+j586hx6suRbzZu7eOHilasIDw3RbKcCnkrVlF2Q6PCynSouNGJOnTolcyw8PBwAZJxnXwe4XC6yqjIBOb+mhmYGuJR3EfOhHYFO/+sGmE7yVxsAwMGjN86fP4ZhIcFq3fPWMGoKASf5KyFyUXEur3DZ7vz580hKSnreNsG4ceMwZ84clTp5mekMK7uDxw/gmXM5DBXYht3j5ODJ0ydwcXaRW64qhBCcuZQO78iJbZ7n5BeGw8dOYtqk8Wr3KfZVtG56itJSaTestgz8BU/+Vqk/mSv6448/YuPGjSgtLcWVK6LcF3Q6XWtp0TqbzrCyu3D/Dxi6KjaLFHjw8MPRPVgfu0Gj/Z48cw72/dv/jJY23XA7Kx0TuVy1bTHu378PAMg8fxTAUaXrizexOorMVV2wYAEWLFiA/fv3Y8qUKTAzM0NRUZFkh0vvYBiBe+22zh4K8x7l4Z6wbT8+Gp2Gm6V/QSgUqr2pJYYQgoz7j+A9vO1QXGI8Alk4dPQE5kyf2v7JbTBhwgQAgIuLC2J3Hket9YtNFLGbVWv3LMLnYZpzAz58ZxY8PT2V6k+hmnj06BHWr1+PhIQELF++HH379tVPazsuD4yBAwCuag8tyvLTqb0gbnzQ0Pb8tNTuKU6cO46JoyZppN9fDv0KNp2J/HvSmydF+Q+l/rYkP+MWgv19lRaqltja2mLBAtHzQMTVuzjLkZ2/t3bPsq7IwobVy5TWzkAbAn327FncuHEDAPDzzz8jNDRULwWa4TtQMofWNgKBALfKb4Bm1f7DloGVAc5mndGIQOfl5WHWzLbziX/3Wazc4z9uTcCDBw/UEmoxM0eGIDnxFohV2zu7AU4mKgkz0IZA02g01NXVoUuXLmhoaHjpdpA0hS7toE+lnESZXQkM0TEb3+y626ipqZEYO6lKfb0oOc/iDd/AyV1aQ3I5zagofopujs5gMKVjCxb9k4fvPv9YUl9dwoYEwDMpBQ+gWKCFdeUYP8FX5T4UCvTq1asxaNAg2NnZoaqqClu2KI6d8LKgjD00AInRS25urk7soX8+exCk53MPkUYBuBVcNBWLjKyaijmgMWgwcRQJFftxE4T0Rmzcvh4zxoqszVS1h7558yYAwMm9FwwMDSFokUKNwTSGf3gUTMzM8ay8BLVVlS8qtli2U9Ue2srKCmVlZZKdZgeDemSX5IFu/OI+tbSTtq/LhZ2VDwQCAQwMDJCXlyfzhVLKHrolfD6flJaWEqFQSHg8Xoftk3WNqvbQYnRpD+08vTvxT+xPPBbL2vmauBoT/8T+xD+xP6EZ0jRmD21paUkAkI0HzxBru+4y5Z/tOUx+vlVIxs37QO7nycjIUNkemhBCtmzZIvtZewWS7jHbFNpBi22cR4wYoZQ9tMLk9bt378a2bdvA5/MhEAhgbGzc4ZQRukacxFyePfTLoqGv3riK7Q+2wNTVGEaWhlIaumDPU7gtdIapm7GUhgYBmAUm+G7hHjCZTJU19Hc//IS9u3di48EzMhoaALq7uMvV0EX5D/HdZ7G4cuUK/Pz8NKKhBQIBYjbuRU33wRByGlGauOTFKoeAh7fdmzF9/GgMHDiwXQ1dUFAAd3d3qeT1Cqcce/fuxdWrVxEXF4c5c+Zg/fr1ik59aVDKHroFbQVa15Q99LErR2E54MWXzdDMAIZmJpL3Jo5MiTADgKmrqIxryUNVQyXGBI+VlClrD21p82Kn1sVT8W6dtZ0DrO0cZI4/LSpCSEiISvbQgOjLJt6BBgAfr0tIpzlIbJ7FqxwG5Q/wwTtz4Or6YkOprYdRpeyhbWxsYG1tjYaGBgQEBKC6ulphw686zdfSUbfrezRfS9daH0/qH6tUj2FlhJsPbqrVt7oe3SUajqrkZiPfGdjBsFFKmFVBoYZ2cHBAUlISmEwmVq9erdcCrYvdwjK2ajlD1K0LAJzna+zy1prbQnx+TbVmg9AM7OmEg0WyKyeuCgRdGRQK9I8//ognT55g4sSJ2LdvHw4dOqR2Zy8r2l6L5vF4KC0rBdeoSaas5SqHPOjGdFRz5TzNK0HtM9G8WNFac7v1a2rU6r81UcNCYPq/n1AHptRxFystCrSBgQHc3d0BAB999JHaHb3MaHst+q+//kLWjzltnlOwR9Y2Wky3xd3V6r9bN1EuCHnr0G0hfij099dsoE5bW1s4Mnmoa34h0IKmegzopb55hU7iQ3M4HLz11luorKxETEwM3n333TaP6xJdWNtVVoo0pNtCZ5g4SmslIVcITiUPTFsj0BnSjzTiFRB1Df5tnrsxObn3gru3Tztny9JbA7uErXG2MsX9khfvTeqeIDI0Ru12dZKSQlEKZFVTI2sy1zf7xEmwT55AzabNWn0oBEQrGaZuJlIv895msAnpCvPeZjJlrYVfVUIDB6tV38BAuehFHcHZWnp6YW/E1UgwIxkNPW/ePIWG3T/9pJojaUZGBubPny+VAtnPz0/h8fbYvHkz1q1bp9JYZKAB/PIK0CASbm1oabFTsKoY0tX7IZW3vNURGuu0E5EUAFztLIHbLx4MHbqqd43EyGjomJgYzJ07FzU1NfD19cUHH3yAoKAgtfIvK0qBrGpqZE15fTdfSwe/pBTg80AzNUU7BnAqo27AQmMDzeRvZDcqZ5NRmtf2vF8dAn28IWh4sTnUzVwzv0YyAh0eHo7w8HDU1dUhNjYW/v7+WLhwIUpLVV86UpQCWdXUyJry+uZmZgGcZhg6OsHI3R2m45TITqME6qb1MGeqZnnWmpL7tzt8LiEEXZnam5F69vKAGe/FVNFWWwItRiAQ4JdffkFeXh72798vCQ6jCopSIHd2amSG70AYurqB4TsQFu+/99LGhe7CUE/Dv2in4wJa9Og+/H00my6iJRYWFjA1EOV3ASHoaqaeZ4wYhVL6888/Y+vWrThw4AA8PT1x5IjqAVCmTZuGGTNmYN++fZg3bx7y8vKQkpKCxYsXSx3XdWrkzgqhqyxWxlZq1Rf79fHZdci+/ifMu3SVlCkyHy3IuAIHPy+1+m0Pc6bol0vY3AB3ZyUcaNugzZ3CESNGwNPTE4GBgWppaGNjY5kUyC9DauRXJQ2Fl6t6N1vs1/ftN9uVqvfD87+qGtu3hylTtHpCOA3o6eKskTYVSunKlStRXV2NW7duwdzcHJs2bcLhw4c10unLgnjLm33yJNgnTkJQXgYDO3sYevYSuWZpWND5TQKV6vmosHbcErFfn5eXF345dQ5OA1/E3BBvnrTcdOFyODAqz8OoqEiJZaE2MDUSCbSxoAnd7TUTf1yhQF+/fh2pqamIiIjAzJkzsWPHDo10+FLBMELT2Yvgl5YBbDbA58PQ3R2C4mIwBg7QmG2H2KySd48Ptons9rfCegq2w5WlpV9f5v1HsJHzBWm56fLg1jW8M2c2HBxkLe80iRlDJH5M0qzy0mJrFAo0j8dDdXU1aDQampqa1Jpy6AplPVZ6cXmgm5rhUXU12I2NAMMIho0NYAYMhgsBuvsO1Ig9tDiYTsHJp8BJ5T9XamoqTExMNBLB30DIRdbVi+hiJRKglgZLYnvo4ntZKPHriZKSErUj+Le2hxZjaWkJU+ZzmWp8JnMdO2IPrZTHyoULF4iXlxextLQkAwYMIMnJySr6k2gfVT1Wmq5eI1UrV5FBFl105rFiN8KGeMf3JE5T7WXKmN0ZxDu+J/GO70lgINuuNiP4bzx4Rq7HiroR/AmR77EyduxYsvhfa3XnsZKTk4P+/fujoqICtra2Gg0LpWlU9Vjpy+GCm5mFe7ez0PCkEIROA93CAqS2Dt6xH8IxOlojGvrevXuYPXs2psRMQq7vHTC6Gkk8VlrS2qew/h4bTw+V4MCBA/D29tZYjpU9B4/AYYAoK4F4Dr3x4BlY2tii8GEunOkNCBkSCED9HCttaeifjxzD2lXLETx2Jr6NXyZVrnGPlbVr16KsrAyzZ8/GzJkz1d7t0gXKeqw8W7UanCtX4FRcAoh3Co2NYWhlDUbaRSA6WqMR/N8IG4HcalGQmdYeK60xdTUB8x9TPEUJvL29pfpRN4J/sH8emh1EblctsbZzQMm9v/HehwtkppiqRvAHZD1WJHU8Rdac7k7dFbatMY+Vo0eP4tSpUxAIBJg0aRJmzpypsOFXFgIIyssBLhd49gyobwCeVYu8ndVz8pCLb39fdKnq2uHz7U3VMxtVxPjR0cjPlG+IZcU00Nnz0uhRo5CRkYEftmjOva/NraPCwkIUFBSAzWZ3SuxmbWM6fhzQxRKg01+8unSBQbduojINY2RkBBfTjrsYOZo5anwMgMhuxpzGlzleVVKEgV66u88Hbpbgi2tsHLhZ0v7JHaTN1Mg9e/bE/Pnz8eWXX2osxtrLhHFwEEwjh6P59xQIuVzQTE1gaGcPZvgwrW20uFp64AHa957n1/PRx0l7O3XuDraor5c2uy3NzcLij7UXl7oley49wtfnH6CJJ0RumWiOvHCY+l8mhVIaEhKCI0eOYOTIkXopzMALazshnw8wmQABDHv21GqcO38Pf/DrZbVja4yLTBE6OExr43gzOgoFt/+SOmZtYqize30muwRNPFF+8SaeEGeyNaOlFY5eXqwKfYN94iR42dlAYyNQUQHS1ATeo0dajXM3kjUKxkVm7Z7Xw9hVq6mSzczMYEp/sXPJrq9DT2ftzNnlMdrHASZGIvEzMaJjtI9mNnEUTjny8vJgZWUFR0dH0Gg00Gg0/PPPPxrp9GWBX1YKUlMDcDiAoSHA5ULYUA/2CdHuhzamHebm5nA2dkEB8to8z93SQ+N9t8bGlIlyItKSpY/uYcScFVrvU4x4enEmuwSjfRw0Mt0A2hDogoICjXTwSkAIwOMBNBroTGPwHxdoNaSBu6V7mwItaBKgn1PH4jirQ9iQwdh1/DwAgEkTqO1ZoywLh/XUmCCLUTjluHPnDsLDw+Hj44PNmzfj5EkV9mxfUsSBZWg0mmj9WQyPBwNHB4mNtLbo69gPgmbFhkoGT40wmvWm1voX4+3VB+xKkeOGBVM9JwRl2HPpESbsvIw9lx5pvG2FAh0bG4uDBw/C1tYW7733Hr744guVOyksLERERARCQkKk8hQCol0+JycnsFgsREZGqtyHMnAzs8DLzQX/n3ygpWOugQEM7OxhnbBJq+ak0WEjQStSvNbbHY4KN3M0CY1Gg/lzewpHO91lIDuTXYLMp7XYdyVf44KtUKCFQiGcnZ1Bo9FgZWWllk2sONf31atXkZKSIrWVKc73nZaW1uHk9ep6fTN8B0L47BloJiaAgQGahEJkc7lo4vHQfD1d697fjo6O6CZQbC7paKG79B8mz004ffppzzulNaN9HODrbAkajYbMp7UaW+EA2hBoNzc3xMfHo7q6Gtu3b1crx4qiXN+Aavm+1c31bRwcBIv33wNj4EAY9umDhwBGV5bjoUAAurGJyNdQy3Q3U/xU72iunQ0VeZgbi1yfOsO0wcvBAr7Olhpb4QDaEOg9e/bA3t4eISEhoNPp2Lt3r8qdKMr1DaiW71sTXt/GwUGwTtiE7ieOgR0kml4IHR3A8B2ok/QUdmbyNbSQK4STlWa8NzqCg4YM65VBPOV41sDF8Q/CNPpgqHAix+Fw8Oabb2Ls2LFISkpCaWmpQoOY1nzzzTdS3i1r1qxBVVUVjh8/juPHj0udGxQUBBMTE9DpdEm+7/bmj0wmE0ymZryEAaA8LAQ4fhS1Sz+B9dtva6zdtlDkJ8it4sF3kGZDb7U5jq6618xijaxJzSxGoYaeOXMmcnJysGzZMjAYDKXCdMXGxuLy5cuSl7W1tdxc38DLke9bbBqrSxPZLkxLEDlhbo3YDPRwUm4KpQ6aVAwdZeGwnhrXzGIUCnRdXR1GjRqFiooKLFu2DAKBav5wgHSubxaLhaKiImRnZ2Pr1q1YunQpEhISMHToUEyZMqVT8n2bPH/gNTZrfwdPU3Sz7gZBk1DmOIPH0Ol8Vt/MGhROOQQCAVatWoWhQ4fi/PnzakVOkpfr28nJSWJ03tn5vg2eP6RqI4abIkyYJiACWQ1NE9DVzt6qDEKh7JfqVUbh13Pfvn1wcXHB6tWr8ezZsw6vQLyK1D33/qjXYVB3Lo8LmpyrT2hErV9DZeHz2zeUepVQKNBOTk4QCAT49NNP8ezZM7i5uelwWLrl2XMblcIc7cVyk+mzthp0Y9nLLzDidzjGH4UsCgV6zpw5qK2txZgxY1BSUoK5c+fqclw6hfdUFGy8+bFqeVBUoYZdDbqh7OXnGDejtEy9FBTKwDMQPRQ+a9RMyITORuEcuqqqSpIKOTo6GiwWS1dj0imEEPCfPAEA8J+oFsVUFWo5NYCcqbKBFR05edno662bnTuPfr44lJwGa2d3nfSnbWRUxJMnT/DkyRPY29vjxIkTqKysREpKitai53Q2hU+ewKqsHABgXFKiVvB0Zahky99AMjAzwKMSzRvtyKO8vhkNAjqsXPtAaMjAvbI6lNerly2gs5HR0C2nFtu2bcO2bdsA6HaNVpdkXbwId57IQ8WloQG309MRNmKE1vstayyXe5xGo6G8UbNp1BRRxeaCTwBDAzoauAIwDQ1QxebCzkK3QTM1iYxAi5fQBAIBKisrYWtrq9PlLF1TW1AAm+fLDTZ0OgqzswEtCzSHw0EFR/E8uaRRc8Y6bWFjykA9hw8agC7GhhAS0bFXGRmBFq9snDlzBjY2NigvL8ebb76Jr776Su3A3S8jgsoXAVvooIFX3r4tibKIw9neunULAHA7JwvlTeUwKjCUmzToQf4D3Lx5E7m5uRofS0vsLIxfaW0sDxmBXrt2Lezs7CQXUygUYsOGDVizZg2+/PJLnQ9Q2/CfSkf04beK8KMJxOFslTEfCNgfIPlfW+Fs9REZgb548SIuX74seU+n07F27VoEBb288ZNVhc1mg/5Q2hVK8OABBAKBRqdZLcPZmpqa4qsDXyLTUpTuuGXyenHWK6FAiPEGUzBt7HSthrPVR2QEWtG2qz5ON84lJWFIRSVarj77FT7FhaNHETV1qsb6aRnOFgCE/+NLktOLEad8k1CLNkNvUchHRqAtLS1x9+5d9O37Yh00Ozv7lYht1x6tA//d+OUQjPl8POSLVjnEfy/v/Qk2Laz+NKklCSEoZbe/iqGrlQ59Q0agt2zZgsmTJ2PEiBHw8PBAQUEBTp8+jaNHj3bG+DRGXl4eevfuLXN8c4v/Y2ue23KcS8amc8lS5z148EAjQl1ZWYl6o9p2M8hVNVeq3dfriIxAe3h44Nq1azh9+jT++ecf9OvXD2vWrJEKU/sqItbMBw4cQK9evXD0o4/AyhdNNpoJwVMBH84GhjB+vt5+oU9vTNv2NXJzczF79myZkK6qkv84HxzzZhijbTvkOo72kl6KKa9vxuNqNggANytTvVjxkLv1bWxsjMmTJ2usk7q6Onh7e8PT0xMGBgZSzrCq5PvmcDjgcF7YHiizu+ft7Y2cEyfwTmERzFo8LwS0EjDnR//g/sWL8I6I6HDbHaGw9AkMzdt/4GwiTeByuVo1Ja1ic1HP4Uv+1weB1ol1d1ue3ark+1bHSVYgEKAp+RzM2jFst6LTUXVS1o5bXRqbGkEzav+yCw0EUl9abWBjyoAF0xAGNKCBw3/lt70BHQl0W57dGRkZCAsLk8r33R7qOMk+e/YM3co69sBlUlbWoS+YMjAYDECOYX9raISu9TjNdhbGCHARBQ0vrWvGvTLd2LFoE51EthZ7dgcEBGDEiBEYM2aMxBFWlXzfqjjJitMllJWVIY/JgEkju906eQwjGObnS9VXFwdbR/DzBDAwbXvawRAydZaIlMMXgjz/+6qjFYFu7fW9atUqREREyPXsVjXft7KIY/UpZdddWQ7Mni2pHxoaqvY4ern1guHvRkA7gZEsmZZaNQgrr29GFZsLOg1gPrfLdrXWXrRTXaEVgY6NjUVsbKzk/bJly2BoaIioqChkZGRg+fLlkjJxvu9+/fohMzMTGzdu1MaQJB43Bw4cwMMTJxCYerHdOjdHRcMjOhqzZ8/WmMeOo6MjLLld0YC2VzHsTLQXL6O8vhl3y+thbEBHs0AIGzMmDOk09LHTjjLRJTqZcixduhSzZ8/GF198gTlz5sDGxgbZ2dk6zfctDm7Tp08f1FXXwKcDqwcFNbXw8vKSqq8uNBoNjmbOeNCOQDuaay8cWBWbKxFmBwumXljZidGJQDs5Ocl4dvv4+HRKvu8zP/yAiTl3RPlU2mHwrb9x/uefNT4Gd0t3PMAdheW8Wj78ewUoLFcXsfB6PP9bxVbdo/9lQ7+CMnQA29NnYNXBWBQOhID52381PoZh/VngPVPsbW1R0gWjh4/WeL9i7CyMYWPKQBWbi8fVbPCFRG+E+uXPd6whHj434WxsaER2K0MreTuFYgzYzzceOhBzr6NEDo2E9e/dUG8tP2yCp0Ufra9wVLG54AsJCABDOo2acrxqnP3PfwAA/6pVLfZG5ukziIqO1shYDA0N4dXVGzdwVaZMyBNioIP2Y9uJNbS+bHmLeW0Eelg4C67X0jHAiAGTVlr4IZ+H2JpqfNPVCr0MZc1kGw3osAjU7Jw2pFcYrpVekjlu+JiJWbGzNdqXPFp6q+SW16GkngMHC+Yrv9Lx2syh5674F2zemgoXJhM+DIbUSyzEvQyNZMosjY1RNX0aRsyapdHxTH1zKro8lY1A6snorZPo/S0pqeeALyAoqX/1Y3O8NgJNp9Px4Z49KPh0KS5Zyw9l2xJCCC5064aqz1Zj0bZtGt/kMDY2Rqj9UPQocwMA9Chzg2ehF4b3jdJoPx3BwYIJQwMaHCx0H4lU09AIIVrIaq1b6urqYGlpidra2g7tNBbk5uLYv1Yg8q8bsKbTkc3lYnRlOc7Y2sGHwUApgMshwXhry1dwdOl4KmNVYLPZuH//vsQ9i6LjyLvvr42Gbolbnz74+NhR/L3wXdw1ltZKWeZmeLTkY3x8+FetCzMAmJqaYtCgQZQwa4jXUqAB0RRk7rp4cFetxOPnu4Z5JsYwXhePacs/1dvAOvrOayvQYsYsWIAnz1O4VbDCETl9eiePiEIdXnuBBoDQOaK8KhExMZ07EAq1oQQagIOzKOtU165dO3cgFGpDCTSFXqETgd66dStYLBZYLBYsLS0lMd6AzkmNTKHHEB2Sk5NDZsyYIXXs2rVr5IsvvlCqnebmZlJbWyt5FRYWEgCktrZWpXE1NjaSjIwM0tjYqFJ9is6htrZW5r7rdMqxceNGrFu3TupYZ6RGbg21Fqw/6EygGxoaUFVVJRN9qLNSI1PoJzpxkl2zZg3q6+sxadIkmXM1kRqZPN+911U6CYqXA/H9Ji2tN3Q131m8eDG5c+eOzPGlS5eS5ORkIhAIyBtvvEEqKyuVbls8h6Zer+ersLBQIgs6s4fOz8+H8/P1XgASJ1l5DrTK4ujoiMLCQlhYWKi0ZV1XV4cePXqgsLBQa2EUXsa+O7t/dfsmhKC+vh6Ojo6SY3phbacuylrr6Uvfnd2/NvqmNlYo9ApKoCn0CkqgIVo1iYuLUzpe3qved2f3r42+qTk0hV5BaWgKvYISaAq9ghJoCr2CEmgKvYISaAq9ghJoCr2CEmgKvUKvBJrP52Px4sUYMWIEoqKicO3aNSxZsgQAwGKxdDKGtLQ0xMfHq1Q3NTUVFRUVSExMRGZmplrjGDZsGLhcUcznmTNn4l//+hcA4PHjx5jVgTh98fHxSEtLU7pfefegpSmxttGr6KOpqalwdHTEd999h7q6OowYMQLp6emdPawOs3//fsTHxyNGA+EUBg8ejOzsbPj7+4PD4eD+8/jYf/31l0aSHylC3j0wNjbGW2+9pbU+W6JXGrpr1644c+YMrly5AnNzc1y6dElGMy9fvhzh4eGYPHky2Gw24uPj8c477yAsLEySzOj06dMICAhAQEAA/v77b9y5cwcsFgshISGSnOf9+vVDWFgYMjIyZMaRlpaGCRMmICoqClFRUeByubh8+TIiIyMRHByML7/8EgCwd+9eBAYGIiQkBI8fP0ZycjKWLVsm0Y4sFguEEHC5XIwaNQrFxcWIjo7G0KFD8e2337Z5LUJCQpCRkYH8/Hy4u7vDxMQEjY2NuHHjhmTsoaGhuHLlCgghmD9/PqKiohAWFoaSkhJJO1u3bsX69etVvgfffPMNMjMzsWzZMsTExGDSpElYt24dDh48iCFDhmDo0KHIzs5GYmIipk+fjuHDh2P682A/Z8+exeDBgzFp0iSEhIR0bACqmeu/vBw+fJiEhYURBwcHsmfPHhIeHk4IISQ8PJzcunWLLFy4kBBCyKFDh8iOHTtIXFwc2bVrFyGEkP79+xNCCImMjCQNDQ3k0aNHJDExkUycOJEUFhYSHo9Hhg8fToRCIenRowcRCAQy/aemppLw8HAybdo0QgghCxcuJFeuXCGJiYmkpqaG8Hg8EhQURAghJCwsjPD5fJKenk5Onz5N5s6dS/Lz80lcXBxJTU0lK1euJOnp6eTMmTNky5YtJDY2lly9epUQQsiECRPadIYoKioiixYtIvv37yfHjx8nW7ZsISkpKWTcuHEkIiKC1NfXk4aGBhIdHU3Ky8vJzz//TAghZNu2bSQpKYnExcWRTz75hLz77rsauwdz584lKSkpRCAQkCFDhhAej0ceP35Mxo0bR/bt20dWrFhBCCFkxIgRpKioiERGRpL6+nqSn59Pevbs2aG+9WrKcffuXYwYMQJTp06VaLOWjq8PHjyQaD4ejwcWiwUjIyP06dMHAKScC8zMzODh4QEPDw9s3boVs5/nK6yoqEBlZSXc3NxAbyNXi7hNBwcHcLlcODo6YtGiRbC0tASHw0F5eTl69OgBAwMDDBkyBABk5ppTpkzBb7/9hurqaqxcuRLvv/8+Vq5cCRqNhurqajx58kShQ4SjoyMqKytx5coVbNy4EQ4ODkhNTQWTyUR6ejrGjBkDQJRZt0uXLrh8+TJ+//13FBcXY8aMGQCA9PR0uLu7a/QeeHh4oKKiAu7u7jA0NISLiwtqamrkXjNAlMfS3Ny8wzGz9WrK8fvvv2PXrl0AAHt7e9ja2kIofJEd1cPDAxMmTEBaWhr+/e9/Y9iwYQAg4+XC5XLR1NSEkpISLFu2DJ6enjh69CguXLiA6dOnw9rauk1hltfmmjVrsH//fnzxxRdobGyEtbU1Hj9+DKFQiJs3b2Lr1q2g0WhS/nH+/v64d+8eSkpK4ObmBg8PD+zevRtpaWmIiYlpN3dily5dUFVVBVtbWwwaNAinT5/GgAED4Ofnh99//x1nz57FjBkzcObMGXTr1g379u1D3759JWPYtGkT6HQ6Ll++3PaFb4G8e2Bqaippk06nw9bWFv/88w/4fD4eP34sEfjW10wgEKCpqQlFRUWorKzsUP96paEXL16MxYsXIyAgAGZmZli8eDG+++47SXlAQAB++eUXhIeHg8fj4ZdffsG1a9dk2lm1ahWGDx8OOp2OnTt3gk6nY/LkyaitrcXUqVNhYNB2WmN5jB49GoGBgbCxsZGkgl6wYAFCQ0PBZDKRlJQEIyMjvP/++wgKCpLU6927N6ysRAHaV65ciQULFqCmpgaDBw+WHFeEn58fcnJyAIjyuhgZGSE0NBTBwcF44403UFdXh08++QT+/v6Ij4/H+fPn0bVrV9jb20va2LRpE6ZPn44///yz3S8xIP8e3LhxAwsXLpScY2BggA8++ABhYWEghCAxMRHXr1+XaWvFihVgsVhwdXXtcJ5IynyU4qVl165dWLRoEUpLSzFr1iyZXJfy0CsNTaFf0Gg0BAQEQCAQYPfu3R2rQ2loCn1Crx4KKVRDnd3Nl61vSqAp9ApqDv0aUFxcjHnz5oHNZmPatGk4deoUrK2tcefOHckq0IULF3D58mUwmUwcO3YMDIZmUiXrvG+lt4EoXjla7zD26tWLVFdXk6tXr5J3332XpKamkilTphBCCPnwww9Jenr6K9s3paFfAx48eCC1w1heXo6uXbtK7cj17dsXgGi3tKmp6ZXtm5pDvwZ0ZIdRW2nsdN03paFfA5TdYQSAlJQUVFRUYObMma9U39Q6NIVeQU05KPQKSqAp9ApKoCn0CkqgKfQKSqAp9ApKoCn0CkqgKfQKSqCVJC0tDd27dweLxUJERASGDBmiUiCVgoICSYgFsdu+PE6cONGhZKQAEBMTo1JwGACorq7Gf//7XwDAkiVLUFpaKvc8cVlOTo5c97XOhhJoFRg5ciTS0tKQmpqKs2fPYsWKFWq1d+jQIYVlX3/9NRobG9VqvyNkZWXh1KlTAIBt27ahe/fucs8Tl/3222/Izc3V+riUhRJoNampqYGZmRkAUZrnyMhIxMfHIysrS6LF582bBx6Ph4aGBowZMwbh4eFYvXq1pA2xfcPVq1cRFBSEwMBALFiwAGfOnEFmZqYkdNfy5csRFhaGsLAwiSY+fPgw/Pz8EB0djbt378qMLysrC1FRUYiMjERQUBAePHgAoVCI9957D0OGDIG/vz+uXLmChIQEJCcnIykpCSwWC3l5efDy8pIYEC1ZsgQnT54Ei8VCeno6EhMTkZCQgO3bt2Pu3LkARHkDfX190dzcrK3L3T7qGQe+fqSmphJ7e3sSHh5OIiIiyJtvvkmuX79OCCHE0NCQ1NbWEkIIGTJkCMnLyyOEELJx40by/fffkx07dpANGzYQQgg5fvy4JACLq6srIYSQgQMHkocPHxJCCElISCAFBQUkPDyc5Ofnk9OnT5O3336bEEJIbW0tGTBgACGEkN69e5Pa2lrC4/HI4MGDSWpqqtR4f/75Z3L//n1CCCFbtmwhGzZsIEeOHCFz584lhBCSm5tLtm7dSlJTUyXHxH0uXbqUHD16lPB4POLj40N4PJ6kLC4ujuzbt48IBALi7e1NGhsbyR9//EEWL16s2QuuJJRxkgqMHDkSiYmJMscdHBwkCSTv37+PBQsWAACampoQFRWFqqoqTJgwAQAQGhqKr7/+Wqp+TU0NevbsCQAy05icnBxcv35dMu9ubm5GWVkZrK2tJX2KA9a0xMXFBXFxcTAxMUFpaSmGDBmChw8fIjg4GIAoTELv3r3lzr3nz5+Pzz//HIaGhoiOjoahoay40Ol0TJkyBceOHUNKSgo+/PBDBVdNN1BTDg3SMm6Fl5cXDh06JPGZY7FY8PT0lDxIyYuJZ2Njg8ePHwMAPv30U1y/fh10Oh1CoRB9+vSRzN3PnTuHt956C127dkVFRQWqq6shFArlRiz96KOPEB8fj3379sHDwwOEEHh5eeHGjRsAgMLCQsyYMUPST0v69u2Lqqoq7NmzB/PmzZP5rOLz582bh4MHD+LRo0cYPHiw6hdQA1AaWkvs2LED06dPB4/Hg5mZGZKSkhAeHo63334bw4YNg6enp0ydnTt3SkKO9e3bF4GBgQgODsb06dNx8eJF/PHHHwgPD0ddXR0WLVoEJpOJ3bt3IzIyEra2tjAyMpJpc8aMGRg3bhzs7OxgZ2cHHo+H+Ph4nDt3DkOHDoVAIMDWrVvRo0cP3Lx5E3v27JGqP2vWLOzdu1dihC9m0KBBWLFiBfr374/AwEBwOBxMnDhRg1dQNSjzUQqNMGLECBw8eLDDMei0BTXloFCLhoYG+Pn54Y033uh0YQYoDU2hZ1AamkKvoASaQq+gBJpCr/h/uyQUWYk/LxAAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 184.96x184.96 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAALQAAADLCAYAAADHubkmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA6d0lEQVR4nO2deVxUVf/HPzPDsC+yCMgOiuAu+6qAKGqauaCiqWBqmj35mNYTagVuZYtmZmr+HgMf08xKy3IpTHBDTFEUEBEVENlkH/YZZs7vj2FGhpmB2dHpvl+veSn33Hu+Z+79zrln+S40QggBBYWWQO/rBlBQqBJKoSm0CkqhKbQKSqEptApKoSWQlpYGW1tbhIeHIyIiAj4+Pli1apVSdSYnJyMxMVEl7YuLi8Po0aMRERGB8PBwREZGory8XO56EhMTkZycjKysLGzbtk3iOe3t7Th06BAAYNu2bcjKylKm6WqHUmgpTJo0CWlpaUhNTcWNGzeQlZWF7Ozsvm6WkJ07dyI1NRVpaWmIjo7Gl19+qXBdo0ePRnx8vMSy8vJyHDhwAAAQHx+P0aNHKyxHE1AKLQNNTU2or6+HmZkZ6uvrER0djaioKHh7e+PEiRMAAA8PDyxcuBB+fn5YtmwZACA3NxcBAQEYP348jh49Kqzv888/h7+/PwIDA7Fz504A/F53+fLlGD9+PKZOnYrPPvsMkZGRGDt2LNrb23tsX11dHUxMTJCcnIzQ0FAEBwcjMzMTX3zxBYKDgxESEoJjx44BAC5cuAAfHx+MHz8eFy5cAMB/I8XFxQEAEhIS4OfnBy8vL5w4cULYK3/22WeIi4tDWloa6uvrMXXqVISFhSEsLAy3b9+Weg9++eUXBAYGIjg4GKtXr1bJ8+gRQiFGamoqsbGxIWFhYWTw4MHEx8eHHDt2jBBCyPXr18nJkycJIYRkZmaS8ePHE0II0dXVJZWVlYTL5RJ3d3dSXl5OZs2aRS5fvkwIIWT16tUkISGB3Llzh4wdO5Z0dHSQjo4OMmbMGHLnzh0SGxtLtm/fTgghZNq0aWTv3r3C/1+9elWkfbGxsWTUqFEkLCyMREREkBUrVpDGxkaSlJRE5syZQwghJCcnh0RERBAul0vYbDbx8/MjdXV1JCAggBQVFRFCCJk+fTpJSkoiqampJDY2lly/fl14TXV1NdmwYQMpLCwkYWFhQrmpqalk7dq1ZP/+/YQQQnJzc4mfn5/UexAdHU1+/fVXQggh33zzDWGz2ap/YF3QUf9P5sVk0qRJSE5OxuPHjzFx4kS4uroCAOzs7PDVV1/h559/Bo1GA4fDAQDY2trC2tpaeE5bWxtyc3Ph6+sLAAgJCUFOTg7y8vIQFBQEBoMBAAgKCkJubi4AwMvLCwBgbm4OT09P4f/b2trE2rdz506Eh4eLHffw8ADAfzs8fPgQ48aNAwA0NzejqKgIT58+hbOzs7BNXXnw4AECAgJAp9NhaWmJLVu2oKioSEzG3bt3sXTpUgDA0KFDUVlZKfUefP755/joo4/w+eefw9vbG0TN+3jUkKMXnJycsHfvXsTExKC1tRWff/45goKCkJycjMmTJwsfEI1GE7vW3d0dV69eBQBkZmYCADw9PZGRkQEejwcul4srV65g0KBBUuuQFzqd/0gHDx4MLy8v4Txg3rx5cHV1hZmZGR4+fCjSJgEeHh7IzMwEIQRNTU2YPHky6HQ6eDyeyHmenp64cuUKAL5ym5mZSW3/gQMHkJCQgIsXL+LBgwe4fv260t+xJ6geWgbCw8MRHh6OrVu3YsqUKVi5ciW+++47ODg4oLa2Vup1e/bsQUxMDPT19WFiYoJRo0Zh5MiReOmllxAcHIyOjg5Mnz5d2IurktGjRyMgIAChoaFoaWnByy+/DDMzMxw8eBAxMTEwNTWFnp6eyDVeXl4YM2YMQkJCwOVysWHDBlhbW4PFYmHTpk3C8zZs2IDFixfj4MGDYLPZ+L//+z+p7fDx8cH06dNhbGwMGxsbtXzXrtCIut8BFBQahBpyUGgVlEJTaBWUQlNoFZRCU2gVlEJTaBVasWzH4/FQVlYGExMTlazlUrwYEELQ2NgIOzs74fq7Vih0WVkZHB0d+7oZFH1ESUkJHBwcAGiJQpuYmADgfzFTU9M+bg2FpmCxWHB0dBQ+f0BLFFowzDA1NaUU+h9I12EmNSmk0Cr+cQrd1RtF8Dl58qSwfOfOnXJ5luTn5yMkJAQhISHYuHGjWPmDBw8QHh6OsLAwREVF4enTpyCE4M0334S/vz+ioqKExkIUKkCtxqkaoqGhgQAgDQ0NvZ4rsP3tDpvNJnFxcWTQoEEkISFBZtlTp04l165dIzwej0yePJncuXNHpDwmJoakpKQQQgjZt28fiY+PJydPniQzZswgPB6P3Lt3j0yZMkVmedpEc3MzyczMJM3NzQpdL+m5a2UP7evri4aGBtTV1UFHRwdVVVVoaGgQs//tCpvNxuzZs7FhwwaJ5T/88INIrx4eHo7S0lLk5ubC398fNBoNEydOxPnz50Wu27VrFyIiIgAAXC4Xenp6yMvLw4QJE0Cj0eDh4YH79++r7su/QNy7dw8+Pj64deuWyurUiklhd6ZMmYJz586Bx+Nh6NChSE1NFR4HgLNnzwqN4/v3748ff/wRRkZGeOmll5CcnCyxzrlz52Lu3Lk9yjUxMUFZWZnIsf79+wMArl+/jt27d+PSpUu4ceMGvv76ayxbtgy3bt1CcXGxEt/2xefBo8IeOxt50EqFfuWVV7B3717QaDRs3rwZZ86cQWtrK+Lj41FZWSn0RpGHH374AXv37hU5dvjwYZEZdmNjo9DYvSu//fYb1q9fj19//RWWlpaYOHEirl69isjISAQFBSEgIECh76ktFBQ/UVldWqnQ3t7euHv3LgwNDTFt2jR8/PHHYDKZGDJkiNBdSF6k9dCenp64ceMGfHx88Mcff2Dz5s0i5X/88Qc2btyIc+fOwcbGBgCQl5cHFxcXXLhwAVlZWXjyRHUP9EWisbERAFBR26iyOrVSoQG+9wWHwwGNRoObmxucnJzUImf79u1YtmwZ2traEBUVBR8fH1RVVWH16tU4fPgw3n77bdDpdOGPYcKECVizZg3Wr1+Pb775Bqampvj222/V0rbnnZt3+L6UlQ3iPpMKo6oZa18izyoHxfPD4n+vIwDI6JkrCI/Hk/v6f8wqB8WLwaPqJgBABcxw67ZqgvhQCk3RJ5Q8KUVhY+eE2sQWP/xxSSX1UgpN0Sfs/+E3sPu58P+g0fB3Ua1YuARFoBS6j7lZcQOHcpNxs+JGXzdFY/B4PFx+UAXQn6nfI6Yzfv79D6XrphS6j8mtyUELpwW5NTl93RSNcfj4b3jAdBU5RjOyxE+XbitdN6XQAFj79qNy6jSw9u3XuOxhlsNhyDTEMMvhGpfdV5xIvwuGkbnY8ZssI9xUcnKotevQ8tB04AB4ZWXgVlbCdMXrGpXtbesLb1v1RhN6nrjy9w3cbukH6ImXcSzcsP/nFOwbNULh+jXSQ7NYLNjb2wuDcwtob2/HK6+8gpCQkB7DSXWnvb0dLBZL5KNw2/btB6+qSuHrKeTj4O8XwTWXvsmVUdqK+vp6hevXiELfvXsXK1asQFpaGv766y/h8WPHjmHSpEm4dOkSfvzxR4lRNiXx8ccfw8zMTPhRxp+w9fffAQ4H0OG/rDQ57HgeJoR5lSz8mV+JvErFOwVZaWlpQWZZS4/n1FkMxX9/ONnjOT2hEYXOzc3F2bNnMXbsWCQlJQmPZ2ZmIjQ0FHQ6HSNGjEBeXp5M9a1btw4NDQ3CT0lJicJtM5g6FUwvL9CtrcErK+MruIZ4HiaEJfWtqGnhoKS+Ve2yjv5yBlUm7j2eQ9PRxd8PKxSWoRGFdnFxwaeffopz587h4MGDqOp8xTc2NsLY2BgAYGRkhKamJpnq09PTE/oPKutHaLriddj8fhLGS5aA6eUFg6lTFa5LXoZZDkdtWw1KWSV90kvnVbLQwuHCQIcOx34Gapf3d8ET0PWNez0vv5qNlpaee3JpaEShAwMDERQUBF1dXQQGBqKwsBAAYGxsjObmZgD8gNx96eAqUGxNTgq9bX3hYOIIcwPLPumlS+pb0crhobWDh4c1zWofejyqlk1JawwdkXrpqkIyNKLQH374IVJSUsDj8ZCZmYmBAwcC4Jt5Xrx4EYQQZGVlCaPP/5Ng0pjIq8nFncrbWP3Xmzh697DGZHftlRvbuWoderS0tKCskSvTuTomVsi890ghORpR6DVr1mDbtm0YM2YMoqOjUVZWhh07dmDu3LlISUmBv78/Zs2aBX19fU0057miqOEROFwO7tbkooXTgozydI3JHmJjClMaG8V52WDy2mFpyFTb0CM37x7q9KxkPr+S1XOiJGloZB3a3t5e6AYlYMQI/lqjIIvUPxVCAwgAG0NrGDINETggWKPyU74/gC8//QgnT55E1Msvq03O9ex70DGzlfn8ykbFbKSpjZU+Roemg6fNFQiwC8K7Aes0Lp8wmACA4hL1es1U1jeDxhB3T5NGbTNbITnU1jf447ubN28qPLNWhtKmJ7AwsMTd6pw+WZNmc/kWbpXVNWqVU9Mkn4LWtnAUsr6jFBrAn6dPw8fHB5cuqcYmVx68+vmgqagZxjDtkzVpdgdfadhc9abaqWuRT6FZNGNUVMi/Hk0pNICivHsAgLLO5URNMoI2Gt+/8SOc69w0bqRUV1cHhr4hAKC9Q7YVCIVlyanQbD1z3H8o//OQqtA//fSTMKmkttNczd/oqe0WU0OTPMl/goXD4jRqqPTwUSGMLfmJMgVDD3XBapVPlxjGFrhfJP+4XqpC3759G4GBgVi1apVKI9s8j5BW/tpra119n7WhrUOFns8yUvCoEGaW/EA4gqGHOiCEoIkt3xuApsNELatZbllSVzk2b96MzZs348KFC9i2bRvy8/OxdOlSLF26VOvWi3mCyWCr5ieFAlo5mpfNamoBs3NtmNB1wGazoaurq3I5WVlZaKipQjseiBzn1JSI/NudEpr8u5ZSFbqlpQXHjx/HoUOHwOPx8M4776CjowMvv/wyUlJS5Bb0PENa+D00r1X9BjrSaO4DheZ0mQfqm/RDaWmpMKe5qigoKIC3t3eP59T8tl3i8T0AVi+Ogbt7zwZNXZGq0N7e3pg9ezb27t0LNzc34fG7d+/KXPmLAum0JyF9sGzX0dEBAGCx1W++2Z2uE0Fjc0sUP36icoUWREeyCFsIXVcfkTLSwUZHQyV0zGxA0xF9M3BqSlDz23bh9bIiVaFfe+01/Oc//xH+vWXLFrz//vv49NNP5RLwIiBQZNIk/5hNWUrLSgEArPZ6jcvuOm42MTNHSXm52mTRLRyhZztIvMBhqErliCn0d999hwMHDuD27ds4e/YsAH4Y2MbGRrz//vsqFf68wOvsBbiNmu8ls/PvAABYtAbU19ejX79+GpPN5vLA7Pw/ncFAa7tiu3OyQGMwez9JBYgp9Jw5cxAaGopPPvkE69atAyEEdDod1tbWGmmQMmRlZQntqwHA3Nwcrq6uaGtrkzhU8vb2BofDQU11NQCgrPIpbt68CYBvw21hYYGqqioxBwITExO4u7uDy+Xi9m1xT+URI0aAyWTi4cOHaGhoECmzt7eHjY0N6urqUFhYiFt5/BWkxvZGHDtxDK8v5puv3rp1C4SIbnYMGTIEBgYGKC4uRk2N6M6ejY0N7O3t0djYiIKCApEyJpMptJ3Jzs4WLsc+LiqCPuvZvKGqplb4/QVYWlrC2dkZra2tYg4YNBoNXl5eAPgBKFu7zUG6Dl86Gp6iveLZpJCuZwSm+QAQbgfYVUXdbyFI5y7h48ePxcocHR3Rv39/1NbWipWJxbY7ePAgIYSQxMREsnHjRpHP84ogxln3z6uvvkoIIaSgoEBiOSGE3M/PJ846OmJlhw4dIoQQsnv3brGyqKioHuU+ffqUEELIyy+/LFa2fft2Qgghx44dEysb4DhA+J10dXXFynNycgghhCxZskSsLD4+nhDCz1DQvcze3l5Yr729vcQ2bz18moybNEXs+JIlSwghhOTk5IiV6erqCuv18vISKz927BjJzMyUKM9gkD9xjv+dOKw6LLHcev42AoAEBgaKle3evZsQQsj+/fvFYtuJ9dCCkK8uLi7i2q8gJSUlWLRoEdrb2zFt2jTEx8cD4DvPDhkyBO7u7mAwGCL+hopw4cIFsR4aABwcHJCZmSnxmvxr1/C2sTFW19cj0cQMTp9uwyh/f+H3nzNnDoKCgkSuEaQRMzIyklivYNjwxRdfiOVrsbe3BwCMHz8emZmZWLwxDndOZsPldQcEmj4L+p2RkSHWQwsm5x988AFWrlwpUiZ4bj4+PmJtYjKfve7PnDkDDoeD4uJinM97Ah4B9ry/CgDgNyYCn23dJHKtpaWlUHb3ervGxj58+LDEHlrgzGExeRV0bZ4tLtD1jIT/2sbtRHcEPfR7770nFjlW4EM6YcIEsevEFHrixIkAgJMnT+K1117D5MmThVk6FeXLL79EYmIiwsLCEBkZiTfffBMmJiZC59kPPvhArvra29vR3v7MXlbg9T169GiJXi/6+vpSl47q7hdgcOcM21dPF/X19SLn9u/fXxiFvzsMBqPHJSmBI4MkzM3NQafTwbHifw8DOz20k2erLIJXuSScnZ3h7OwssczExKTHNgmGHqWVTzEyJAJVpc+GU/rGplKvNTAw6LHeIUOGSDwuUGjd/k4SJ4U0ho7E44LhiZOTk1S5FhYWYsekampCQgLOnz8PX19fxMfHK5UHZP369QgNDQXAn2DqdHpYS3Oe7Q1Ven1zuoyP6aCBUyw+ZlMX17Oug2P2bIfwaetT4TKeuql4WgNj034ixzigq9HcQTNmQ1KljBw5Etu3b8fZs2fx6NEj+Pj4SDu1VywsLMBgMLB//354e3vDwIDvFSHNebY3VOn1zS0SNYDpKCpSuC55uV2QBR3TZy9JlhELRcWakd/MFv/hmNk6If8FT2DUo3HSjBkzEB0djcjISJSWliol6LvvvsMvv/yCbdu2CY9Jc57tDVV5fZcUF8OqW8Ie/UePlAp0Ig91bbUi41CaOUH2fdXESe6NZrZ4TzzAZSBudUbVVzldvqc6karQ6enp2Lp1Ky5evIjly5cr5ZF97do1HDx4ED///LOIrYA051lNcenI9xjdLvpg/ViNOH/0qEbkN7SLLunpmDBQWKaYc6g8cLlctHSI2z/rMHVR26Sm3VINKbTYpPB///sfFi1aBDMzM/z0008iZR9++KFCQrZs2YKKigpMnjwZALBx40ZkZmZizZo1WLBgATZt2oRFixYJZ9Saoun632B2u9FGdDpqL18BVqxQu3xWu+hGDo1GEzumDnJy78LcXvIWd7OcVnHPGzIv29GU+IX99ttvYsfCwsIAQMx5VlM8fvQIlnckv94Nbt5ETU2N2n9gDex6oNttrW+rU6tMALhxOxsDhkjOC9jCBTgcjshS34uE2JBDsGx37tw5xMbGIjY2FosWLcKPP/6o8capk5Sv98CvRbJ1XVB9A0599ZVa5fN4PNS0Vosdl3RM1dS3tIPOYEgss3IejMybL679u1gP/d///hdbt25FRUUFrly5AgCg0+kqy/T5PNDW1gbepYugS3nr6NJoaE1NBfeDD8CQ8uCV5d79PLBMGoBub/iKpnIQQpR6I/YGq60DdlLKrB1dcCsrFYEB/iqRJdhskWbzLA3B+fI6LosptMCI/+DBg4iOjoaRkRFKS0uFO1zawC9f7cbYklKRlAjdCcovwO/ffotXli1TSxvOXDkDnQF0cLp5GVUZVCK/4B48B0veqFCWkpIS6JhJH0rRaDSw2lS3Fl3UuQwqzea5N27fvi3cw5AFqeajDx8+xObNm7Ft2za8++67GDp0qFZY23G5XNSd+h1Gvex+WjAYuHz8OMjSpWrpLfOe3gXNTrxehiMdP//1MzYMVs+9/utSOpyH9Oy3yGLz0NHRIdwAUwbBXMzy5bVgWsq+ASawh+bJuSEjtcVnzpzB9evXAQBHjhxBSEiIVij0qaQkBOcX9Ng7C/DKzkHqiRMYN3OmStvQ2NiIvEbJjhI0Og13KrNUKq8rVQ3NcNCVED6/C7buI3DpylVEhI1RWp5gE41pKcUeuhceV8oXL0TqU6XRaEIbiaamJjFDmReVshO/wFJG2xRHAtz//nuVtyH5xLdocZHuifGAfh9381XvGcThcNDQ3vuynJWdA7LvP1SJTFlDJEujrkk+52GpPfT69evh7e0Na2tr1NTUYPt2xcZAmqQ3e+i8mzeBGzeQzeP/OEd0bvI86bSfeNDxbOzooKMDczodhn9fx68//QTHLm5oytpDX3x4AcSGoKW8Fa1lfOOk1rJ20HRpMLDTB5x4+OSbj/HvBW+LXKuMPbS+vj7+SDmH8qcNqDv9LJ5gVRl/8lVa+AANtdVoaqgHAJTnZcHSkClMJaKoPXR5RSUA8UmhrPbQXW3UBchlD92Vjo4OUlFRQXg8HuFwOHJYKGsWRe2hS+wcSImdA/GUYA/9ZT9zUmLnQDabmqnUHvqtf79FvD4aRtxWOorbCDvrE5/k4cQneTihMWgqs4e2sbGR2E5ZP/fv31fYHnr7LnF7cshhD93f3lkue2gaIZLHEvv27cPOnTvR0dEBLpcLfX19mVNGaBoWiwUzMzOJ9tBde+j/W7gQQ8ue+c0N0uFvHlxpb8PWRhbeNTaBY+dEyIbOgAmdjgYeDxedHLH0wH+F1ynTQ+8+vgt33DPR0cwFu4qN1rJ2FO1/ApfXHWDoos/voQE03m/GTIM5mP3SHOG1ivbQBQUFiImJwSsLXoPTMB9wuc+GHR2cdvB4BPaug9Da3CTsoUEICq+n4czJE8jMzMSQIUMU6qF37vk/bHr/PbFJoSw9dOX/1iBkZhx2bXhLpEzQQxcVFcHV1RUNDQ1C0wypQ44DBw4gPT0dCQkJWLRoETZv3izt1OeGnuyhW1pasK8Xj/XPmqSMa2ur8TqdjtGjR4scltcemsVi4T6PH3ZMx4gBHaNnsZgN7PSEygwAJoON8OjRA4n1K2oPPSRoHEYGhUltb3fYdZUA+MMTRe2hOzp/PNImhb3ZQxuamMtlDy1VoS0tLWFhYYGmpib4+fmhrk79W7Lq5HY6P5D4rn4WGNRtOaqNEDzhdsCBoQP9bkt0Dzo6sKq+Fjk3b4optLx8+/MBtLg2gi7jUlQBIx85eTkYPkQ18e5M+okrQE8Y9ld+78HEyFCp63UY8i2ZSlXoAQMG4NChQ9DT08P69etfeIWu7zR/HaSjI5wMdsVPUibILlQqaSddUFCAP/8+izZb0cj0XSeF3SGEYM/hr7Fny16lZEsZVfaKtaOLUnIBwNHORqnrTfTkWwuXevZ///tfPH78GDNmzEBSUhKOasikUl2QdsVSHAhoUyK5Z0FBAQYPHtzjOUX7JQcmvId9eDt2jVzRg7rzREFbdlVsKLk6S0+yKQvmhvIZSUlVaAaDIXRDf+utt6Sd9uKgQPDsrhCu4maVgug/TrEDYOQq+grmsXlor+ZAz4oJuq7oUEQwYaysrFRKobPv3lP4WkDxHh6A0nEQR7i7yHW+RlJStLe3Y86cOaiurkZcXByWddpHSDuurRi5GsLQRTwpj3HPnTdy7mXLZc/QnaZW5SKbKuPBIzAuYleKb9T0FgoMAHxGyTd/0IhCC1IgL1++HJMmTcLChQuhr68v9XhvSPP6VifqtH7rjYKSB72f1AMcJaPzl5Urntn13j3+26H2jGLmuPLapIsp9OLFi6U+vG+//VahRmVmZmLJkiUiKZC9vLykHu+Njz/+GBs3bpSvETQlvY77UKGb2uULWNgdQYiv0kL5fhiC88srKxWWPX36dLDZbHxx6haarEV7W4EBkkTDJR4XC1zZcg+1xBQ6Li4OAD+WRnh4OEJCQpCZmYmLFy/K9026IC0FsqKpkdetW4c1a9YI/2axWL2HMlAytghR9gehBMrGjhaEOhMElJGXEgnhuGTFysoKK1euxJWS7bgCycZJktaodZ/mYdXShXLLE1NogWvU5s2bsWoV/wb4+Pjg2LFjclcuQFoKZEVTI+vp6UFPr+dlNjHoyvWwNCWvV4ZWrnJxq80t+UHNV27ZBXtX2S3eSgsfYM/7q9BfBXENRziY4/JjLmh02RwmBhlz4Ogg/zq41DE0l8vF999/D19fX6SnpytlGytIgTxs2DBkZWVh69atPR5XB+zOACpdDZBkQXA+W0MBYCTRylFOofV0+Utf9q6D4DpkhNzXqyIH+8yosdj/ySnw+sv2g/KwVUymVC09cuQIduzYge+++w7u7u5K+RTOnTsX8+bNQ1JSEhYvXoyCggKkpKRg5cqVIsfVmeqiunNjaFW9YhtE1TXy2eV2RTBplbR50hOC85vblItbzVDy7WKughC/7oMGwp7RCFkcsXhtzRjtNUAhOT3uFEZFRcHd3R3+/v5K9dD6+vpiKZA1nRrZqrOX2dXPXGiUJAsPOjhYVV8HKyUeam4uP3iLtM2T3qivVW6XVpeh3Ph/gI3yQw4ajQYXKyOUyLAdoN9QhElhsQrJkaql8fHxqKurw82bN2FsbIyPPvpIqXF0X6PbuUoxSIcpceu7N5jd4w3IgSDYo8vrDjCwk33sL9hYYZoqF1JAn6nc6qyqwjk4mBsCMrzobHTZwnAa8iL1m167dg2pqamIiIjA/Pnz8ZWa3frVDeEoNwYmco69uyKYwBrY6UncWOkNuo5ynuceg9x6P6kHVLUGP3CABXgVbaAzex5a2pnJf48ESFVoDoeDuro60Gg0tLa2qsRhUt305LFSWvUUgOikUNBTP+zgoIUnuvkg8FhpEHhOVFQIPSfktYcWBJNsLWsH05wJppmO0B66K0KPFQAtxc+8Wdpq23Dz5k2F7aFBnr3nHxfkgdttgmvr5AoDI2PUPi1HQ82zuCBd160V9VgxNzdHZWUlSktLYW6gA86Dq4C5vTA+NCDBm8WwClwuFwwGAwUFBWKJgxTyWDl//jzx9PQkZmZmZOTIkeTs2bMq8i9RPcp6rHgzxaPlCzxW3jAyVqnHikOMLfFJHt67x4qO6jxWrK2tCcCP0m9hbStW/v7+Y+TIzRIybfGbEr9PZmamwh4rhBCyfft2iR4rtnE7pXrJCLxQoqKiVOOxkpOTg+HDh6OqqgpWVlZ9uvXbG7J4rKx+5RV88+efIpNCWXroy21tmFdbjTdefhlLO6Pxy9tD5+XlYcGCBXB53QGmw4xl76FL+WPoiJhwfP7udoV76OzsbMTFxWHJhm0wNDaRzWMFQO7VNFw4dVwpj5WuPTQAvPlJEooMB4GuZwReezMqkleL7BRym+uQGOmA5UsX99pDy+Wx8uGHH6KyshILFizA/PnzYWZmJu3U54aePFbsOycZkiaFA3tY9TDr3GG0s7UV85yQ12PFwE4PTDP+Le/usdIdQ2cDfn8EwM3JTUSOvB4rgqHSga3xUuX1hImJicIeKwD/xyaY5Hl6XEQ5h78WLfBK6bpTaFBxB6/GzBZGrOpp61suj5Xjx4+jtrYWR44cwcyZM2FjY4MjR45Irfx5h2ag+ERDFdcrg4O54hkKAL49BQCkZ+djUMA4kTLBbmD3XUQel4v2wiwsmjdXKdPV7tj3MwB6iGtvo89TqvPscaZXUlKCoqIitLS0aDx2s6ohnYl+ciSkXOjNBQsALJyUM1QHAB5HMZvsIK+g3k/qASsrKyxduhTN+5NgLWWnsPsu4oOsv7FscSyclEj3IYlhLgPALWkCQ99YYrl9P+U213pMjTxw4EAsWbIEn376qdKJg/oa0jm2/k+DYpsUIwMDFZYtsAluud0OOlP0PvZm4A9AatIieRnpMQg5ZU9gZefQ67lcVrXKlRkAIkMDYfB7Mtj6wySWO1oYSTwuK1IVOjg4GLt371aq8ueJ5W+8gQffJmF6SysMxHph/m6gtF3EGw72SkVfFdgEl/xeBvwu//WCNHLKMiYkCH99nSyTQpsbqCc+tLW1Nez0OCiSUMZta8JIN2lxUWVDqkLn5+ejubkZRkbK/WKeF2xsbBDm7wf/q9ekniNtF7F8uHJe14Ix7Jlbp5FnfUekrGtcju67iLwOHqKN5qtsDKujowMT3d7ftI31dXBzUMyWQhacLY1QJMGjzaC+EFFhS5WqW6pCFxQUwNzcHHZ2dqDRaKDRaHj0SP35P9SJzmAPkPQMuZYgOwiB/lDlQtsKxrDcIx141JQPhp64UknaReQ9oOGNRW8oJbs7pno6vcaffnL3FubHqjZAZVecLQ2Bp+LH7Q14sLKyUqpuqT/XoqIisNlsFBUVobCw8IVXZgAY9co03JdzLnBLj4mQuXNVIj964mwwi2W3I3Gmu8DeTrVxuUcP9UDVk+Iez2Hy2CLr+arG19MVHU3iu3wulsrF8AB6UOjc3FyEhYVhxIgR+Pjjj3Hy5EmlhfU1o/z98cDVRa5rqgZ7wNlNOVsIAZaWlnDVlX21yMPcUyVyu+Lr442nhT17gZvqq9fMISI0EKbN4paHjubKL41KVehVq1bh8OHDsLKywooVK7Bp0yZpp/ZKSUkJIiIiEBwcLJKnEODv8gkiXEZGRiosQxZoNBp0R4+S+XxCiFzny4KHxRCZwgJ0VHMxbvR4lcoG+JFIjXSkDzfaWppha6HeTTQTExPYG4reA25THXyGSM7MJQ9SFZrH48HBwQE0Gg3m5uZKzbQFub7T09ORkpIispUpyPedlpYmc/L69vZ2sFgskY+seL7yCgplHEPfYTIRtGCBzHXLwoywmeCW9X6eZbU1Iseo5wdu2kM0ouK7WQgPDVaL3K44dVueM2kqQUSI4kujAqQqtIuLCxITE1FXV4cvv/xSqRwr0nJ9A4rl+1Ym13dgZCRyB8n22n8ybCg8RsjvstQTo0eMhkNL7+31NB+itoRFbo4D0FArJdtWC0thW2R5cLYyBLq8qRyMVePqJVWh9+/fDxsbGwQHB4NOp+PAgQMKC5GW6xtQLN+3Mrm+aTQa9AMDwOvltc8mBMbByu3QSZM/1ErypoKAjqYOhLgrHlimN8YEB6Hs3h2JZSZqHj8LCBk1BNymZ5tczpaqWR6WqtDt7e2YMmUK1q1bh6amJlRUyB5sZNeuXQgNDRV+zpw5IzHXN6BYvm9lc32Pf+MNXOslKuYVc3NMXrlSrnplZYLXRHRUS3c4MHlijpmTZ6lFNsAfw+oQtthxHo8HMwP5vXkUISTAF2bszs6LEDhZKL/CAfSg0PPnz0dOTg7Wrl0LXV1ducJ0rVq1CpcvXxZ+LCwsJOb6Bvom37eDszNqfXx6PKfV30+iNZcqiBwTCatqW6nlg00Hq9VhGACMmeLDmYqih/Aa0fPbQ1Xo6+ujvxG/DR1NtQgZrZo0dlIVmsViYfLkyaiqqsLatWtFbGjlpWuu7/DwcJSWliI7Oxs7duzAmjVrsG3bNowZMwbR0dEay/ftOmM6KohkY6EHdBpGvjpfbbLpdDo8LCQvyfE6eBhuM1JtsgWYGuiC1y2AZe2Thxg1QjWxqGXBxpT/ozVtr0agb+8Rs2Shx7gc69atw5gxY3Du3Dmw2eKvKFmRlOvb3t5e6PndF/m+x8+ejf1ff42pD8WHOPeGDsVyNS8hjrb3QnqTeDQqWrEO5qxQzUZOT4wa5olLj0W/uwEDGnW1szLmb/VbGNBgaKjmIUdSUhKcnJywfv161NbWyrwC8aJAp9OhFxwsNjlkEwKTMaFq99CZOXGWxF1DJzipfHdQEqNGjEBNiWhEUGM5g4srywBz/lKwuZGcUbB6QKpC29vbg8vl4p133kFtba0wI6g2EbF0KW52CymWYWqCiStWqF22paUlnJguYscHmvcSW1dF6OrqwoD2bMjB7ehAvx48aNTBYFf+8qWZnuos+6Qq9KJFi9DQ0ICpU6eivLwcsbGKBf54nnEZNAhVw0QnI80jhittICMrbv1EJ8DcFi5GOqp2Z7InjLss0VWVFiPAW3OyAWCALX+920SFpqpSFbqmpgbvv/8+Jk6ciI0bNwqdHLUN3REjQDqd93ggMBiluYc6wnEkuG3PJtvMJ3qYOm6qxuSbGxmA2xnWobm6Ah69pM1QNQJXK0Nd1W0giSn048eP8fjxY9jY2ODXX39FdXU1UlJSVOpX9jwxZPJklHWOowsYDARGR2tM9kvhU8B4+qx3smM4qG2pUBKBvl6ofFwEADBkMjTulSSQZ6zCHlpsFtB1aLFz507s3LkTQN9GsFcnPiEhOGnVH6ipRrWtLQZ6eGhMtqWlJazo1gD4YQeczCR7cquLwe7uaK37HwDV9pLyokzcwO6IKbRgCY3L5aK6uhpWVlZqsyl4HmAwGKDbDwDy86Bjb6/xH+4AI75nCOERuPRT3tpMHmg0Ggw6N1hMDFS30iAvlubK23AIEHvHcLlcvP322xg6dChmzJgBDw8P/Pvf/wZHgre0tsAYwFcqhp363I6kYW3Ej+zZUcNFhO+4Xs5WPXo6fBUYrGT8O2UwN1OjQn/44YewtrZGfn4+0tPTcf/+fVhaWmLDhg0qE/q8YdE5GbIdpplt364Mc+VvLhk1mmDEMNVa9smCfqdCq8PDW1YMVBjzRGzIceHCBVy+fFn4N51Ox4cffohAJdz4n3eWvPMOQqZOhaen6j1EemOYB/9HZKFv2ScBMY0N+NvP3W1sNIkq7VbE7qC0L8Zkqset/XnA0NCwxzBX6kTgVW+ub94n8m2sVRPzQxnkzpfTA2JDDjMzM9y9e1fkWHZ29gsR2+5FxlRXdeNIeXB26ruhhgBVdpZiPfT27dsxa9YsREVFwc3NDUVFRTh16hSOHz+uMqEU4pjo9Y1CBwUFITMzs0+GWwJUubIk1kO7ubnh6tWr8Pf3B4vFwrBhw5Ceno5Bg2RPB0YhP6b6ffMGFAy3VGXt1tdInIXo6+tj1izVeUywWCwMGTIE7u7uYDAYIs6wiuT77ovUyOrGxkr9fnz/CNQUVF+Eq1evkk2bNkks+9///kf27NlDuFwumTBhAmltbe21voSEhB6jvr9INDc3k4yMDNLY2NjXTdE4zc3NJDMzkzQ3Nyt0vSCDQtfnrpHN+548uzMzMxEaGiqS77s3lHGSfd4wNDQE05mBE8U/4WbFjb5ujkZRx3BHIwrdk2e3Ivm+lXWSfd7IrclBC6cFuTU5fd2UFx61KHR3r++2tjapnt2K5vtWJax9+1E5dRpY+/ZrXDYADLMcDkOmIYZZas6fT1tRi0J39/o+f/68VM9uQb5vQgiysrLgoUFrNwGtv/8Ozq1baP1dgeDNFM8VGhlySPLsFnh9z507FykpKfD398esWbPU7r4vCYOpU8H08oLBVM0Z13eFGnKoDqlp3V4kBGnduqb3epG4WXEDuTU5GGY5HN62vn3dnBcGSc/9+U8P+w/A29aXUmQV8WJnAqKg6Aal0BRaBaXQFFoFpdAUWgWl0BRaBaXQFFoFpdAUWgWl0BRaBaXQFFoFpdAUWgWl0BRaBaXQFFqFRhR6x44dCA8PR3h4OMzMzHDz5k1hmSZTI1NoPxo1H83NzcXWrVtx5MgR4bGMjAykpKTggw8+kLkeSV7fjo6OL6z5KIViSDIf1eiQY+vWrdi4caPIMU2nRqbQbjSm0E1NTaipqRHLBKDp1MgU2o1aDPx37dqFY8eOCf/esGEDGhsbMXPmTLFzAwMDYWBgADqdLnSg7d+/5wCCenp6IgH+BKMmbQg4QyE7guctMmpWOlqIjKxcuZLk5uaKHV+zZg05e/Ys4XK5ZPz48aS6ulruuktKSiQGnqE+/4xPSUmJUBc05oJVWFgIBwcH4d/Z2dlISUnBmjVrsGDBAmzatAmLFi1SKDWynZ0dSkpKYGJiolDgP8GksqSkROOTyr6U3dfylZVNCEFjYyPs7OyEx7TCSVZZ+tLJtq8dfLXtu1MbKxRaBaXQFFoFpdDgr5okJCSoNDXCiyC7r+WrQzY1hqbQKqgemkKroBSaQqugFJpCq6AUmkKroBSaQqugFJpCq6AUmkKr0CqF7ujowMqVKxEVFYUJEybg6tWrWL16NQAgPDxcI21IS0tDYmKiQtempqaiqqoKycnJyMrKUqodY8eOBZvNBgDMnz8f//nPfwAAxcXFePXVV3u9PjExEWlpaXLLlfQMupoSqxutCniempoKOzs77NmzBywWC1FRUcjIyOjrZsnMwYMHkZiYiLi4OKXr8vX1RXZ2Nnx8fNDe3o579+4BAP7++2+EhIQoXb80JD0DfX19zJkzR20yu6JVPXS/fv1w+vRpXLlyBcbGxrh48aJYz/zuu+8iLCwMs2bNQktLCxITE/Haa68hNDQU7777LgDg1KlT8PPzg5+fH27duoXc3FyEh4cjODhYmPN82LBhCA0NRWZmplg70tLSMH36dEyYMAETJkwAm83G5cuXERkZiaCgIHz66acAgAMHDsDf3x/BwcEoLi7G2bNnsXbtWmHvGB4eDkII2Gw2Jk+ejLKyMkycOBFjxozB7t27e7wXwcHByMzMRGFhIVxdXWFgYIDm5mZcv35d2PaQkBBcuXIFhBAsWbIEEyZMQGhoKMrLy4X17NixA5s3b1b4GezatQtZWVlYu3Yt4uLiMHPmTGzcuBGHDx9GQEAAxowZg+zsbCQnJyMmJgbjxo1DTEwMAODMmTPw9fXFzJkzERwcLFsDFDPXf345duwYCQ0NJQMGDCD79+8nYWFhhBBCwsLCyM2bN8nrr79OCCHk6NGj5KuvviIJCQlk7969hBBChg8fTgghJDIykjQ1NZGHDx+S5ORkMmPGDFJSUkI4HA4ZN24c4fF4xNHRkXC5XDH5qampJCwsjMydO5cQQsjrr79Orly5QpKTk0l9fT3hcDgkMDCQEEJIaGgo6ejoIBkZGeTUqVMkNjaWFBYWkoSEBJKamkri4+NJRkYGOX36NNm+fTtZtWoVSU9PJ4QQMn369B6dIUpLS8ny5cvJwYMHyS+//EK2b99OUlJSyLRp00hERARpbGwkTU1NZOLEieTp06fkyJEjhBBCdu7cSQ4dOkQSEhLI22+/TZYtW6ayZxAbG0tSUlIIl8slAQEBhMPhkOLiYjJt2jSSlJRE3nvvPUIIIVFRUaS0tJRERkaSxsZGUlhYSAYOHCiTbK0acty9exdRUVGYPXu2sDfrmqX0/v37wp6Pw+EgPDwcTCZTmEquq3OBkZER3Nzc4Obmhh07dmDBggUAgKqqKlRXV8PFxQV0uvQXnKDOAQMGgM1mw87ODsuXL4eZmRna29vx9OlTODo6gsFgICAgAADExprR0dH46aefUFdXh/j4eLzxxhuIj48HjUZDXV0dHj9+LNUhws7ODtXV1bhy5Qq2bt2KAQMGIDU1FXp6esjIyMDUzoxftbW1MDU1xeXLl/Hnn3+irKwM8+bNA8D3yHd1dVXpM3Bzc0NVVRVcXV2ho6MDJycn1NfXS7xnAD+PpbGxca9ueQK0asjx559/Yu/evQAAGxsbWFlZgcfjCcvd3Nwwffp0pKWl4ZNPPsHYsWMBQMzLhc1mo7W1FeXl5Vi7di3c3d1x/PhxnD9/HjExMbCwsOhRmSXVuWHDBhw8eBCbNm1Cc3MzLCwsUFxcDB6Phxs3bmDHjh2g0Wgi/nE+Pj7Iy8tDeXk5XFxc4Obmhn379iEtLQ1xcXFwcXHpsQ2mpqaoqamBlZUVvL29cerUKYwcORJeXl74888/cebMGcybNw+nT59G//79kZSUhKFDhwrb8NFHH4FOp+Py5cs93/guSHoGhoaGwjrpdDqsrKzw6NEjdHR0oLi4WKjw3e8Zl8tFa2srSktLUV1dLZN8reqhV65ciZUrV8LPzw9GRkZYuXIl9uzZIyz38/PD999/j7CwMHA4HHz//fe4evWqWD3r1q3DuHHjQKfT8fXXX4NOp2PWrFloaGjA7NmzwWAw5G7bSy+9BH9/f1haWgpTQS9duhQhISHQ09PDoUOHwGQy8cYbbyAwMFB43eDBg2Fubg4AiI+Px9KlS1FfXw9fX1/hcWl4eXkhJ4ef+1BHRwdMJhMhISEICgrC+PHjwWKx8Pbbb8PHxweJiYk4d+4c+vXrBxsbG2EdH330EWJiYnDp0qVef8SA5Gdw/fp1vP7668JzGAwG3nzzTYSGhoIQguTkZFy7dk2srvfeew/h4eFwdnaGgYFBr7IBynyU4jlm7969WL58OSoqKvDqq68iNTW112u0qoem0C5oNBr8/PzA5XKxb98+2a6hemgKbUKrJoUUiqHM7ubzJptSaAqtghpD/wMoKyvD4sWL0dLSgrlz5+K3336DhYUFcnNzhatA58+fx+XLl6Gnp4cTJ05AV1f3xZQt9zYQxQtH9x3GQYMGkbq6OpKenk6WLVtGUlNTSXR0NCGEkH/9618kIyPjhZVN9dD/AO7fvy+yw/j06VP069dPZEdu6NChAPi7pa2trS+sbGoM/Q9Alh1GRWICPo+yqR76H4C8O4wAkJKSgqqqKsyfP/+Fkk2tQ1NoFdSQg0KroBSaQqugFJpCq6AUmkKroBSaQqugFJpCq6AUmkKroBRaTtLS0mBra4vw8HBEREQgICBAoUAqRUVFwhALArd9Sfz6668yJSMFgLi4OIWCwwBAXV0dfv75ZwDA6tWrUVFRIfE8QVlOTo5E97W+hlJoBZg0aRLS0tKQmpqKM2fO4L333lOqvqNHj0ot++KLL9Dc3KxU/bJw+/Zt/PbbbwCAnTt3wtbWVuJ5grKffvoJ+fn5am+XvFAKrST19fUwMjICwE/zHBkZicTERNy+fVvYiy9evBgcDgdNTU2YOnUqwsLCsH79emEdAvuG9PR0BAYGwt/fH0uXLsXp06eRlZUlDN317rvvIjQ0FKGhocKe+NixY/Dy8sLEiRNx9+5dsfbdvn0bEyZMQGRkJAIDA3H//n3weDysWLECAQEB8PHxwZUrV7Bt2zacPXsWhw4dQnh4OAoKCuDp6Sk0IFq9ejVOnjyJ8PBwZGRkIDk5Gdu2bcOXX36J2NhYAPy8gaNHj0ZbW5u6bnfvKGcc+M8jNTWV2NjYkLCwMBIREUGmTJlCrl27RgghREdHhzQ0NBBCCAkICCAFBQWEEEK2bt1KvvnmG/LVV1+RLVu2EEII+eWXX4QBWJydnQkhhIwaNYo8ePCAEELItm3bSFFREQkLCyOFhYXk1KlTZOHChYQQQhoaGsjIkSMJIYQMHjyYNDQ0EA6HQ3x9fUlqaqpIe48cOULu3btHCCFk+/btZMuWLeTHH38ksbGxhBBC8vPzyY4dO0hqaqrwmEDmmjVryPHjxwmHwyEjRowgHA5HWJaQkECSkpIIl8slQ4YMIc3NzeSvv/4iK1euVO0NlxPKOEkBJk2ahOTkZLHjAwYMECaQvHfvHpYuXQoAaG1txYQJE1BTU4Pp06cDAEJCQvDFF1+IXF9fX4+BAwcCgNgwJicnB9euXROOu9va2lBZWQkLCwuhTEHAmq44OTkhISEBBgYGqKioQEBAAB48eICgoCAA/DAJgwcPljj2XrJkCT744APo6Ohg4sSJ0NERVxc6nY7o6GicOHECKSkp+Ne//iXlrmkGasihQrrGrfD09MTRo0eFPnPh4eFwd3cXTqQkxcSztLREcXExAOCdd97BtWvXQKfTwePx4OHhIRy7//HHH5gzZw769euHqqoq1NXVgcfjSYxY+tZbbyExMRFJSUlwc3MDIQSenp64fv06AKCkpATz5s0TyunK0KFDUVNTg/3792Px4sVi31Vw/uLFi3H48GE8fPgQvr6+it9AFUD10Griq6++QkxMDDgcDoyMjHDo0CGEhYVh4cKFGDt2LNzd3cWu+frrr4Uhx4YOHQp/f38EBQUhJiYGFy5cwF9//YWwsDCwWCwsX74cenp62LdvHyIjI2FlZQUmkylW57x58zBt2jRYW1vD2toaHA4HiYmJ+OOPPzBmzBhwuVzs2LEDjo6OuHHjBvbv3y9y/auvvooDBw4IjfAFeHt747333sPw4cPh7++P9vZ2zJgxQ4V3UDEo81EKlRAVFYXDhw/LHINOXVBDDgqlaGpqgpeXF8aPH9/nygxQPTSFlkH10BRaBaXQFFoFpdAUWsX/A3g3aTuD4Tr6AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 184.96x184.96 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "for batch, df in predictions.groupby(\"data_batch_name\"):\n",
    "    df = df.reset_index(level=\"data_batch_name\")\n",
    "    # Get the corresponding measurement info\n",
    "    batch_measured = activity_df.loc[df.index]\n",
    "    fig, ax = plot_report_accuracy(batch_measured, df[\"pred\"], batch)\n",
    "    plot_utils.save_fig(fig, os.path.join(figures_dir, f\"predictedObservedClasses_{batch}\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "78130c86-60d3-4535-a4a4-a0fe3bc7ff41",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.9.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}