ReinforcementLearningDemo.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "f5fb88d5-81de-42e7-898b-3738191ecbd9",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2207de6b-303c-451c-bbce-eedd57208e0e",
   "metadata": {},
   "source": [
    "# 1. Rescorla-Wagner model of classical conditioning \n",
    "The Rescorla-Wagner model of associative learning states that animals learn by storing associations between stimuli and rewards, and that these associations are updated using prediction errors. We can implement this model in a few simple lines of code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "id": "762ccebb-62b1-4916-85df-217473e48cd1",
   "metadata": {},
   "outputs": [],
   "source": [
    "def rescorla_wagner_update(w, stimuli, reward):\n",
    "    expected_reward = np.dot(w, stimuli)\n",
    "    prediction_error = reward - expected_reward\n",
    "    new_w = w + learning_rate * prediction_error * stimuli\n",
    "    return new_w"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2f19ac12-4de0-4f6c-9aa9-c2ba1513f7ba",
   "metadata": {},
   "source": [
    "Let's now simulate a simple conditioning + extinction paradigm whereby one stimulus is paired with reward for 100 trials, followed by 100 trials of extinction (the stimulus is shown but no reward is received). \n",
    "\n",
    "### Exercise 1.1\n",
    "Can you reason about why the learning and extinction curves look like a logarithmic function? "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 100,
   "id": "c3c93cec-8227-4ea8-863d-92b2c37b680f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x16878fd30>"
      ]
     },
     "execution_count": 100,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGwCAYAAABVdURTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA+50lEQVR4nO3de3wU9b3/8fcmkpCAWQgJScAgkXLRcpVLjJSLNRUVqRZbOWoNRcWqqGi0hfBTUNtDoh4FK6ititqjVvSU2laUFoOJUlKgXLxUpYLBWEwCEdmEEAjszu8PzDabbDabZHdnd/b1fDxizOxM9jtMdva93+93PmMzDMMQAACARcSY3QAAAIBAItwAAABLIdwAAABLIdwAAABLIdwAAABLIdwAAABLIdwAAABLOcXsBoSay+XSl19+qVNPPVU2m83s5gAAAD8YhqG6ujr169dPMTG++2aiLtx8+eWXyszMNLsZAACgE7744guddtppPteJunBz6qmnSjr5j5OUlGRyawAAgD9qa2uVmZnpfh/3JerCTdNQVFJSEuEGAIAI48+UEiYUAwAASyHcAAAASyHcAAAAS4m6OTcAAIQTp9Op48ePm92MsBAXF9fuZd7+INwAAGACwzBUVVWlQ4cOmd2UsBETE6OsrCzFxcV16fcQbgAAMEFTsOnbt68SExOjvrBsU5HdyspKDRgwoEv/HoQbAABCzOl0uoNNnz59zG5O2EhNTdWXX36pEydOqFu3bp3+PUwoBgAgxJrm2CQmJprckvDSNBzldDq79HsINwAAmCTah6JaCtS/B+EGAABYiqnh5p133tGMGTPUr18/2Ww2vfbaa+1uU1JSorPPPlvx8fH61re+peeeey7o7QQAAJHD1AnF9fX1GjVqlK699lrNnDmz3fXLy8s1ffp03XjjjXrxxRdVXFys66+/XhkZGZo2bVoIWgxEpkpHg8pr6pWV0kOSVF5Trx5xsapvdLq/+3oslOsEYvsMe0JI/30BhBdTw81FF12kiy66yO/1n3zySWVlZenhhx+WJJ155pnauHGjli1bRrhBRGsKH8EIDK+/X6lnNpbLZUhNo9mGlzb4eiyU63R1+xibVDhzhGaNH+DjNwAItqNHj+rGG2/Utm3b9PHHH+uSSy7xa4QmECLqUvCysjLl5uZ6LJs2bZpuv/32Nrc5duyYjh075v65trY2WM0D3PzpKWn6/sE+hx548xO5fLybBypUdPaxUK7T1e1dhrRozYeaPCSVHhxEjebnnHD5u3c6nUpISNBtt92m3//+9yF97ogKN1VVVUpLS/NYlpaWptraWjU0NCghofUBLSws1H333ReqJiIKtOxlaRlg/O0p6YhAhYpo4TQM7a05EjYneSCYVm+tUMGaD+Qygt9z+frrr+vHP/6xvvrqK8XGxmrnzp0aM2aMFixYoKKiIknS9ddfr6NHj+qFF17QE088IUn629/+FtJKzBEVbjqjoKBA+fn57p9ra2uVmZlpYosQCdrqeWkeXJr4CjAEDnPE2mwamEL9EFhfpaPBHWyk4PdcTpo0SXV1ddqxY4fGjRun0tJSpaSkqKSkxL1OaWmpFixYEPDn7oiICjfp6emqrq72WFZdXa2kpCSvvTaSFB8fr/j4+FA0DxHI21yXjva8RFKAsX3zH8NLo309Fsp1urp9rM2mpTOH02uDqFBeU99qSDuYPZd2u12jR49WSUmJxo0bp5KSEt1xxx267777dPjwYTkcDu3evVtTpkwJ+HN3RESFm5ycHL3xxhsey9avX6+cnByTWoRI1BRo/JnrEi7Bpatv+DGSrp+cpTkTsyRJe2uOKDEuRkcaXe7vTT0d3h4L5TqB2J5gg2iRldJDMTZ5nMeC3XM5ZcoUlZSU6M4779S7776rwsJCvfLKK9q4caMOHjyofv36afDgwUF7fn+YGm4OHz6s3bt3u38uLy/Xzp07lZycrAEDBqigoED79u3Tb3/7W0nSjTfeqBUrVujnP/+5rr32Wm3YsEGvvPKK1q5da9YuIAI0753xNqwUTP6EklibTT+/cKhGntYrqIGh+Ru+rzd/f4JBKNYJ1HMAVpZhT1DhzBFatOZDOQ0jJD2XU6dO1apVq/Tee++pW7duGjZsmKZOnaqSkhJ9/fXXpvfaSCaHm3/84x8677zz3D83zY2ZPXu2nnvuOVVWVqqiosL9eFZWltauXas77rhDjz76qE477TQ9/fTTXAYOt5ZzZVZtLA96mPEWYPzpKelMTwNv+ABamjV+gCYPSdXemiMh6blsmnezbNkyd5CZOnWqioqK9PXXX+vOO+8M6vP7w9RwM3XqVBk+PtJ6qz48depU7dixI4itQiSqdDR4BBmbAjOk5Cu4TB+R0WaPSUd6SgCgqzLsCSE7z/Tu3VsjR47Uiy++qBUrVkiSJk+erCuuuELHjx/36Ln56KOP1NjYqIMHD6qurk47d+6UJI0ePTqobYyoOTdAS02h5ul3yz3CTFeDja+el7Y+GRFgAESLKVOmaOfOnZo6daokKTk5WWeddZaqq6s1dOhQ93oXX3yxPv/8c/fPY8aMkSSfHRuBYDOC/Qxhpra2Vna7XQ6HQ0lJSWY3B53QfEJw0RufdDrIxNikBRcO8zrXhUmpAILp6NGjKi8vV1ZWlrp37252c8KGr3+Xjrx/03ODiNFy6KmzmvfKEGAAwHoINwh7bQ09+dJ8roy3OTKEGgCwLsINwlZnQ83cFnNlCDMAEF0INwgrnZ1P09ZQE6EGAKIP4QZho/nN3/zRfEIwvTMAgCaEG4SF9774WgvXfOCzkm+T5kNPBBoAQEuEG5iqI/NqCDUAAH8QbmCa1VsrtPD3H7Qbarh0GwDQEYQbmMI9DNXG48ynAQB0VozZDUD0qHQ0aNOeGv36nT26bOWmNufXxEj6w83n6qdTBilnUB+CDQBEoJKSEl166aXKyMhQjx49NHr0aL344osheW56bhAS/l4JFWOTCmeO0KjM3qFpGAAgKDZt2qSRI0dqwYIFSktL0+uvv668vDzZ7XZdcsklQX1uwg2Czp8roZgsDABd4NgnHdwjJQ+S7P2D9jSvv/66fvzjH+urr75SbGysdu7cqTFjxmjBggUqKiqSJF1//fU6evSoXnjhBY9t58+fr7/+9a9as2ZN0MMNw1IIqtVbK3wOQUkn/whfm3euFl18FsEGADpq+2+l5cOl52ec/L79t0F7qkmTJqmurk47duyQJJWWliolJUUlJSXudUpLS913C2/J4XAoOTk5aO1rQrhB0LQ3aViSYm02FV7OMBQAdIpjn/Tn+ZLhOvmz4ZL+fPvJ5UFgt9s1evRod5gpKSnRHXfcoR07dujw4cPat2+fdu/erSlTprTa9pVXXtHWrVs1Z86coLStOcINgsJXj02MTSq4aJh+N/ccbVx4nmaNHxD6BgKAFRzc859g08RwSgc/C9pTTpkyRSUlJTIMQ++++65mzpypM888Uxs3blRpaan69eunwYMHe2zz9ttva86cOXrqqaf07W9/O2hta8KcGwScrx6bpiuh6KlBsDTdnywrpQfDnLC+5EGSLcYz4NhipeQzgvaUU6dO1apVq/Tee++pW7duGjZsmKZOnaqSkhJ9/fXXrXptSktLNWPGDC1btkx5eXlBa1dz9NwgoNrrsWEICsG0emuFJhZt0FVPbdbEog1avbXC7CYBwWXvL8149GSgkU5+n7E8qJOKm+bdLFu2zB1kmsJNSUmJx3ybkpISTZ8+XQ888IBuuOGGoLWpJXpuEDCVjgYV0GMDkzT9/TWVG3AZ0qI1H2rykFR6cGBtZ+dJg84/ORSVfEZQg40k9e7dWyNHjtSLL76oFStWSJImT56sK664QsePH3cHnrfffluXXHKJ5s+fr8svv1xVVVWSpLi4uKBPKqbnBgFR6WjQ6+9/6bWODT02CIXymvpWf39Ow9DemiPmNAgIJXt/KWtS0INNkylTpsjpdLp7aZKTk3XWWWcpPT1dQ4cOlSQ9//zzOnLkiAoLC5WRkeH+mjlzZtDbZzMMf+7DbB21tbWy2+1yOBxKSkoyuzmW4KtAX4ykP8yjxwbBV+lo0MSiDR5/h7E2mzYuPI+eG4Sdo0ePqry8XFlZWerevbvZzQkbvv5dOvL+Tc8NuqRp8rC3YMNl3gilDHuCCmeOUKzNJunk39/SmcMJNkAUYs4NOs3XXb3vmX6mLh6ZwRsLQmrW+AGaPCRVe2uOcMNVIIoRbtApviYPx9psBBuYJsOewN8eEOUYlkKnbPv86zYnDzMUAAAwEz036LCm4aiWuNwbADomyq7paVeg/j3ouUGHtDUcxeXeAOC/bt26SZKOHKFUQXONjY2SpNjY2C79Hnpu4DdftWx+9V9jdMmofqFvFABEoNjYWPXq1Uv79++XJCUmJsr2zZV+0crlcunAgQNKTEzUKad0LZ4QbuAXX7VsYm02jR1Ijw0AdER6erokuQMOpJiYGA0YMKDLQY9wg3a1LGvfHLVEAKBzbDabMjIy1LdvXx0/ftzs5oSFuLg4xcR0fcYM4Qbt8lbWXqKWDQAEQmxsbJfnmMATE4rRrh5xsWrZQ0gtGwBAuCLcwKfVWyv0g8c3yWhxvx6GogAA4YphKbTJ21ybGElrbs7hkm8AQNii5wZt8laF2CXpSKPLlPYAAOAPwg28Wr21Qre+tKPV8libTQNTEk1oEQAA/iHcoBVfVYiZawMACHfMuUErbV36TRViAEAkoOcGrbR16TdViAEAkYBwAw9c+g0AiHQMS8GNS78BAFZAzw3cvM214dJvRLJKR4M27alRpaPB7KYACCF6buCWldJDMTZ5BBwu/Uakan4n+xibVDhzhGaNH2B2swCEAD03kHTyE255Tb0WXDRMsd/MJmauDSJVyyFWlyEtWvMhPThAlKDnBq0+4S64cJhGntZLA1MSCTaISN6GWJ2Gob01R/ibBqIAPTdRztsn3AfX7SLYIKI1DbE2xxArED0IN1HO1ydcIFJl2BNUOHMEQ6xAlGJYKso1FexrWdeGT7iIdLPGD9DkIanaW3OEnkggyhBuoljTXBsK9sGqMuwJ/C0DUYhwE6Uo2AcAsCrm3EQpCvYBAKyKcBOluJoEAGBVhJsodt13stwBh7k2AACrYM5NFGpetM8m6YbJWZozMYtgAwCwBHpuokzLicSGpGfe3WtmkwAACCjTw83KlSs1cOBAde/eXdnZ2dqyZYvP9ZcvX66hQ4cqISFBmZmZuuOOO3T06NEQtTbyUbQPAGB1poab1atXKz8/X0uWLNH27ds1atQoTZs2Tfv37/e6/ksvvaSFCxdqyZIl+vjjj/XMM89o9erVWrRoUYhbHrmYSAwAsDpTw80jjzyiuXPnas6cOTrrrLP05JNPKjExUatWrfK6/qZNmzRx4kRdddVVGjhwoC644AJdeeWV7fb24D8oSw8AsDrTJhQ3NjZq27ZtKigocC+LiYlRbm6uysrKvG5z7rnn6oUXXtCWLVs0YcIEffbZZ3rjjTd0zTXXtPk8x44d07Fjx9w/19bWBm4nIlClo0GZyYlac3OOjjS6KEsPALAc08JNTU2NnE6n0tLSPJanpaXpk08+8brNVVddpZqaGn3nO9+RYRg6ceKEbrzxRp/DUoWFhbrvvvsC2vZI1fwqqRibVDhzhHIG9TG7WQAABJTpE4o7oqSkREuXLtXjjz+u7du3a82aNVq7dq1+8YtftLlNQUGBHA6H++uLL74IYYvDR8urpFyGtGjNh6p0NJjbMAAAAsy0npuUlBTFxsaqurraY3l1dbXS09O9bnPPPffommuu0fXXXy9JGjFihOrr63XDDTfo//2//6eYmNZZLT4+XvHx8YHfgQjj6yophqUAAFZiWs9NXFycxo4dq+LiYvcyl8ul4uJi5eTkeN3myJEjrQJMbGysJMkwDG+b4BtcJQUAiBamDkvl5+frqaee0vPPP6+PP/5YN910k+rr6zVnzhxJUl5enseE4xkzZuiJJ57Qyy+/rPLycq1fv1733HOPZsyY4Q458I6rpBDNKh0N2rSnhmFYIEqYevuFWbNm6cCBA1q8eLGqqqo0evRorVu3zj3JuKKiwqOn5u6775bNZtPdd9+tffv2KTU1VTNmzNB///d/m7ULEaHS0aDymnpNHpKqjQvP096aI1wlhajhbSL9rPEDzG4WgCCyGVE2nlNbWyu73S6Hw6GkpCSzmxN0nNgRzSodDZpYtMFjvlmszaaNC88j3AMRpiPv3xF1tRQ6hiukEO243QgQnQg3FsaJHdGOifRAdCLcWBgndkQ7JtID0cnUCcUIrqYT+6I1H8ppGJzYEZVmjR+gyUNSmUgPRBHCjYVxHyngpAx7An/7QBQh3FgU95ECAEQr5txYEFdJAQCiGeHGgrhKCgAQzQg3FsRVUgCAaEa4sSAufwUARDMmFFsUl78CAKIV4caCmm6UmZXSgyukAABRh3BjMdwoEwAQ7ZhzYyFcAg4AAOHGUrgEHAAAwo2lcAk4AACEG0vhEnAAAJhQbDlcAg4AiHaEGwviDsiAd83LJPAaAayLcGMhnLiBtlEmAYgehBuL4MQNtK2tMgmTh6TyQQCwICYUWwD1bQDfKJMARBfCjQVw4gZ8o0wCEF0INxbAiRvwjTIJQHRhzo0FNJ24F635UE7D4MQNeEGZBCB6EG4sghM30D7KJADRgXBjIZy4AQBgzk3Eq3Q0aNOeGq6MAgDgG/TcRDBq2wAA0Bo9NxGK2jYAAHhHuIlQ1LYBAMA7wk2EorYNAADeEW4iFEXJAADwjgnFEYzaNgAAtEa4iXDUtgEAwBPDUgAAwFIINwAAwFIINxGKysRA5/H6AayNOTcRiMrEQOfx+gGsj56bCENlYqDzeP0A0YFwE2GoTAx0Hq8fIDoQbiIMlYmBzuP1A0QHwk2EoTIx0Hm8foDoYDMMw2h/Neuora2V3W6Xw+FQUlKS2c3ptEpHA5WJgU7i9QNEno68f3O1VISiMjHQebx+AGtjWAoAAFgK4SbCUHwMAADfGJaKIBQfAwCgffTcRAiKjwEA4B/CTYSg+BgAAP4h3EQIio8BAOAfwk2EoPgYAAD+YUJxBJk1foAmD0ml+BgAAD4QbiIMxccAAPCNYSkAAGAphBsAUY3CmID1MCwFIGpRGBOwJtN7blauXKmBAweqe/fuys7O1pYtW3yuf+jQIc2bN08ZGRmKj4/XkCFD9MYbb4SotebgkyUQeBTGBKzL1J6b1atXKz8/X08++aSys7O1fPlyTZs2Tbt27VLfvn1brd/Y2Kjvfe976tu3r/7v//5P/fv31+eff65evXqFvvEhwidLIDh8FcZk0j4Q2UwNN4888ojmzp2rOXPmSJKefPJJrV27VqtWrdLChQtbrb9q1SodPHhQmzZtUrdu3SRJAwcO9Pkcx44d07Fjx9w/19bWBm4HgqytT5aTh6Ry8gW6qKkwZvOAQ2FMwBpMG5ZqbGzUtm3blJub+5/GxMQoNzdXZWVlXrf505/+pJycHM2bN09paWkaPny4li5dKqfT2ebzFBYWym63u78yMzMDvi/Bwi0XgOChMCZgXab13NTU1MjpdCotLc1jeVpamj755BOv23z22WfasGGDrr76ar3xxhvavXu3br75Zh0/flxLlizxuk1BQYHy8/PdP9fW1kZMwOGTJRBcFMYErCmirpZyuVzq27evfvOb3yg2NlZjx47Vvn379NBDD7UZbuLj4xUfHx/ilgZG0yfLRWs+lNMw+GQJBAGFMQHrMS3cpKSkKDY2VtXV1R7Lq6urlZ6e7nWbjIwMdevWTbGxse5lZ555pqqqqtTY2Ki4uLigttkMfLIEAKBjTJtzExcXp7Fjx6q4uNi9zOVyqbi4WDk5OV63mThxonbv3i2Xy+Ve9q9//UsZGRmWDDZNMuwJyhnUh2ADAIAfTK1zk5+fr6eeekrPP/+8Pv74Y910002qr693Xz2Vl5engoIC9/o33XSTDh48qPnz5+tf//qX1q5dq6VLl2revHlm7QIAAAgzps65mTVrlg4cOKDFixerqqpKo0eP1rp169yTjCsqKhQT85/8lZmZqb/85S+64447NHLkSPXv31/z58/XggULzNoFAAAQZmyGYRjtr2YdtbW1stvtcjgcSkpKMrs5AADADx15/zb99gsAAACBRLgJY9xTCgCAjouoOjfRhHtKAQDQOfTchCHuVgyEHj2lgHXQcxOGuFsxEFr0lALWQs9NGGq6p1Rz3FMKCA56SgHrIdyEIe5WDISOr55SAJGJYakwxT2lgNBo6iltHnDoKQUiGz03YYx7SgHBR08pYD303ACIevSUAtZCuAEAnezBIdQA1sCwFAAAsBTCDQAAsBTCDQAAsBTCTRiiDDwAAJ3HhOIwQxl4AAC6hp6bMEIZeAAAuo5wE0YoAw8AQNcRbsIIN8wEzMecNyDyEW7CCGXgAXOt3lqhiUUbdNVTmzWxaINWb60wu0kAOsFmGIbR/mrWUVtbK7vdLofDoaSkJLOb41Wlo4Ey8ECIVToaNLFoQ6sbaG5ceB6vQyAMdOT9m6ulwhBl4IHQ8zXnjdcjEFk6NSyVl5enZ599Vnv27Al0ewDAFMx5A6yjU+EmLi5OhYWFGjx4sDIzM/XjH/9YTz/9tD799NNAtw8AQoI5b4B1dGnOzb59+/TOO++otLRUpaWl+te//qWMjAz9+9//DmQbAyoS5twAMA9z3oDwFLI5N71791afPn3Uu3dv9erVS6eccopSU1O78isBwFTMeQMiX6eGpRYtWqRzzz1Xffr00cKFC3X06FEtXLhQVVVV2rFjR6DbCAAA4LdODUvFxMQoNTVVd9xxh2bOnKkhQ4YEo21BEc7DUpWOBpXX1CsrpQefHAEAaCbow1I7duxQaWmpSkpK9PDDDysuLk5TpkzR1KlTNXXq1IgKO+GCG2YCABAYASni995772nZsmV68cUX5XK55HQ6A9G2oAjHnhuKhwEA4FvQe24Mw9COHTtUUlKikpISbdy4UbW1tRo5cqSmTJnSqUZHM4qHAQAQOJ0KN8nJyTp8+LBGjRqlKVOmaO7cuZo0aZJ69eoV4OZFh6biYS17bigeBgBAx3Uq3LzwwguaNGlS2AzrRLqm4mGL1nwop2FQPAwIE0zyByITN84MIxQPA8IHk/yB8NKR9+9O1blBcGTYE5QzqA/BBjBZpaPBHWykk0PGi9Z8qEpHg7kNA+AXwg0AtOBrkj+A8Ee4AYAWuEM4ENkINwDQAncIByJbl26cCQBWNWv8AE0eksokfyACEW4AoA3cIRyITAxLAQAASyHcAAAASyHcmKzS0aBNe2qonwEAQIAw58ZEVEAFACDw6LkxCRVQAQAIDsKNSaiACkQOho+ByMKwlEmaKqA2DzhUQAXCD8PHQOSh58YkVEAFwh/Dx0BkoufGRFRABcKbr+FjXq9A+CLcmIwKqED4YvgYiEwMSwFAGxg+BiITPTcA4APDx0DkIdwAQDsYPgYiC8NSAADAUgg3AADAUgg3AADAUsIi3KxcuVIDBw5U9+7dlZ2drS1btvi13csvvyybzabLLrssuA0EAAARw/Rws3r1auXn52vJkiXavn27Ro0apWnTpmn//v0+t9u7d6/uuusuTZo0KUQtBRDtuMcUEBlshmEY7a8WPNnZ2Ro/frxWrFghSXK5XMrMzNStt96qhQsXet3G6XRq8uTJuvbaa/Xuu+/q0KFDeu2117yue+zYMR07dsz9c21trTIzM+VwOJSUlBTw/QFgTdxjCjBXbW2t7Ha7X+/fpvbcNDY2atu2bcrNzXUvi4mJUW5ursrKytrc7v7771ffvn113XXXtfschYWFstvt7q/MzMyAtL2r+AQIRA7uMQVEFlPr3NTU1MjpdCotLc1jeVpamj755BOv22zcuFHPPPOMdu7c6ddzFBQUKD8/3/1zU8+NmfgECEQW7jEFRJaIKuJXV1ena665Rk899ZRSUlL82iY+Pl7x8fFBbpn/2voEOHlIKidJIExxjykgspgablJSUhQbG6vq6mqP5dXV1UpPT2+1/p49e7R3717NmDHDvczlckmSTjnlFO3atUuDBg0KbqO7iE+AQORpusfUojUfymkY3GMKCHOmhpu4uDiNHTtWxcXF7su5XS6XiouLdcstt7Raf9iwYfrggw88lt19992qq6vTo48+avpwkz/4BAhEJu4xBUQO04el8vPzNXv2bI0bN04TJkzQ8uXLVV9frzlz5kiS8vLy1L9/fxUWFqp79+4aPny4x/a9evWSpFbLwxWfAIHIxT2mgMhgeriZNWuWDhw4oMWLF6uqqkqjR4/WunXr3JOMKyoqFBNjejmegOITIAAAwWN6nZtQ68h18gAAIDxETJ0bAIhE1KkCwpvpw1IAEEmoUwWEP3puAMBPVCoGIgPhBgD85KtOFYDwQbgBAD811alqjjpVQPgh3ACAn5rqVMXaTiYc6lQB4YkJxQDQAdSpAsIf4QYAOohKxUB4Y1gKAABYCuEmxCj+BQBAcDEsFUIU/wKsp9LRoPKaemWl9GCoCggThJsQaav41+QhqZwQgQjFBxYgPDEsFSIU/wKshWrFQPgi3IQIxb8Aa+EDCxC+CDchQvEvwFr4wAKEL+bchBDFvwDraPrAsmjNh3IaBh9YgDBCuAkxin8B1sEHFiA8EW4AoAv4wAKEH+bcAAAASyHcAEAAUH0cCB8MSwFAF1HMDwgv9NwAQBdQzA8IP4QbAOgCivkB4YdwAwBdQDE/IPwQbgCgC6g+DoQfJhQDQBdRzA8IL4QbAAgAivkB4YNhKQAAYCmEGwAIIIr5AeZjWCpEKh0NKq+pV1ZKD7quAYuimB8QHgg3IcAJD7C+tor5TR6SygcaIMQYlgoyqpcC0YFifkD4INwEGSc8IDpQzA8IH4SbIOOEB0QHivkB4YM5N0HWdMJbtOZDOQ2DEx5gYRTzA8KDzTAMo/3VrKO2tlZ2u10Oh0NJSUkhe95KRwMnPCCKcIUkEFgdef+m5yZEqF4KRA+ukATMxZwbAAggrpAEzEe4AYAA4gpJwHyEGwAIIK6QBMxHuAGAAOKScMB8TCgGgADjknDAXIQbAAgCrpAEzMOwFAAEWaWjQZv21HDFFBAi9NwAQBBR8wYIPXpuACBIqHkDmINwAwBBQs0bwByEGwAIEmreAOYg3ABAkFDzBjAHE4oBIIia17xJjItRfaNTlY4GAg4QRISbIKp0NKi8pl5ZKT04kQFRLMOeoHf+dYCrpoAQIdwECZd/AmjS1lVTk4ek8sEHCALm3AQBl38CaI6rpoDQItwEAScyAM1x1RQQWoSbIOBEBqA5rpoCQos5N0HQdCJbtOZDOQ2DExkArpoCQigsem5WrlypgQMHqnv37srOztaWLVvaXPepp57SpEmT1Lt3b/Xu3Vu5ubk+1zfLrPEDtHHhefrd3HO0ceF5TCYGoAx7gioO1usHj2/SVU9t1sSiDVq9tcLsZgGWY3q4Wb16tfLz87VkyRJt375do0aN0rRp07R//36v65eUlOjKK6/U22+/rbKyMmVmZuqCCy7Qvn37Qtzy9mXYE5QzqA+fzABI4mIDIFRshmEY7a8WPNnZ2Ro/frxWrFghSXK5XMrMzNStt96qhQsXtru90+lU7969tWLFCuXl5bW7fm1trex2uxwOh5KSkrrcfgDw16Y9Nbrqqc2tlv9u7jnKGdTHhBYBkaMj79+m9tw0NjZq27Ztys3NdS+LiYlRbm6uysrK/PodR44c0fHjx5WcnOz18WPHjqm2ttbjCwDMwMUGQGiYGm5qamrkdDqVlpbmsTwtLU1VVVV+/Y4FCxaoX79+HgGpucLCQtntdvdXZmZml9sNAJ3BVVNAaET01VJFRUV6+eWXVVJSou7du3tdp6CgQPn5+e6fa2trCTgATMNVU0DwmRpuUlJSFBsbq+rqao/l1dXVSk9P97nt//zP/6ioqEhvvfWWRo4c2eZ68fHxio+PD0h7ASAQuNcUEFymDkvFxcVp7NixKi4udi9zuVwqLi5WTk5Om9s9+OCD+sUvfqF169Zp3LhxoWgqAAQMV00BwWX6sFR+fr5mz56tcePGacKECVq+fLnq6+s1Z84cSVJeXp769++vwsJCSdIDDzygxYsX66WXXtLAgQPdc3N69uypnj17mrYfAOAvX7doYXgK6DrTw82sWbN04MABLV68WFVVVRo9erTWrVvnnmRcUVGhmJj/dDA98cQTamxs1A9/+EOP37NkyRLde++9oWw6AHRK01VTzQMOV00BgWN6nZtQo84NgHCwemuFxy1afn7hUI04za6slB703gBedOT92/SeGwCIRs2vmnp/3yE98OYnTC4GAsT02y8AQLTKsCdoYEqiO9hITC4GAoFwAwAm8jW5GEDnEG4AwETebskQIykxjtMz0Fm8egDARC1vySBJLkk/eHyTVm+tMK9hQAQj3ARBpaNBm/bUMGYOwC+zxg/Qmptz1CzfMPcG6AKulgqw1VsrKKkOoMPqG51qWZiDwn5A59BzE0CUVAfQWd7m3lDYD+gcwk0AcdUDgM5qOfcmRtJ13xloapuASEW4CSA+eQHoilnjB2jjwvN0w6QzJJv0m3fLNbFoAxOLgQ4i3ARQy09esTabls4czng5gA55euNnDG8DXcCE4gBrXlJ9YEoiwQZAh3DHcKDrCDdBkGFP4CQEoFO83TE8RtJX9cdU6Wjg3AL4gWEpAAgjLYe3bZIMSbe8tIP5N4CfCDcAEGaaJhavuHKMbLaT4UZi/g3gL8INAIShDHuCknvGeZ1/s23v1+Y0CogQhBsACFPeyktI0m0v72B4CvCBcAMAYapp/k3LEzXDU4BvhBsACGOzxg/Qr64a02o5w1NA2wg3ABDmxp7em+EpoAMINwAQ5hieAjqGcAMAEcDX8BQ35wU8EW4AIEJ4G56KkZQYx6kcaI5XBABEiJbViyXJJekHj29i7g3QDOEGACLIrPEDtObmHDXLN3IZUsHvP9Dr73/J/BtAhBsAiDj1jU4ZLSoXu8T9p4AmhBsAiDBtVS6WuIIKkAg3ABBxvM29aY4Cf4h2p5jdAABAx80aP0CTh6Rq296vddvLO1rdYPO2l3eovvGEZo0fYE4DARPRcwMAESrDnqBLRvVrs8Bfwe8/0Htf0IOD6EO4CaBKR4M27alhrBtASLVV4M8l6TIuE0cUYlgqQFZvrVDBmg/kMqQYm1Q4cwTdwQBCpqnAX8vhKeObCcaTh6Qqw55gTuOAEKPnJgAqHQ3uYCNxtQKA0Gvr/lMSE4wRfQg3AVBeU9/q0xL3ewEQarPGD9Af5p0rbxdRcQdxRBPCTQB4qzkRa7NpYEqiOQ0CELVGZfZWkY8JxlQxRjQg3ARAy5oTsTabls4czvg2AFP4mmB8y0s7dG7hBi1d+xEhB5ZlM4yWRbytrba2Vna7XQ6HQ0lJSQH93ZWOBu2tOaKBKYkEGwCmqnQ0aGLRhlZD5s1x8QMiSUfev+m5CaAMe4JyBvUh2AAwXXtVjCVq4cC6CDcAYFGzxg/QxoXnacWVY9q+F5WohQPrIdwAgIX5qmLcxKAHBxbDnBsAiBKVjgY9u3Gvnnr3M3k78dts0sILh2nEaXZlpfRgiB1hpSPv34QbAIgy733xtS57fJN8nf2ZbIxww4RiAECb2qqF0xyTjRHJCDcAEIV8VTNu4pJ06cpN1MRBxCHcAECUaurB8XW5uCT95t1yTSzaoF+X7tGmPTUEHYQ95twAQJRrKkD6/r5DeuCNT+RqZ32bpLmTsjTnO1lMOkbIMKHYB8INALTNn8nGTQg5CCUmFAMAOsWfycZNDJ0csuJeVQg39NwAAFppqonz9MbPfN6fqrkYm7Tgmzo5PeJiVd/opF4OAoZhKR8INwDgv47Ox2mJoSsECuHGB8INAHROexWOfWkeciSpvKaeXh10COHGB8INAHRNV0OOdHK+TlPgmT4ygyEstItw4wPhBgACozPzcnyhdwe+EG58INwAQGA1n5fz4Ju75Ozi20p7vTvSyeDDpOXoQrjxgXADAMHTFHQS42K09v2qTg1d+dI8+DRfRgCyPsKND4QbAAidlkNXNimgYcebzgaglt8JROEl4sLNypUr9dBDD6mqqkqjRo3SY489pgkTJrS5/quvvqp77rlHe/fu1eDBg/XAAw/o4osv9uu5CDcAEHpNPToDUxIlKaBzdTrDWwBqyVvdHm8BSPIekvwNUP6sQ8iKsHCzevVq5eXl6cknn1R2draWL1+uV199Vbt27VLfvn1brb9p0yZNnjxZhYWFuuSSS/TSSy/pgQce0Pbt2zV8+PB2ny+o4caxTzq4R+rWQzpe/5/vyYNOPu7tsUCtE4rnoB2R245IaivtiJp2HGg8RftrDqpbQk9t3fVvvfL+V0q0HdNeV7oMSVkxVTrsilfPmGPu7+WudKmNx8xYZ+836wz0so6vxzqyzudGum6e+i2N7vmV4hJOVWNDXavvqaefJUk68PlHHXosmOuknn6W0k775ngHQESFm+zsbI0fP14rVqyQJLlcLmVmZurWW2/VwoULW60/a9Ys1dfX6/XXX3cvO+ecczR69Gg9+eST7T5f0MLN9t9Kf54vGd5KXPnzGaGr64TiOWhH5LbDn3VoB+0wtx1NE4iNb9axyZBhSDab3N9P9vTYFGNr/Vi4rRPt7XAaNm0bea8mXH67j78H/0XMvaUaGxu1bds25ebmupfFxMQoNzdXZWVlXrcpKyvzWF+Spk2b1ub6x44dU21trcdXwDn2+Qg20smXbHsZsqvrhOI5aEfktsOfdWgH7TC3HTb3d0O2b9axfbOw6XuMTYqxeX8s3NaJ9nbE2gyd/f59qv73HoWaqeGmpqZGTqdTaWlpHsvT0tJUVVXldZuqqqoOrV9YWCi73e7+yszMDEzjmzu4x0ewAQAgOp1ic6nm809C/ryWvyt4QUGBHA6H++uLL74I/JMkD5Jslv+nBACgQ04YMUo5fVjIn9fUd+SUlBTFxsaqurraY3l1dbXS09O9bpOent6h9ePj45WUlOTxFXD2/tKMRyVbbBsr2P7TZ9emrq4TiuegHZHbDn/WoR20g3a0XKdp/k/TQFrr777W6er2J7lkk8v4Zp1vFrb87jTaXsfXY8Fc54QRo+0jlwR0UrG/Tgn5MzYTFxensWPHqri4WJdddpmkkxOKi4uLdcstt3jdJicnR8XFxbr99tvdy9avX6+cnJwQtNiHs/OkQedLBz+TuiVKx4/853vyGSfX8fZYoNYJxXPQjshtRyS1lXbQjjBqh+2bdWzfLGv1vdnv6dBjHVgnJvkMVdcdVc3nn6hbQg8db6hv9b2pd8TbOr4eC+Y6KacP0wQTgo0k86+WWr16tWbPnq1f//rXmjBhgpYvX65XXnlFn3zyidLS0pSXl6f+/fursLBQ0slLwadMmaKioiJNnz5dL7/8spYuXRoel4IDAICg6Mj7t6k9N9LJS7sPHDigxYsXq6qqSqNHj9a6devck4YrKioUE/Of0bNzzz1XL730ku6++24tWrRIgwcP1muvveZXsAEAANZnes9NqNFzAwBA5ImYOjcAAACBRrgBAACWQrgBAACWQrgBAACWQrgBAACWQrgBAACWQrgBAACWQrgBAACWQrgBAACWYvrtF0KtqSBzbW2tyS0BAAD+anrf9ufGClEXburq6iRJmZmZJrcEAAB0VF1dnex2u891ou7eUi6XS19++aVOPfVU2Wy2gP7u2tpaZWZm6osvvrDkfausvn8S+2gFVt8/iX20AqvvnxT4fTQMQ3V1derXr5/HDbW9ibqem5iYGJ122mlBfY6kpCTL/rFK1t8/iX20Aqvvn8Q+WoHV908K7D6212PThAnFAADAUgg3AADAUgg3ARQfH68lS5YoPj7e7KYEhdX3T2IfrcDq+yexj1Zg9f2TzN3HqJtQDAAArI2eGwAAYCmEGwAAYCmEGwAAYCmEGwAAYCmEmwBZuXKlBg4cqO7duys7O1tbtmwxu0mdVlhYqPHjx+vUU09V3759ddlll2nXrl0e60ydOlU2m83j68YbbzSpxR1z7733tmr7sGHD3I8fPXpU8+bNU58+fdSzZ09dfvnlqq6uNrHFHTdw4MBW+2iz2TRv3jxJkXn83nnnHc2YMUP9+vWTzWbTa6+95vG4YRhavHixMjIylJCQoNzcXH366ace6xw8eFBXX321kpKS1KtXL1133XU6fPhwCPeibb727/jx41qwYIFGjBihHj16qF+/fsrLy9OXX37p8Tu8HfeioqIQ70nb2juGP/nJT1q1/8ILL/RYJ5yPodT+Pnp7XdpsNj300EPudcL5OPrz/uDPObSiokLTp09XYmKi+vbtq5/97Gc6ceJEwNpJuAmA1atXKz8/X0uWLNH27ds1atQoTZs2Tfv37ze7aZ1SWlqqefPm6e9//7vWr1+v48eP64ILLlB9fb3HenPnzlVlZaX768EHHzSpxR337W9/26PtGzdudD92xx136M9//rNeffVVlZaW6ssvv9TMmTNNbG3Hbd261WP/1q9fL0n60Y9+5F4n0o5ffX29Ro0apZUrV3p9/MEHH9SvfvUrPfnkk9q8ebN69OihadOm6ejRo+51rr76av3zn//U+vXr9frrr+udd97RDTfcEKpd8MnX/h05ckTbt2/XPffco+3bt2vNmjXatWuXvv/977da9/777/c4rrfeemsomu+X9o6hJF144YUe7f/d737n8Xg4H0Op/X1svm+VlZVatWqVbDabLr/8co/1wvU4+vP+0N451Ol0avr06WpsbNSmTZv0/PPP67nnntPixYsD11ADXTZhwgRj3rx57p+dTqfRr18/o7Cw0MRWBc7+/fsNSUZpaal72ZQpU4z58+eb16guWLJkiTFq1Civjx06dMjo1q2b8eqrr7qXffzxx4Yko6ysLEQtDLz58+cbgwYNMlwul2EYkX38DMMwJBl/+MMf3D+7XC4jPT3deOihh9zLDh06ZMTHxxu/+93vDMMwjI8++siQZGzdutW9zptvvmnYbDZj3759IWu7P1runzdbtmwxJBmff/65e9npp59uLFu2LLiNCxBv+zh79mzj0ksvbXObSDqGhuHfcbz00kuN7373ux7LIuk4tnx/8Occ+sYbbxgxMTFGVVWVe50nnnjCSEpKMo4dOxaQdtFz00WNjY3atm2bcnNz3ctiYmKUm5ursrIyE1sWOA6HQ5KUnJzssfzFF19USkqKhg8froKCAh05csSM5nXKp59+qn79+umMM87Q1VdfrYqKCknStm3bdPz4cY/jOWzYMA0YMCBij2djY6NeeOEFXXvttR43i43k49dSeXm5qqqqPI6b3W5Xdna2+7iVlZWpV69eGjdunHud3NxcxcTEaPPmzSFvc1c5HA7ZbDb16tXLY3lRUZH69OmjMWPG6KGHHgpoV38olJSUqG/fvho6dKhuuukmffXVV+7HrHYMq6urtXbtWl133XWtHouU49jy/cGfc2hZWZlGjBihtLQ09zrTpk1TbW2t/vnPfwakXVF348xAq6mpkdPp9DhIkpSWlqZPPvnEpFYFjsvl0u23366JEydq+PDh7uVXXXWVTj/9dPXr10/vv/++FixYoF27dmnNmjUmttY/2dnZeu655zR06FBVVlbqvvvu06RJk/Thhx+qqqpKcXFxrd4w0tLSVFVVZU6Du+i1117ToUOH9JOf/MS9LJKPnzdNx8bb67DpsaqqKvXt29fj8VNOOUXJyckRd2yPHj2qBQsW6Morr/S4IeFtt92ms88+W8nJydq0aZMKCgpUWVmpRx55xMTW+u/CCy/UzJkzlZWVpT179mjRokW66KKLVFZWptjYWEsdQ0l6/vnndeqpp7Ya9o6U4+jt/cGfc2hVVZXX12rTY4FAuIFP8+bN04cffugxJ0WSxxj3iBEjlJGRofPPP1979uzRoEGDQt3MDrnooovc/z9y5EhlZ2fr9NNP1yuvvKKEhAQTWxYczzzzjC666CL169fPvSySj1+0O378uK644goZhqEnnnjC47H8/Hz3/48cOVJxcXH66U9/qsLCwogo8/9f//Vf7v8fMWKERo4cqUGDBqmkpETnn3++iS0LjlWrVunqq69W9+7dPZZHynFs6/0hHDAs1UUpKSmKjY1tNRO8urpa6enpJrUqMG655Ra9/vrrevvtt3Xaaaf5XDc7O1uStHv37lA0LaB69eqlIUOGaPfu3UpPT1djY6MOHTrksU6kHs/PP/9cb731lq6//nqf60Xy8ZPkPja+Xofp6emtJvmfOHFCBw8ejJhj2xRsPv/8c61fv96j18ab7OxsnThxQnv37g1NAwPsjDPOUEpKivvv0grHsMm7776rXbt2tfvalMLzOLb1/uDPOTQ9Pd3ra7XpsUAg3HRRXFycxo4dq+LiYvcyl8ul4uJi5eTkmNiyzjMMQ7fccov+8Ic/aMOGDcrKymp3m507d0qSMjIygty6wDt8+LD27NmjjIwMjR07Vt26dfM4nrt27VJFRUVEHs9nn31Wffv21fTp032uF8nHT5KysrKUnp7ucdxqa2u1efNm93HLycnRoUOHtG3bNvc6GzZskMvlcoe7cNYUbD799FO99dZb6tOnT7vb7Ny5UzExMa2GciLFv//9b3311Vfuv8tIP4bNPfPMMxo7dqxGjRrV7rrhdBzbe3/w5xyak5OjDz74wCOoNoX1s846K2ANRRe9/PLLRnx8vPHcc88ZH330kXHDDTcYvXr18pgJHkluuukmw263GyUlJUZlZaX768iRI4ZhGMbu3buN+++/3/jHP/5hlJeXG3/84x+NM844w5g8ebLJLffPnXfeaZSUlBjl5eXG3/72NyM3N9dISUkx9u/fbxiGYdx4443GgAEDjA0bNhj/+Mc/jJycHCMnJ8fkVnec0+k0BgwYYCxYsMBjeaQev7q6OmPHjh3Gjh07DEnGI488YuzYscN9tVBRUZHRq1cv449//KPx/vvvG5deeqmRlZVlNDQ0uH/HhRdeaIwZM8bYvHmzsXHjRmPw4MHGlVdeadYuefC1f42Njcb3v/9947TTTjN27tzp8bpsurpk06ZNxrJly4ydO3cae/bsMV544QUjNTXVyMvLM3nP/sPXPtbV1Rl33XWXUVZWZpSXlxtvvfWWcfbZZxuDBw82jh496v4d4XwMDaP9v1PDMAyHw2EkJiYaTzzxRKvtw/04tvf+YBjtn0NPnDhhDB8+3LjggguMnTt3GuvWrTNSU1ONgoKCgLWTcBMgjz32mDFgwAAjLi7OmDBhgvH3v//d7CZ1miSvX88++6xhGIZRUVFhTJ482UhOTjbi4+ONb33rW8bPfvYzw+FwmNtwP82aNcvIyMgw4uLijP79+xuzZs0ydu/e7X68oaHBuPnmm43evXsbiYmJxg9+8AOjsrLSxBZ3zl/+8hdDkrFr1y6P5ZF6/N5++22vf5ezZ882DOPk5eD33HOPkZaWZsTHxxvnn39+q33/6quvjCuvvNLo2bOnkZSUZMyZM8eoq6szYW9a87V/5eXlbb4u3377bcMwDGPbtm1Gdna2Ybfbje7duxtnnnmmsXTpUo9gYDZf+3jkyBHjggsuMFJTU41u3boZp59+ujF37txWHxLD+RgaRvt/p4ZhGL/+9a+NhIQE49ChQ622D/fj2N77g2H4dw7du3evcdFFFxkJCQlGSkqKceeddxrHjx8PWDtt3zQWAADAEphzAwAALIVwAwAALIVwAwAALIVwAwAALIVwAwAALIVwAwAALIVwAwAALIVwAwAALIVwAyDg7r33Xo0ePbpD29hsNr322mtBaU9XTZ06VbfffrvZzQDgJ8INAL905A3+rrvu8rhxHgCE0ilmNwCAdRiGIafTqZ49e6pnz55mNyesOZ1O2Ww2xcTwGRMINF5VANr1k5/8RKWlpXr00Udls9lks9m0d+9elZSUyGaz6c0339TYsWMVHx+vjRs3thqW2rp1q773ve8pJSVFdrtdU6ZM0fbt2zvUhqlTp+q2227Tz3/+cyUnJys9PV333nuv+/G9e/fKZrNp586d7mWHDh2SzWZTSUmJJLnb+5e//EVjxoxRQkKCvvvd72r//v168803deaZZyopKUlXXXWVjhw54vH8J06c0C233CK73a6UlBTdc889an5rvmPHjumuu+5S//791aNHD2VnZ7ufV5Kee+459erVS3/605901llnKT4+XhUVFR36NwDgH8INgHY9+uijysnJ0dy5c1VZWanKykplZma6H1+4cKGKior08ccfa+TIka22r6ur0+zZs7Vx40b9/e9/1+DBg3XxxRerrq6uQ+14/vnn1aNHD23evFkPPvig7r//fq1fv77D+3PvvfdqxYoV2rRpk7744gtdccUVWr58uV566SWtXbtWf/3rX/XYY4+1eu5TTjlFW7Zs0aOPPqpHHnlETz/9tPvxW265RWVlZXr55Zf1/vvv60c/+pEuvPBCffrpp+51jhw5ogceeEBPP/20/vnPf6pv374dbjuA9jEsBaBddrtdcXFxSkxMVHp6eqvH77//fn3ve99rc/vvfve7Hj//5je/Ua9evVRaWqpLLrnE73aMHDlSS5YskSQNHjxYK1asUHFxsc/n9uaXv/ylJk6cKEm67rrrVFBQoD179uiMM86QJP3whz/U22+/rQULFri3yczM1LJly2Sz2TR06FB98MEHWrZsmebOnauKigo9++yzqqioUL9+/SSdnHe0bt06Pfvss1q6dKkk6fjx43r88cc1atSoDrUXQMfQcwOgy8aNG+fz8erqas2dO1eDBw+W3W5XUlKSDh8+3OFhmZa9QhkZGdq/f3+H29v896SlpSkxMdEdbJqWtfy955xzjmw2m/vnnJwcffrpp3I6nfrggw/kdDo1ZMgQ93yjnj17qrS0VHv27HFvExcX57VnC0Bg0XMDoMt69Ojh8/HZs2frq6++0qOPPqrTTz9d8fHxysnJUWNjY4eep1u3bh4/22w2uVwuSXJPzG0+D+b48ePt/h6bzebz9/rj8OHDio2N1bZt2xQbG+vxWPOJ1QkJCR4BCUBwEG4A+CUuLk5Op7NT2/7tb3/T448/rosvvliS9MUXX6impiaQzVNqaqokqbKyUmPGjJEkj8nFXbV582aPn5vmDsXGxmrMmDFyOp3av3+/Jk2aFLDnBNA5hBsAfhk4cKA2b96svXv3qmfPnkpOTvZ728GDB+t///d/NW7cONXW1upnP/uZEhISAtq+hIQEnXPOOSoqKlJWVpb279+vu+++O2C/v6KiQvn5+frpT3+q7du367HHHtPDDz8sSRoyZIiuvvpq5eXl6eGHH9aYMWN04MABFRcXa+TIkZo+fXrA2gGgfcy5AeCXu+66S7GxsTrrrLOUmpraofkyzzzzjL7++mudffbZuuaaa3TbbbcF5UqhVatW6cSJExo7dqxuv/12/fKXvwzY787Ly1NDQ4MmTJigefPmaf78+brhhhvcjz/77LPKy8vTnXfeqaFDh+qyyy7T1q1bNWDAgIC1AYB/bEbzAWoAAIAIR88NAACwFMINAACwFMINAACwFMINAACwFMINAACwFMINAACwFMINAACwFMINAACwFMINAACwFMINAACwFMINAACwlP8PmpCwRwKM96UAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# initialize some parameters \n",
    "learning_rate = .08\n",
    "n_stimuli = 2\n",
    "\n",
    "# intialize the weights to zeros (no prior reward expectation)\n",
    "w = np.zeros(n_stimuli)\n",
    "\n",
    "ws = []\n",
    "for trial in range(200):\n",
    "    stimuli = np.array([1, 0])\n",
    "    if trial < 100:\n",
    "        reward = 1\n",
    "    else:\n",
    "        reward = 0\n",
    "\n",
    "    w = rescorla_wagner_update(w, stimuli, reward)\n",
    "\n",
    "    ws.append(w)\n",
    "\n",
    "plt.plot(ws, '.')\n",
    "plt.xlabel('trial number')\n",
    "plt.ylabel('w')\n",
    "plt.legend(['w1', 'w2'])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f30037e2-beac-4e6f-a8aa-2d4f5e3186e3",
   "metadata": {},
   "source": [
    "## Forward blocking, overshadowing and overexpectation\n",
    "\n",
    "In forward blocking, a first cue A is first paired with reward until asymptote, after which A is paired with reward in compound with B. Subsequently we measure conditioned responding to cue B \\\n",
    "A —> +  \\\n",
    "AB —> + \\\n",
    "B —> ? \n",
    "\n",
    "### Exercise 1.2\n",
    "What do you expect the conditioned responding to B to be like, given the Rescorla-Wagner learning rule? \n",
    "\n",
    "Let's simulate this below: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "id": "f5bef417-a147-48d4-a474-546305096cd4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# intialize the weights to zeros (no prior reward expectation)\n",
    "w = np.zeros(n_stimuli)\n",
    "\n",
    "ws = []\n",
    "for trial in range(200):\n",
    "    reward = 1\n",
    "    if trial < 100:\n",
    "        stimuli = np.array([1, 0])\n",
    "    else:\n",
    "        stimuli = np.array([1, 1])\n",
    "\n",
    "    w = rescorla_wagner_update(w, stimuli, reward)\n",
    "\n",
    "    ws.append(w)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 102,
   "id": "9fad2c20-8c66-4624-ae45-7b0e2df9895e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "conditioned response to B: 0.00011960593412453592\n"
     ]
    }
   ],
   "source": [
    "print(f'conditioned response to B: {np.dot(np.array([0,1]), w)}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 103,
   "id": "a3a024df-d1d0-47f6-9f1d-b803e5d2bb23",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x16887b370>"
      ]
     },
     "execution_count": 103,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGwCAYAAABVdURTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA5C0lEQVR4nO3df3wU9Z3H8ffsSkICZkMIJAGDiZQfUgggP2L0+GFNRUSqF1tzag2litWiImgL4RTU9gjqKXgFpVVRe2rF9tC2qLQYSZSSAgWi2CoVGgzFJBCRDSSBwO7cHzQrS0J+kN2d3cnr+XhsaWa+k/lMJtl9+53vfMcwTdMUAACATTisLgAAACCQCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWzrG6gFDzer36/PPPde6558owDKvLAQAAbWCapg4fPqw+ffrI4Wi5b6bThZvPP/9cqampVpcBAADOwt69e3Xeeee12KbThZtzzz1X0skfTlxcnMXVAACAtqipqVFqaqrvc7wlnS7cNF6KiouLI9wAABBh2jKkhAHFAADAVgg3AADAVgg3AADAVjrdmJu28ng8On78uNVlhIWoqKhWb7sDACBcEG5OY5qmKisrdejQIatLCRsOh0Pp6emKioqyuhQAAFpFuDlNY7Dp3bu3YmNjO/1Ef42THlZUVKhfv36d/ucBAAh/hJtTeDweX7Dp2bOn1eWEjV69eunzzz/XiRMn1KVLF6vLAQCgRQykOEXjGJvY2FiLKwkvjZejPB6PxZUAANA6wk0zuPTij58HACCSEG4AAICtWBpu3nvvPU2dOlV9+vSRYRh64403Wt2mqKhIF110kaKjo/W1r31NL7zwQtDrBAAAkcPSAcW1tbUaPny4vv/97ysnJ6fV9mVlZZoyZYpuv/12vfzyyyosLNStt96qlJQUTZo0KQQVA5Gvwl2vsupadYtyqrbB4/s3PbGbJDW7LpRtqIM6IqGOSKrVqjrSE7spxRXTyjtScFgabiZPnqzJkye3uf2KFSuUnp6uxx9/XJJ04YUXasOGDVqyZAnhpo2OHj2q22+/XVu3btXHH3+sq6++uk09ZgiuMwWOQL/x7Njn1iNvfyKv2bSGxpFVzawKaRvqoI5IqKMtbTp7HQ5DKsgZptwx/VrYIjgi6lbwkpISZWdn+y2bNGmS7rnnnjNuc+zYMR07dsz3dU1NTbDKa6LxA8vK9Ho6j8ejmJgY3X333fq///s/q8uxrVPPvXT2gaNRoN54WtKW7ULRhjqoIxLqaEubzl6H15Tmr/5I4wf2CvlnYESFm8rKSiUlJfktS0pKUk1Njerr6xUT0/SHV1BQoIceeihUJfqs2lKu/NU75DWDn17XrFmj7373u/riiy/kdDpVWlqqkSNHau7cuVq8eLEk6dZbb9XRo0f10ksv6emnn5Yk/elPf2Im5rNwei/L6QFmzYcVem5DmbxmxwNHo0C98QBAKHlMU3uq6wg3gZafn685c+b4vq6pqVFqampQ91nhrvcFGyn46XXcuHE6fPiwtm/frtGjR6u4uFiJiYkqKirytSkuLtbcuXMDvm+7OlPPy6nBpVFLAYbAAaAzcxqG0hJDP3dcRIWb5ORkVVVV+S2rqqpSXFxcs702khQdHa3o6OhQlOdTVl3b5BJDMNOry+XSiBEjVFRUpNGjR6uoqEizZ8/WQw89pCNHjsjtdmvXrl2aMGFCwPdtF6f2xrS358UuAcb41/+YrV0eC3Ib6qCOSKgjkmq1qg6nYWhRzlBLhmVEVLjJysrSW2+95bds3bp1ysrKsqii5qUndpPDkF/ACXZ6nTBhgoqKinTvvffq/fffV0FBgV577TVt2LBBBw8eVJ8+fTRgwICg7T+SnN4rs3JDWZPemEbhElwC9cbjNAz9+MpByjgvXrFRDtU1eH3/Nv5+7qmua7IulG2ogzoioY5IqtWqOtISYzvn3VJHjhzRrl27fF+XlZWptLRUCQkJ6tevn/Lz87Vv3z798pe/lCTdfvvtWrZsmX784x/r+9//vt5991299tprevPNN606hGaluGJUkDNM81d/JI9phiS9Tpw4UStXrtQHH3ygLl26aPDgwZo4caKKior05Zdfdupem5Z6ZYIdXjoaOILxxtPa72Fbfk9D0YY6qCMS6mhLG+oIPUvDzV/+8hdddtllvq8bx8ZMmzZNL7zwgioqKlReXu5bn56erjfffFOzZ8/Wk08+qfPOO0/PPvtsWN4Gnjumn8YP7BWy9No47mbJkiW+IDNx4kQtXrxYX375pe69996g7j9cWNEr01yAcUi6dXy6pl+aLqnjgaNRJLypAIDVLA03EydOlNnCf9I2N/vwxIkTtX379iBWFTgprpiQfdD06NFDGRkZevnll7Vs2TJJ0vjx43X99dfr+PHjfj03f/vb39TQ0KCDBw/q8OHDKi0tlSSNGDEiJLUGQ4W73i/IBKpXpqXgMmVYyhl7TE4PLAQOAAidiBpzg5ZNmDBBpaWlmjhxoiQpISFBQ4YMUVVVlQYNGuRrd9VVV+mzzz7zfT1y5EhJajFohqvGUPPs+2V+YaajR9JSz8uZeloIMAAQHgwzEj/ROqCmpkYul0tut1txcXF+644ePaqysjKlp6era9euFlUYfsLt59J46WnHPrcWv/XJWQeZU3tlmuuNIawAQPho6fP7dPTcIGKcfunpbDXXK0OYAQD7INwg7J3p0lNL2tMrQ6gBAHsh3CBsnW2omUGvDAB0aoQbhJWzHU9z6qUmemUAoHMj3CBsnPqw0bZwGNLcKwcr47x4emcAAD6EG4SFD/Z+qXmrd7Q4k2+jUy89EWgAAKcj3MBS7RlXQ6gBALQF4QaWWbWlXPP+b0eroeZM42kAAGgO4QaW8F2GOsN6xtMAAM6Ww+oCEFpFRUW65pprlJKSom7dumnEiBF6+eWXQ7LvCne9Nu6u1s/f261rl2884/gah6TXf3iJfjChv7L69yTYAADahZ6bTmbjxo3KyMjQ3LlzlZSUpDVr1igvL08ul0tXX3110Pbb1juhHIZUkDNMw1N7BK0WAIC9EW6Cyb1POrhbSugvufoGbTdr1qzRd7/7XX3xxRdyOp0qLS3VyJEjNXfuXC1evFiSdOutt+ro0aN66aWX/LadNWuW/vjHP2r16tVBCzdtuROKwcIAgEDhslSwbPultHSo9OLUk/9u+2XQdjVu3DgdPnxY27dvlyQVFxcrMTFRRUVFvjbFxcW+p4Wfzu12KyEhISi1rdpS3uIlKOnkL+EbMy/R/KuGEGwAAB1GuAkG9z7p97Mk03vya9Mr/f6ek8uDwOVyacSIEb4wU1RUpNmzZ2v79u06cuSI9u3bp127dmnChAlNtn3ttde0ZcsWTZ8+PeB1tTZoWJKchqGC67gMBQAIHMJNMBzc/VWwaWR6pIP/CNouJ0yYoKKiIpmmqffff185OTm68MILtWHDBhUXF6tPnz4aMGCA3zbr16/X9OnT9cwzz+jrX/96QOtpqcfGYUj5kwfrVzMu1oZ5lyl3TL+A7hsA0Lkx5iYYEvpLhsM/4BhOKeGCoO1y4sSJWrlypT744AN16dJFgwcP1sSJE1VUVKQvv/yySa9NcXGxpk6dqiVLligvLy+gtbTUY9N4JxQ9NQCAYKHnJhhcfaWpT54MNNLJf6cuDeqg4sZxN0uWLPEFmcZwU1RU5DfepqioSFOmTNEjjzyi2267LaB1tNZjwyUoAECw0XMTLBflSf0vP3kpKuGCoAYbSerRo4cyMjL08ssva9myZZKk8ePH6/rrr9fx48d9gWf9+vW6+uqrNWvWLF133XWqrKyUJEVFRXV4UHGFu1759NgAACxGz00wufpK6eOCHmwaTZgwQR6Px9dLk5CQoCFDhig5OVmDBg2SJL344ouqq6tTQUGBUlJSfK+cnJwO7bvCXa81H37e7Dw29NgAAELJMM22PIfZPmpqauRyueR2uxUXF+e37ujRoyorK1N6erq6du1qUYXhp7WfS0sT9DkkvT6THhsAQMe09Pl9Oi5LoUNamqDPaRhalDOUYAMACCnCDc5aS0/1fmDKhboqI4VJ+QAAIceYG5yVlgYPOw2DYAMAsAzhBmdl62dfnnHw8KKcoQQbAIBluCzVjE42xrpVp/88Gi9HnY7bvQEA4YCem1N06dJFklRXV2dxJeGloaFBkuR0Os94OYrbvQEA4YKem1M4nU7Fx8dr//79kqTY2FgZhmFxVdbyer06cOCAYmNjdaD2uN7cUdHs5aj/+Y+Runp4n9AXCADAaQg3p0lOTpYkX8CB5HA4VPqlUz9esb7ZYOM0DI1Ko8cGABAeCDenMQxDKSkp6t27t44fP251OWHh4FGvfvzCmYMNA4gBAOGEcHMGTqdTTqfT6jLCwp591c0GG+ayAQCEIwYUo1Xdopw6fegRc9kAAMIV4QYtWrWlXP/+1Ea/xytwKQoAEM64LIUzarzt+9RLUg5Jq3+YxS3fAICwRc8Nzqi5WYi9kuoavJbUAwBAWxBu0KxVW8p11yvbmyx3GobSEmMtqAgAgLYh3KCJlmYhZqwNACDcMeYGTZRV1zILMQAgYtFzgybOdOs3sxADACIB4QZ+uPUbABDpuCwFH279BgDYAT038GlurA23fgMAIg3hBj7pid3kaGasDbd+AwAiCeEGkk5ekiqrrtXcyYPl/NdoYsbaAAAiEWNuoFVbyn1jbRyGNPfKwco4L15pibEEGwBAxKHnppM7fRCx15QeXbuTYAMAiFiEm06uuUHEHtPUnuo6awoCAKCDCDed3Jkm7GMQMQAgUhFuOjEm7AMA2BEDijspJuwDANgVPTedFBP2AQDsinDTSTFhHwDArgg3ndgt/5buCziMtQEA2AVjbjqhUyftMyTdNj5d0y9NJ9gAAGyBnptO5vSBxKak597fY2VJAAAElOXhZvny5UpLS1PXrl2VmZmpzZs3t9h+6dKlGjRokGJiYpSamqrZs2fr6NGjIao28jFpHwDA7iwNN6tWrdKcOXO0cOFCbdu2TcOHD9ekSZO0f//+Ztu/8sormjdvnhYuXKiPP/5Yzz33nFatWqX58+eHuPLIxUBiAIDdWRpunnjiCc2YMUPTp0/XkCFDtGLFCsXGxmrlypXNtt+4caMuvfRS3XjjjUpLS9MVV1yhG264odXeHnwlxRWjgpxhPPkbAGBblg0obmho0NatW5Wfn+9b5nA4lJ2drZKSkma3ueSSS/TSSy9p8+bNGjt2rP7xj3/orbfe0s0333zG/Rw7dkzHjh3zfV1TUxO4g4hAFe56pSbEavUPs1TX4OUBmQAA27Es3FRXV8vj8SgpKclveVJSkj755JNmt7nxxhtVXV2tf/u3f5Npmjpx4oRuv/32Fi9LFRQU6KGHHgpo7ZHq1LukHIZUkDNMWf17Wl0WAAABZfmA4vYoKirSokWL9NRTT2nbtm1avXq13nzzTf3kJz854zb5+flyu92+1969e0NYcfg4/S4prynNX/2RKtz11hYGAECAWdZzk5iYKKfTqaqqKr/lVVVVSk5ObnabBx54QDfffLNuvfVWSdKwYcNUW1ur2267Tf/5n/8ph6NpVouOjlZ0dHTgDyDCtHSXFJelAAB2YlnPTVRUlEaNGqXCwkLfMq/Xq8LCQmVlZTW7TV1dXZMA43Q6JUmmaTa3Cf6Fu6QAAJ2FpZel5syZo2eeeUYvvviiPv74Y91xxx2qra3V9OnTJUl5eXl+A46nTp2qp59+Wq+++qrKysq0bt06PfDAA5o6daov5KB53CUFAOgsLH38Qm5urg4cOKAFCxaosrJSI0aM0Nq1a32DjMvLy/16au6//34ZhqH7779f+/btU69evTR16lT913/9l1WHEBEq3PUqq67V+IG9tGHeZdpTXcddUgAA2zLMTnY9p6amRi6XS263W3FxcVaXE3TN3SGVO6af1WUBANAu7fn8jqi7pdA+3CEFAOiMCDc2xnOkAACdEeHGxrhDCgDQGRFubIw7pAAAnZGld0shuHiOFACgMyLc2BTPkQIAdFZclrIh7pICAHRmhBsb4i4pAEBnRrixIe6SAgB0ZoQbG+IuKQBAZ8aAYpvKHdNP4wf24jlSAIBOh3BjQ40PykxP7MYdUgCATodwYzM8KBMA0Nkx5sZGuAUcAADCja1wCzgAAIQbW+EWcAAACDe2wi3gAAAwoNh2uAUcANDZEW5sKMUVQ6gBAHRaXJaykQp3vTburubuKABAp0bPjU0wvw0AACfRc2MDzG8DAMBXCDc2wPw2AAB8hXBjA8xvAwDAVwg3NsD8NgAAfIUBxTbB/DYAAJxEuLER5rcBAIDLUhGPuW0AAPBHz00EY24bAACaoucmQjG3DQAAzSPcRCjmtgEAoHmEmwjF3DYAADSPcBOhmNsGAIDmMaA4gjG3DQAATRFuIhxz2wAA4I/LUgAAwFYINwAAwFYINxGKmYkBAGgeY24iEDMTAwBwZvTcRBhmJgYAoGWEmwjDzMQAALSMcBNhmJkYAICWEW4iDDMTAwDQMgYURyBmJgYA4MwINxGKmYkBAGgel6UAAICtEG4iDJP3AQDQMi5LRRAm7wMAoHX03EQIJu8DAKBtCDcRgsn7AABoG8JNhGDyPgAA2oZwEyGYvA8AgLZhQHEEYfI+AABaR7iJMEzeBwBAy7gsBQAAbIVwAwAAbIVwAwAAbMXycLN8+XKlpaWpa9euyszM1ObNm1tsf+jQIc2cOVMpKSmKjo7WwIED9dZbb4WoWmvwyAUAANrO0gHFq1at0pw5c7RixQplZmZq6dKlmjRpknbu3KnevXs3ad/Q0KBvfvOb6t27t37zm9+ob9+++uyzzxQfHx/64kOERy4AANA+hmmaZuvNgiMzM1NjxozRsmXLJEler1epqam66667NG/evCbtV6xYoccee0yffPKJunTp0qZ9HDt2TMeOHfN9XVNTo9TUVLndbsXFxQXmQIKkwl2vSxe/6zczsdMwtGHeZdwxBQDoVGpqauRyudr0+W3ZZamGhgZt3bpV2dnZXxXjcCg7O1slJSXNbvO73/1OWVlZmjlzppKSkjR06FAtWrRIHo/njPspKCiQy+XyvVJTUwN+LMHCIxcAAGg/y8JNdXW1PB6PkpKS/JYnJSWpsrKy2W3+8Y9/6De/+Y08Ho/eeustPfDAA3r88cf105/+9Iz7yc/Pl9vt9r327t0b0OMIJh65AABA+0XUJH5er1e9e/fWL37xCzmdTo0aNUr79u3TY489poULFza7TXR0tKKjo0NcaWA0PnJh/uqP5DFNHrkAAEAbWBZuEhMT5XQ6VVVV5be8qqpKycnJzW6TkpKiLl26yOl0+pZdeOGFqqysVENDg6KiooJasxV45AIAAO1j2WWpqKgojRo1SoWFhb5lXq9XhYWFysrKanabSy+9VLt27ZLX6/Ut+/vf/66UlBRbBptGKa4YZfXvSbABAKANLJ3nZs6cOXrmmWf04osv6uOPP9Ydd9yh2tpaTZ8+XZKUl5en/Px8X/s77rhDBw8e1KxZs/T3v/9db775phYtWqSZM2dadQgAACDMWDrmJjc3VwcOHNCCBQtUWVmpESNGaO3atb5BxuXl5XI4vspfqamp+sMf/qDZs2crIyNDffv21axZszR37lyrDgEAAIQZS+e5sUJ77pMHAADhISLmuQEAAAgGwk0Y45lSAAC0X0TNc9OZ8EwpAADODj03YajCXe8LNpLkNaX5qz+iBwcAgDYg3IQhnikFAMDZI9yEIZ4pBQDA2SPchKHGZ0o5jZMJh2dKAQDQdgwoDlM8UwoAgLNDuAljKa4YQg0AAO3EZSkAAGArhBsAAGArhBsAAGArhBsAAGArhBsAAGArhJswxAMzAQA4e9wKHmZ4YCYAAB1Dz00Y4YGZAAB0HOEmjPDATAAAOo5wE0Z4YCYAAB1HuAkjPDATAICOY0BxmOGBmQAAdAzhJgzxwEwAAM7eWV2WysvL0/PPP6/du3cHuh4AAIAOOatwExUVpYKCAg0YMECpqan67ne/q2effVaffvppoOsDAABoF8M0TbP1Zs3bt2+f3nvvPRUXF6u4uFh///vflZKSon/+85+BrDGgampq5HK55Ha7FRcXZ3U5AACgDdrz+d2hu6V69Oihnj17qkePHoqPj9c555yjXr16deRbAgAAdMhZhZv58+frkksuUc+ePTVv3jwdPXpU8+bNU2VlpbZv3x7oGgEAANrsrC5LORwO9erVS7Nnz1ZOTo4GDhwYjNqCIpwvS1W461VWXav0xG7cLQUAwCna8/l9VreCb9++XcXFxSoqKtLjjz+uqKgoTZgwQRMnTtTEiRMjKuyECx6YCQBAYHRoQHGjDz74QEuWLNHLL78sr9crj8cTiNqCIhx7birc9bp08bt+z5VyGoY2zLuMHhwAABSCnhvTNLV9+3YVFRWpqKhIGzZsUE1NjTIyMjRhwoSzKroza+mBmYQbAADa56zCTUJCgo4cOaLhw4drwoQJmjFjhsaNG6f4+PgAl9c5ND4w8/SeGx6YCQBA+51VuHnppZc0bty4sLmsE+kaH5g5f/VH8pgmD8wEAKADAjLmJpKE45ibRhXueh6YCQBAM4I+5gbBwQMzAQDouA7NUAwAABBuCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcWq3DXa+PualW4660uBQAAW2ASPwut2lKu/NU75DUlhyEV5AxT7ph+VpcFAEBEo+fGIhXuel+wkU4+NHP+6o/owQEAoIMINxYpq671ewq4JHlMU3uq66wpCAAAmyDcWCQ9sZschv8yp2EoLTHWmoIAALAJwo1FUlwxKsgZJqdxMuE4DUOLcoby4EwAADqIAcUWyh3TT+MH9tKe6jqlJcYSbAAACADCjcVSXDGEGgAAAojLUgAAwFYINwAAwFYINwAAwFYINwAAwFYINwAAwFYINwAAwFbCItwsX75caWlp6tq1qzIzM7V58+Y2bffqq6/KMAxde+21wS0QAABEDMvDzapVqzRnzhwtXLhQ27Zt0/DhwzVp0iTt37+/xe327Nmj++67T+PGjQtRpQAAIBJYHm6eeOIJzZgxQ9OnT9eQIUO0YsUKxcbGauXKlWfcxuPx6KabbtJDDz2kCy64IITVAgCAcGdpuGloaNDWrVuVnZ3tW+ZwOJSdna2SkpIzbvfwww+rd+/euuWWW1rdx7Fjx1RTU+P3CgcV7npt3F2tCne91aUAAGArlj5+obq6Wh6PR0lJSX7Lk5KS9MknnzS7zYYNG/Tcc8+ptLS0TfsoKCjQQw891NFSA2rVlnLlr94hryk5DKkgZ5hyx/SzuiwAAGzB8stS7XH48GHdfPPNeuaZZ5SYmNimbfLz8+V2u32vvXv3BrnKllW4633BRpK8pjR/9Uf04AAAECCW9twkJibK6XSqqqrKb3lVVZWSk5ObtN+9e7f27NmjqVOn+pZ5vV5J0jnnnKOdO3eqf//+fttER0crOjo6CNWfnbLqWl+waeQxTe2pruMBmgAABIClPTdRUVEaNWqUCgsLfcu8Xq8KCwuVlZXVpP3gwYO1Y8cOlZaW+l7f+ta3dNlll6m0tFSpqamhLP+spCd2k8PwX+Y0DKUlxlpTEAAANmNpz40kzZkzR9OmTdPo0aM1duxYLV26VLW1tZo+fbokKS8vT3379lVBQYG6du2qoUOH+m0fHx8vSU2Wh6sUV4wKcoZp/uqP5DFNOQ1Di3KG0msDAECAWB5ucnNzdeDAAS1YsECVlZUaMWKE1q5d6xtkXF5eLocjooYGtSp3TD+NH9hLe6rrlJYYS7ABACCADNM0zdab2UdNTY1cLpfcbrfi4uKsLgcAALRBez6/7dUlAgAAOj3CDQAAsBXCDQAAsBXCDQAAsBXCDQAAsBXCDQAAsBXCDQAAsBXCDQAAsBXCDQAAsBXCTYhVuOu1cXe1Ktz1VpcCAIAtWf5sqc5k1ZZy5a/eIa8pOQypIGeYcsf0s7osAABshZ6bEKlw1/uCjSR5TWn+6o/owQEAIMAINyFSVl3rCzaNPKapPdV11hQEAIBNEW5CJD2xmxyG/zKnYSgtMdaaggAAsCnCTYikuGJUkDNMTuNkwnEahhblDFWKK8biygAAsBcGFIdQ7ph+Gj+wl/ZU1yktMZZgAwBAEBBuQizFFUOoAQAgiLgsBQAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwAwAAbIVwEyIV7npt3F2tCne91aUAAGBrPBU8BFZtKVf+6h3ympLDkApyhil3TD+rywIAwJbouQmyCne9L9hIkteU5q/+iB4cAACChHATZGXVtb5g08hjmtpTXWdNQQAA2BzhJsjSE7vJYfgvcxqG0hJjrSkIAACbI9wEWYorRgU5w+Q0TiYcp2FoUc5QpbhiLK4MAAB7YkBxCOSO6afxA3tpT3Wd0hJjCTYAAAQR4SZEUlwxhBoAAEKAy1IAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDcAAMBWCDdBVOGu18bd1apw11tdCgAAncY5VhdgV6u2lCt/9Q55TclhSAU5w5Q7pp/VZQEAYHv03ARBhbveF2wkyWtK81d/RA8OAAAhQLgJgrLqWl+waeQxTe2prrOmIAAAOhHCTRCkJ3aTw/Bf5jQMpSXGWlMQAACdCOEmCFJcMSrIGSancTLhOA1Di3KGKsUVY3FlAADYX1iEm+XLlystLU1du3ZVZmamNm/efMa2zzzzjMaNG6cePXqoR48eys7ObrG9VXLH9NOGeZfpVzMu1oZ5lzGYGACAELE83KxatUpz5szRwoULtW3bNg0fPlyTJk3S/v37m21fVFSkG264QevXr1dJSYlSU1N1xRVXaN++fSGuvHUprhhl9e9Jjw0AACFkmKZptt4seDIzMzVmzBgtW7ZMkuT1epWamqq77rpL8+bNa3V7j8ejHj16aNmyZcrLy2u1fU1NjVwul9xut+Li4jpcPwAACL72fH5b2nPT0NCgrVu3Kjs727fM4XAoOztbJSUlbfoedXV1On78uBISEppdf+zYMdXU1Pi9AACAfVkabqqrq+XxeJSUlOS3PCkpSZWVlW36HnPnzlWfPn38AtKpCgoK5HK5fK/U1NQO1w0AAMKX5WNuOmLx4sV69dVX9frrr6tr167NtsnPz5fb7fa99u7dG+IqAQBAKFn6+IXExEQ5nU5VVVX5La+qqlJycnKL2/73f/+3Fi9erHfeeUcZGRlnbBcdHa3o6OiA1AsAAMKfpT03UVFRGjVqlAoLC33LvF6vCgsLlZWVdcbtHn30Uf3kJz/R2rVrNXr06FCUCgAAIoTlD86cM2eOpk2bptGjR2vs2LFaunSpamtrNX36dElSXl6e+vbtq4KCAknSI488ogULFuiVV15RWlqab2xO9+7d1b17d8uOAwAAhAfLw01ubq4OHDigBQsWqLKyUiNGjNDatWt9g4zLy8vlcHzVwfT000+roaFB3/72t/2+z8KFC/Xggw+GsnQAABCGLJ/nJtSY5wYAgMgTMfPcAAAABBrhBgAA2ArhBgAA2ArhBgAA2ArhBgAA2ArhJggq3PXauLtaFe56q0sBAKDTsXyeG7tZtaVc+at3yGtKDkMqyBmm3DH9rC4LAIBOg56bAKpw1/uCjSR5TWn+6o/owQEAIIQINwFUVl3rCzaNPKapPdV11hQEAEAnRLgJoPTEbnIY/suchqG0xFhrCgIAoBMi3ARQiitGBTnD5DROJhynYWhRzlCluGIsrgwAgM6DAcUBljumn8YP7KU91XVKS4wl2AAAEGKEmyBIccUQagAAsAiXpQAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgAAgK0QbgKowl2vjburVeGut7oUAAA6rXOsLsAuVm0pV/7qHfKaksOQCnKGKXdMP6vLAgCg06HnJgAq3PW+YCNJXlOav/ojenAAALAA4SYAyqprfcGmkcc0tae6zpqCAADoxAg3AZCe2E0Ow3+Z0zCUlhhrTUEAAHRihJsASHHFqCBnmJzGyYTjNAwtyhmqFFeMxZUBAND5MKA4QHLH9NP4gb20p7pOaYmxBBsAACxCuAmgFFcMoQYAAItxWQoAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANgK4QYAANhKWISb5cuXKy0tTV27dlVmZqY2b97cYvtf//rXGjx4sLp27aphw4bprbfeClGlAAAg3J1jdQGrVq3SnDlztGLFCmVmZmrp0qWaNGmSdu7cqd69ezdpv3HjRt1www0qKCjQ1VdfrVdeeUXXXnuttm3bpqFDh1pwBKdw75MO7pa6dJOO1371b0L/k+ubWxeoNqHYB3VEbh2RVCt1UEck1BFJtVpVR0J/ydVXVjBM0zQt2fO/ZGZmasyYMVq2bJkkyev1KjU1VXfddZfmzZvXpH1ubq5qa2u1Zs0a37KLL75YI0aM0IoVK1rdX01NjVwul9xut+Li4gJ3INt+Kf1+lmR6m1lp/Ovfln7UHW0Tin1QR+TW0ZY21EEd1BHYNp28DsMhTX1SuiivhfZt157Pb0svSzU0NGjr1q3Kzs72LXM4HMrOzlZJSUmz25SUlPi1l6RJkyadsf2xY8dUU1Pj9wo4974Wgo108hehtQzZ0Tah2Ad1RG4dbWlDHdRBHYFt08nrML3S7+85+RkZYpaGm+rqank8HiUlJfktT0pKUmVlZbPbVFZWtqt9QUGBXC6X75WamhqY4k91cHcLwQYAgE7K9EgH/xHy3YbFgOJgys/Pl9vt9r327t0b+J0k9D/Z/QYAAL5iOKWEC0K+W0s/kRMTE+V0OlVVVeW3vKqqSsnJyc1uk5yc3K720dHRiouL83sFnKvvyeuKhvMMDQzJMM6wLlBtQrEP6ojcOtrShjqogzoC26aT12E4palLLRlUbOndUlFRURo1apQKCwt17bXXSjo5oLiwsFB33nlns9tkZWWpsLBQ99xzj2/ZunXrlJWVFYKKW3BRntT/8pPdb11ipeN1X/3bmFqbWxeoNqHYB3VEbh2RVCt1UEck1BFJtVpVR8IFnfduqVWrVmnatGn6+c9/rrFjx2rp0qV67bXX9MknnygpKUl5eXnq27evCgoKJJ28FXzChAlavHixpkyZoldffVWLFi1q863gQbtbCgAABE17Pr8tn+cmNzdXBw4c0IIFC1RZWakRI0Zo7dq1vkHD5eXlcji+unp2ySWX6JVXXtH999+v+fPna8CAAXrjjTesn+MGAACEBct7bkKNnhsAACJPxMxzAwAAEGiEGwAAYCuEGwAAYCuEGwAAYCuEGwAAYCuEGwAAYCuEGwAAYCuEGwAAYCuEGwAAYCuWP34h1BonZK6pqbG4EgAA0FaNn9ttebBCpws3hw8fliSlpqZaXAkAAGivw4cPy+Vytdim0z1byuv16vPPP9e5554rwzAC+r1ramqUmpqqvXv32vK5VXY/PoljtAO7H5/EMdqB3Y9PCvwxmqapw4cPq0+fPn4P1G5Op+u5cTgcOu+884K6j7i4ONv+skr2Pz6JY7QDux+fxDHagd2PTwrsMbbWY9OIAcUAAMBWCDcAAMBWCDcBFB0drYULFyo6OtrqUoLC7scncYx2YPfjkzhGO7D78UnWHmOnG1AMAADsjZ4bAABgK4QbAABgK4QbAABgK4QbAABgK4SbAFm+fLnS0tLUtWtXZWZmavPmzVaXdNYKCgo0ZswYnXvuuerdu7euvfZa7dy506/NxIkTZRiG3+v222+3qOL2efDBB5vUPnjwYN/6o0ePaubMmerZs6e6d++u6667TlVVVRZW3H5paWlNjtEwDM2cOVNSZJ6/9957T1OnTlWfPn1kGIbeeOMNv/WmaWrBggVKSUlRTEyMsrOz9emnn/q1OXjwoG666SbFxcUpPj5et9xyi44cORLCozizlo7v+PHjmjt3roYNG6Zu3bqpT58+ysvL0+eff+73PZo774sXLw7xkZxZa+fwe9/7XpP6r7zySr824XwOpdaPsbm/S8Mw9Nhjj/nahPN5bMvnQ1veQ8vLyzVlyhTFxsaqd+/e+tGPfqQTJ04ErE7CTQCsWrVKc+bM0cKFC7Vt2zYNHz5ckyZN0v79+60u7awUFxdr5syZ+vOf/6x169bp+PHjuuKKK1RbW+vXbsaMGaqoqPC9Hn30UYsqbr+vf/3rfrVv2LDBt2727Nn6/e9/r1//+tcqLi7W559/rpycHAurbb8tW7b4Hd+6deskSd/5znd8bSLt/NXW1mr48OFavnx5s+sfffRR/c///I9WrFihTZs2qVu3bpo0aZKOHj3qa3PTTTfpr3/9q9atW6c1a9bovffe02233RaqQ2hRS8dXV1enbdu26YEHHtC2bdu0evVq7dy5U9/61reatH344Yf9zutdd90VivLbpLVzKElXXnmlX/2/+tWv/NaH8zmUWj/GU4+toqJCK1eulGEYuu666/zahet5bMvnQ2vvoR6PR1OmTFFDQ4M2btyoF198US+88IIWLFgQuEJNdNjYsWPNmTNn+r72eDxmnz59zIKCAgurCpz9+/ebkszi4mLfsgkTJpizZs2yrqgOWLhwoTl8+PBm1x06dMjs0qWL+etf/9q37OOPPzYlmSUlJSGqMPBmzZpl9u/f3/R6vaZpRvb5M03TlGS+/vrrvq+9Xq+ZnJxsPvbYY75lhw4dMqOjo81f/epXpmma5t/+9jdTkrllyxZfm7fffts0DMPct29fyGpvi9OPrzmbN282JZmfffaZb9n5559vLlmyJLjFBUhzxzht2jTzmmuuOeM2kXQOTbNt5/Gaa64xv/GNb/gti6TzePrnQ1veQ9966y3T4XCYlZWVvjZPP/20GRcXZx47diwgddFz00ENDQ3aunWrsrOzfcscDoeys7NVUlJiYWWB43a7JUkJCQl+y19++WUlJiZq6NChys/PV11dnRXlnZVPP/1Uffr00QUXXKCbbrpJ5eXlkqStW7fq+PHjfudz8ODB6tevX8Sez4aGBr300kv6/ve/7/ew2Eg+f6crKytTZWWl33lzuVzKzMz0nbeSkhLFx8dr9OjRvjbZ2dlyOBzatGlTyGvuKLfbLcMwFB8f77d88eLF6tmzp0aOHKnHHnssoF39oVBUVKTevXtr0KBBuuOOO/TFF1/41tntHFZVVenNN9/ULbfc0mRdpJzH0z8f2vIeWlJSomHDhikpKcnXZtKkSaqpqdFf//rXgNTV6R6cGWjV1dXyeDx+J0mSkpKS9Mknn1hUVeB4vV7dc889uvTSSzV06FDf8htvvFHnn3+++vTpow8//FBz587Vzp07tXr1agurbZvMzEy98MILGjRokCoqKvTQQw9p3Lhx+uijj1RZWamoqKgmHxhJSUmqrKy0puAOeuONN3To0CF973vf8y2L5PPXnMZz09zfYeO6yspK9e7d22/9Oeeco4SEhIg7t0ePHtXcuXN1ww03+D2Q8O6779ZFF12khIQEbdy4Ufn5+aqoqNATTzxhYbVtd+WVVyonJ0fp6enavXu35s+fr8mTJ6ukpEROp9NW51CSXnzxRZ177rlNLntHynls7vOhLe+hlZWVzf6tNq4LBMINWjRz5kx99NFHfmNSJPld4x42bJhSUlJ0+eWXa/fu3erfv3+oy2yXyZMn+/5/RkaGMjMzdf755+u1115TTEyMhZUFx3PPPafJkyerT58+vmWRfP46u+PHj+v666+XaZp6+umn/dbNmTPH9/8zMjIUFRWlH/zgByooKIiIaf7/4z/+w/f/hw0bpoyMDPXv319FRUW6/PLLLawsOFauXKmbbrpJXbt29VseKefxTJ8P4YDLUh2UmJgop9PZZCR4VVWVkpOTLaoqMO68806tWbNG69ev13nnnddi28zMTEnSrl27QlFaQMXHx2vgwIHatWuXkpOT1dDQoEOHDvm1idTz+dlnn+mdd97Rrbfe2mK7SD5/knznpqW/w+Tk5CaD/E+cOKGDBw9GzLltDDafffaZ1q1b59dr05zMzEydOHFCe/bsCU2BAXbBBRcoMTHR93tph3PY6P3339fOnTtb/duUwvM8nunzoS3vocnJyc3+rTauCwTCTQdFRUVp1KhRKiws9C3zer0qLCxUVlaWhZWdPdM0deedd+r111/Xu+++q/T09Fa3KS0tlSSlpKQEubrAO3LkiHbv3q2UlBSNGjVKXbp08TufO3fuVHl5eUSez+eff169e/fWlClTWmwXyedPktLT05WcnOx33mpqarRp0ybfecvKytKhQ4e0detWX5t3331XXq/XF+7CWWOw+fTTT/XOO++oZ8+erW5TWloqh8PR5FJOpPjnP/+pL774wvd7Genn8FTPPfecRo0apeHDh7faNpzOY2ufD215D83KytKOHTv8gmpjWB8yZEjACkUHvfrqq2Z0dLT5wgsvmH/729/M2267zYyPj/cbCR5J7rjjDtPlcplFRUVmRUWF71VXV2eapmnu2rXLfPjhh82//OUvZllZmfnb3/7WvOCCC8zx48dbXHnb3HvvvWZRUZFZVlZm/ulPfzKzs7PNxMREc//+/aZpmubtt99u9uvXz3z33XfNv/zlL2ZWVpaZlZVlcdXt5/F4zH79+plz5871Wx6p5+/w4cPm9u3bze3bt5uSzCeeeMLcvn27726hxYsXm/Hx8eZvf/tb88MPPzSvueYaMz093ayvr/d9jyuvvNIcOXKkuWnTJnPDhg3mgAEDzBtuuMGqQ/LT0vE1NDSY3/rWt8zzzjvPLC0t9fu7bLy7ZOPGjeaSJUvM0tJSc/fu3eZLL71k9urVy8zLy7P4yL7S0jEePnzYvO+++8ySkhKzrKzMfOedd8yLLrrIHDBggHn06FHf9wjnc2iarf+emqZput1uMzY21nz66aebbB/u57G1zwfTbP099MSJE+bQoUPNK664wiwtLTXXrl1r9urVy8zPzw9YnYSbAPnZz35m9uvXz4yKijLHjh1r/vnPf7a6pLMmqdnX888/b5qmaZaXl5vjx483ExISzOjoaPNrX/ua+aMf/ch0u93WFt5Gubm5ZkpKihkVFWX27dvXzM3NNXft2uVbX19fb/7whz80e/ToYcbGxpr//u//blZUVFhY8dn5wx/+YEoyd+7c6bc8Us/f+vXrm/29nDZtmmmaJ28Hf+CBB8ykpCQzOjravPzyy5sc+xdffGHecMMNZvfu3c24uDhz+vTp5uHDhy04mqZaOr6ysrIz/l2uX7/eNE3T3Lp1q5mZmWm6XC6za9eu5oUXXmguWrTILxhYraVjrKurM6+44gqzV69eZpcuXczzzz/fnDFjRpP/SAznc2iarf+emqZp/vznPzdjYmLMQ4cONdk+3M9ja58Pptm299A9e/aYkydPNmNiYszExETz3nvvNY8fPx6wOo1/FQsAAGALjLkBAAC2QrgBAAC2QrgBAAC2QrgBAAC2QrgBAAC2QrgBAAC2QrgBAAC2QrgBAAC2QrgBEHAPPvigRowY0a5tDMPQG2+8EZR6OmrixIm65557rC4DQBsRbgC0SXs+4O+77z6/B+cBQCidY3UBAOzDNE15PB51795d3bt3t7qcsObxeGQYhhwO/hsTCDT+qgC06nvf+56Ki4v15JNPyjAMGYahPXv2qKioSIZh6O2339aoUaMUHR2tDRs2NLkstWXLFn3zm99UYmKiXC6XJkyYoG3btrWrhokTJ+ruu+/Wj3/8YyUkJCg5OVkPPvigb/2ePXtkGIZKS0t9yw4dOiTDMFRUVCRJvnr/8Ic/aOTIkYqJidE3vvEN7d+/X2+//bYuvPBCxcXF6cYbb1RdXZ3f/k+cOKE777xTLpdLiYmJeuCBB3Tqo/mOHTum++67T3379lW3bt2UmZnp268kvfDCC4qPj9fvfvc7DRkyRNHR0SovL2/XzwBA2xBuALTqySefVFZWlmbMmKGKigpVVFQoNTXVt37evHlavHixPv74Y2VkZDTZ/vDhw5o2bZo2bNigP//5zxowYICuuuoqHT58uF11vPjii+rWrZs2bdqkRx99VA8//LDWrVvX7uN58MEHtWzZMm3cuFF79+7V9ddfr6VLl+qVV17Rm2++qT/+8Y/62c9+1mTf55xzjjZv3qwnn3xSTzzxhJ599lnf+jvvvFMlJSV69dVX9eGHH+o73/mOrrzySn366ae+NnV1dXrkkUf07LPP6q9//at69+7d7toBtI7LUgBa5XK5FBUVpdjYWCUnJzdZ//DDD+ub3/zmGbf/xje+4ff1L37xC8XHx6u4uFhXX311m+vIyMjQwoULJUkDBgzQsmXLVFhY2OK+m/PTn/5Ul156qSTplltuUX5+vnbv3q0LLrhAkvTtb39b69ev19y5c33bpKamasmSJTIMQ4MGDdKOHTu0ZMkSzZgxQ+Xl5Xr++edVXl6uPn36SDo57mjt2rV6/vnntWjRIknS8ePH9dRTT2n48OHtqhdA+9BzA6DDRo8e3eL6qqoqzZgxQwMGDJDL5VJcXJyOHDnS7ssyp/cKpaSkaP/+/e2u99Tvk5SUpNjYWF+waVx2+ve9+OKLZRiG7+usrCx9+umn8ng82rFjhzwejwYOHOgbb9S9e3cVFxdr9+7dvm2ioqKa7dkCEFj03ADosG7durW4ftq0afriiy/05JNP6vzzz1d0dLSysrLU0NDQrv106dLF72vDMOT1eiXJNzD31HEwx48fb/X7GIbR4vdtiyNHjsjpdGrr1q1yOp1+604dWB0TE+MXkAAEB+EGQJtERUXJ4/Gc1bZ/+tOf9NRTT+mqq66SJO3du1fV1dWBLE+9evWSJFVUVGjkyJGS5De4uKM2bdrk93Xj2CGn06mRI0fK4/Fo//79GjduXMD2CeDsEG4AtElaWpo2bdqkPXv2qHv37kpISGjztgMGDND//u//avTo0aqpqdGPfvQjxcTEBLS+mJgYXXzxxVq8eLHS09O1f/9+3X///QH7/uXl5ZozZ45+8IMfaNu2bfrZz36mxx9/XJI0cOBA3XTTTcrLy9Pjjz+ukSNH6sCBAyosLFRGRoamTJkSsDoAtI4xNwDa5L777pPT6dSQIUPUq1evdo2Xee655/Tll1/qoosu0s0336y77747KHcKrVy5UidOnNCoUaN0zz336Kc//WnAvndeXp7q6+s1duxYzZw5U7NmzdJtt93mW//8888rLy9P9957rwYNGqRrr71WW7ZsUb9+/QJWA4C2McxTL1ADAABEOHpuAACArRBuAACArRBuAACArRBuAACArRBuAACArRBuAACArRBuAACArRBuAACArRBuAACArRBuAACArRBuAACArfw//QhAtfiNz9oAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(ws, '.')\n",
    "plt.xlabel('trial number')\n",
    "plt.ylabel('w')\n",
    "plt.legend(['w1', 'w2'])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fb5fd6ef-a890-4d1b-a3cd-58cd250121db",
   "metadata": {},
   "source": [
    "### Exercise 1.3\n",
    "\n",
    "In overshadowing, two stimuli are presented together. What should we expect conditioned responding to B to be, after this paradigm? \\\n",
    "AB—> + \\\n",
    "B —> ?\n",
    "\n",
    "In overexpectation, we pair A and B independently with reward, until asymptote, before the compound phase. What happens to the A and B associations during this compound learning phase? \\\n",
    "A—>+ \\\n",
    "B—>+ \\\n",
    "AB—>+ \\\n",
    "B—>?\n",
    "\n",
    "Can you show what the RW model predicts in these cases? "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "id": "bd8052ac-edbd-4195-8671-8e15d3775de5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# write your code here "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9e808983-0388-40b1-a108-0631fb990e47",
   "metadata": {},
   "source": [
    "# 2. Kalman Filter model of classical conditioning\n",
    "\n",
    "The Rescorla-Wagner model gives the maximum-likelihood estimate of the weights, but does not give an estimate of the uncertainty the animal might have about these weights. A Bayesian generalization of the RW model leads to the Kalman Filter. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 105,
   "id": "c54f620b-fb1d-4ef1-8471-b7f625afe524",
   "metadata": {},
   "outputs": [],
   "source": [
    "def kalman_filter_update(mean_w, covariance, stimuli, reward):\n",
    "\n",
    "    # Prediction step;\n",
    "    a_priori_covariance = covariance + transition_noise_variance * np.eye(len(mean_w))\n",
    "\n",
    "    # Compute statistics of interest;\n",
    "    expected_reward = np.dot(stimuli, mean_w)\n",
    "    prediction_error = reward - expected_reward  # the \"prediction error\". Actually expected value of the prediction error\n",
    "    residual_cov = np.dot(stimuli, np.matmul(a_priori_covariance, stimuli)) + observation_noise_variance\n",
    "\n",
    "    # Correction step;\n",
    "    kalman_gain = np.matmul(a_priori_covariance, stimuli) * residual_cov**-1\n",
    "    mean_w = mean_w + kalman_gain * prediction_error  # weight update\n",
    "    covariance = a_priori_covariance - np.outer(kalman_gain, residual_cov*kalman_gain)  # covariance update \n",
    "\n",
    "    return mean_w, covariance"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a756398f-c917-47a1-a391-553d046c9fee",
   "metadata": {},
   "source": [
    "Let's simulate a very simple paradigm whereby a single stimulus is presented many times without any reward for 100 trials, after which the stimulus is paired with a reward for another 100 trials. \n",
    "\n",
    "### Exercise 2.1\n",
    "What happens to the variance estimate during the first 100 trials? What does this mean for subsequent conditioning? How does this differ from the original Rescorla-Wagner model?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 111,
   "id": "3f19ddb9-c623-46aa-9174-3f82e0d6f520",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'w')"
      ]
     },
     "execution_count": 111,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAGwCAYAAAC5ACFFAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABQMUlEQVR4nO3de1hUdeI/8PfMAMNNZkAuA4qKihfyGiTSRd0kId3v5je3tCwva7i1khlaSk9qaYVmmWVubrve+qbbbdNuG0UYlEloKGupsemjosmASjBc5DZzfn98fjMwAsNF4HCG9+t55mHmnM8cPscRzpvP7agkSZJARERERE1Sy10BIiIiou6MYYmIiIjIAYYlIiIiIgcYloiIiIgcYFgiIiIicoBhiYiIiMgBhiUiIiIiB1zkroAzsFgsuHjxInr16gWVSiV3dYiIiKgVJElCWVkZQkJCoFY3337EsNQBLl68iNDQULmrQURERO1w/vx59O3bt9n9DEsdoFevXgDEP7aPj4/MtSEiIqLWMJlMCA0NtV3Hm8Ow1AGsXW8+Pj4MS0RERArT0hAaDvAmIiIicoBhiYiIiMgBhiUiIiIiBzhmiYiIyImYzWbU1tbKXY1uwdXVFRqN5rqPw7BERETkBCRJgtFoRElJidxV6Vb0ej0MBsN1rYPIsEREROQErEEpMDAQnp6ePX6RZEmSUFlZiaKiIgBAcHBwu4/FsERERKRwZrPZFpR69+4td3W6DQ8PDwBAUVERAgMD290lxwHeRERECmcdo+Tp6SlzTbof67/J9YzjYlgiIiJyEj29660pHfFvoqiw9M033+B//ud/EBISApVKhX379rX4noyMDNx4443QarUYPHgwdu7c2ajMli1bMGDAALi7uyM6OhqHDh3q+MoTERGRIikqLFVUVGD06NHYsmVLq8qfOXMG06ZNw+9+9zvk5uZiyZIleOihh/DFF1/Yyrz77rtISkrC6tWrceTIEYwePRpxcXG2AWFERETUs6kkSZLkrkR7qFQq7N27F9OnT2+2zPLly/HZZ5/hp59+sm2bNWsWSkpKkJqaCgCIjo7GTTfdhNdffx0AYLFYEBoaikcffRQrVqxo8rjV1dWorq62vbbeiK+0tJT3hiMioi5XVVWFM2fOICwsDO7u7nJXp1tx9G9jMpmg0+lavH4rqmWprbKyshAbG2u3LS4uDllZWQCAmpoa5OTk2JVRq9WIjY21lWlKSkoKdDqd7REaGto5J0BEREQAROiZN28eRo4cCRcXF4eNJR3NqcOS0WhEUFCQ3bagoCCYTCZcvXoVly9fhtlsbrKM0Whs9rjJyckoLS21Pc6fP98p9QeA6mpAmW1/REREHcdsNsPDwwOLFy9u1BDS2bjOUjtotVpotdou+V4nTgD9+wN+fl3y7YiIyElIElBZKc/39vQEWjMJ7dNPP8UDDzyAK1euQKPRIDc3F2PHjsXy5cuxbt06AMBDDz2EqqoqvP3223jjjTcAAN99912XrlTu1GHJYDCgsLDQblthYSF8fHzg4eEBjUYDjUbTZBmDwdCVVW2W2QzU1cldCyIiUprKSsDbW57vXV4OeHm1XO62225DWVkZjh49iqioKGRmZsLf3x8ZGRm2MpmZmVi+fHnnVbYVnLobLiYmBunp6Xbb0tLSEBMTAwBwc3NDZGSkXRmLxYL09HRbGblZLOyGIyIi56TT6TBmzBhbOMrIyMDjjz+Oo0ePory8HL/++itOnTqFiRMnylpPRbUslZeX49SpU7bXZ86cQW5uLvz8/NCvXz8kJyfj119/xVtvvQUAePjhh/H666/jySefxJ/+9Cfs378f7733Hj777DPbMZKSkjB37lxERUVh3Lhx2LRpEyoqKjB//vwuP7+mSBLDEhERtZ2np2jhket7t9bEiRORkZGBpUuX4ttvv0VKSgree+89HDhwAMXFxQgJCUF4eHjnVbYVFBWWfvjhB/zud7+zvU5KSgIAzJ07Fzt37kRBQQHy8/Nt+8PCwvDZZ5/h8ccfx6uvvoq+ffviH//4B+Li4mxlZs6ciUuXLmHVqlUwGo0YM2YMUlNTGw36lpPFIncNiIhIaVSq1nWFyW3SpEnYvn07/vOf/8DV1RXDhg3DpEmTkJGRgd9++032ViVAYWFp0qRJcLQsVFOrc0+aNAlHjx51eNzExEQkJiZeb/U6BbvhiIjImVnHLb3yyiu2YDRp0iSsW7cOv/32G5YuXSpzDZ18zJIzYDccERE5M19fX4waNQq7d+/GpEmTAAATJkzAkSNH8N///teuZenEiRPIzc1FcXExSktLkZubi9zc3E6vo6JalnoiSWI3HBERObeJEyciNzfXFpb8/PwQERGBwsJCDB061FZu6tSpOHfunO312LFjAcBhr1NHYMtSN8eWJSIicnabNm2CJEkYNmyYbVtubi4KCgrsyp09exaSJDV6dDaGpW6OLUtERETyYljq5jjAm4iISF4MSwrAsERERCQfhiUFMJvlrgERESlBV4zfUZqO+DdhWFIAhiUiInLE1dUVAFAp151zuzHrv4n136g9uHSAAjAsERGRIxqNBnq9HkVFRQAAT09PqFQqmWslL0mSUFlZiaKiIuj1emg0mnYfi2FJATgbjoiIWmIwGADAFphI0Ov1tn+b9mJYUgCGJSIiaolKpUJwcDACAwNRW1srd3W6BVdX1+tqUbJiWFIAdsMREVFraTSaDgkIVI8DvBWALUtERETyYVhSALYsERERyYdhSQHYskRERCQfhiUFYMsSERGRfBiWFIAtS0RERPJhWFIASeL94YiIiOTCsKQAFgvDEhERkVwYlhRAktgVR0REJBeGJQVgNxwREZF8GJYUgGGJiIhIPgxLCsBuOCIiIvkwLCkAW5aIiIjkw7CkAGxZIiIikg/DkgJw6QAiIiL5MCwpALvhiIiI5MOwpBDshiMiIpIHw5ICsBuOiIhIPgxLCsAB3kRERPJhWFIAjlkiIiKSj+LC0pYtWzBgwAC4u7sjOjoahw4darbspEmToFKpGj2mTZtmKzNv3rxG++Pj47viVFqNYYmIiEg+LnJXoC3effddJCUlYevWrYiOjsamTZsQFxeHvLw8BAYGNir/4Ycfoqamxvb6ypUrGD16NO655x67cvHx8dixY4fttVar7byTaAd2wxEREclHUS1LGzduREJCAubPn4+IiAhs3boVnp6e2L59e5Pl/fz8YDAYbI+0tDR4eno2CktardaunK+vb1ecTqtxgDcREZF8FBOWampqkJOTg9jYWNs2tVqN2NhYZGVlteoY27Ztw6xZs+Dl5WW3PSMjA4GBgRg6dCgeeeQRXLlyxeFxqqurYTKZ7B6djWGJiIhIHooJS5cvX4bZbEZQUJDd9qCgIBiNxhbff+jQIfz000946KGH7LbHx8fjrbfeQnp6OtavX4/MzEzceeedMJvNzR4rJSUFOp3O9ggNDW3fSbUSu+GIiIjko6gxS9dj27ZtGDlyJMaNG2e3fdasWbbnI0eOxKhRozBo0CBkZGRg8uTJTR4rOTkZSUlJttcmk6lTAxO74YiIiOSjmJYlf39/aDQaFBYW2m0vLCyEwWBw+N6Kigq88847WLBgQYvfZ+DAgfD398epU6eaLaPVauHj42P36ExsWSIiIpKPYsKSm5sbIiMjkZ6ebttmsViQnp6OmJgYh+99//33UV1djQceeKDF73PhwgVcuXIFwcHB113njsKlA4iIiOSjmLAEAElJSfj73/+OXbt24eTJk3jkkUdQUVGB+fPnAwDmzJmD5OTkRu/btm0bpk+fjt69e9ttLy8vxxNPPIHvv/8eZ8+eRXp6Ou666y4MHjwYcXFxXXJOrcGwREREJB9FjVmaOXMmLl26hFWrVsFoNGLMmDFITU21DfrOz8+HWm2f//Ly8nDgwAF8+eWXjY6n0Whw7Ngx7Nq1CyUlJQgJCcGUKVOwdu3abrXWErvhiIiI5KOSJLZZXC+TyQSdTofS0tIOH7+UmQlcuADExAADB3booYmIiHq01l6/FdUN11OxZYmIiEg+DEsKwbBEREQkD4YlhXCwRiYRERF1IoYlhWBYIiIikgfDkgKoVOyGIyIikgvDkkKwZYmIiEgeDEsKoFazZYmIiEguDEsKoFKxZYmIiEguDEsKwDFLRERE8mFYUgCGJSIiIvkwLCkAu+GIiIjkw7CkAGxZIiIikg/DkgJYwxJveUxERNT1GJYUQKUSQYlhiYiIqOsxLCmAWi2CErviiIiIuh7DkgKwZYmIiEg+DEsKwZYlIiIieTAsKYD1didsWSIiIup6DEsKoFKJrwxLREREXY9hSQGsY5bYDUdERNT1GJYUgOssERERyYdhSQHYskRERCQfhiUFUP//T4ktS0RERF2PYUkh2A1HREQkD4YlBWA3HBERkXwYlhTAersTtiwRERF1PYYlBeDtToiIiOTDsKQQ1uUDiIiIqGsxLCkEW5aIiIjkwbCkIGxZIiKinkaSgNOn5a2Di7zfntqCLUtEROTszp4FDhwADh0CcnKAH38EysqAc+eAfv3kqZPiWpa2bNmCAQMGwN3dHdHR0Th06FCzZXfu3AmVSmX3cHd3tysjSRJWrVqF4OBgeHh4IDY2Fr/88ktnn0a7MCwREZEzuXwZ+Phj4KmngClTgMBAICwMePBBYPNm4OBBEZRcXIC8PPnqqaiWpXfffRdJSUnYunUroqOjsWnTJsTFxSEvLw+BgYFNvsfHxwd5Df6FVSqV3f4XX3wRr732Gnbt2oWwsDCsXLkScXFxOHHiRKNgJScO8CYiIiWrrBStRVlZwOHDwNGjohXpWmo10L8/EB4ODBwIDB4MjBwJjB/f5VW2UVRY2rhxIxISEjB//nwAwNatW/HZZ59h+/btWLFiRZPvUalUMBgMTe6TJAmbNm3C008/jbvuugsA8NZbbyEoKAj79u3DrFmzOudE2oEDvImISCksFtF9duCACEY5OaJlqLa2cdmgIGDIEBGMhgwBIiIAX19ArxcPT0/Ay0t8lYtiwlJNTQ1ycnKQnJxs26ZWqxEbG4usrKxm31deXo7+/fvDYrHgxhtvxAsvvIAbbrgBAHDmzBkYjUbExsbayut0OkRHRyMrK6vZsFRdXY3q6mrba5PJdL2n1ypsWSIiou6otBT49lsRjrKzRatRaWnjcj4+IhANGiRajCIiAINBhCJf3/pQ5OlZf1/U7kAxYeny5cswm80ICgqy2x4UFISff/65yfcMHToU27dvx6hRo1BaWoqXXnoJN998M44fP46+ffvCaDTajnHtMa37mpKSkoJnn332Os+o7diyREREcrNYRCvRN9+IMUWHDgH//W/jP+jd3EQgsj6GDgUGDBCBqXdvwNu7Phy5dPM00s2rd31iYmIQExNje33zzTdj+PDh+Nvf/oa1a9e2+7jJyclISkqyvTaZTAgNDb2uurZGXV2nfwsiIiI7lZUiFH37LfD996JL7cqVxuX8/YHhw8VYo2HDxMMajHQ6EYy8vACttuvP4XopJiz5+/tDo9GgsLDQbnthYWGzY5Ku5erqirFjx+LUqVMAYHtfYWEhgoOD7Y45ZsyYZo+j1Wqh7eJPW6MBGvT8ERERdYrCQmD/fiAzU4SjEycajzXSaERr0ZAh9eOM+vUT4cjfX4Qib2/Aw6N7dae1l2LCkpubGyIjI5Geno7p06cDACwWC9LT05GYmNiqY5jNZvz444+YOnUqACAsLAwGgwHp6em2cGQymZCdnY1HHnmkM06j3VxcGJaIiKjjnTkjwpG1W+306cbDPnQ6+1aj4cMBP7/6ViNrl5qbmzzn0NkUE5YAICkpCXPnzkVUVBTGjRuHTZs2oaKiwjY7bs6cOejTpw9SUlIAAGvWrMH48eMxePBglJSUYMOGDTh37hweeughAGKm3JIlS/Dcc88hPDzctnRASEiILZB1F2xZIiKi62Udb/TVV/XdaufPNy7Xty9www31rUbWsUYBASIYWcORM7QatYaiwtLMmTNx6dIlrFq1CkajEWPGjEFqaqptgHZ+fj7UDT653377DQkJCTAajfD19UVkZCQOHjyIiIgIW5knn3wSFRUVWLhwIUpKSnDrrbciNTW1W62xBNSHJUkSay4RERG1xGIBjh0DvvyyflXsa0azQKUS0/YjIsQg7BEjgD59RIuRv399OPLwkOccugOVJHGO1fUymUzQ6XQoLS2Fj49Phx47M1MM7HZxEX3GEycCrq4d+i2IiMhJSBLw888iHH39tehWu3TJvoxGI1qMhg8XXWojR4oWo969xfT9Xr1EOHLWLrWGWnv9VlTLUk+m0QBVVSI4MSwREREgwtGZMyIc7d8PfPcdcPGifRlXVxGMIiLE1xtuEKGod28x7sjbWwSk7j59X078p1EIjUYEpdrant0USkTU0128CHzxRX04OnPGfr9aLVqMRowQwWjECBGO/PzEw9pyxHDUevynUgiNBjCbudYSEVFPU1EhBmR/8YXoWsvLs5+tplKJFbFHjhStRw271RiOOgb/6RRCoxED9RiWiIicm9kM5OYC//43kJYm7q1WVWVfpn//+nA0YgQQElLftebjw261jsZ/SoVhWCIicj4XLwKffy5aj779Frj2jlu+vsDYsSIgjRolpvbrdEBgYH046gkDsuXCsKQwDEtERMpXXQ1kZIiAtH8/cPy4/b3VXF1Fi9GoUeIxdKgIR0FBIhz5+HD8aldiWFIYhiUiImU6fx74+GPRvfbtt0BZmf3+/v2BMWPqW4969xZjj/z8RDjy9uY6e3JhWFIQtZqreBMRKUVdHZCVJQLSl18CP/1k33rUq1d919ro0eLeanq9fesRxx11D/wYFMS61hIREXVPV64An34KfPaZmLl2+bL9/kGDgBtvFOFoxAgRjhq2Hnl6ylJtagHDkoLwZrpERN2LJAGnTgEffCBC0uHDYj08K61WhKMxY8TX0FAxWDswUIxB8vERfwhT98awpCC8mS4RkfzMZuCHH4B//Uu0IJ04Yb8/OBi46SYx7mj0aNFqFBgoxiDpdGw9UiKGJQWxruJtvVccERF1jepqID0d2LsXSE0FLlyw3x8RIQLSjTeKrjadDjAYRCuSjw9vU6V0vOQqiIsLcPUqwxIRUVcoKwM++gjYt09M7//tt/p9rq5icHZkpHj07StakIKCxDgkb28xKYecAy+5CsJbnhARda6SEhGQPvhABKTKyvp93t7AuHGi9WjsWBGM/P3FAG29nt1rzoxhSUEadsMREVHHKC4WrUcffABkZtoHJH9/4OabRUC64QbRemQwiPFHer0YwE3Oj2FJQdiyRETUMYqLxQDtf/0L+OYbMcTBKjAQuOUWICpKBCS93j4gcRhEz8OPXEGsK7c2nJZKREStU1IiwtG774qA1HB2cVCQCEg33SQGa+v19Ten1es5vb+nY1hSILYsERG1ztWrYnr/nj1iFe2Kivp9BgNw660iIA0fLkJRcLDoatPpGJCoHsOSAjEsERE1r7ZW3KT27bdFULpypX5fYCAwYYIYqG29OW2fPiIg6fWcwUZNY1hSGJWK3XBERNeyWICcHOD//k8M1j5/vn5fr14iIMXEiDFIvr6ii826SCRbkKglDEsKw/vDERHVO3MG2LkTeO894Oef67drtWIWW0yMmObv5ycCkr8/B2lT2/G/i8Lw/nBE1NOVlYlp/m+9BXz7rZglDIg/JqOixEDtyEjR5RYUJL76+XEVbWo/hiWFcXERa4BIUv3sOCIiZ2c2AwcOADt2iEUjS0rq9w0bBvzud2KgtrX1yDpQ291dtiqTE2FYUhitVoSlqirAw0Pu2hARdR5JAs6dE91s//wn8N//1u/z8wMmTxZdbeHh9d1svXuLlbaJOhLDksK4uYm/qBiWiMhZXb0KvP8+sGuXWA/JOgPYxUWEowkTgDFjREDq21e0JOl0nMlGnYdhSWFcXERzNAd5E5EzkSTgp5+AN98Ui0ZeulS/b8gQ4PbbgfHjxdpIQUH1K2pzHBJ1BYYlhWJYIiJnUFkpZrL94x9AVpZYAgAQLUWxsWKw9uDBIhj16SNakby85K0z9TwMSwrk4gKYTHLXgoiofSQJOH4c+NvfGrcijR0rQlJUlAhIwcGiJcnXl91sJB+GJQXSakVY4ow4IlKSqirgww+BrVuB776zb0WaMgWYOBEYMECEJOtYJM5mo+6AYUmB3NzEL52aGhGciIi6s/PnxViknTuBCxfqt1tbkSIjRTDq00e0Iun1/EOQuheGJQXSaoHychGYGJaIqDsym4GDB4HXXwc++UTMcAPEeKP4eLEuUr9+QEAAW5Go+2NYUiBXV3F/uKoq0XxNRNRdVFYCe/aIrracnPrtYWHAtGliRltAgGhFMhg4FomUQXH/Rbds2YIBAwbA3d0d0dHROHToULNl//73v+O2226Dr68vfH19ERsb26j8vHnzoFKp7B7x8fGdfRrXTZLq/1IjIpLb+fPAihXAwIFAQoIISmq1WBNp7Vpg40bg3nvF7LZbbxU3tO3dm0GJlEFRLUvvvvsukpKSsHXrVkRHR2PTpk2Ii4tDXl4eAgMDG5XPyMjAfffdh5tvvhnu7u5Yv349pkyZguPHj6NPnz62cvHx8dixY4fttVYBfVsuLqIrjohILhYLcOyYCEIffFD/B5xOB0ydCkyaVN+CZJ32zxvYkhKpJEmS5K5Ea0VHR+Omm27C66+/DgCwWCwIDQ3Fo48+ihUrVrT4frPZDF9fX7z++uuYM2cOANGyVFJSgn379rW7XiaTCTqdDqWlpfDx8Wn3cZqSmSlWr9Xr7bdfvgz06iVWsyUi6kp1dUBamghJ+/fXz2oLCwPuugsYN65+wHZIiAhPHLBN3VFrr9+Kyfg1NTXIyclBcnKybZtarUZsbCyysrJadYzKykrU1tbCz8/PbntGRgYCAwPh6+uL22+/Hc899xx69+7d7HGqq6tRXV1te22SYdEjrVb8FVdbyxVsiahrXL0K7N4NbN4sWpSsoqKA3/8eGDlSdK317w8EBgKenvLVlagjKSYsXb58GWazGUFBQXbbg4KC8PPPP7fqGMuXL0dISAhiY2Nt2+Lj43H33XcjLCwMp0+fxlNPPYU777wTWVlZ0Gg0TR4nJSUFzz77bPtPpgO4u4vWpYqKxq1OREQd6bffgL/+VQzatk7912jEtP/4eDFOKTCwfnYb/4AjZ6OYsHS91q1bh3feeQcZGRlwbzA/ddasWbbnI0eOxKhRozBo0CBkZGRg8uTJTR4rOTkZSUlJttcmkwmhoaGdV/kmWGfEmUwMS0TUOQoKgFdfFbciuXJFbOvVS7QiTZ4sVtcODhZT/zlYm5yZYsKSv78/NBoNCgsL7bYXFhbCYDA4fO9LL72EdevW4auvvsKoUaMclh04cCD8/f1x6tSpZsOSVqvtFoPA3dxE61K/fnLXhIicydmzwMsvA2+9VX9rpYAAYMYMMZPN318EpD59+Mca9QyKCUtubm6IjIxEeno6pk+fDkAM8E5PT0diYmKz73vxxRfx/PPP44svvkBUVFSL3+fChQu4cuUKgoODO6rqncbLCyguFit5u7nJXRsiUjJJAvLygPXrxf3arDPb+vYF/vhHsT6Sv7/44yw4GPD2lre+RF1JMWEJAJKSkjB37lxERUVh3Lhx2LRpEyoqKjB//nwAwJw5c9CnTx+kpKQAANavX49Vq1Zhz549GDBgAIxGIwDA29sb3t7eKC8vx7PPPosZM2bAYDDg9OnTePLJJzF48GDExcXJdp6t5ekpmslNJvFLjIiorSQJ+PFH4PnngX37xB9fADBokAhJkZGiVal/f7EEgIeHrNUlkoWiwtLMmTNx6dIlrFq1CkajEWPGjEFqaqpt0Hd+fj7UDTrN33jjDdTU1OCPf/yj3XFWr16NZ555BhqNBseOHcOuXbtQUlKCkJAQTJkyBWvXru0W3Wwt0WjElF2GJSJqK0kCfvpJLBi5d69YDgAQi0X+8Y9iZpu/v1gOICiIt1aink1R6yx1V3Kss2R16ZLYN358h35bInJSkgQcPy5akv71LzFRBABuvBG45x5gyBAxs23AABGSOLONnJnTrbNETfPyAkpLxf2YuKYJETVHkoCffwZeeAF477367rYxY4CZM4Fhw+xbkrjSNlE9/jgonIeHmBFXWsqwRERN++UXEZL++U/Aup7uyJHArFnA8OH2LUkMSUSN8cdC4VQq0Ux+8aIYfMlbChCR1fnzwHPPiSUAqqrEthtuECEpIkKEJLYkEbWMPx5OQK8HCgvFQG+dTu7aEJHciouBl14CtmypXydp2DDgvvuAESPE7LawMPEHFkMSUcv4Y+IEPDzEQO+iIoYlop6sogJ4803gxReB/79SCvr1A+bOBUaPFmOSBg4UIYkDt4laj2HJSfTqJZrc+/fnApVEPU1NjRi0vWaNGJ8EiNuPPPggEBNTH5JCQvj7gag9GJachI+PGLd0+bL4hUhEzs9sBr76Cnj6aeCHH8Q2Ly8xu23y5Prutr59xc23iah9GJachEYjmtXPnRODNjkOgch5WVfdXrEC+Pxzsc3FBfjf/wWmTRMDtvv1Ey3NXl7y1pXIGfCS6kT8/cXtTwoKgNBQuWtDRJ2hsFB0t23bVr8MwJQpYtXt4GDRijRgAG9wS9SRGJaciIuLGOx9+rRofmezO5HzqKoCtm4FUlLEZA4AGDUK+NOfRDgyGMS4JH9/LiFC1NEYlpyMn58Y6H3hAjB4sNy1IaLrZTYDX3whutx+/FFsCwoCHnoIGDtWhKPBg0VY0mjkrSuRs2JYcjJqNeDrC5w6Jb727i13jYiovf77X+DJJ4GPPxbjlNzdgfvvB+64Q/xsDxokut04w42oczEsOSEfHzFu6eRJIDJSdM0RkXKUl4tFJTduBMrKxLa4ODHLLTBQDN4OCwO8veWtJ1FPwbDkpIKCRFdcXp64BxSb54m6P7NZzG5btkz87AJAeDjw5z+LVqTgYPHV31/eehL1NAxLTkqtFoHpzBnRdD9kiNhGRN3T6dOiy23vXtHl5uUlBm/fdlv9uKSQEP7hQyQHhiUnptWKX7I//yxmxwwZwlkyRN1NVRWwebO44a31Pm7x8fVdbgMGiC43dqcTyYdhycl5eoq/Uk+eFF8HD+aClUTdgSQB2dnA4sXA4cNi28CBwF/+Un9rkvBwMcOViOTFy2YP4OUlWpROngSuXgWGD+caTERyKi0Fnn0W2LJF3NfNzQ2YP1/cosTPT4SkPn3Y5UbUXTAs9RCenmIdlrNnxZ3Jhw4VC1cSUdexrpn02GNieQ9AzFhduFAM3h4wQLQq8RYlRN0Lw1IP4uYm1mQpKgIOHRKzagYMYCsTUVcoLASeeALYvRuwWMQSH488Atx0k5iMER4uxihxXCFR98Ow1MOo1aKFqbwcOHFCrMdknZLs6ip37Yicj9kMfPSRaE26cEFsu+MO4IEHROvuwIHijxatVtZqEpEDDEs9lLe3aOovLgZycsQ4ibAw8Rcuf2kTdYzCQrFm0p49ojXJ3x9YskSMG7QO4OYq+0TdH8NSD6ZSiV/Uej3w228iNPn4iK66oCBAp2OXAFF7mM3Ap5+KmW75+WLb1KniViXWNZP69WNrLpFSMCwRNBrxC9xiEeu8nDghBp/6+Ym/fn19gV69GJyIWuPyZTE26a23xM9U796iNemGG8TP09Ch4g8UIlIOhiWyUavFL3G9HqiuBkpKAKNRdMvp9WLwqU4nWp/YVUdkz2IB0tPFzLazZ8W2KVOABx8UY5MGDwb69+c6Z0RKxB9bapJWK8IRIFYYLisTs+hUKrEMgU4nWqO8vMTD05O3U6Geq7wcWLMG2LQJqK0VrbLXtib5+spdSyJqL4YlapG7e/3yAmazWNjyyhXg4kURnrRasd/XV7Q6eXiIh7s7W6DI+R07Bjz0UP0q3LfcIm58ax2bNGAAxyYRKR3DErWJRiNm0nl7i9eSJLrsqqrEtOjaWrHd1VWs6+TuLgKUt7cITtaHm5t4sDWKlKq2Fvj734EVK0TLq1YLLFoE3HyzWJ5j6FARmIhI+RiW6LqoVPYtT1a1tfUhymQC6upEsFKpxJgNV1fx8PAQXXieniI8WbdbH9ayHFxO3cnFi+Iebh99JF4PGQIsXSrWKxs0SDzYqkrkPBiWqFNYw861JEkEp5oaEahKS8XsIbNZ7ANEMNJoRFBycRHPG7ZKabX1QcparuHXhg+2XFFHkiQxiHvBgvolAe6/H7jrLjGIe9gw0arEcE/kXBiWqEupVM0HKStJEuGprq7+UV4uglVdnZh1ZC1npVbXh6OGQUmjabq1yhqs1Gr79177vKmHSsUQ1hNVVwMbNgBr14qw7+8PLF9ef9ugIUN4TzciZ9WusDRnzhz87ne/w4QJEzBo0KCOrpNDW7ZswYYNG2A0GjF69Ghs3rwZ48aNa7b8+++/j5UrV+Ls2bMIDw/H+vXrMXXqVNt+SZKwevVq/P3vf0dJSQluueUWvPHGGwgPD++K06EmWLvq2jLF2mIRActsFs+tr2tqxIB06zbro2ErVsPn1q8NQ1HD5w23NQxl1vpa695wn0rV+P3Xbrv2uaNtrXlQx7pwQSwJ8Pnn4nV0NJCYKALTkCFiSQCNRt46ElHnaVdYcnNzQ0pKChYsWIA+ffpg4sSJmDRpEiZOnNipIePdd99FUlIStm7diujoaGzatAlxcXHIy8tDoHWeewMHDx7Efffdh5SUFPz+97/Hnj17MH36dBw5cgQjRowAALz44ot47bXXsGvXLoSFhWHlypWIi4vDiRMn4M47zCqGNcBc76wjSaoPUw2D1bVfrV2J1tcN91mfW4/XMJBZNRXSGp5Lw9DTmmDU1PZrw13Df6drQ1hzobCp79GW563d35lf27qv4XNJAjIygD/9SXS7qdWiC27KFNHdNmwYB3ET9QQqSbr213jr/frrr/jmm2+QmZmJzMxM/Pe//0VwcDAuWO8W2cGio6Nx00034fXXXwcAWCwWhIaG4tFHH8WKFSsalZ85cyYqKirw6aef2raNHz8eY8aMwdatWyFJEkJCQrB06VIsW7YMAFBaWoqgoCDs3LkTs2bNalW9TCYTdDodSktL4ePj0wFnKvz6K3DggGgh6cDDUjdlDVrXPqz7mirTcFtryjX3tannDcs0DA/W1y09t3K07dqg6ChoXVvG0ftaOmZrjgUAR48Cr70muuD8/YHk5Pput6FDG09sICJlae31+7rGLPn6+qJ3797w9fWFXq+Hi4sLAgICrueQzaqpqUFOTg6Sk5Nt29RqNWJjY5GVldXke7KyspCUlGS3LS4uDvv27QMAnDlzBkajEbGxsbb9Op0O0dHRyMrKajYsVVdXo7q62vbaZDK197QcmjwZyMvrlEMTURtER4tlAQIDRWtSaCjHrRH1JO0KS0899RQyMjJw9OhRDB8+HBMnTsSKFSswYcIE+HbSMrWXL1+G2WxGUFCQ3fagoCD8/PPPTb7HaDQ2Wd5oNNr2W7c1V6YpKSkpePbZZ9t8Dm1lXY8I4DgUIjm4uwP33gvccYfodouIEPd6I6KepV1had26dQgICMDq1atx9913Y8iQIR1dr24tOTnZrsXKZDIhNDS0w7/Pf/4DZGaK8TG88SZ1d466BRuWaa47sKltjo7T0nsdlXdUtuF2i0V0wQ0YIFqUPDyaLk9Ezq1dYeno0aPIzMxERkYGXn75Zbi5udkGeU+aNKlTwpO/vz80Gg0KCwvtthcWFsJgMDT5HoPB4LC89WthYSGCg4PtyowZM6bZumi1Wmi54hy1U1ODwFvzuqmxRe3Z1rAeVteOSbJua+51U8+vHQvUngHXbS3T1HtaGkRufW7tRmtuFqF1f+/enO1G1NO1KyyNHj0ao0ePxuLFiwEA//nPf/DKK69g0aJFsFgsMJvNHVpJQMzAi4yMRHp6OqZPnw5ADPBOT09HYmJik++JiYlBeno6lixZYtuWlpaGmJgYAEBYWBgMBgPS09Nt4chkMiE7OxuPPPJIh58DdU8NZ71ZF8dsaUacNcg0fN5UeGjua1sfQONZbE2t/dTctmuftzSjrrWvW7O9NdvaU8aqtaGKiKi92hWWJEnC0aNHkZGRgYyMDBw4cAAmkwmjRo3CxIkTO7qONklJSZg7dy6ioqIwbtw4bNq0CRUVFZg/fz4Asf5Tnz59kJKSAgB47LHHMHHiRLz88suYNm0a3nnnHfzwww948803AQAqlQpLlizBc889h/DwcNvSASEhIbZARt2XNdxYA4t1wUrrGksNv1qDTlOtJUDzQaNhyGi4flJTz1uafu9oX1Ovm9vGEEBE1LXaFZb8/PxQXl6O0aNHY+LEiUhISMBtt90GfScPrJk5cyYuXbqEVatWwWg0YsyYMUhNTbUN0M7Pz4e6wRSVm2++GXv27MHTTz+Np556CuHh4di3b59tjSUAePLJJ1FRUYGFCxeipKQEt956K1JTU7nGkkysK3Y3XMHbuthkXV19OWvoabgoZMNQ4+FRfzNf68rd194GxdGq3de+ZkAhIuq52rXO0meffYbbbrutQ9cUUrLOWmcJcK4B3haLuB/ctQ9rNxbQ+J5wbm71N+pteE+4hmWaujccERFRSzp1naVp06a1u2Lk/GprxQyimpr6G+YConWmYWuPTifupeXhYd8C1PDBtWyIiEhuvJEuXRfrvdeqqsRzSRIhR6sVrUGBgUCvXvWv3dzqW4jYtUVERErAsERtUlUFVFSIgGQNRh4eQFAQ4Osrnlsfbm4MREREpHwMS+SQ2SzCUVmZeO7uDnh7A/36ifvVeXmJhwv/JxERkZPiJY4akSQRjqy3vPP2Fovy+fuLgOTtzRYjIiLqORiWyKamBigpEV1tvXoBAweKMUd6vRhnRERE1BMxLBGqq4HLl8XzgACgb1/xlUtNERERMSz1aHV1wKVLotutTx8gNFR0tXG6PhERUT2GpR6quBgoLweCg0V3W0AAxyERERE1hWGph6mtBYxGMUg7MlK0KHHFayIiouYxLPUglZWi261/fyA8XMxsIyIiIscYlnqIkhLR7RYRAQwezHWRiIiIWouXzB7AuhzAmDFiMUmOTSIiImo9zntycqWlovtt5EjR/cagRERE1DYMS06sokJ0vY0cKVqUiIiIqO0YlpxUXZ1YaHLIENGiRERERO3DsOSkjEaxLMCgQex6IyIiuh4MS06ouFisozR8OODqKndtiIiIlI1hycnU1QFlZWJ5AK6jREREdP0YlpzM5ctAUBAQEiJ3TYiIiJwDw5ITqa4WLUuDBrH7jYiIqKMwLDmRK1eAvn2BwEC5a0JEROQ8GJacRF0dYLEAoaGAmp8qERFRh+Fl1UmUlAC+voCfn9w1ISIici4MS06iokKs0s0b5BIREXUshiUnUF4OeHkB/v5y14SIiMj5MCw5gdJSwGAQC1ESERFRx2JYUjiLBZAkEZaIiIio4zEsKVxlJeDpCeh0cteEiIjIOTEsKVxFBaDXA+7ucteEiIjIOTEsKVx1NRAQIHctiIiInBfDkoLV1YmlAnjDXCIios6jmLBUXFyM2bNnw8fHB3q9HgsWLEB5ebnD8o8++iiGDh0KDw8P9OvXD4sXL0ZpaaldOZVK1ejxzjvvdPbpdIiKCjFeqVcvuWtCRETkvBSzhOHs2bNRUFCAtLQ01NbWYv78+Vi4cCH27NnTZPmLFy/i4sWLeOmllxAREYFz587h4YcfxsWLF/HBBx/Yld2xYwfi4+Ntr/V6fWeeSoexLkTJm+YSERF1HkWEpZMnTyI1NRWHDx9GVFQUAGDz5s2YOnUqXnrpJYSEhDR6z4gRI/Cvf/3L9nrQoEF4/vnn8cADD6Curg4uDZa61uv1MLRh7n11dTWqq6ttr00mU3tO67rV1QG9e8vyrYmIiHoMRXTDZWVlQa/X24ISAMTGxkKtViM7O7vVxyktLYWPj49dUAKARYsWwd/fH+PGjcP27dshSZLD46SkpECn09keoaGhbTuhDlBdDbi5cbwSERFRZ1NEWDIajQgMDLTb5uLiAj8/PxiNxlYd4/Lly1i7di0WLlxot33NmjV47733kJaWhhkzZuAvf/kLNm/e7PBYycnJKC0ttT3Onz/fthPqANXVYrkAT88u/9ZEREQ9iqzdcCtWrMD69esdljl58uR1fx+TyYRp06YhIiICzzzzjN2+lStX2p6PHTsWFRUV2LBhAxYvXtzs8bRaLbRa7XXX63rU1IiFKHnjXCIios4l66V26dKlmDdvnsMyAwcOhMFgQFFRkd32uro6FBcXtzjWqKysDPHx8ejVqxf27t0L1xZGQ0dHR2Pt2rWorq6WPRA5UlUFyND7R0RE1OPIGpYCAgIQ0IoVFWNiYlBSUoKcnBxERkYCAPbv3w+LxYLo6Ohm32cymRAXFwetVouPP/4Y7q1Y5jo3Nxe+vr7dOigB4n5wvHEuERFR51NEJ87w4cMRHx+PhIQEbN26FbW1tUhMTMSsWbNsM+F+/fVXTJ48GW+99RbGjRsHk8mEKVOmoLKyEm+//TZMJpNt1lpAQAA0Gg0++eQTFBYWYvz48XB3d0daWhpeeOEFLFu2TM7TbZF1/DlvcUJERNT5FBGWAGD37t1ITEzE5MmToVarMWPGDLz22mu2/bW1tcjLy0NlZSUA4MiRI7aZcoMHD7Y71pkzZzBgwAC4urpiy5YtePzxxyFJEgYPHoyNGzciISGh606sHWprxdpKDEtERESdTyW1NE+eWmQymaDT6WxLE3SkzEyxnlLDdTLLy0VgmjBBLB9AREREbdfa67cilg4ge1VVYskABiUiIqLOx7CkQDU1XIySiIioqzAsKZDZzJlwREREXYVhSaE4uJuIiKhrMCwpTG2tWLWbYYmIiKhrMCwpTHU1oNUCHh5y14SIiKhnYFhSmJoaEZY4E46IiKhrMCwpTF2daFVSqeSuCRERUc/AsKQwZjPHKxEREXUlhiWFYVgiIiLqWgxLCiNJ4r5wRERE1DUYlhTIRTG3PyYiIlI+hiUFYlgiIiLqOgxLCiJJYhYcwxIREVHXYVhSELMZ0GgYloiIiLoSw5KC1NWJsMQB3kRERF2HYUlB2LJERETU9RiWFIRhiYiIqOsxLCmI2Sy64DQauWtCRETUczAsKYjZLG6iS0RERF2HYUlB6up4qxMiIqKuxrCkIHV1bFkiIiLqagxLCiJJgJub3LUgIiLqWRiWFIYz4YiIiLoWw5LCcEFKIiKirsWwpDBsWSIiIupaDEsKYTYDajXDEhERUVdjWFIIs1kEJYYlIiKirsWwpBC81QkREZE8GJYUoq6OLUtERERyYFhSCLYsERERyUMxYam4uBizZ8+Gj48P9Ho9FixYgPLycofvmTRpElQqld3j4YcftiuTn5+PadOmwdPTE4GBgXjiiSdQV1fXmafSLnV1YkFKlUrumhAREfUsimmnmD17NgoKCpCWloba2lrMnz8fCxcuxJ49exy+LyEhAWvWrLG99vT0tD03m82YNm0aDAYDDh48iIKCAsyZMweurq544YUXOu1c2sNs5n3hiIiI5KCIsHTy5Emkpqbi8OHDiIqKAgBs3rwZU6dOxUsvvYSQkJBm3+vp6QmDwdDkvi+//BInTpzAV199haCgIIwZMwZr167F8uXL8cwzz8CtG91bxGzmfeGIiIjkoIhuuKysLOj1eltQAoDY2Fio1WpkZ2c7fO/u3bvh7++PESNGIDk5GZWVlXbHHTlyJIKCgmzb4uLiYDKZcPz48WaPWV1dDZPJZPfobBYLW5aIiIjkoIiWJaPRiMDAQLttLi4u8PPzg9FobPZ9999/P/r374+QkBAcO3YMy5cvR15eHj788EPbcRsGJQC2146Om5KSgmeffba9p9MukiQWpSQiIqKuJWtYWrFiBdavX++wzMmTJ9t9/IULF9qejxw5EsHBwZg8eTJOnz6NQYMGtfu4ycnJSEpKsr02mUwIDQ1t9/Fai2GJiIio68kalpYuXYp58+Y5LDNw4EAYDAYUFRXZba+rq0NxcXGz45GaEh0dDQA4deoUBg0aBIPBgEOHDtmVKSwsBACHx9VqtdDKMICIM+GIiIi6nqxhKSAgAAEBAS2Wi4mJQUlJCXJychAZGQkA2L9/PywWiy0AtUZubi4AIDg42Hbc559/HkVFRbZuvrS0NPj4+CAiIqKNZ9P5GJaIiIi6niI6doYPH474+HgkJCTg0KFD+O6775CYmIhZs2bZZsL9+uuvGDZsmK2l6PTp01i7di1ycnJw9uxZfPzxx5gzZw4mTJiAUaNGAQCmTJmCiIgIPPjgg/jPf/6DL774Ak8//TQWLVokS8tRS9gNR0RE1PUUc/ndvXs3hg0bhsmTJ2Pq1Km49dZb8eabb9r219bWIi8vzzbbzc3NDV999RWmTJmCYcOGYenSpZgxYwY++eQT23s0Gg0+/fRTaDQaxMTE4IEHHsCcOXPs1mXqLlQqtiwRERHJQSVJkiR3JZTOZDJBp9OhtLQUPj4+HXrszEyxendZGTB+PNCGIVpERETkQGuv34ppWSK2LBEREcmBYUlBGJaIiIi6HsOSgnCANxERUdfj5VcBJIkDvImIiOTCsKQAFosISmxZIiIi6nq8/CoAW5aIiIjkw7CkANab6DIsERERdT2GJQWwroTFbjgiIqKux8uvArBliYiISD4MSwrAMUtERETyYVhSAGtYYjccERFR1+PlVwHYskRERCQfhiUF4DpLRERE8uHlVyHYskRERCQPhiUFsLYsMSwRERF1PYYlBeAAbyIiIvnw8qsAXGeJiIhIPgxLCsDZcERERPJhWFIAa8sSERERdT1eghXAYgE0GrlrQURE1DMxLCkEwxIREZE8GJYUgN1wRERE8uElWAEsFoYlIiIiufASrACSxG44IiIiuTAsKQRbloiIiOTBS7ACcMwSERGRfHgJVggXF7lrQERE1DMxLCkEW5aIiIjkwUuwQnCANxERkTwYlhSA94UjIiKSD8OSAjAsERERyYdhSQFUKo5ZIiIikotiLsHFxcWYPXs2fHx8oNfrsWDBApSXlzdb/uzZs1CpVE0+3n//fVu5pva/8847XXFKraZWs2WJiIhILoqZkD579mwUFBQgLS0NtbW1mD9/PhYuXIg9e/Y0WT40NBQFBQV22958801s2LABd955p932HTt2ID4+3vZar9d3eP2vF1uWiIiI5KGIsHTy5Emkpqbi8OHDiIqKAgBs3rwZU6dOxUsvvYSQkJBG79FoNDAYDHbb9u7di3vvvRfe3t522/V6faOyjlRXV6O6utr22mQyteV02oxjloiIiOSjiPaKrKws6PV6W1ACgNjYWKjVamRnZ7fqGDk5OcjNzcWCBQsa7Vu0aBH8/f0xbtw4bN++HZIkOTxWSkoKdDqd7REaGtq2E2ojdsMRERHJRxFhyWg0IjAw0G6bi4sL/Pz8YDQaW3WMbdu2Yfjw4bj55pvttq9Zswbvvfce0tLSMGPGDPzlL3/B5s2bHR4rOTkZpaWltsf58+fbdkJtxAHeRERE8pG1G27FihVYv369wzInT5687u9z9epV7NmzBytXrmy0r+G2sWPHoqKiAhs2bMDixYubPZ5Wq4VWq73uerUWu+GIiIjkI2tYWrp0KebNm+ewzMCBA2EwGFBUVGS3va6uDsXFxa0aa/TBBx+gsrISc+bMabFsdHQ01q5di+rq6i4NRI6wZYmIiEg+soalgIAABAQEtFguJiYGJSUlyMnJQWRkJABg//79sFgsiI6ObvH927Ztwx/+8IdWfa/c3Fz4+vp2m6AEsGWJiIhIToqYDTd8+HDEx8cjISEBW7duRW1tLRITEzFr1izbTLhff/0VkydPxltvvYVx48bZ3nvq1Cl88803+Pe//93ouJ988gkKCwsxfvx4uLu7Iy0tDS+88AKWLVvWZefWGgxLRERE8lFEWAKA3bt3IzExEZMnT4ZarcaMGTPw2muv2fbX1tYiLy8PlZWVdu/bvn07+vbtiylTpjQ6pqurK7Zs2YLHH38ckiRh8ODB2LhxIxISEjr9fNqC3XBERETyUUktzZOnFplMJuh0OpSWlsLHx6dDj52ZCZSVAbfeCnTDtTKJiIgUq7XXb7ZXKAC74YiIiOTDsKQQ7IYjIiKSBy/BCsAVvImIiOTDsKQAHOBNREQkH16CFYBjloiIiOTDsKQADEtERETyYVhSAHbDERERyYeXYAXgAG8iIiL5MCwpBFuWiIiI5MFLsAJwzBIREZF8GJYUgN1wRERE8mFYUgC2LBEREcmHYUkBNBq5a0BERNRzMSwpAAd3ExERyYeXYQVgyxIREZF8GJYUgC1LRERE8uFlWAEYloiIiOTDy7ACsBuOiIhIPgxLCsCWJSIiIvnwMqwAbFkiIiKSD8OSAjAsERERyYdhSQHYDUdERCQfXoYVgC1LRERE8mFYUgDeF46IiEg+DEvdnFrNbjgiIiI58TLczalUbFkiIiKSE8NSN8ewREREJC+GpW5OpWI3HBERkZx4Ge7m2LJEREQkL4albo4DvImIiOTFy7ACsGWJiIhIPooJS88//zxuvvlmeHp6Qq/Xt+o9kiRh1apVCA4OhoeHB2JjY/HLL7/YlSkuLsbs2bPh4+MDvV6PBQsWoLy8vBPOoH3UaoYlIiIiOSkmLNXU1OCee+7BI4880ur3vPjii3jttdewdetWZGdnw8vLC3FxcaiqqrKVmT17No4fP460tDR8+umn+Oabb7Bw4cLOOIV24QBvIiIieakkSZLkrkRb7Ny5E0uWLEFJSYnDcpIkISQkBEuXLsWyZcsAAKWlpQgKCsLOnTsxa9YsnDx5EhERETh8+DCioqIAAKmpqZg6dSouXLiAkJCQVtXJZDJBp9OhtLQUPj4+13V+1zp2DOjfH9DpOvSwREREPV5rr99O22Zx5swZGI1GxMbG2rbpdDpER0cjKysLAJCVlQW9Xm8LSgAQGxsLtVqN7OzsZo9dXV0Nk8lk9+gsI0cyKBEREcnJacOS0WgEAAQFBdltDwoKsu0zGo0IDAy02+/i4gI/Pz9bmaakpKRAp9PZHqGhoR1c+3ocr0RERCQvWcPSihUroFKpHD5+/vlnOavYpOTkZJSWltoe58+fl7tKRERE1Elc5PzmS5cuxbx58xyWGThwYLuObTAYAACFhYUIDg62bS8sLMSYMWNsZYqKiuzeV1dXh+LiYtv7m6LVaqHVattVLyIiIlIWWcNSQEAAAgICOuXYYWFhMBgMSE9Pt4Ujk8mE7Oxs24y6mJgYlJSUICcnB5GRkQCA/fv3w2KxIDo6ulPqRURERMqimDFL+fn5yM3NRX5+PsxmM3Jzc5Gbm2u3JtKwYcOwd+9eAIBKpcKSJUvw3HPP4eOPP8aPP/6IOXPmICQkBNOnTwcADB8+HPHx8UhISMChQ4fw3XffITExEbNmzWr1TDgiIiJybrK2LLXFqlWrsGvXLtvrsWPHAgC+/vprTJo0CQCQl5eH0tJSW5knn3wSFRUVWLhwIUpKSnDrrbciNTUV7u7utjK7d+9GYmIiJk+eDLVajRkzZuC1117rmpMiIiKibk9x6yx1R525zhIRERF1jh6/zhIRERFRR2BYIiIiInKAYYmIiIjIAYYlIiIiIgcYloiIiIgcYFgiIiIicoBhiYiIiMgBxSxK2Z1Zl6oymUwy14SIiIhay3rdbmnJSYalDlBWVgYACA0NlbkmRERE1FZlZWXQ6XTN7ucK3h3AYrHg4sWL6NWrF1QqVYcd12QyITQ0FOfPn3falcGd/Ryd/fwAnqMzcPbzA3iOzqAzzk+SJJSVlSEkJARqdfMjk9iy1AHUajX69u3bacf38fFxyv/4DTn7OTr7+QE8R2fg7OcH8BydQUefn6MWJSsO8CYiIiJygGGJiIiIyAGGpW5Mq9Vi9erV0Gq1clel0zj7OTr7+QE8R2fg7OcH8BydgZznxwHeRERERA6wZYmIiIjIAYYlIiIiIgcYloiIiIgcYFgiIiIicoBhqRvbsmULBgwYAHd3d0RHR+PQoUNyV6ldUlJScNNNN6FXr14IDAzE9OnTkZeXZ1dm0qRJUKlUdo+HH35Yphq33TPPPNOo/sOGDbPtr6qqwqJFi9C7d294e3tjxowZKCwslLHGbTdgwIBG56hSqbBo0SIAyvsMv/nmG/zP//wPQkJCoFKpsG/fPrv9kiRh1apVCA4OhoeHB2JjY/HLL7/YlSkuLsbs2bPh4+MDvV6PBQsWoLy8vAvPwjFH51hbW4vly5dj5MiR8PLyQkhICObMmYOLFy/aHaOpz33dunVdfCZNa+kznDdvXqO6x8fH25VR8mcIoMmfSZVKhQ0bNtjKdOfPsDXXh9b8/szPz8e0adPg6emJwMBAPPHEE6irq+uwejIsdVPvvvsukpKSsHr1ahw5cgSjR49GXFwcioqK5K5am2VmZmLRokX4/vvvkZaWhtraWkyZMgUVFRV25RISElBQUGB7vPjiizLVuH1uuOEGu/ofOHDAtu/xxx/HJ598gvfffx+ZmZm4ePEi7r77bhlr23aHDx+2O7+0tDQAwD333GMro6TPsKKiAqNHj8aWLVua3P/iiy/itddew9atW5GdnQ0vLy/ExcWhqqrKVmb27Nk4fvw40tLS8Omnn+Kbb77BwoULu+oUWuToHCsrK3HkyBGsXLkSR44cwYcffoi8vDz84Q9/aFR2zZo1dp/ro48+2hXVb1FLnyEAxMfH29X9n//8p91+JX+GAOzOraCgANu3b4dKpcKMGTPsynXXz7A114eWfn+azWZMmzYNNTU1OHjwIHbt2oWdO3di1apVHVdRibqlcePGSYsWLbK9NpvNUkhIiJSSkiJjrTpGUVGRBEDKzMy0bZs4caL02GOPyVep67R69Wpp9OjRTe4rKSmRXF1dpffff9+27eTJkxIAKSsrq4tq2PEee+wxadCgQZLFYpEkSdmfIQBp7969ttcWi0UyGAzShg0bbNtKSkokrVYr/fOf/5QkSZJOnDghAZAOHz5sK/P5559LKpVK+vXXX7us7q117Tk25dChQxIA6dy5c7Zt/fv3l1555ZXOrVwHaOr85s6dK911113NvscZP8O77rpLuv322+22KeUzlKTG14fW/P7897//LanVasloNNrKvPHGG5KPj49UXV3dIfViy1I3VFNTg5ycHMTGxtq2qdVqxMbGIisrS8aadYzS0lIAgJ+fn9323bt3w9/fHyNGjEBycjIqKyvlqF67/fLLLwgJCcHAgQMxe/Zs5OfnAwBycnJQW1tr93kOGzYM/fr1U+znWVNTg7fffht/+tOf7G4erfTP0OrMmTMwGo12n5lOp0N0dLTtM8vKyoJer0dUVJStTGxsLNRqNbKzs7u8zh2htLQUKpUKer3ebvu6devQu3dvjB07Fhs2bOjQ7o3OlpGRgcDAQAwdOhSPPPIIrly5YtvnbJ9hYWEhPvvsMyxYsKDRPqV8htdeH1rz+zMrKwsjR45EUFCQrUxcXBxMJhOOHz/eIfXijXS7ocuXL8NsNtt98AAQFBSEn3/+WaZadQyLxYIlS5bglltuwYgRI2zb77//fvTv3x8hISE4duwYli9fjry8PHz44Ycy1rb1oqOjsXPnTgwdOhQFBQV49tlncdttt+Gnn36C0WiEm5tbowtQUFAQjEajPBW+Tvv27UNJSQnmzZtn26b0z7Ah6+fS1M+gdZ/RaERgYKDdfhcXF/j5+Snyc62qqsLy5ctx33332d2kdPHixbjxxhvh5+eHgwcPIjk5GQUFBdi4caOMtW2d+Ph43H333QgLC8Pp06fx1FNP4c4770RWVhY0Go3TfYa7du1Cr169GnXxK+UzbOr60Jrfn0ajscmfVeu+jsCwRF1q0aJF+Omnn+zG8wCwGyMwcuRIBAcHY/LkyTh9+jQGDRrU1dVsszvvvNP2fNSoUYiOjkb//v3x3nvvwcPDQ8aadY5t27bhzjvvREhIiG2b0j/Dnqy2thb33nsvJEnCG2+8YbcvKSnJ9nzUqFFwc3PDn//8Z6SkpHT722rMmjXL9nzkyJEYNWoUBg0ahIyMDEyePFnGmnWO7du3Y/bs2XB3d7fbrpTPsLnrQ3fAbrhuyN/fHxqNptFo/8LCQhgMBplqdf0SExPx6aef4uuvv0bfvn0dlo2OjgYAnDp1qiuq1uH0ej2GDBmCU6dOwWAwoKamBiUlJXZllPp5njt3Dl999RUeeughh+WU/BlaPxdHP4MGg6HRhIu6ujoUFxcr6nO1BqVz584hLS3NrlWpKdHR0airq8PZs2e7poIdaODAgfD397f9n3SWzxAAvv32W+Tl5bX4cwl0z8+wuetDa35/GgyGJn9Wrfs6AsNSN+Tm5obIyEikp6fbtlksFqSnpyMmJkbGmrWPJElITEzE3r17sX//foSFhbX4ntzcXABAcHBwJ9euc5SXl+P06dMIDg5GZGQkXF1d7T7PvLw85OfnK/Lz3LFjBwIDAzFt2jSH5ZT8GYaFhcFgMNh9ZiaTCdnZ2bbPLCYmBiUlJcjJybGV2b9/PywWiy0odnfWoPTLL7/gq6++Qu/evVt8T25uLtRqdaPuKyW4cOECrly5Yvs/6QyfodW2bdsQGRmJ0aNHt1i2O32GLV0fWvP7MyYmBj/++KNd8LUG/4iIiA6rKHVD77zzjqTVaqWdO3dKJ06ckBYuXCjp9Xq70f5K8cgjj0g6nU7KyMiQCgoKbI/KykpJkiTp1KlT0po1a6QffvhBOnPmjPTRRx9JAwcOlCZMmCBzzVtv6dKlUkZGhnTmzBnpu+++k2JjYyV/f3+pqKhIkiRJevjhh6V+/fpJ+/fvl3744QcpJiZGiomJkbnWbWc2m6V+/fpJy5cvt9uuxM+wrKxMOnr0qHT06FEJgLRx40bp6NGjtplg69atk/R6vfTRRx9Jx44dk+666y4pLCxMunr1qu0Y8fHx0tixY6Xs7GzpwIEDUnh4uHTffffJdUqNODrHmpoa6Q9/+IPUt29fKTc31+5n0zqD6ODBg9Irr7wi5ebmSqdPn5befvttKSAgQJozZ47MZyY4Or+ysjJp2bJlUlZWlnTmzBnpq6++km688UYpPDxcqqqqsh1DyZ+hVWlpqeTp6Sm98cYbjd7f3T/Dlq4PktTy78+6ujppxIgR0pQpU6Tc3FwpNTVVCggIkJKTkzusngxL3djmzZulfv36SW5ubtK4ceOk77//Xu4qtQuAJh87duyQJEmS8vPzpQkTJkh+fn6SVquVBg8eLD3xxBNSaWmpvBVvg5kzZ0rBwcGSm5ub1KdPH2nmzJnSqVOnbPuvXr0q/eUvf5F8fX0lT09P6X//93+lgoICGWvcPl988YUEQMrLy7PbrsTP8Ouvv27y/+XcuXMlSRLLB6xcuVIKCgqStFqtNHny5EbnfeXKFem+++6TvL29JR8fH2n+/PlSWVmZDGfTNEfneObMmWZ/Nr/++mtJkiQpJydHio6OlnQ6neTu7i4NHz5ceuGFF+zChpwcnV9lZaU0ZcoUKSAgQHJ1dZX69+8vJSQkNPqDU8mfodXf/vY3ycPDQyopKWn0/u7+GbZ0fZCk1v3+PHv2rHTnnXdKHh4ekr+/v7R06VKptra2w+qp+v+VJSIiIqImcMwSERERkQMMS0REREQOMCwREREROcCwREREROQAwxIRERGRAwxLRERERA4wLBERERE5wLBERERE5ADDEhEpwjPPPIMxY8a06T0qlQr79u3rlPpcr0mTJmHJkiVyV4OIWoFhiYhk05bAsGzZMrubaRIRdRUXuStAROSIJEkwm83w9vaGt7e33NXp1sxmM1QqFdRq/h1M1JH4E0VEspg3bx4yMzPx6quvQqVSQaVS4ezZs8jIyIBKpcLnn3+OyMhIaLVaHDhwoFE33OHDh3HHHXfA398fOp0OEydOxJEjR9pUh0mTJmHx4sV48skn4efnB4PBgGeeeca2/+zZs1CpVMjNzbVtKykpgUqlQkZGBgDY6vvFF19g7Nix8PDwwO23346ioiJ8/vnnGD58OHx8fHD//fejsrLS7vvX1dUhMTEROp0O/v7+WLlyJRrerrO6uhrLli1Dnz594OXlhejoaNv3BYCdO3dCr9fj448/RkREBLRaLfLz89v0b0BELWNYIiJZvPrqq4iJiUFCQgIKCgpQUFCA0NBQ2/4VK1Zg3bp1OHnyJEaNGtXo/WVlZZg7dy4OHDiA77//HuHh4Zg6dSrKysraVI9du3bBy8sL2dnZePHFF7FmzRqkpaW1+XyeeeYZvP766zh48CDOnz+Pe++9F5s2bcKePXvw2Wef4csvv8TmzZsbfW8XFxccOnQIr776KjZu3Ih//OMftv2JiYnIysrCO++8g2PHjuGee+5BfHw8fvnlF1uZyspKrF+/Hv/4xz9w/PhxBAYGtrnuROQYu+GISBY6nQ5ubm7w9PSEwWBotH/NmjW44447mn3/7bffbvf6zTffhF6vR2ZmJn7/+9+3uh6jRo3C6tWrAQDh4eF4/fXXkZ6e7vB7N+W5557DLbfcAgBYsGABkpOTcfr0aQwcOBAA8Mc//hFff/01li9fbntPaGgoXnnlFahUKgwdOhQ//vgjXnnlFSQkJCA/Px87duxAfn4+QkJCAIhxW6mpqdixYwdeeOEFAEBtbS3++te/YvTo0W2qLxG1HluWiKhbioqKcri/sLAQCQkJCA8Ph06ng4+PD8rLy9vcDXVtq1VwcDCKioraXN+GxwkKCoKnp6ctKFm3XXvc8ePHQ6VS2V7HxMTgl19+gdlsxo8//giz2YwhQ4bYxmt5e3sjMzMTp0+ftr3Hzc2tyZY3Iuo4bFkiom7Jy8vL4f65c+fiypUrePXVV9G/f39otVrExMSgpqamTd/H1dXV7rVKpYLFYgEA20DphuOIamtrWzyOSqVyeNzWKC8vh0ajQU5ODjQajd2+hgPdPTw87AIXEXU8hiUiko2bmxvMZnO73vvdd9/hr3/9K6ZOnQoAOH/+PC5fvtyR1UNAQAAAoKCgAGPHjgUAu8He1ys7O9vutXXslUajwdixY2E2m1FUVITbbrutw74nEbUdwxIRyWbAgAHIzs7G2bNn4e3tDT8/v1a/Nzw8HP/3f/+HqKgomEwmPPHEE/Dw8OjQ+nl4eGD8+PFYt24dwsLCUFRUhKeffrrDjp+fn4+kpCT8+c9/xpEjR7B582a8/PLLAIAhQ4Zg9uzZmDNnDl5++WWMHTsWly5dQnp6OkaNGoVp06Z1WD2IyDGOWSIi2SxbtgwajQYREREICAho03ijbdu24bfffsONN96IBx98EIsXL+6UmWDbt29HXV0dIiMjsWTJEjz33HMdduw5c+bg6tWrGDduHBYtWoTHHnsMCxcutO3fsWMH5syZg6VLl2Lo0KGYPn06Dh8+jH79+nVYHYioZSqpYWc8EREREdlhyxIRERGRAwxLRERERA4wLBERERE5wLBERERE5ADDEhEREZEDDEtEREREDjAsERERETnAsERERETkAMMSERERkQMMS0REREQOMCwREREROfD/AKCI+pqYfXOZAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# initialize some parameters \n",
    "n_stimuli = 2\n",
    "transition_noise_variance = 1e-5\n",
    "observation_noise_variance = 1.\n",
    "prior_covariance = 100 * np.eye(n_stimuli)\n",
    "trials = np.arange(200)\n",
    "\n",
    "# intialize the weights to zeros (no prior reward expectation)\n",
    "w = np.zeros(n_stimuli)\n",
    "\n",
    "ws = []\n",
    "covs = []\n",
    "covariance = prior_covariance \n",
    "for trial in trials:\n",
    "    stimuli = np.array([1, 0])\n",
    "    if trial < 100:\n",
    "        reward = 0\n",
    "    else:\n",
    "        reward = 1\n",
    "\n",
    "    w, covariance = kalman_filter_update(w, covariance, stimuli, reward)\n",
    "\n",
    "    ws.append(w)\n",
    "    covs.append(covariance)\n",
    "\n",
    "# now plot the uncertainty as shaded errorbars\n",
    "w1_mean = np.array(ws)[:, 0]\n",
    "w1_var = np.array(covs)[:, 0, 0]\n",
    "\n",
    "plt.plot(trials, w1_mean, color='blue')\n",
    "plt.legend(['w1', 'w2'])\n",
    "\n",
    "# Add shaded uncertainty area around w1\n",
    "plt.fill_between(trials, w1_mean - w1_var, w1_mean + w1_var, color='blue', alpha=0.2)\n",
    "\n",
    "plt.xlabel('trial number')\n",
    "plt.ylabel('w')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2e980c4e-96f7-4a55-92ba-e15fb310d03c",
   "metadata": {},
   "source": [
    "### Exercise 2.2\n",
    "In part 1, we simulated the forward blocking phenomenon using Rescorla-Wagner. Animals also show a 'backward blocking' phenomenon, where blocking occurs even if the compound phase is shown before A is paired with reward alone, as below:\n",
    "\n",
    "AB —> + \\\n",
    "A —> +  \\\n",
    "B —> ? \n",
    "\n",
    "Simulate this paradigm using both the RW and KF models, and explain what you see. Why do or don't these two models show backward blocking?\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7dbee652-767d-405f-9264-432cecd5cf66",
   "metadata": {},
   "outputs": [],
   "source": [
    "# write your code here "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ff8ce5d0-cd16-47d4-b6c3-de59d3e6d5ff",
   "metadata": {},
   "source": [
    "# 3. Temporal difference (TD) learning\n",
    "\n",
    "With the TD learning model, we can model reinforcement learning over time within the trial. This allows us to model what the TD prediction error looks like at different times within the trial, across different stages of learning. Below is an implementation of the TD update. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 183,
   "id": "a908c6fe-2913-4ca0-be5b-70db43ff4a6b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def temporal_difference_update(w, current_stimulus, next_stimulus, reward):\n",
    "    expected_reward =  np.dot(w, current_stimulus) - gamma * np.dot(w, next_stimulus) \n",
    "    prediction_error = reward - expected_reward\n",
    "    new_w = w + learning_rate * prediction_error * current_stimulus\n",
    "    return new_w, prediction_error"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f18d3a95-afbc-494e-bf06-9dfb5fd8349f",
   "metadata": {},
   "source": [
    "Let's use this model to simulate a simple conditioning paradigm where a single stimulus is reliably followed by reward later in the same trial. We can show what the prediction error looks like at different points in the trial, early during learning as well as late. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 189,
   "id": "bfd23828-b960-4529-8e60-f7512cfb5143",
   "metadata": {},
   "outputs": [],
   "source": [
    "# First, we need to represent the stimuli, as well as a representation of the stimuli over time. We'll use the \n",
    "# complete serial compound representation (see slides).\n",
    "\n",
    "def construct_stimulus(trial_length, onset, duration):\n",
    "    \"\"\"\n",
    "    Args:\n",
    "        - 'trial_length': length of trial in time steps \n",
    "        - 'onset': onset of each stimulus (time step). This is a list with length n_stimuli.\n",
    "        - 'duration': duration of each stimulus (as above) \n",
    "    \"\"\"\n",
    "    n_stimuli = len(onset)\n",
    "    stimulus = np.zeros((trial_length, n_stimuli))\n",
    "    for i in range(n_stimuli):\n",
    "        stimulus[onset[i]:onset[i]+duration[i], i] = 1\n",
    "    return stimulus\n",
    "\n",
    "\n",
    "def construct_csc(stimulus):\n",
    "    x = np.array([]).reshape(stimulus.shape[0], 0)\n",
    "    for d in range(stimulus.shape[1]):\n",
    "        diag_matrix = np.diag(stimulus[:, d])\n",
    "        x = np.hstack((x, diag_matrix))\n",
    "    return x\n",
    "\n",
    "\n",
    "# set some parameters\n",
    "learning_rate = .1 \n",
    "gamma = .99\n",
    "n_trials = 1000\n",
    "n_time_steps = 20\n",
    "\n",
    "stim = construct_stimulus(n_time_steps, [5], [20])\n",
    "x = construct_csc(stim)\n",
    "reward = construct_stimulus(n_time_steps, [15], [1])\n",
    "\n",
    "# initialize weight to zero\n",
    "d = x.shape[0]\n",
    "w = np.zeros(d)\n",
    "\n",
    "ws = []\n",
    "errors = []\n",
    "\n",
    "# learn for\n",
    "prediction_errors = []\n",
    "for trial in range(n_trials):\n",
    "    trial_errors = []\n",
    "    for time_step in range(n_time_steps-1):\n",
    "        w, delta = temporal_difference_update(w, x[time_step], x[time_step+1], reward[time_step+1])\n",
    "        trial_errors.append(delta)\n",
    "    prediction_errors.append(trial_errors)\n",
    "\n",
    "\n",
    "prediction_errors = np.array(prediction_errors)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 198,
   "id": "cda65eb0-cd80-4957-8ab6-132a3237fd6a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAn4AAAHVCAYAAABv4/bQAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABm1UlEQVR4nO3de1xUdf4/8NeZAWYAYUCRywgK4l0RDJOwdTUlwfVbUltq31qVLWvddDO2UtuEzP3lpbZHF0lb8/qrTXO3bH9dMCWxG+omkmLqQuGdAUGZ4SK3mfP7Q2d0YLgMzMwZZl7Px2MeMed8zmfe53QY33xuRxBFUQQRERERuTyZ1AEQERERkWMw8SMiIiJyE0z8iIiIiNwEEz8iIiIiN8HEj4iIiMhNMPEjIiIichNM/IiIiIjcBBM/IiIiIjfBxI+IiIjITTDxIyKXJAgCXnzxRZvUlZubC0EQkJuba5P6ACAyMhLz5s2zWX32MmnSJEyaNEnqMIjIRpj4EZFDbN26FYIgtPk6ePCg1CESEbk8D6kDICL38tJLLyEqKqrV9kGDBkkQTef8+te/xrVr1+Dl5SV1KA735ZdfSh0CEdkQEz8icqhp06Zh7NixdqnbYDCgsbERSqXSpvXKZDKb1ymFrlwfd0x2iVwZu3qJyOm8+uqrGD9+PPr06QNvb2/Ex8fjn//8Z6tygiBg4cKFeP/99zFy5EgoFApkZ2e3Krd//34IgoCPP/641b5//OMfEAQBeXl5bcZjaYzfpEmTMGrUKPz000+466674OPjg379+mHt2rVdO2kAVVVVWLx4MSIiIqBQKDBo0CCsWbMGBoPBrFx3r4+x2/27775Deno6+vbtC19fX9x33324fPmyWR0tx/gZr8WHH36I//N//g/Cw8OhVCoxZcoUFBcXt4ohKysLAwcOhLe3N8aNG4dvvvmG4waJJMQWPyJyKK1Wi4qKCrNtgiCgT58+pvdvvPEG7r33Xjz88MNobGzEjh078OCDD+LTTz/F9OnTzY796quv8OGHH2LhwoUICgpCZGRkq8+cNGkSIiIi8P777+O+++4z2/f+++8jOjoaiYmJVp/L1atXkZKSgvvvvx8zZ87EP//5TyxZsgQxMTGYNm2aVXXV1dVh4sSJuHjxIp544gn0798f33//PZYtW4bS0lK8/vrrprLdvT4FBQUAgEWLFiEwMBCZmZk4c+YMXn/9dSxcuBA7d+7sMN7Vq1dDJpPhmWeegVarxdq1a/Hwww/j0KFDpjLr16/HwoULMWHCBDz99NM4c+YMUlNTERgYiPDwcKuuDxHZiEhE5ABbtmwRAVh8KRQKs7J1dXVm7xsbG8VRo0aJkydPNtsOQJTJZOKJEydafR4AMTMz0/R+2bJlokKhEKuqqkzbysvLRQ8PD7Nyluzfv18EIO7fv9+0beLEiSIAcfv27aZtDQ0NYmhoqPjb3/623fpEURQHDBggzp071/R+5cqVoq+vr/jf//7XrNzSpUtFuVwunjt3zrStu9fH+P8iKSlJNBgMpu1PP/20KJfLza7RxIkTxYkTJ7a6FsOHDxcbGhpM29944w0RgHj8+HHTtejTp494++23i01NTaZyW7duFQGY1UlEjsOuXiJyqKysLOzdu9fs9cUXX5iV8fb2Nv189epVaLVaTJgwAfn5+a3qmzhxIkaMGNHh586ZMwcNDQ1mXaI7d+5Ec3MzHnnkkS6dS69evcyO9fLywrhx4/DLL79YXdeuXbswYcIEBAYGoqKiwvRKSkqCXq/H119/bSprq+vz+OOPQxAE0/sJEyZAr9fj7NmzHcablpZmNv5vwoQJAGA69x9++AGVlZWYP38+PDxudi49/PDDCAwM7LB+IrIPdvUSkUONGzeuw8kdn376Kf7617+ioKAADQ0Npu23JilGlmYIWzJs2DDcfvvteP/99/Hoo48CuN7Ne8cdd3R5RnF4eHirmAIDA3Hs2DGr6yoqKsKxY8fQt29fi/vLy8tNP9vq+vTv39/svTEhu3r1aofxdnSsMXlseW09PDwsdscTkWMw8SMip/LNN9/g3nvvxa9//Wu8/fbbCAsLg6enJ7Zs2YJ//OMfrcrf2vrVkTlz5uCpp57ChQsX0NDQgIMHD2LdunVdjlUul1vcLoqi1XUZDAbcfffdeO655yzuHzJkCADbXp/uxG/Lcycix2HiR0RO5V//+heUSiX27NkDhUJh2r5ly5Zu1z179mykp6fjgw8+wLVr1+Dp6YlZs2Z1u15biI6ORk1NDZKSktotZ8/rY0sDBgwAABQXF+Ouu+4ybW9ubsaZM2cwevRoqUIjcmsc40dETkUul0MQBOj1etO2M2fOYPfu3d2uOygoCNOmTcN7772H999/HykpKQgKCup2vbYwc+ZM5OXlYc+ePa32VVVVobm5GYB9r48tjR07Fn369MHGjRtNsQPXu9c705VMRPbBFj8icqgvvvgCp06darV9/PjxGDhwIKZPn47XXnsNKSkp+N///V+Ul5cjKysLgwYN6tLYuZbmzJmDBx54AACwcuXKbtdnK88++yz+/e9/43/+538wb948xMfHo7a2FsePH8c///lPnDlzBkFBQXa/Prbi5eWFF198EYsWLcLkyZMxc+ZMnDlzBlu3bkV0dLTF8YhEZH9M/IjIoTIyMixu37JlCwYOHIjJkydj06ZNWL16NRYvXoyoqCisWbMGZ86csUlic8899yAwMBAGgwH33ntvt+uzFR8fHxw4cAAvv/wydu3ahe3bt8Pf3x9DhgzBihUroFKpAMDu18eWFi5cCFEU8be//Q3PPPMMYmNj8e9//xt/+tOfXOJJKEQ9kSByJC4RuZHm5mao1Wrcc8892LRpk9ThuB2DwYC+ffvi/vvvx8aNG6UOh8jtcIwfEbmV3bt34/Lly5gzZ47Uobi8+vr6VrN8t2/fjitXrvCRbUQSYYsfEbmFQ4cO4dixY1i5ciWCgoIsLnZMtpWbm4unn34aDz74IPr06YP8/Hxs2rQJw4cPx5EjR8wWgCYix+AYPyJyC+vXr8d7772HuLg4bN26Vepw3EJkZCQiIiLw5ptv4sqVK+jduzfmzJmD1atXM+kjkghb/IiIiIjcBMf4EREREbkJdvVaYDAYcOnSJfj5+XGtKSIiInJ6oiiiuroaarUaMlnb7XpM/Cy4dOkSIiIipA6DiIiIyCrnz59HeHh4m/uZ+Fng5+cH4PrF8/f3lzgaIiIiovbpdDpERESYcpi2MPGzwNi96+/vz8SPiIiIeoyOhqhxcgcRERGRm5A08fv6669xzz33QK1WQxAE7N69u8NjcnNzcdttt0GhUGDQoEEW1+PKyspCZGQklEolEhIScPjwYdsHT0RERNTDSJr41dbWIjY2FllZWZ0qX1JSgunTp+Ouu+5CQUEBFi9ejMceewx79uwxldm5cyfS09ORmZmJ/Px8xMbGIjk5GeXl5fY6DSIiIqIewWkWcBYEAR9//DFSU1PbLLNkyRJ89tlnKCwsNG2bPXs2qqqqkJ2dDQBISEjA7bffjnXr1gG4vjRLREQEFi1ahKVLl3YqFp1OB5VKBa1WyzF+REQkmcZmA46cvYpxUb0hl3F5MWpbZ3OXHjXGLy8vD0lJSWbbkpOTkZeXBwBobGzEkSNHzMrIZDIkJSWZyljS0NAAnU5n9iIiIpLa37/+GQ9tPIgt35VIHQq5iB6V+Gk0GoSEhJhtCwkJgU6nw7Vr11BRUQG9Xm+xjEajabPeVatWQaVSmV5cw4+IiJzBd8WVAIC8nysljoRcRY9K/Oxl2bJl0Gq1ptf58+elDomIiNycKIo4cUkLADhxiT1RZBs9ah2/0NBQlJWVmW0rKyuDv78/vL29IZfLIZfLLZYJDQ1ts16FQgGFQmGXmImIiLriwtVr0NU3AwA0unpU1jSgTy/+W0Xd06Na/BITE5GTk2O2be/evUhMTAQAeHl5IT4+3qyMwWBATk6OqQwREVFP0LKVj61+ZAuSJn41NTUoKChAQUEBgOvLtRQUFODcuXMArnfBzpkzx1T+D3/4A3755Rc899xzOHXqFN5++218+OGHePrpp01l0tPTsXHjRmzbtg0nT57EggULUFtbi7S0NIeeGxERUXf8dKOb14iJH9mCpF29P/zwA+666y7T+/T0dADA3LlzsXXrVpSWlpqSQACIiorCZ599hqeffhpvvPEGwsPD8e677yI5OdlUZtasWbh8+TIyMjKg0WgQFxeH7OzsVhM+iIiInFnhjURPrVLikrYehS0SQaKucJp1/JwJ1/EjIiKpJby8D2W6BvxpymC8mVOEqCBf7H9mktRhkZNyyXX8iIiI3EFFTQPKdA0QBODB+HAAQElFLWoamiWOjHo6Jn5EREROxjieL6qPLyJ6+yDUXwkAOFnKcX7UPUz8iIiInIxx/b4R6utddiNv/PfERY7zo+5h4kdERORkjC1+I9WqG//1N9tO1FVM/IiIiJzMT6bE73rCN+JGAsjEj7qLiR8REZETqa5vQklFLYCbid+oftf/+9+yajQ06yWLjXo+Jn5ERERO5GRpNQAgTKU0PaKtX4A3VN6eaDaIKCqrkTI86uGY+BERETkR48QOY2sfAAiCcMs4P07woK5j4kdEROREjOP4jOP6jDjBg2yBiR8REZETOdFiYofRSE7wIBtg4kdEROQkGpr1KCq7PsavdeJ3/f3JUh30Bj5tlbqGiR8REZGT+K+mBs0GESpvT/QL8DbbN7BvLyg9Zahr1Jtm/RJZi4kfERGRkzBO3BjVzx+CIJjtk8sEDA/jBA/qHiZ+RERETqLlEztaMnb3/sRxftRFTPyIiIichKWlXG7FCR7UXUz8iIiInIDeIJoWb2478bvZ1SuKnOBB1mPiR0RE5ARKKmpxrUkPb085ooJ6WSwzJMQPcpmAq3VNKNXWOzhCcgVM/IiIiJyAsZt3WNj15M4Spaccg4N73SjP7l6yHhM/IiIiJ2BM5Ea1MbHDyDjOr/AiZ/aS9Zj4EREROYGOJnYY8dFt1B1M/IiIiCQmimKHS7kY3VzShS1+ZD0mfkRERBK7pK1HVV0TPGQChoRanthhNOJG4ndJW4+rtY2OCI9cCBM/IiIiiZ24MV5vUHAvKDzk7Zb1U3piQB+f68exu5esxMSPiIhIYp3t5jW6dT0/Imsw8SMiIpLYzcSv/YkdRnyCB3UVEz8iIiKJGVvuRvWzrsWvkC1+ZCWnSPyysrIQGRkJpVKJhIQEHD58uM2ykyZNgiAIrV7Tp083lZk3b16r/SkpKY44FSIiIqtcqW00PYVjeJhfp44xtviVVNSitqHZbrGR65E88du5cyfS09ORmZmJ/Px8xMbGIjk5GeXl5RbLf/TRRygtLTW9CgsLIZfL8eCDD5qVS0lJMSv3wQcfOOJ0iIiIrGJs7Yvs4wM/pWenjunrp0CwnwKiCJzSsLuXOk/yxO+1117D/PnzkZaWhhEjRmDDhg3w8fHB5s2bLZbv3bs3QkNDTa+9e/fCx8enVeKnUCjMygUGBrYZQ0NDA3Q6ndmLiIjIEayd2GHEhZypKyRN/BobG3HkyBEkJSWZtslkMiQlJSEvL69TdWzatAmzZ8+Gr6+v2fbc3FwEBwdj6NChWLBgASorK9usY9WqVVCpVKZXRERE106IiIjISsbEbUQnJ3YYmSZ4XGTiR50naeJXUVEBvV6PkJAQs+0hISHQaDQdHn/48GEUFhbiscceM9uekpKC7du3IycnB2vWrMGBAwcwbdo06PV6i/UsW7YMWq3W9Dp//nzXT4qIiMgKnX1UW0umFr9STvCgzvOQOoDu2LRpE2JiYjBu3Diz7bNnzzb9HBMTg9GjRyM6Ohq5ubmYMmVKq3oUCgUUCoXd4yUiIrpVbUMzSipqAVjf1WucAXxaU43GZgO8PCQfvUU9gKR3SVBQEORyOcrKysy2l5WVITQ0tN1ja2trsWPHDjz66KMdfs7AgQMRFBSE4uLibsVLRERkSydLdRBFIMRfgb5+1jVAhAd6w1/pgSa9iKLyajtFSK5G0sTPy8sL8fHxyMnJMW0zGAzIyclBYmJiu8fu2rULDQ0NeOSRRzr8nAsXLqCyshJhYWHdjpmIiMhWujqxAwAEQTCNC+QED+osyduF09PTsXHjRmzbtg0nT57EggULUFtbi7S0NADAnDlzsGzZslbHbdq0CampqejTp4/Z9pqaGjz77LM4ePAgzpw5g5ycHMyYMQODBg1CcnKyQ86JiIioM7o6vs/ImDD+xMSPOknyMX6zZs3C5cuXkZGRAY1Gg7i4OGRnZ5smfJw7dw4ymXl+evr0aXz77bf48ssvW9Unl8tx7NgxbNu2DVVVVVCr1Zg6dSpWrlzJcXxERORUrH1UW0t8Zi9ZSxBFUZQ6CGej0+mgUqmg1Wrh79+1X0YiIqL2NDYbMDIzG016Ed88dxcievtYXcdpTTWSX/8avl5yHH8xGTKZYIdIqSfobO4ieVcvERGROyoqr0aTXoS/0gPhgd5dqiO6ry8UHjLUNupx9kqdjSMkV8TEj4iISALGhZdHqlUQhK611HnIZRgWdr11p/Aiu3upY0z8iIiIJNDdiR1GfHQbWYOJHxERkQRMEzv62SrxY4sfdYyJHxERkYMZDCJOlnZ9Db9b3bqkC+drUkeY+BERETnYmcpa1DbqofCQYWCQb7fqGhbqB7lMQGVtI8p0DTaKkFwVEz8iIiIHM3bzDgvzh4e8e/8UKz3liO7re6NedvdS+5j4EREROVh3F25uydjdywke1BEmfkRERA5mbJkb1c3xfUbGBJJLulBHmPgRERE5kCiKbPEjyTDxIyIiciCNrh5XahshlwkYGupnkzpH3EggL1ZdQ1Vdo03qJNfExI+IiMiBjE/sGNS3F5SecpvUqfL2RETv6499+4mtftQOJn5EREQOZOtuXqORYezupY4x8SMiInIg48SOEbZO/PgED+oEJn5EREQOZGyRG9XPNjN6jYz1FbLFj9rBxI+IiMhBrtY24mLVNQD2a/H75XINrjXqbVo3uQ6rE7+mpiZ4eHigsLDQHvEQERG5rJ9uPJ+3f28f+Cs9bVp3sL8SQb0UMIjASQ1b/cgyqxM/T09P9O/fH3o9/5ogIiKyhnH8na0ndhjdHOfHxI8s61JX71/+8hc8//zzuHLliq3jISIicln2mtFrZKz3J07woDZ4dOWgdevWobi4GGq1GgMGDICvr6/Z/vz8fJsER0RE5EpuJn62ndhhxCd4UEe6lPilpqbaOAwiIiLXVtfYjF8u1wCwf4vfKU01mvQGeMo5h5PMdSnxy8zMtHUcRERELu1kaTUMItDXT4Fgf6VdPqN/bx/4KTxQ3dCM4vIaDA+zT4JJPVeXEj+jI0eO4OTJkwCAkSNHYsyYMTYJioiIyNX8ZOeJHQAgkwkYrvbH4ZIrOHFJx8SPWulS4ldeXo7Zs2cjNzcXAQEBAICqqircdddd2LFjB/r27WvLGImIiHo8e0/sMBppSvy0eCA+3K6fRT1Plzr/Fy1ahOrqapw4cQJXrlzBlStXUFhYCJ1Ohz/96U+2jpGIiKjHs/fEDiNO8KD2dKnFLzs7G/v27cPw4cNN20aMGIGsrCxMnTrVZsERERG5gia9Aac11QAc0+IHACcv6WAwiJDJBLt+HvUsXWrxMxgM8PRsveK4p6cnDAaD1fVlZWUhMjISSqUSCQkJOHz4cJtlt27dCkEQzF5KpfkgWVEUkZGRgbCwMHh7eyMpKQlFRUVWx0VERGQLxeU1aNQb4KfwQESgj10/a1BwL3h5yFDd0IzzV+vs+lnU83Qp8Zs8eTKeeuopXLp0ybTt4sWLePrppzFlyhSr6tq5cyfS09ORmZmJ/Px8xMbGIjk5GeXl5W0e4+/vj9LSUtPr7NmzZvvXrl2LN998Exs2bMChQ4fg6+uL5ORk1NfXW3eiRERENmDsdh2h9rd7C5ynXIZhoX5mn0tk1KXEb926ddDpdIiMjER0dDSio6MRFRUFnU6Ht956y6q6XnvtNcyfPx9paWkYMWIENmzYAB8fH2zevLnNYwRBQGhoqOkVEhJi2ieKIl5//XW88MILmDFjBkaPHo3t27fj0qVL2L17t8X6GhoaoNPpzF5ERES2UnjROKPXvuP7jIzdvcbPJTLqUuIXERGB/Px8fPbZZ1i8eDEWL16Mzz//HPn5+QgP7/wMosbGRhw5cgRJSUk3A5LJkJSUhLy8vDaPq6mpwYABAxAREYEZM2bgxIkTpn0lJSXQaDRmdapUKiQkJLRZ56pVq6BSqUyviIiITp8DERFRR35y0IxeoxGc4EFtsHpyR1NTE7y9vVFQUIC7774bd999d5c/vKKiAnq93qzFDgBCQkJw6tQpi8cMHToUmzdvxujRo6HVavHqq69i/PjxOHHiBMLDw6HRaEx1tKzTuK+lZcuWIT093fRep9Mx+SMiIpswGET8VHoj8evnmMTPmGAy8aOWrE78PD090b9/f+j1envE06HExEQkJiaa3o8fPx7Dhw/HO++8g5UrV3apToVCAYVCYasQiYiITM5dqUNNQzO8PGSI7tvLIZ85PNQfMgGoqGlAua7ebk8KoZ6nS129f/nLX/D888/jypUr3frwoKAgyOVylJWVmW0vKytDaGhop+rw9PTEmDFjUFxcDACm47pTJxERka0YW92Ghfo57Nm53l5yDLyRZLLVj27V5ckdX3/9NdRqNYYOHYrbbrvN7NVZXl5eiI+PR05OjmmbwWBATk6OWatee/R6PY4fP46wsDAAQFRUFEJDQ83q1Ol0OHToUKfrJCIispUTDnhUmyU3u3s5wYNu6tICzqmpqTYLID09HXPnzsXYsWMxbtw4vP7666itrUVaWhoAYM6cOejXrx9WrVoFAHjppZdwxx13YNCgQaiqqsIrr7yCs2fP4rHHHgNwfcbv4sWL8de//hWDBw9GVFQUli9fDrVabdO4iYiIOqPQQU/saGmUWoVPCi6h8CJb/OgmqxO/5uZmCIKA3//+91bN4G3LrFmzcPnyZWRkZECj0SAuLg7Z2dmmyRnnzp2DTHazYfLq1auYP38+NBoNAgMDER8fj++//x4jRowwlXnuuedQW1uLxx9/HFVVVfjVr36F7OzsVgs9ExER2ZMoivhJ6ha/Urb40U2CKIqitQf5+fnh+PHjiIyMtENI0tPpdFCpVNBqtfD3d+wvKhERuY4yXT0SXs6BTABOrEiBt5fcYZ9dVdeIuJf2AgB+zJwKlXfrJ26R6+hs7tLlJ3ccOHCgy8ERERG5A+P4uui+vRya9AFAgI8X+gV4A7i5jiBRl8b4TZs2DUuXLsXx48cRHx8PX19fs/333nuvTYIjIiLqyU5cdOzCzS2NVPvjYtU1nLikRWJ0H0liIOfSpcTvj3/8I4Drj1trSRAEydb4IyIiciYnJJrYYTRSrcKXP5WxxY9MupT4GQwGW8dBRETkcowTKxz1xI6WRvXjEzzInFVj/H7zm99Aq705O2j16tWoqqoyva+srDSbXUtEROSutHVNOH/lGgBgZJh0LX4AUHy5BvVN7I0jKxO/PXv2oKGhwfT+5ZdfNnt6R3NzM06fPm276IiIiHooY2tfeKA3VD7SzKgN8Vegj68X9AYRpzTVksRAzsWqxK/lyi9dWAmGiIjILfx0SdqJHcD1cfcj+AQPuoVjHhpIRETkZqSe2GFk/HyO8yPAysRPEAQIgtBqGxEREZmT6hm9Ld18Zi8TP7JyVq8oipg3bx4UCgUAoL6+Hn/4wx9M6/jdOv6PiIjIXdU36fHz5VoAztDidz3xO1WqQ7PeAA85O/vcmVWJ39y5c83eP/LII63KzJkzp3sRERER9XCnNNXQG0QE9fJCiL9C0lgi+/jC10uO2kY9fqmoxZAQP0njIWlZlfht2bLFXnEQERG5jMKL17t5R6hVkg+JksmuT/D4z5mrKLyoZeLn5tjeS0REZGMnnGBG7604wYOMmPgRERHZ2E9OMrHDiEu6kBETPyIiIhtq1htMiyVLPbHDyJiA/nRJxzV43RwTPyIiIhv6+XItGpoN6KXwwIDePlKHAwAYHOwHT7kAXX0zLly9JnU4JCEmfkRERDZk7E4dEeYPmcw51rr18pBhaOj1SR3s7nVvTPyIiIhsqPDi9QkUI5xkfJ/RyLDr3c7G+Mg9MfEjIiKyIWd5YkdLI/txggcx8SMiIrIZURTxU6lzPKO3JT66jQAmfkRERDZz/so1VNc3w0suw+CQXlKHY2ZYqD8EASivbsDlaj5i1V0x8SMiIrIRYzfqkNBe8HSyZ+L6KjwQFeQLgN297sy57koiIqIezPTEjjDn6uY14hM8iIkfERGRjRhb0kb1c66JHUajblnImdwTEz8iIiIbKbxkXMrFuVv8CtnV67acIvHLyspCZGQklEolEhIScPjw4TbLbty4ERMmTEBgYCACAwORlJTUqvy8efMgCILZKyUlxd6nQUREbqy8uh6XqxsgCMDwMD+pw7HIOLP3bGUddPVNEkdDUpA88du5cyfS09ORmZmJ/Px8xMbGIjk5GeXl5RbL5+bm4qGHHsL+/fuRl5eHiIgITJ06FRcvXjQrl5KSgtLSUtPrgw8+cMTpEBGRmzKOmxsY5AsfLw+Jo7Es0NcLapUSAHCS3b1uSfLE77XXXsP8+fORlpaGESNGYMOGDfDx8cHmzZstln///ffxxz/+EXFxcRg2bBjeffddGAwG5OTkmJVTKBQIDQ01vQIDAx1xOkRE5KaM4+acbf2+lkZwgodbkzTxa2xsxJEjR5CUlGTaJpPJkJSUhLy8vE7VUVdXh6amJvTu3dtse25uLoKDgzF06FAsWLAAlZWVbdbR0NAAnU5n9iIiIrKGsz6xoyUu5OzeJE38KioqoNfrERISYrY9JCQEGo2mU3UsWbIEarXaLHlMSUnB9u3bkZOTgzVr1uDAgQOYNm0a9Hq9xTpWrVoFlUplekVERHT9pIiIyC2d6CEtfjcTP07wcEfOOQihk1avXo0dO3YgNzcXSqXStH327Nmmn2NiYjB69GhER0cjNzcXU6ZMaVXPsmXLkJ6ebnqv0+mY/BERUafp6ptwtrIOgPO3+I3qdz0xLS6vQX2THkpPucQRkSNJ2uIXFBQEuVyOsrIys+1lZWUIDQ1t99hXX30Vq1evxpdffonRo0e3W3bgwIEICgpCcXGxxf0KhQL+/v5mLyIios4yju/rF+CNQF8viaNpX5hKiUAfTzQbRPy3rFrqcMjBJE38vLy8EB8fbzYxwzhRIzExsc3j1q5di5UrVyI7Oxtjx47t8HMuXLiAyspKhIWF2SRuIiKiW50wrd/n/A0HgiDwCR5uTPJZvenp6di4cSO2bduGkydPYsGCBaitrUVaWhoAYM6cOVi2bJmp/Jo1a7B8+XJs3rwZkZGR0Gg00Gg0qKmpAQDU1NTg2WefxcGDB3HmzBnk5ORgxowZGDRoEJKTkyU5RyIicm09ZWKHEcf5uS/Jx/jNmjULly9fRkZGBjQaDeLi4pCdnW2a8HHu3DnIZDfz0/Xr16OxsREPPPCAWT2ZmZl48cUXIZfLcezYMWzbtg1VVVVQq9WYOnUqVq5cCYVC4dBzIyIi99BTlnIxGsGZvW5LEEVRlDoIZ6PT6aBSqaDVajnej4iI2lXfpMfIzD3QG0R8v3Qy1AHeUofUoeLyGiS9dgDennIUrkiGXCZIHRJ1U2dzF8m7eomIiHqy/5ZVQ28Q0dvXC2EqZccHOIGoIF/4eMlxrUmPkooaqcMhB2LiR0RE1A2FF43dvP4QhJ7RciaXCRgedr1VyBg/uQcmfkRERN1gnCDRE2b03ooTPNwTEz8iIqJu6ClP7GiJj25zT0z8iIiIukhvEHFKc7Ortye5dS0/zvN0H0z8iIiIuuiXyzWobzLAx0uOqD6+UodjlcEhveAhE6C91oSLVdekDocchIkfERFRFxm7SYeH+UPWw5ZEUXjIMTjEDwC7e90JEz8iIqIuMk6MGNXDunmNRnGcn9th4kdERNRFN5dy6VkTO4xMEzwucmavu2DiR0RE1AWiKPbYpVyMRva7OcGD3AMTPyIioi64cPUadPXN8JQLGHJjrFxPMzzMH4IAaHT1qKxpkDoccgAmfkRERF1gbCUbHOwHL4+e+c9pL4UHIm/MRmarn3vomXcqERGRxH660c3b09bva2kEJ3i4FSZ+REREXWBMlEb165kTO4xGmRZy5gQPd8DEj4iIqAtuPqqtZ7f4GeP/iS1+boGJHxERkZUqahqg0dVDEK5PkOjJjInfLxW1qGloljgasjcmfkRERFYytvZF9fGFr8JD4mi6p08vBUL9lQCAk6Vs9XN1TPyIiIis1NPX72uJCzm7DyZ+REREVro5vq9nT+wwGsmZvW6DiR8REZGVfnKRiR1GI9R8goe7YOJHRERkher6JpRU1AJwncRvVL/r51FUXo2GZr3E0ZA9MfEjIiKywsnSagBAmEqJPr0UEkdjG/0CvKHy9kSTXkRRWY3U4ZAdMfEjIiKywgkXeWLHrQRBuGWcHyd4uDImfkRERFYwjoMb4SITO4w4wcM9MPEjIiKygqs8saOlkZzg4RaY+BEREXVSQ7MeRWXXx/i5XuJ3/XxOluqgN4gSR0P24hSJX1ZWFiIjI6FUKpGQkIDDhw+3W37Xrl0YNmwYlEolYmJi8Pnnn5vtF0URGRkZCAsLg7e3N5KSklBUVGTPUyAiIjdQVFaDZoOIAB9P9AvwljocmxrYtxeUnjLUNepxprJW6nDITiRP/Hbu3In09HRkZmYiPz8fsbGxSE5ORnl5ucXy33//PR566CE8+uijOHr0KFJTU5GamorCwkJTmbVr1+LNN9/Ehg0bcOjQIfj6+iI5ORn19fWOOi0iInJBt07sEARB4mhsSy4TTM8dZnev6xJEUZS0PTchIQG333471q1bBwAwGAyIiIjAokWLsHTp0lblZ82ahdraWnz66aembXfccQfi4uKwYcMGiKIItVqNP//5z3jmmWcAAFqtFiEhIdi6dStmz57dYUw6nQ4qlQparRb+/vZpyj/0S6VpHSgiIuoZPjteim+KKvD4rwfi+d8Mlzocm3th93G8d/AckoYHI2l4iNThuJSxkb0xKLiX3ervbO4i6ZOlGxsbceTIESxbtsy0TSaTISkpCXl5eRaPycvLQ3p6utm25ORk7N69GwBQUlICjUaDpKQk036VSoWEhATk5eVZTPwaGhrQ0NBgeq/T2f8vnX8euYBdRy7Y/XOIiMj2XG18n9GoGxM89p0sx76TlnveqGtW3R9j18SvsyRN/CoqKqDX6xESYv5XRUhICE6dOmXxGI1GY7G8RqMx7Tdua6tMS6tWrcKKFSu6dA5dNULtj6Q6/jVFRNTThPgrkDwyVOow7OKeWDUKzlehoqZR6lBcjrOMCZU08XMWy5YtM2tF1Ol0iIiIsOtnpt0ZhbQ7o+z6GURERNbwVXhg9W9HSx0G2ZGkkzuCgoIgl8tRVlZmtr2srAyhoZb/mgoNDW23vPG/1tSpUCjg7+9v9iIiIiJyNZImfl5eXoiPj0dOTo5pm8FgQE5ODhITEy0ek5iYaFYeAPbu3WsqHxUVhdDQULMyOp0Ohw4darNOIiIiIncgeVdveno65s6di7Fjx2LcuHF4/fXXUVtbi7S0NADAnDlz0K9fP6xatQoA8NRTT2HixIn429/+hunTp2PHjh344Ycf8Pe//x3A9ecNLl68GH/9618xePBgREVFYfny5VCr1UhNTZXqNImIiIgkJ3niN2vWLFy+fBkZGRnQaDSIi4tDdna2aXLGuXPnIJPdbJgcP348/vGPf+CFF17A888/j8GDB2P37t0YNWqUqcxzzz2H2tpaPP7446iqqsKvfvUrZGdnQ6lUOvz8iIiIiJyF5Ov4OSOtVouAgACcP3+e4/2IiIjI6RknplZVVUGlUrVZTvIWP2dUXX39OYz2ntlLREREZEvV1dXtJn5s8bPAYDDg0qVL8PPzs9sjeYyZOVsVW+O1aRuvjWW8Lm3jtbGM16VtvDaWOft1EUUR1dXVUKvVZkPkWmKLnwUymQzh4eEO+SwuH9M2Xpu28dpYxuvSNl4by3hd2sZrY5kzX5f2WvqMJF3OhYiIiIgch4kfERERkZtg4icRhUKBzMxMKBQKqUNxOrw2beO1sYzXpW28NpbxurSN18YyV7kunNxBRERE5CbY4kdERETkJpj4EREREbkJJn5EREREboKJHxEREZGbYOJHRERE5CaY+NlRVlYWIiMjoVQqkZCQgMOHD7dbfteuXRg2bBiUSiViYmLw+eefOyhSx1m1ahVuv/12+Pn5ITg4GKmpqTh9+nS7x2zduhWCIJi9lEqlgyJ2nBdffLHVeQ4bNqzdY9zhnomMjGx1XQRBwJNPPmmxvCvfL19//TXuueceqNVqCIKA3bt3m+0XRREZGRkICwuDt7c3kpKSUFRU1GG91n5XOZv2rktTUxOWLFmCmJgY+Pr6Qq1WY86cObh06VK7dXbl99EZdXTPzJs3r9V5pqSkdFivK98zACx+5wiCgFdeeaXNOnvKPcPEz0527tyJ9PR0ZGZmIj8/H7GxsUhOTkZ5ebnF8t9//z0eeughPProozh69ChSU1ORmpqKwsJCB0duXwcOHMCTTz6JgwcPYu/evWhqasLUqVNRW1vb7nH+/v4oLS01vc6ePeugiB1r5MiRZuf57bfftlnWXe6Z//znP2bXZO/evQCABx98sM1jXPV+qa2tRWxsLLKysizuX7t2Ld58801s2LABhw4dgq+vL5KTk1FfX99mndZ+Vzmj9q5LXV0d8vPzsXz5cuTn5+Ojjz7C6dOnce+993ZYrzW/j86qo3sGAFJSUszO84MPPmi3Tle/ZwCYXY/S0lJs3rwZgiDgt7/9bbv19oh7RiS7GDdunPjkk0+a3uv1elGtVourVq2yWH7mzJni9OnTzbYlJCSITzzxhF3jlFp5ebkIQDxw4ECbZbZs2SKqVCrHBSWRzMxMMTY2ttPl3fWeeeqpp8To6GjRYDBY3O8u9wsA8eOPPza9NxgMYmhoqPjKK6+YtlVVVYkKhUL84IMP2qzH2u8qZ9fyulhy+PBhEYB49uzZNstY+/vYE1i6NnPnzhVnzJhhVT3ueM/MmDFDnDx5crtleso9wxY/O2hsbMSRI0eQlJRk2iaTyZCUlIS8vDyLx+Tl5ZmVB4Dk5OQ2y7sKrVYLAOjdu3e75WpqajBgwABERERgxowZOHHihCPCc7iioiKo1WoMHDgQDz/8MM6dO9dmWXe8ZxobG/Hee+/h97//PQRBaLOcu9wvtyopKYFGozG7J1QqFRISEtq8J7ryXeUKtFotBEFAQEBAu+Ws+X3syXJzcxEcHIyhQ4diwYIFqKysbLOsO94zZWVl+Oyzz/Doo492WLYn3DNM/OygoqICer0eISEhZttDQkKg0WgsHqPRaKwq7woMBgMWL16MO++8E6NGjWqz3NChQ7F582Z88skneO+992AwGDB+/HhcuHDBgdHaX0JCArZu3Yrs7GysX78eJSUlmDBhAqqrqy2Wd8d7Zvfu3aiqqsK8efPaLOMu90tLxv/v1twTXfmu6unq6+uxZMkSPPTQQ/D392+znLW/jz1VSkoKtm/fjpycHKxZswYHDhzAtGnToNfrLZZ3x3tm27Zt8PPzw/33399uuZ5yz3hIHQC5ryeffBKFhYUdjoFITExEYmKi6f348eMxfPhwvPPOO1i5cqW9w3SYadOmmX4ePXo0EhISMGDAAHz44Yed+kvTHWzatAnTpk2DWq1us4y73C9kvaamJsycOROiKGL9+vXtlnWX38fZs2ebfo6JicHo0aMRHR2N3NxcTJkyRcLInMfmzZvx8MMPdzhJrKfcM2zxs4OgoCDI5XKUlZWZbS8rK0NoaKjFY0JDQ60q39MtXLgQn376Kfbv34/w8HCrjvX09MSYMWNQXFxsp+icQ0BAAIYMGdLmebrbPXP27Fns27cPjz32mFXHucv9Yvz/bs090ZXvqp7KmPSdPXsWe/fubbe1z5KOfh9dxcCBAxEUFNTmebrTPQMA33zzDU6fPm319w7gvPcMEz878PLyQnx8PHJyckzbDAYDcnJyzFoibpWYmGhWHgD27t3bZvmeShRFLFy4EB9//DG++uorREVFWV2HXq/H8ePHERYWZocInUdNTQ1+/vnnNs/TXe4Zoy1btiA4OBjTp0+36jh3uV+ioqIQGhpqdk/odDocOnSozXuiK99VPZEx6SsqKsK+ffvQp08fq+vo6PfRVVy4cAGVlZVtnqe73DNGmzZtQnx8PGJjY60+1mnvGalnl7iqHTt2iAqFQty6dav4008/iY8//rgYEBAgajQaURRF8Xe/+524dOlSU/nvvvtO9PDwEF999VXx5MmTYmZmpujp6SkeP35cqlOwiwULFogqlUrMzc0VS0tLTa+6ujpTmZbXZsWKFeKePXvEn3/+WTxy5Ig4e/ZsUalUiidOnJDiFOzmz3/+s5ibmyuWlJSI3333nZiUlCQGBQWJ5eXloii67z0jitdnDfbv319csmRJq33udL9UV1eLR48eFY8ePSoCEF977TXx6NGjptmpq1evFgMCAsRPPvlEPHbsmDhjxgwxKipKvHbtmqmOyZMni2+99ZbpfUffVT1Be9elsbFRvPfee8Xw8HCxoKDA7HunoaHBVEfL69LR72NP0d61qa6uFp955hkxLy9PLCkpEfft2yfedttt4uDBg8X6+npTHe52zxhptVrRx8dHXL9+vcU6euo9w8TPjt566y2xf//+opeXlzhu3Djx4MGDpn0TJ04U586da1b+ww8/FIcMGSJ6eXmJI0eOFD/77DMHR2x/ACy+tmzZYirT8tosXrzYdB1DQkLE3/zmN2J+fr7jg7ezWbNmiWFhYaKXl5fYr18/cdasWWJxcbFpv7veM6Ioinv27BEBiKdPn261z53ul/3791v8/TGev8FgEJcvXy6GhISICoVCnDJlSqtrNmDAADEzM9NsW3vfVT1Be9elpKSkze+d/fv3m+poeV06+n3sKdq7NnV1deLUqVPFvn37ip6enuKAAQPE+fPnt0rg3O2eMXrnnXdEb29vsaqqymIdPfWeEURRFO3apEhEREREToFj/IiIiIjcBBM/IiIiIjfBxI+IiIjITTDxIyIiInITTPyIiIiI3AQTPyIiIiI3wcSPiIiIyE0w8SMiIiJyE0z8iIiIiNwEEz8iIiIiN8HEj4iIiMhNMPEjIiIichNM/IiIiIjcBBM/IiIiIjfBxI+IiIjITTDxIyIiInITTPyIiIiI3AQTPyKiNgiCgBdffNFm9c2bNw+RkZE2q89eXnzxRQiCIHUYRGQHTPyIyOG2bt0KQRDwww8/dLuuuro6vPjii8jNze1+YERELs5D6gCIiLqjrq4OK1asAABMmjTJpnVfu3YNHh7u9zX5wgsvYOnSpVKHQUR2wBY/IqI2KJVKl0j86urqrCrv4eEBpVJpp2iISEpM/IjIKTU2NiIjIwPx8fFQqVTw9fXFhAkTsH//flOZM2fOoG/fvgCAFStWQBCEVuPyTp06hQceeAC9e/eGUqnE2LFj8e9//7tTMbSsyzj2rbi4GPPmzUNAQABUKhXS0tKsTq6MDAYDXn/9dYwcORJKpRIhISF44okncPXqVbNyn3zyCaZPnw61Wg2FQoHo6GisXLkSer3erNykSZMwatQoHDlyBL/+9a/h4+OD559/HmfOnIEgCHj11Vfx97//HdHR0VAoFLj99tvxn//8x6wOS2P8BEHAwoULsXv3bowaNQoKhQIjR45EdnZ2q3PKzc3F2LFjoVQqER0djXfeeYfjBomcRM//U5aIXJJOp8O7776Lhx56CPPnz0d1dTU2bdqE5ORkHD58GHFxcejbty/Wr1+PBQsW4L777sP9998PABg9ejQA4MSJE7jzzjvRr18/LF26FL6+vvjwww+RmpqKf/3rX7jvvvu6FNvMmTMRFRWFVatWIT8/H++++y6Cg4OxZs0aq+t64oknsHXrVqSlpeFPf/oTSkpKsG7dOhw9ehTfffcdPD09AVwfF9mrVy+kp6ejV69e+Oqrr5CRkQGdTodXXnnFrM7KykpMmzYNs2fPxiOPPIKQkBDTvn/84x+orq7GE088AUEQsHbtWtx///345ZdfTJ/Vlm+//RYfffQR/vjHP8LPzw9vvvkmfvvb3+LcuXPo06cPAODo0aNISUlBWFgYVqxYAb1ej5deesmUoBORxEQiIgfbsmWLCED8z3/+02aZ5uZmsaGhwWzb1atXxZCQEPH3v/+9advly5dFAGJmZmarOqZMmSLGxMSI9fX1pm0Gg0EcP368OHjw4A7jbFlvZmamCMDs80VRFO+77z6xT58+HdY3d+5cccCAAab333zzjQhAfP/9983KZWdnt9peV1fXqr4nnnhC9PHxMTu/iRMnigDEDRs2mJUtKSkRAYh9+vQRr1y5Ytr+ySefiADE//f//l+r87wVANHLy0ssLi42bfvxxx9FAOJbb71l2nbPPfeIPj4+4sWLF03bioqKRA8Pj1Z1EpHjsauXiJySXC6Hl5cXgOvdoVeuXEFzczPGjh2L/Pz8Do+/cuUKvvrqK8ycORPV1dWoqKhARUUFKisrkZycjKKiIly8eLFLsf3hD38wez9hwgRUVlZCp9NZVc+uXbugUqlw9913m+KrqKhAfHw8evXqZdat7e3tbfrZeD4TJkxAXV0dTp06ZVavQqFAWlqaxc+cNWsWAgMDzWIHgF9++aXDeJOSkhAdHW16P3r0aPj7+5uO1ev12LdvH1JTU6FWq03lBg0ahGnTpnVYPxHZH7t6ichpbdu2DX/7299w6tQpNDU1mbZHRUV1eGxxcTFEUcTy5cuxfPlyi2XKy8vRr18/q+Pq37+/2XtjInX16lX4+/t3up6ioiJotVoEBwe3GZ/RiRMn8MILL+Crr75qlWBqtVqz9/369TMlzdbE3pGWxxqPNx5bXl6Oa9euYdCgQa3KWdpGRI7HxI+InNJ7772HefPmITU1Fc8++yyCg4Mhl8uxatUq/Pzzzx0ebzAYAADPPPMMkpOTLZbpajIil8stbhdF0ap6DAYDgoOD8f7771vcbxwXV1VVhYkTJ8Lf3x8vvfQSoqOjoVQqkZ+fjyVLlpjO1ejW1kFbxm6r8yYi6TDxIyKn9M9//hMDBw7ERx99ZDYbNDMz06xcWzNFBw4cCADw9PREUlKS/QLthujoaOzbtw933nlnu8labm4uKisr8dFHH+HXv/61aXtJSYkjwuy04OBgKJVKFBcXt9pnaRsROR7H+BGRUzK2Lt3amnTo0CHk5eWZlfPx8QFwvVXsVsHBwZg0aRLeeecdlJaWtqr/8uXLNo7YejNnzoRer8fKlStb7Wtubjadk6Vr0djYiLffftshcXaWXC5HUlISdu/ejUuXLpm2FxcX44svvpAwMiIyYosfEUlm8+bNFteBe+qpp/A///M/+Oijj3Dfffdh+vTpKCkpwYYNGzBixAjU1NSYynp7e2PEiBHYuXMnhgwZgt69e2PUqFEYNWoUsrKy8Ktf/QoxMTGYP38+Bg4ciLKyMuTl5eHChQv48ccfHXm6rUycOBFPPPEEVq1ahYKCAkydOhWenp4oKirCrl278MYbb+CBBx7A+PHjERgYiLlz5+JPf/oTBEHA//2//9cpu1hffPFFfPnll7jzzjuxYMEC6PV6rFu3DqNGjUJBQYHU4RG5PSZ+RCSZ9evXW9w+b948zJs3DxqNBu+88w727NmDESNG4L333sOuXbtaPZf33XffxaJFi/D000+jsbERmZmZGDVqFEaMGIEffvgBK1aswNatW1FZWYng4GCMGTMGGRkZDjjDjm3YsAHx8fF455138Pzzz8PDwwORkZF45JFHcOeddwIA+vTpg08//RR//vOf8cILLyAwMBCPPPIIpkyZ0ub4RanEx8fjiy++wDPPPIPly5cjIiICL730Ek6ePNlq9jEROZ4gOuOfjERE5FJSU1Nx4sQJFBUVSR0KkVvjGD8iIrKpa9eumb0vKirC559/jkmTJkkTEBGZsMWPiIhsKiwsDPPmzcPAgQNx9uxZrF+/Hg0NDTh69CgGDx4sdXhEbo1j/IiIyKZSUlLwwQcfQKPRQKFQIDExES+//DKTPiInwBY/IiIiIjfBMX5EREREboJdvRYYDAZcunQJfn5+bT4VgIiIiMhZiKKI6upqqNVqyGRtt+sx8bPg0qVLiIiIkDoMIiIiIqucP38e4eHhbe5n4meBn58fgOsXz9/fX+JoiIiIiNqn0+kQERFhymHawsTPAmP3rr+/PxM/IiIi6jE6GqLGyR1EREREbkLSxO/rr7/GPffcA7VaDUEQsHv37g6Pyc3NxW233QaFQoFBgwZh69atrcpkZWUhMjISSqUSCQkJOHz4sO2DJyIiIuphJE38amtrERsbi6ysrE6VLykpwfTp03HXXXehoKAAixcvxmOPPYY9e/aYyuzcuRPp6enIzMxEfn4+YmNjkZycjPLycnudBhEREVGP4DQLOAuCgI8//hipqaltllmyZAk+++wzFBYWmrbNnj0bVVVVyM7OBgAkJCTg9ttvx7p16wBcX5olIiICixYtwtKlSy3W29DQgIaGBtN74wBJrVbLMX5ERETk9HQ6HVQqVYe5S48a45eXl4ekpCSzbcnJycjLywMANDY24siRI2ZlZDIZkpKSTGUsWbVqFVQqlenFpVyIiIjIFfWoxE+j0SAkJMRsW0hICHQ6Ha5du4aKigro9XqLZTQaTZv1Llu2DFqt1vQ6f/68XeInIiIikhKXcwGgUCigUCikDoOIiIjIrnpU4hcaGoqysjKzbWVlZfD394e3tzfkcjnkcrnFMqGhoY4MlYiIiMjp9Kiu3sTEROTk5Jht27t3LxITEwEAXl5eiI+PNytjMBiQk5NjKkNERETkriRN/GpqalBQUICCggIA15drKSgowLlz5wBcH3s3Z84cU/k//OEP+OWXX/Dcc8/h1KlTePvtt/Hhhx/i6aefNpVJT0/Hxo0bsW3bNpw8eRILFixAbW0t0tLSHHpuRERERM5G0q7eH374AXfddZfpfXp6OgBg7ty52Lp1K0pLS01JIABERUXhs88+w9NPP4033ngD4eHhePfdd5GcnGwqM2vWLFy+fBkZGRnQaDSIi4tDdnZ2qwkfRERERO7GadbxcyadXQuHiIiIyBm45Dp+RERERNR1TPyIiIiI3AQTPyIiIiI3wcSPiIiIyE0w8SMiIiJyE0z8iLrov2XVePebX6A3cGI8ERH1DD3qkW1EzuT5j47jh7NXEdRLgdQx/aQOh4iIqENs8SPqgsZmA45d0AIAjpy9KnE0REREncPEj6gLTml0aNQbAAA/XqiSNhgiIqJOYuJH1AUF56tMP58s1aG+SS9dMERERJ3ExI+oC25N/Jr0In4q1UkXDBERUScx8SPqgh9vJH5+Cg+z90RERM6MiR+RlbTXmvDz5VoAwG/jwwGYtwASERE5KyZ+RFY6fmM2b//ePpg8LBgAW/yIiKhnYOJHZCXjLN7YiACMDlcBAM5U1qGqrlHCqIiIiDrGxI/ISkfPVQEAYsNVCPDxQlSQLwB29xIRkfNj4kdkBVEUTQnemP4BAIC4iOv//fG8VpqgiIiIOomJH5EVLmnrUVHTALlMwEj19W7e2BvdvQXn+QQPIiJybkz8iKxgnMQxLNQPSk85gOtj/QDgxwtaiKIoUWREREQdY+JHZAVj4mfs3gWAEWp/eMoFXKltxIWr16QJjIiIqBOY+BFZ4eiNxC/2lsRP4SHHiDB/s/1ERETOiIkfUSc16w2mNfxubfEDbunuZeJHREROjIkfUScVX67BtSY9eik8EN23l9m+OCZ+RETUAzDxI+qkghvr98X0U0EuE8z2GVv8jl/UoklvcHBkREREncPEj6iTbn1iR0tRfXzhp/RAQ7MBpzXVjg2MiIiok5j4EXVSwXnL4/sAQCYTbnb33kgQiYiInI1TJH5ZWVmIjIyEUqlEQkICDh8+3GbZSZMmQRCEVq/p06ebysybN6/V/pSUFEecCrmousZmnNboAFhO/AAgNvz6dmOXMBERkbPxkDqAnTt3Ij09HRs2bEBCQgJef/11JCcn4/Tp0wgODm5V/qOPPkJjY6PpfWVlJWJjY/Hggw+alUtJScGWLVtM7xUKhf1Oglxe4UUdDCIQ4q9AqEppsUwsW/yIiMjJSZ74vfbaa5g/fz7S0tIAABs2bMBnn32GzZs3Y+nSpa3K9+7d2+z9jh074OPj0yrxUygUCA0N7VQMDQ0NaGhoML3X6XTWnga5OEsLN7cUG3H90W1F5TWoaWhGL4Xkv15ERERmJO3qbWxsxJEjR5CUlGTaJpPJkJSUhLy8vE7VsWnTJsyePRu+vr5m23NzcxEcHIyhQ4diwYIFqKysbLOOVatWQaVSmV4RERFdOyFyWQUWFm5uKdhPiX4B3hBF4Bhb/YiIyAlJmvhVVFRAr9cjJCTEbHtISAg0Gk2Hxx8+fBiFhYV47LHHzLanpKRg+/btyMnJwZo1a3DgwAFMmzYNer3eYj3Lli2DVqs1vc6fP9/1kyKXZEz84m6M42uLsdXvxxsTQYiIiJxJj+6L2rRpE2JiYjBu3Diz7bNnzzb9HBMTg9GjRyM6Ohq5ubmYMmVKq3oUCgXHAFKbLlc34GLVNQgCEBOuardsXEQAPj+u4ULORETklCRt8QsKCoJcLkdZWZnZ9rKysg7H59XW1mLHjh149NFHO/ycgQMHIigoCMXFxd2Kl9yTMYkb1LcX/JSe7ZY1zexl4kdERE5I0sTPy8sL8fHxyMnJMW0zGAzIyclBYmJiu8fu2rULDQ0NeOSRRzr8nAsXLqCyshJhYWHdjpncT3sLN7c0qp8KMgHQ6Oqh0dbbNzAiIiIrSb6OX3p6OjZu3Iht27bh5MmTWLBgAWpra02zfOfMmYNly5a1Om7Tpk1ITU1Fnz59zLbX1NTg2WefxcGDB3HmzBnk5ORgxowZGDRoEJKTkx1yTuRaCjoxo9fIV+GBISF+ALisCxEROR/Jx/jNmjULly9fRkZGBjQaDeLi4pCdnW2a8HHu3DnIZOb56enTp/Htt9/iyy+/bFWfXC7HsWPHsG3bNlRVVUGtVmPq1KlYuXIlx/GR1QwGsVNLudwqLiIApzTVKDhfheSRnVtSiIiIyBEEURRFqYNwNjqdDiqVClqtFv7+/lKHQxL65XINJv/tABQeMhSuSIanvONG8h2Hz2HpR8cxProP/jH/DgdESURE7q6zuYvkXb1EzszYzTuqn6pTSR9wcyzgsQta6A38u4qIiJwHEz+idhi7eWM7WL/vVoODe8HbU46ahmb8crnGPoERERF1ARM/onYUXLi+EHNc/4BOH+Mhl5nW++OyLkRE5EyY+BG1oaFZj5OXrj+3uaMndrRknAjCxI+IiJwJEz+iNpwsrUaj3oDevl6I6O1t1bHGrmEu6UJERM6EiR9RG26O71NBEASrjjV2DZ8qrUZ9k+VnRBMRETkaEz+iNhi7aTvzxI6W1Colgnop0GwQceKS1raBERERdRETP6I2/NiNxE8QBMRFGCd4MPEjIiLnwMSPyAJtXRN+qagFYP3EDiPjBI8fOcGDiIicBBM/IguMkzIG9PFBoK9Xl+qI5cxeIiJyMkz8iCzoysLNLY2+cey5K3W4UtvY/aCIiIi6iYkfkQXGFr+4LozvM1J5e2JgX1+z+oiIiKTExI+oBVEUuzWj91bG8YEF56q6VQ8REZEtMPEjauFi1TVU1DTCQyZgpNq/W3UZE0e2+BERkTNg4kfUwo83ll8ZHuYPpae8W3XdOrNXFMXuhkZERNQtTPyIWig4fxUAEHtjHb7uGBbmBy+5DFfrmnDuSl236yMiIuoOqxO/pqYmeHh4oLCw0B7xEEnO2OLXnRm9RgoPOYbf6C7msi5ERCQ1qxM/T09P9O/fH3o9nz9KrqdZb8Dxi9cTvzE3nrfbXWNM3b18ggcREUmrS129f/nLX/D888/jypUrto6HSFL/LavBtSY9/BQeGBjUyyZ1xpoe3XbVJvURERF1lUdXDlq3bh2Ki4uhVqsxYMAA+Pr6mu3Pz8+3SXBEjmacfTs6QgWZTLBJnXERgQCAwks6NOkN8JRzaC0REUmjS4lfamqqjcMgcg7G9fZsMb7PKLKPD/yVHtDVN+NUaTViwrs/aYSIiKgrupT4ZWZm2joOIqdgbPHr7sLNtxIEAbERAfimqAIFF6qY+BERkWS6lPgZHTlyBCdPngQAjBw5EmPGjLFJUERSqG1oxn/LqgHcnJBhK2NuJH4/nq/C7+4YYNO6iYiIOqtLiV95eTlmz56N3NxcBAQEAACqqqpw1113YceOHejbt68tYyRyiOMXtTCIQJhKiWB/pU3rNrYgckkXIiKSUpdGmS9atAjV1dU4ceIErly5gitXrqCwsBA6nQ5/+tOfbB0jkUP8aHw+rw3H9xkZE7+fL9dAV99k8/qJiIg6o0uJX3Z2Nt5++20MHz7ctG3EiBHIysrCF198YXV9WVlZiIyMhFKpREJCAg4fPtxm2a1bt0IQBLOXUmneOiOKIjIyMhAWFgZvb28kJSWhqKjI6rjIvRjH98XZaP2+WwX1UiA80BuiCBRe4Hp+REQkjS4lfgaDAZ6enq22e3p6wmAwWFXXzp07kZ6ejszMTOTn5yM2NhbJyckoLy9v8xh/f3+UlpaaXmfPnjXbv3btWrz55pvYsGEDDh06BF9fXyQnJ6O+vt6q2Mi92GNG762MrX5H2d1LREQS6VLiN3nyZDz11FO4dOmSadvFixfx9NNPY8qUKVbV9dprr2H+/PlIS0vDiBEjsGHDBvj4+GDz5s1tHiMIAkJDQ02vkJAQ0z5RFPH666/jhRdewIwZMzB69Ghs374dly5dwu7du60+V3IP5bp6XNLWQxBgt1m3cTcSyh+Z+BERkUS6lPitW7cOOp0OkZGRiI6ORnR0NKKioqDT6fDWW291up7GxkYcOXIESUlJNwOSyZCUlIS8vLw2j6upqcGAAQMQERGBGTNm4MSJE6Z9JSUl0Gg0ZnWqVCokJCS0WWdDQwN0Op3Zi9zLjze6X4cE+6GXoluT3dtk7EI2dikTERE5Wpf+hYuIiEB+fj727duHU6dOAQCGDx9ulmx1RkVFBfR6vVmLHQCEhISY6m1p6NCh2Lx5M0aPHg2tVotXX30V48ePx4kTJxAeHg6NRmOqo2Wdxn0trVq1CitWrLAqdnItxsepGR+vZg8j1f6QywSU6RpQqr2GMJW33T6LiIjIEqsTv6amJnh7e6OgoAB333037r77bnvE1abExEQkJiaa3o8fPx7Dhw/HO++8g5UrV3apzmXLliE9Pd30XqfTISIiotuxUs/x4/nrLX62XLi5JR8vDwwJ8cPJUh1+PF/FxI+IiBzO6q5eT09P9O/fH3q9vtsfHhQUBLlcjrKyMrPtZWVlCA0N7XQ8Y8aMQXFxMQCYjrOmToVCAX9/f7MXuQ+DQbw5o9eOid+t9Rec58xeIiJyvC6N8fvLX/6C559/HleuXOnWh3t5eSE+Ph45OTmmbQaDATk5OWateu3R6/U4fvw4wsLCAABRUVEIDQ01q1On0+HQoUOdrpPcyy8Vtaiub4bSU4YhIX52/ay4G13Jxq5lIiIiR+rSGL9169ahuLgYarUaAwYMgK+vr9n+/Pz8TteVnp6OuXPnYuzYsRg3bhxef/111NbWIi0tDQAwZ84c9OvXD6tWrQIAvPTSS7jjjjswaNAgVFVV4ZVXXsHZs2fx2GOPAbg+43fx4sX461//isGDByMqKgrLly+HWq1GampqV06XXJxxlu0otQqe8i79LdRpxq7k4xe00BtEyGWCXT+PiIjoVl1K/GyZQM2aNQuXL19GRkYGNBoN4uLikJ2dbZqcce7cOchkN/8xvnr1KubPnw+NRoPAwEDEx8fj+++/x4gRI0xlnnvuOdTW1uLxxx9HVVUVfvWrXyE7O7vVQs9EABzWzQsAg4P94OMlR22jHj9frrF7CyMREdGtBFEURWsOaG5uxssvv4zf//73CA8Pt1dcktLpdFCpVNBqtRzv5wbuXfctjl3Q4q2HxuCeWLXdP2/WO3k4VHIFa387GjNv5yQiIiLqvs7mLlb3a3l4eOCVV15Bc3NztwIkcgb1TXqcLL2+bqMjWvxu/ZwCrudHREQO1uUndxw4cMDWsRA53MlSHZr0Ivr4eiE80DHLqxgTPz7Bg4iIHK1LY/ymTZuGpUuX4vjx44iPj281uePee++1SXBE9lZwI/mKjQiAIDhmooVxgscpTTWuNerh7SV3yOcSERF1KfH74x//COD6c3ZbEgTBJmv8ETmCsdUt9sZzdB0hTKVEXz8FLlc34MQlLcZG9nbYZxMRkXvrUlevwWBo88Wkj3oS4zN6jc/RdQRBEG5ZyLnKYZ9LRERkVeL3m9/8BlrtzScOrF69GlVVVab3lZWVZsuqEDmzqrpGlFTUAgBiw+33jF5LmPgREZEUrEr89uzZg4aGBtP7l19+2ezpHc3NzTh9+rTtoiOyI2NrX1SQLwJ8vBz62aYJHpzZS0REDmRV4tdyyT8rlwAkcioF56oAOL61DwBibnzm+SvXUFnT0EFpIiIi27Dv86mInJixtS3WQev33cpf6Ynovr5mcRAREdmbVYmfIAitlrxw1BIYRLYkiqJpRq+jFm5uKS4iEABQcF7bQUkiIiLbsGo5F1EUMW/ePCgUCgBAfX09/vCHP5jW8bt1/B+RM7tw9RoqaxvhKRcwPEyax/LFRajwr/wLnOBBREQOY1XiN3fuXLP3jzzySKsyc+bM6V5ERA5gTLaGh/lD6SnNAsqxtzzBQxRFtp4TEZHdWZX4bdmyxV5xEDmU1N28ADAs1B9eHjJorzXhbGUdIoN8Oz6IiIioGzi5g9xSgQRP7GjJy0OGkWp/s3iIiIjsiYkfuZ0mvQGFl65PqJBiRu+tjIknEz8iInIEJn7kdv5bVo36JgP8lB4YKHH36pgbj4rjki5EROQITPzI7dzazSuTSTuhwtjid+KSDo3NBkljISIi18fEj9yOcWJHbITjn9jR0oA+Pgjw8URjswGnNDqpwyEiIhfHxI/czo83Fkw2LqAsJUEQTK1+P3KcHxER2RkTP3IrNQ3N+G95NQBpntFriXGCyVEmfkREZGdM/MitHL+ghSgCapUSwf5KqcMBcP0JHgBb/IiIyP6Y+JFbMc6ejbsxm9YZGLt6f75cC119k7TBEBGRS2PiR26l4FwVAGkXbm6pTy8FInp7AwCO3Rh/SEREZA9M/MitGFv8pF64uSXTBA+u50dERHbExI/cRpmuHqXaesgEIKafc0zsMDI+M5hP8CAiInti4kduw5hUDQnxg6/CQ9pgWrg18RNFUdpgiIjIZTlF4peVlYXIyEgolUokJCTg8OHDbZbduHEjJkyYgMDAQAQGBiIpKalV+Xnz5kEQBLNXSkqKvU+DnNyPtzyxw9mMVKsglwm4XN2AUm291OEQEZGLkjzx27lzJ9LT05GZmYn8/HzExsYiOTkZ5eXlFsvn5ubioYcewv79+5GXl4eIiAhMnToVFy9eNCuXkpKC0tJS0+uDDz5wxOmQEzO2+DnTjF4jby85hoX6AWB3LxER2Y/kid9rr72G+fPnIy0tDSNGjMCGDRvg4+ODzZs3Wyz//vvv449//CPi4uIwbNgwvPvuuzAYDMjJyTErp1AoEBoaanoFBkr/lAaSjsEg4tiF6zNmnbHFD7g54YTr+RERkb1Imvg1NjbiyJEjSEpKMm2TyWRISkpCXl5ep+qoq6tDU1MTevfubbY9NzcXwcHBGDp0KBYsWIDKyso262hoaIBOpzN7kWv5paIGNQ3N8PaUY0hIL6nDsYgTPIiIyN4kTfwqKiqg1+sREhJitj0kJAQajaZTdSxZsgRqtdoseUxJScH27duRk5ODNWvW4MCBA5g2bRr0er3FOlatWgWVSmV6RUREdP2kyCkdvbF+X0w/FTzkkjd0W2RM/I5f1EJv4AQPIiKyPeea2mil1atXY8eOHcjNzYVSefPxW7Nnzzb9HBMTg9GjRyM6Ohq5ubmYMmVKq3qWLVuG9PR003udTsfkz8XcXL/PuZZxuVV0317w9ZKjtlGPovJqDAv1lzokIiJyMZI2fQQFBUEul6OsrMxse1lZGUJDQ9s99tVXX8Xq1avx5ZdfYvTo0e2WHThwIIKCglBcXGxxv0KhgL+/v9mLXMuPN56IERfhvGM95TIBo40LObO7l4iI7EDSxM/Lywvx8fFmEzOMEzUSExPbPG7t2rVYuXIlsrOzMXbs2A4/58KFC6isrERYWJhN4qaepb5Jj5Ol18dtOnOLH3BzggfH+RERkT1IPtgpPT0dGzduxLZt23Dy5EksWLAAtbW1SEtLAwDMmTMHy5YtM5Vfs2YNli9fjs2bNyMyMhIajQYajQY1NTUAgJqaGjz77LM4ePAgzpw5g5ycHMyYMQODBg1CcnKyJOdI0jpxSYdmg4igXl7oF+AtdTjtiruRmBbwmb1ERGQHko/xmzVrFi5fvoyMjAxoNBrExcUhOzvbNOHj3LlzkMlu5qfr169HY2MjHnjgAbN6MjMz8eKLL0Iul+PYsWPYtm0bqqqqoFarMXXqVKxcuRIKhcKh50bOwdhtGhcRAEEQpA2mA8au6P+WVaOusRk+XpL/ihIRkQtxin9VFi5ciIULF1rcl5uba/b+zJkz7dbl7e2NPXv22CgycgUFTvzEjpZCVUqE+CtQpmtA4UUdxkX17vggIiKiTpK8q5fI3m7O6A2QNI7OiuUEDyIishMmfuTSrtY24mxlHYCe0eIH3HykXMGNhJWIiMhWmPiRSzMmTwODfKHy8ZQ2mE6Ku5GgFtxYdJqIiMhWmPiRSzN2l/aUbl4AiAlXQRCAi1XXcLm6QepwiIjIhTDxI5d264zensJP6YlBfa8/T/gYu3uJiMiGmPiRyxJF8eaM3h6U+AFcyJmIiOyDiR+5rPNXruFqXRO85DIMD/OTOhyrMPEjIiJ7YOJHLss4sWO42h8KD7m0wVhpzI3E78fzVRBFUdpgiIjIZTDxI5dlnBUbF+7cz+e1ZGioH7w8ZNDVN6OkolbqcIiIyEUw8SOX1dMWbr6Vp1yGUWp/ADfPg4iIqLuY+JFLatIbUHhRC6Bnzei9lfG5vVzPj4iIbIWJH7mk05pqNDQb4K/0QGQfX6nD6ZLYiOtd1AUXtBJHQkREroKJH7mkW5dxkckEaYPpojE3WvxOXtKhoVkvcTREROQKmPiRSyrogQs3txTR2xuBPp5o1BtwsrRa6nCIiMgFMPEjl2R6VNuN5972RIIgmCam/Mj1/IiIyAaY+JHLqa5vQvHlGgA9c0bvreKY+BERkQ0x8SOXc/yCFqII9AvwRl8/hdThdAuf4EFERLbExI9cjvGJHT15fJ+Rsav6l4paaOuapA2GiIh6PCZ+5HJ+dIGJHUa9fb0woI8PAODYxSppgyEioh6PiR+5nFuXcnEFxlY/LuRMRETdxcSPXIpGW48yXQPkMgGj+vlLHY5NmGb28tFtRETUTUz8yKUYW/uGhPjBx8tD2mBsJM40wUMLURSlDYaIiHo0Jn7kUm4u3KySNhAbGqn2h4dMQEVNAy5WXZM6HCIi6sGY+JFLcYWFm1tSesoxLMwPAPDjeT63l4iIuo6JH7kMvUHE8YvXE6O4/gHSBmNjcRznR0RENsDEj1zGz5drUNPQDB8vOQYH+0kdjk1xZi8REdmCUyR+WVlZiIyMhFKpREJCAg4fPtxu+V27dmHYsGFQKpWIiYnB559/brZfFEVkZGQgLCwM3t7eSEpKQlFRkT1PgZyAcXzfqH4qyGWCtMHYmLHF7/hFLZr1BmmDISKiHkvyxG/nzp1IT09HZmYm8vPzERsbi+TkZJSXl1ss//333+Ohhx7Co48+iqNHjyI1NRWpqakoLCw0lVm7di3efPNNbNiwAYcOHYKvry+Sk5NRX1/vqNMiCRjH941xkfX7bhXdtxd6KTxwrUmPovIaqcMhIqIeShAlXh8iISEBt99+O9atWwcAMBgMiIiIwKJFi7B06dJW5WfNmoXa2lp8+umnpm133HEH4uLisGHDBoiiCLVajT//+c945plnAABarRYhISHYunUrZs+e3arOhoYGNDQ0mN7rdDpERERAq9XC398+a8H937wzyPul0i51u6tDv1xBZW0j3n74NvwmJkzqcGzufzcexPc/V2JM/wCEqZRSh0NERFb433ED8KvBQXarX6fTQaVSdZi7SLrQWWNjI44cOYJly5aZtslkMiQlJSEvL8/iMXl5eUhPTzfblpycjN27dwMASkpKoNFokJSUZNqvUqmQkJCAvLw8i4nfqlWrsGLFChucUecdu6DF58c1Dv1Md+AhE3Bb/0Cpw7CLhKg++P7nShw9V4WjUgdDRERWmTC4r9QhAJA48auoqIBer0dISIjZ9pCQEJw6dcriMRqNxmJ5jUZj2m/c1laZlpYtW2aWTBpb/OzpvjH9MDrcddaacxZDQvwQ6qKtYfN/HYXwQG/UNTZLHQoREVlp7ADnaJRwjUcbdJNCoYBCoXDoZ44fFITxg+zX5Euux8fLA7+ND5c6DCIi6sEkndwRFBQEuVyOsrIys+1lZWUIDQ21eExoaGi75Y3/taZOIiIiIncgaeLn5eWF+Ph45OTkmLYZDAbk5OQgMTHR4jGJiYlm5QFg7969pvJRUVEIDQ01K6PT6XDo0KE26yQiIiJyB5J39aanp2Pu3LkYO3Ysxo0bh9dffx21tbVIS0sDAMyZMwf9+vXDqlWrAABPPfUUJk6ciL/97W+YPn06duzYgR9++AF///vfAQCCIGDx4sX461//isGDByMqKgrLly+HWq1GamqqVKdJREREJDnJE79Zs2bh8uXLyMjIgEajQVxcHLKzs02TM86dOweZ7GbD5Pjx4/GPf/wDL7zwAp5//nkMHjwYu3fvxqhRo0xlnnvuOdTW1uLxxx9HVVUVfvWrXyE7OxtKpWsO+iciIiLqDMnX8XNGWq0WAQEBOH/+vN3W8SMiIiKyFeOKJFVVVVCp2l41RPIWP2dUXV0NAHZf0oWIiIjIlqqrq9tN/NjiZ4HBYMClS5fg5+cHQbDPM1+NmTlbFVvjtWkbr41lvC5t47WxjNelbbw2ljn7dRFFEdXV1VCr1WZD5Fpii58FMpkM4eGOWS/N39/fKW8gZ8Br0zZeG8t4XdrGa2MZr0vbeG0sc+br0l5Ln5Gky7kQERERkeMw8SMiIiJyE0z8JKJQKJCZmenwR8X1BLw2beO1sYzXpW28NpbxurSN18YyV7kunNxBRERE5CbY4kdERETkJpj4EREREbkJJn5EREREboKJHxEREZGbYOJnR1lZWYiMjIRSqURCQgIOHz7cbvldu3Zh2LBhUCqViImJweeff+6gSB1n1apVuP322+Hn54fg4GCkpqbi9OnT7R6zdetWCIJg9lIqlQ6K2HFefPHFVuc5bNiwdo9xh3smMjKy1XURBAFPPvmkxfKufL98/fXXuOeee6BWqyEIAnbv3m22XxRFZGRkICwsDN7e3khKSkJRUVGH9Vr7XeVs2rsuTU1NWLJkCWJiYuDr6wu1Wo05c+bg0qVL7dbZld9HZ9TRPTNv3rxW55mSktJhva58zwCw+J0jCAJeeeWVNuvsKfcMEz872blzJ9LT05GZmYn8/HzExsYiOTkZ5eXlFst///33eOihh/Doo4/i6NGjSE1NRWpqKgoLCx0cuX0dOHAATz75JA4ePIi9e/eiqakJU6dORW1tbbvH+fv7o7S01PQ6e/asgyJ2rJEjR5qd57fffttmWXe5Z/7zn/+YXZO9e/cCAB588ME2j3HV+6W2thaxsbHIysqyuH/t2rV48803sWHDBhw6dAi+vr5ITk5GfX19m3Va+13ljNq7LnV1dcjPz8fy5cuRn5+Pjz76CKdPn8a9997bYb3W/D46q47uGQBISUkxO88PPvig3Tpd/Z4BYHY9SktLsXnzZgiCgN/+9rft1tsj7hmR7GLcuHHik08+aXqv1+tFtVotrlq1ymL5mTNnitOnTzfblpCQID7xxBN2jVNq5eXlIgDxwIEDbZbZsmWLqFKpHBeURDIzM8XY2NhOl3fXe+app54So6OjRYPBYHG/u9wvAMSPP/7Y9N5gMIihoaHiK6+8YtpWVVUlKhQK8YMPPmizHmu/q5xdy+tiyeHDh0UA4tmzZ9ssY+3vY09g6drMnTtXnDFjhlX1uOM9M2PGDHHy5Mntlukp9wxb/OygsbERR44cQVJSkmmbTCZDUlIS8vLyLB6Tl5dnVh4AkpOT2yzvKrRaLQCgd+/e7ZarqanBgAEDEBERgRkzZuDEiROOCM/hioqKoFarMXDgQDz88MM4d+5cm2Xd8Z5pbGzEe++9h9///vcQBKHNcu5yv9yqpKQEGo3G7J5QqVRISEho857oyneVK9BqtRAEAQEBAe2Ws+b3sSfLzc1FcHAwhg4digULFqCysrLNsu54z5SVleGzzz7Do48+2mHZnnDPMPGzg4qKCuj1eoSEhJhtDwkJgUajsXiMRqOxqrwrMBgMWLx4Me68806MGjWqzXJDhw7F5s2b8cknn+C9996DwWDA+PHjceHCBQdGa38JCQnYunUrsrOzsX79epSUlGDChAmorq62WN4d75ndu3ejqqoK8+bNa7OMu9wvLRn/v1tzT3Tlu6qnq6+vx5IlS/DQQw/B39+/zXLW/j72VCkpKdi+fTtycnKwZs0aHDhwANOmTYNer7dY3h3vmW3btsHPzw/3339/u+V6yj3jIXUA5L6efPJJFBYWdjgGIjExEYmJiab348ePx/Dhw/HOO+9g5cqV9g7TYaZNm2b6efTo0UhISMCAAQPw4YcfduovTXewadMmTJs2DWq1us0y7nK/kPWampowc+ZMiKKI9evXt1vWXX4fZ8+ebfo5JiYGo0ePRnR0NHJzczFlyhQJI3MemzdvxsMPP9zhJLGecs+wxc8OgoKCIJfLUVZWZra9rKwMoaGhFo8JDQ21qnxPt3DhQnz66afYv38/wsPDrTrW09MTY8aMQXFxsZ2icw4BAQEYMmRIm+fpbvfM2bNnsW/fPjz22GNWHecu94vx/7s190RXvqt6KmPSd/bsWezdu7fd1j5LOvp9dBUDBw5EUFBQm+fpTvcMAHzzzTc4ffq01d87gPPeM0z87MDLywvx8fHIyckxbTMYDMjJyTFribhVYmKiWXkA2Lt3b5vleypRFLFw4UJ8/PHH+OqrrxAVFWV1HXq9HsePH0dYWJgdInQeNTU1+Pnnn9s8T3e5Z4y2bNmC4OBgTJ8+3arj3OV+iYqKQmhoqNk9odPpcOjQoTbvia58V/VExqSvqKgI+/btQ58+fayuo6PfR1dx4cIFVFZWtnme7nLPGG3atAnx8fGIjY21+linvWeknl3iqnbs2CEqFApx69at4k8//SQ+/vjjYkBAgKjRaERRFMXf/e534tKlS03lv/vuO9HDw0N89dVXxZMnT4qZmZmip6enePz4calOwS4WLFggqlQqMTc3VywtLTW96urqTGVaXpsVK1aIe/bsEX/++WfxyJEj4uzZs0WlUimeOHFCilOwmz//+c9ibm6uWFJSIn733XdiUlKSGBQUJJaXl4ui6L73jChenzXYv39/ccmSJa32udP9Ul1dLR49elQ8evSoCEB87bXXxKNHj5pmp65evVoMCAgQP/nkE/HYsWPijBkzxKioKPHatWumOiZPniy+9dZbpvcdfVf1BO1dl8bGRvHee+8Vw8PDxYKCArPvnYaGBlMdLa9LR7+PPUV716a6ulp85plnxLy8PLGkpETct2+feNttt4mDBw8W6+vrTXW42z1jpNVqRR8fH3H9+vUW6+ip9wwTPzt66623xP79+4teXl7iuHHjxIMHD5r2TZw4UZw7d65Z+Q8//FAcMmSI6OXlJY4cOVL87LPPHByx/QGw+NqyZYupTMtrs3jxYtN1DAkJEX/zm9+I+fn5jg/ezmbNmiWGhYWJXl5eYr9+/cRZs2aJxcXFpv3ues+Ioiju2bNHBCCePn261T53ul/2799v8ffHeP4Gg0Fcvny5GBISIioUCnHKlCmtrtmAAQPEzMxMs23tfVf1BO1dl5KSkja/d/bv32+qo+V16ej3sado79rU1dWJU6dOFfv27St6enqKAwYMEOfPn98qgXO3e8bonXfeEb29vcWqqiqLdfTUe0YQRVG0a5MiERERETkFjvEjIiIichNM/IiIiIjcBBM/IiIiIjfBxI+IiIjITTDxIyIiInITTPyIiIiI3AQTPyIiIiI3wcSPiIiIyE0w8SMisqF58+YhNTVV6jCIiCzykDoAIqKeQhCEdvdnZmbijTfeAB+IRETOiokfEVEnlZaWmn7euXMnMjIycPr0adO2Xr16oVevXlKERkTUKezqJSLqpNDQUNNLpVJBEASzbb169WrV1Ttp0iQsWrQIixcvRmBgIEJCQrBx40bU1tYiLS0Nfn5+GDRoEL744guzzyosLMS0adPQq1cvhISE4He/+x0qKiocfMZE5GqY+BER2dm2bdsQFBSEw4cPY9GiRViwYAEefPBBjB8/Hvn5+Zg6dSp+97vfoa6uDgBQVVWFyZMnY8yYMfjhhx+QnZ2NsrIyzJw5U+IzIaKejokfEZGdxcbG4oUXXsDgwYOxbNkyKJVKBAUFYf78+Rg8eDAyMjJQWVmJY8eOAQDWrVuHMWPG4OWXX8awYcMwZswYbN68Gfv378d///tfic+GiHoyjvEjIrKz0aNHm36Wy+Xo06cPYmJiTNtCQkIAAOXl5QCAH3/8Efv377c4XvDnn3/GkCFD7BwxEbkqJn5ERHbm6elp9l4QBLNtxtnCBoMBAFBTU4N77rkHa9asaVVXWFiYHSMlIlfHxI+IyMncdttt+Ne//oXIyEh4ePBrmohsh2P8iIiczJNPPokrV67goYcewn/+8x/8/PPP2LNnD9LS0qDX66UOj4h6MCZ+RERORq1W47vvvoNer8fUqVMRExODxYsXIyAgADIZv7aJqOsEkUvMExEREbkF/ulIRERE5CaY+BERERG5CSZ+RERERG6CiR8RERGRm2DiR0REROQmmPgRERERuQkmfkRERERugokfERERkZtg4kdERETkJpj4EREREbkJJn5EREREbuL/A03JXywVupAcAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axes = plt.subplots(2, 1, sharey=True)\n",
    "axes[0].plot(prediction_errors[0, :])\n",
    "axes[0].set_title('Early in learning')\n",
    "axes[1].plot(prediction_errors[-1, :])\n",
    "axes[1].set_title('Late in learning')\n",
    "plt.xlabel('Time')\n",
    "plt.ylabel('Error')\n",
    "axes[0].set_ylabel('Error')\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d390b261-65c8-4292-ba2d-136663865fd8",
   "metadata": {},
   "source": [
    "### Exercise 3.1\n",
    "Can you guess what would happen to this error if the reward was omitted? Simulate this to verify. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "369fab21-a7cd-463e-8e28-d865d33dac31",
   "metadata": {},
   "source": [
    "### Exercise 3.2\n",
    "Show what happens to the prediction error between the first and last trials. Is this a realistic model of what dopamine neurons do? "
   ]
  }
 ],
 "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.10.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}