DifferentiatingPerceptron.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import sklearn.datasets\n",
    "import sklearn.linear_model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "5hIbr52I7Z7U"
   },
   "source": [
    "# Differentiating Perceptron\n",
    "\n",
    "A very short notebook that gives an idea of how the simple neural network model known as perceptron can be trained by computing a closed form gradient expression.\n",
    "The chain rule then allows to generalize gradient calculation to arbitrary big multi layer network based on this idea."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "cellView": "both",
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "apJbCsBHl-2A"
   },
   "outputs": [],
   "source": [
    "# First, we load a toy dataset and play a bit with it\n",
    "data, label = sklearn.datasets.load_digits(n_class=10, return_X_y=True)\n",
    "nb_train = 1200\n",
    "nb_test = 300\n",
    "ids = np.random.choice(data.shape[0],nb_train+nb_test,replace=False) \n",
    "train_dataset = data[ids[:nb_train],:]\n",
    "train_labels = label[ids[:nb_train]]\n",
    "test_dataset = data[ids[nb_train:],:]\n",
    "test_labels = label[ids[nb_train:]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/gnthibault/anaconda3/lib/python3.6/site-packages/sklearn/linear_model/sag.py:326: ConvergenceWarning: The max_iter was reached which means the coef_ did not converge\n",
      "  \"the coef_ did not converge\", ConvergenceWarning)\n"
     ]
    }
   ],
   "source": [
    "#try a simple linear model\n",
    "def trainModel(nSamples,data,label):\n",
    "  # skl is wants data to be: dimension 0: list of samples, dimension 1: features\n",
    "  X = data[:nSamples,:]\n",
    "  y = label[:nSamples]\n",
    "\n",
    "  # run the classifier\n",
    "  model = sklearn.linear_model.LogisticRegression(multi_class='multinomial',solver='saga')\n",
    "  model.fit(X, y)\n",
    "  return model\n",
    "\n",
    "def testModel(model, data, label):\n",
    "  X = data.reshape(-1, np.prod(data.shape[1:3]))\n",
    "  test=label\n",
    "  correct = model.predict(X) == label\n",
    "  return np.count_nonzero(correct)/correct.size\n",
    "    \n",
    "# plot result function of training samples\n",
    "nbSamplesConfig=[50, 100, 500, 1200]\n",
    "\n",
    "scores = np.zeros(len(nbSamplesConfig))\n",
    "for idx, nbSample in enumerate(nbSamplesConfig):\n",
    "  model = trainModel(nbSample,train_dataset,train_labels)\n",
    "  scores[idx]= testModel(model,test_dataset,test_labels)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "colab_type": "text",
    "id": "L8oww1s4JMQx"
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "No handles with labels found to put in legend.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcAAAAEdCAYAAACbhbeCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nOzdd5xcdb3/8dc7m7rpDbIkpIcmSL10SUBULhaK6A0CEkWKgHoV708QxBhRruWKDZAihF5EwQhREE1ApCVI6MRsetlAetu03f38/vh+Z3Mymdk9u9kys/t5Ph7zmJlTv2fOmfnM93u+RWaGc8451950aO0EOOecc63BA6Bzzrl2yQOgc865dskDoHPOuXbJA6Bzzrl2yQOgc865dqnJAqCkCZIsPvbJMX9sYv7JTbTP4XF7Exqx7nRJ05siHa5wJa7L4c2w7XGSJkoqyj+SkhZIure101EfSftJ+ruk9fFcnp5nuWY7H3G7jWoz5r816UiaLGlBS+6zOb64G4Dzckw/P85zriU9ARwDVDTDtscB38VLUprbz4CRwGcJ5/KZPMuNo/nOx+1x341xaXy4AtMcF8ofgHMlKTNBUjfgLOD3zbA/B0jq0tppKERmtsLMXjSzra2dlvaoia7L/YFnzewv8VyuaYJ0dUr+RtXHzJaY2YuN2ZeZvW1mbzdmXde8miMA3gMMA45PTDsj7itnAJR0rqTXJG2RtFLSPZLKspYplXSTpFWSNkqaAgzJs72xkv4maYOkTZKelHRgYw5G0vck/SsWv6yMRTFH51huYEzfYklb4/M9yR8ASQdLejQew2ZJsyVdlZi/QNLkHNs2SRMT7yfGaQfGY9sIPBznfVTSVEkVkiolvSnpCkklObZ7YTy2zZLWSHpG0rGSukhaIemGHOtkihT3a+DnOFDSLZL+HdO1WNL9kgbnWPZsSe/G6+ENSZ/KLkaS1FXSDfH4NkpaLulP2enKVQSaKfqTNF7SO/EamSnp+Kx1/0PSXxPna56kmzLngJDbANge91FnEVlc5jpJX5U0P16fz0j6QNZyDb0O9ovXwSZJiyR9Ic4/L36OGyVNkzQqT7oulFQeP+9/SToxxzL1fqfiOXpO0iclvSppK3XkfBSC0HXxeLfF5+skdYrzx8XPdDhwXl2fcV3nQztulVwq6ceSlgFbgT5pr0vlKAJtwPnMvnbHxXU/JenXCr8rK+M12Sdr3YGSHlD4/Vkj6c64nkkal++zjevmvX4T20597I29zrTj+1bvdZbjGEol/Sh+vtvi89VKFHNL6iHpVzFNWyW9L+lppfiN6ljfAo2wEHiWUAz6jzjt88CjwMbshSVdBNwCPARcBewF/BA4StJhZpZZ5xbgv4DvATOAjwD359jex4E/Eoq+zo2TvwX8Q9IHzWxxA49nMHADsAToHrf5rKTDzeyNuM++wPNAP+A64HVgD+A0oDOwVdKRwHSgHPh63N4Y4IMNTE/SH4HfAj8CauK0kcDfgF8BW4AjgInAQODKzIqSfgpcEdf/blz/aGComT0v6U7gAklXmdmWxD4vBp4xs3cbmNZ+MT1XASsI5/kK4J+S9svsQ9JHgPuAKcA3Yrp/DnQF/p3YXhegJ+HzrojbvxR4QdL+Zra8nvR8CNgX+E5M1/eBxyUNN7O1knoATwIvAxMIxffDgWPj+rcT/oBdQPizV53yczgXmA18jXBt/AT4Y/wMqlJuI9vvgNuAnxI+gzskjSEUCV4JdAJ+Qfi+HJW17jjgcOBqQlD4FvBnSQeb2Wxo8HdqH+CXhM9zHrC6jnTfRSjW/CHwHOGzvZpwDX8O+Beh2HEK4Tv//Tq2leZ8XB23cxFQQjjvQ0lxXdZhd87nL4DH47HuC/w4pvv8xDJ/AA6K6SsHPk34btcpxfULKb+TCc16neU4ho7xGA4gnPs3CL9R34lpvyIuegPwKeDbwBygP3Ac0If6mFmTPAgfsgGjgS8Cawg/WmVAFSFgjYvLnBzXKQHeA6Zlbev4uNxX4/t9CRfGlVnL3RyXm5CYVg78LWu5XsBK4OeJadOB6Q08xhLCn4bZwC8S0yfF9B1ax7rPAouB0jqWWQBMzjHdgImJ9xPjtK/Vk17F9F4dz0eHOH10TO/P6lh3ZFzmvMS0D8b9jm+C66UE2Dtu74zE9OeBNwElph0el8t7vuL2Sglf9K/nuC6HZ33Oa4C+iWlHxOU+l/X+g3XsM3MeOqY8ZiN8QTslpp0Vpx+7G9fB5xPT+hK+b6uAXonpX43LDsvazzZg78S0noSgdU8jv1M1wCEpPosDs48nTr8m+3Mn/Fnc5fNIez4IP/xGCKiqZxv5rsuJgDXyfE5PXrvs+B28K2t7vyYEJMX3H43LfTZruSlx+rg6jqPe67ehx94C19lkYEHi/XlxWydkpfPquL094vs3qeO3rK5Hc928/x3hH/ongXOA5YRcSbZ9CTml+5ITzew5Qk5ybJx0FKEI9eGs9R9Mvon/RkYB90nqmHkAlcALwAkNPRBJJ8ds/SrCCd9O+Je7b2KxjwIzzOzVPNsoJfwjuc/MKhuahjo8mmNfZbFYYyHhItlOyCX1IXzWACcTPs9b823YzOYR/n1dnJh8MeGf4h8ak1hJX1Yo6t5I+CwXxVn7xvklhC/u7y1e2TEtrwDzc2zvs5JekrQ2bm8T0IOdz00+L9jO95LeiM9D4/McYC1wi0IR/d5pj7MefzWz7XXstzH+nHkRj+l94EUzW59YJpNjzz6OFy2RgzOzDeyoONSY79QCM5uVIs2Z9bJroWbej6XpPZa8rjLquy7rsTvn84ms928Qfjf3jO+PJvwJzf6eP5Ji26mu3wYee7NdZ3mcQogDz2dde08RcpuZW1EzgAmSvi3pCOW43ZNPswTAeHCPESL45wk//DU5Fu0Xn3PV0FuemJ+5H/he1jLZ7zM/8L8l/PAnH58gZI1Tk3QYMJVQdHsB4QP/D+A1Qu42oz/hX2o+fQmfdV3LNMZOn1ssF59CONbrgJMI6f1BXCST5sznUF96bgKOU7jXmCn+vdPMtjU0oZK+Erf3NHAmcCQ7LuBMugYQLuz3c2xip3Mt6ZOEYvN3CEVIRxGOdQU7n5t8diqasx2VZLrG9+uAE4FlMd2LFO43fjrFtlPvl1AcVLvfRsquFLItz7Rc+8n+DmWmZe4DNfQ7lba2bb7v/vKs+U1pl7SlvC7rsjvns751y4A1WQEWcp+znaS5fhtx7M15neWyB6E+SfZ193Kcn7n2vkK4RfZFQjB8X6F+QGkd2waa5x5gxt2ECN8BODvPMpkLYFCOeYOAV+LrzIW7J+G+Aon3Savi81WEk5qtoT/cnyb8KzozeRHGe35rE8utpO4TuYZQNFTXMhCKPzonJ0iqK2hn/5sdRchBnWdmtf+sY7BIWhmfBxOKc/OZSii+uJgQ9HtSR66xHuMJxWiZcnskjciRru3s+NFN2pMd/04z2ys3swmJ7XWiCX84Y07m0/Ff5xGE6+rheN/izabaTw4NvQ52R/Z3KDNtaXzd0O9U2rZyye/+3MT0QVnzm1KutKW5LltLBdBXUqesIJjrnO0ixfXbksde33WWyypCyc9n88xfAGChnshVwFWShhGKof+XcG1+q65ENWf7pb8Siix/Y2Zv5VlmNuFfwPjkREnHEiL/9DjpJUIAyf4gxme9n034UD5gZjNzPF5v4DGUEoogar84kk5i1+KNp4AjJR2cayOx2PM5QvOQbnXsbyHh3kjSxxuYXghBJJPeToRi6KSnCZ/nRXVtLObabyHk5C8HnjazuXWtU0/asv/JfiFrf9XATMKXNtmM5nAg+4tZSvhzknQe4T5GkzKzKgtV4L9D+M7sH2dl/rHXdU4bY3evg4Y4Olk8Jqln3NcLcVJTf6cyno3P2d/hzLU6vRHbbMz5qPe6bEUvEq7nM7Kmf6YhG6nj+m3JY6/vOsvlL4Si1I15rr2V2SuY2UIz+z9CcXK9Nf+bLQcYf8zy5fxql5F0LaGc+l5C+f9gQpHdHOCOuNxsSfcDk2Ix3wzCfbdTs7Znki4j1MLqTAjAKwn/NI4FFpnZzxpwGH8B/huYHGtF7kO4iLL/tdxAKIZ7WtJ1hA9/AKEW6CWxSPibhAa8L0j6P0Lx40hChYGvxO08SKhZdQOhdtjBhEocab1D+PH8gaRqwsX99eyFzGxu3Mc34oU4hRDojwTeNbOHEov/lnAT/GBCjngnClWxpwFfMLPJdaTtL8C3JH2bUIRxEuGfWrbvEv5QPCrpVsLnOJFQNJYsRv8LcHriszqCUBSyliYg6ROEPwiPEf6Fdifc4N/Aji9tpm3XFZL+DFSb2cwm2P3uXgcN8R7wlEIzgkztvO7EGpfN8J0ibvdNSQ8AE2MO5XnC/aDvAA9YrGHdQI05H2mvyxZnZk9J+idwq6QBhMpIZxGuB9j5+7CTlNdvSx57nddZHvcRAvLf4m/ma4SSkVGEWp+nm1mlpBcIv2FvEG5XjSV8RnfVm6rG1JzJ9SBRC7SOZcaRqAWamH5uPLithGzvPUBZ1jKlhFqfq+NBTiFULNmpFmhc9hjCD8caQnHSAsKPyjGJZaaTohYo4Ud1PrCZEHhPzrUuodjuVkKxxTZCjc+7gC6JZQ4F/kT4kd5MuGH8rcT8DsC1hCBWSaiEMor8tf92qX0IHELIbVYSguwk4Etk1YSMy15CaLKxNX6u05OfUWK5Jwn3EnLt7+Nx26fU8zl2i+dvBeFL+DghV5erJuDnCDmPrcBbhH/ArwKPZn1W18V0VRL+XBxKVg1K8tcCvTdHGmvTQqgE8FA891tiuqcCRyWWLwFuJNyzrCGrlmCe7V+XNW04u9Zk3q3rINfxkeO7l1kuXh9z4+f9KnBSjrSn/U4914DfjM7xHC4k/FlbGN93yloubS3QnOcj8Rl/qbHXJflrgaY5n9PJXQs0+3dwArteqwPj57yB8LtxN6GZhAEH1/FZpLl+G3TszX2dkVULNE7rGvf/Ljt+p2bEaR3jMj+K21tHqAj3BrEFQX2PTHVb53KK9zsXEaq7fyfH/B8S/o0dZM10MUkaQvj3+wMzq+sfo3NtnqRfE3JG/awIejhS6N/zOTM7t75lW1pzVoJxRUzSQMK/yK8RciQ35Vl0LPDDpgp+8R7pzwj3KVcSion/HyEndHtT7MO5YqHQ0X9vQklIZ0LTgC8DPymG4FfoPAC6fD4O3EnI/Z1vZjmrt5vZcU2832pCTcBfE6o5byL0KPSZfGlwrg3bRKiHMIrQRnA+oceTn7RmotoKLwJ1zjnXLrXLHKCkVwk3l8tbOy3OOVdERgMrzOzQ1k5IU2iXARAY2Lt378GHHHJIfQ3TnXPORbNmzWLdunWtnYwm014DYPkhhxwyePr06a2dDuecKxrjxo3jmWeeaTMlZz6StXPOuXbJA6Bzzrl2yQOgc865dskDoHPOuXbJA6Bzzrl2qcUDoKRTJM2WVC7pyhzzh0n6m6TXJU2P/UBm5lVLmhUfUxLTR8SRwcslPRR7rXfOOefyatEAGIeqvxH4T+AA4GxJB2Qt9lPgbjP7IGEkg+sT8zab2SHx8anE9B8BN5jZaEJv9Rc020E455xrE1q6HeCRhFG85wFIepAwZt7biWUOAL4RX08jjGeVVxw49STCEDoQhiCaSBjmwznn2q0NW7azYGUl81dtYv6KTSxYtYl5KzfxjY/sw9h9BrZ28lpdSwfAwYRx8jKWAEdlLfMacCbwC8I4cD0l9TezVUBXSTMJI4H/r5k9Rugwea2ZVSW2mbOHF0kXEQaJ3LeiwvtVds4Vvy3bq1m4qpL5Kzcyf2V4XrCyknkrN7FyY+4BI2YvX+8BkMLsCeabwK/jMCDPEkZfr47zhpnZUkkjgb9LeoMwCGIqZnYrYXTl6WVlZWObON3OOdcstlfXsHh1ZcjBxZzc/JWbWLCykmXrNpNvTIPOHTswvH8pIwZ0Z/iA7owc0J3h/buz76CeLXsABaqlA+BSYO/E+yFxWi0zW0bIASKpB/BpM1sb5y2Nz/MkTSeMAP57oI+kjjEXuMs2nXOu0NXUGMvWbY6BLRRVLlgZAt3iNZuprskd5Uo6iKH9S2Og68GIAeF5+IBS9urdjQ4d1MJHUjxaOgDOAMZIGkEIUuPZce8OAEkDgNVmVgNcBdwRp/cFKs1sa1zmOODHZmaSpgFnAQ8C5wN/bKkDcs65tMyMFRu2Mj8GtuS9uQWrKtlWVZNzPQkG9+kWc3I7B7ohfbvRqcRbtDVGiwZAM6uSdDnwJFAC3GFmb0maBMw0synAOOB6SUYoAr0srr4/cIukGkLt1f81s0zlmW8BD0q6DngV+G2LHZRzzmVZW7ltpxxc5rFg5SY2bavOu97Anl0Y0b97bZHliPgY1r+Urp1KWvAI2ocWvwdoZlOBqVnTrk28fgR4JMd6zwMH5dnmPEINU+ecaxEbt1btFOAWZHJ0KzextnJ73vV6d+vEiMz9uKwg17NrpxY8AleIlWCcc64gbNlezaLVlTsqnqzYEeRWbMhdwxKgtHPJLhVPRgzszoj+3enb3fvpKBQeAJ1z7dr26hqWrNm8S8WT+Ss31VvDcli/0tocXDLgDezZhdBE2RUyD4DOuTavpsaoWL9lRw4u0ZRg8epKquqoYbl3v247FVWOiDm6vfp0o8RrWBY1D4DOuTbBzFixcWvo+SSrUfiCVZvYmqeGJcBevbsyYmAsqkwEuiF9S+nc0WtYtlUeAJ1zRWVd5XbmrdyYuCe3I9Bt3FqVd70BPbrEpgOJe3MDujOsX3e6dfYalu2RB0DnXMHZtLUq0dvJzvfm1tRRw7JX146MGNijtuLJ8AGljIyNwr2GpcvmAdA51yq2VlWzaFXljmYEiW6+3lufv4Zlt04lu1Q8ybzuW9rJK5+41DwAOueaTVV1DUvXbs7ZKHzZ2s3kqXtC55IODO1fukvFk5EDu7OH17B0TcQDoHNut9TUGMvXb9m5QXhsFL54dSXbq3NHuQ6CYf1Ld6l4MmKA17B0LcMDoHOuXmbGyo3bdm4MXtuH5Sa2bM9fw7Ksd9ecjcL39hqWrpV5AHTO1Vq3eXvOiicLVm5iQ501LDvX5uSSNSyH9/calq5weQB0rp2p3FYV28pt2mV8udWbtuVdr2fXjrv0X5kJeL28hqUrQh4AnWuDtlZVs3h1ZW1j8GSj8OXrt+Rdr1unkhjgdr031697Z6984toUD4DOFanqGmPpms2hUXhtxZMQ6JauyV/DslOJGNpv18FTRw7owZ69vIalaz88ADpXwGpqjPc25KhhuXITi+qpYTm0X2mi4kkpIwb2YET/7uzVpysdfQBV5zwAOtfazIzVm7btPHBqvDe3cFUlm7fnH0B1UK+uu1Q8GTGgO3v360aXjl75xLm6eAB0roWs37I95wjh81ZuYsOW/DUs+3fvnHM0guEDSint7F9h5xrLvz3ONaHN26pD27jsZgSrNrFyYx01LLt0zDkawfAB3endzWtYOtccUgdAhTvjnwROAPoDE81soaSxwBwzW9ZMaXSuoGyrqmHR6sqdejzJNAqvWJe/hmXXTh12aiuXDHT9vYalcy0uVQCU1BeYChwFbAB6AL8CFgIXAquBrzZTGp1rcdU1xrI8fVguWVOZt4Zlxw4KfVjmaBQ+qFdXOnj3Xs4VjLQ5wJ8AewPHATOAZFnO08D/pN2hpFOAXwAlwO1m9r9Z84cBdwADCYH1XDNbIukQ4GagF1AN/MDMHorrTAbGAuviZiaY2ay0aXJu87ZqHp65mH/MWcmCVZtYtKqSbdW5u/eSYEjfbrsUVY4c0J3Bfbp5DUvnikTaAHga8E0ze0FSdtWyRYTgWK+47o3AR4AlwAxJU8zs7cRiPwXuNrO7JJ0EXA+cB1QCnzezOZL2Al6R9KSZrY3r/Y+ZPZLyeJwDQtdf9764kN8+N3+XXlD27NUl52gEe/cr9RqWzrUBaQNgD2BpnnldgbTlOkcC5WY2D0DSg4TgmgyABwDfiK+nAY8BmNm/MwuY2TJJ7xNyiWtxroFWbNjKHf+cz70vLKzt4/LgIb35/DHD2a+sJ8P7d6d7F68j5lxblvYbPhv4KKG4M9tY4I2U2xkMLE68X0K4r5j0GnAmoZj0DKCnpP5mtiqzgKQjgc7A3MR6P5B0LfA34Eoz22VETUkXARcB+1ZUVKRMsmtLlqyp5LZn5/HgjMVsrQpFnMeO6s+l40Zz3Oj+XhHFuXYkbQC8Cfi1pHXA/XFaH0lfAC4nBJWm8s24rwnAs4ScZ21LYEllwD3A+WaWuUlzFbCcEBRvBb4FTMresJndCtwqaXpZWdnYJkyzK3Dl72/k5ulz+eOspVTFGiwn778nl544isOG9m3l1DnnWkOqAGhmt0oaCXyPHYHlr0AN8GMzuy/l/pay8/3CIWQVrcbmFGcCSOoBfDpzn09SL+AJ4GozezGxTiY7t1XSnYQg6hxvLFnHTdPL+ctbyzELXYSddshefHncKPYb1Ku1k+eca0Wpb3KY2ZWSbiZUYNkDWAX8NXM/L6UZwBhJIwiBbzzwueQCkgYAq2Pu7ipCjVAkdQYeJVSQeSRrnTIzq4htFU8H3mxAmlwbY2a8PH81N06fy7P/XgFA55IOnHXEEC4+YSTD+ndv5RQ65wpBvQEwBp4fAfeb2Qzg9sbuzMyqJF0OPEloBnGHmb0laRIw08ymAOOA6yUZoQj0srj6Z4mN8GPxKOxo7nCfpIGEyjizgEsam0ZXvMyM6bNXcOO0cmYuXANAaecSzjlqKF/60Ej27NW1lVPonCsk9QZAM9sm6WJC7mu3mdlUQqP65LRrE68fAXZpzmBm9wL35tnmSU2RNlecqmuMqW9UcNP0ubxTsR6A3t06MeHY4Uw4djh9u3du5RQ65wpR2iLQV4GDCDky5wrCtqoaHn11Cb95Zh7zV24CYI+eXbjwQyM5+6ih9PBmDM65OqT9hbgCeEDSQuAJM8vTEZRzza9yWxUPvryY2/4xr7bvzb37deOSsaP49GFD6NrJG6k75+qXNgD+DugN/BHYLmkFkAyCZmbDmjpxziWtq9zO3S8s4I5/zmdN5XYA9tmzB5eOG80nPljmXZA55xokbQD8GzsHPOdazIoNW/ntc/O598WFbMz02rJ3Hy4bN4qT99/TO5h2zjVK2naAE5o5Hc7tYsmaSm59dh4PJXptOW50fy4bN5pjRnmvLc653eO1BFzBKX9/AzdPn7dTry0fOWBPLh03ikO91xbnXBNpyIC4BwHfJfT92RdYQ+is+vtmlrYvUOfyen3JWm6aNpcn3w69tpR0EGccOphLxo5i30E9Wzt5zrk2Ju2AuP8BPANsBqYQ+t0cRBgh/uOSTjCzV5otla7NMjNemr+aG6eV8485K4HQa8tnjhjCxSeMYmj/0lZOoXOurUqbA7ye0L3Yh81sQ2aipJ6EESKuJ4wW4VwqZsbf332fG6eV869FYUSr7p1LOOfoYXzp+BHs4b22OOeaWdoAeDRwXjL4AZjZBkk/Au5q8pS5Nqm6xnjijQpumlbOu8vD5dSntBNfOHYE5x87jD6l3muLc65lpA2A9TWB8CYSrk5bq6p59F9L+c0zc1mwqhIIvbZcdMJIzj5yqA8+65xrcWl/dV4Cvi3p6awi0O6EsfdezLuma9cqt1XxwMuLue3ZeSxfH3ptGdqvNPTacvhgunT0Xlucc60jbQD8NjAdWCjpcaCCUAnmVKCUMIKDc7XWVW7nrhcWcGei15Z99+zJpSeO4uMHea8tzrnWl7Yh/MuSjgauBT4G9ANW480gXJb3N2wJvba8sJBN26oBOGTvPlx+4mhO2m8P77XFOVcwGjIg7uvAWc2YFlfEFq+u5JZn5/LwzCVsi722HD96AJeeOIpjRnqvLc65wpO2HeBAoK+Z/TvHvH0II7ivbOrEucI3570N3Dx9Ln98bRnVsdeWj31gTy4dN5qD9+7Tyqlzzrn80uYAbyIUeV6cY97Xgf6EEdtdO/Ha4rXcNL2cJ996Dwi9tpx56GAuGTeKffb0Xlucc4UvbQA8Hrgsz7yngF83TXJcITMzXpy3mpumJ3pt6diB/zpiby46YSR79/NeW5xzxSNtAOwLrMszbz0hB+jaKDPjb++8z43Ty3k10WvLuccM44LjR7BHT++1xTlXfNIGwCXAUYRxAbMdRWgW4dqYquoannijgpunz92p15YvHjeC848ZTu/STq2cQueca7y0AfAR4CpJr5nZE5mJkj4OXAncnHaHkk4BfgGUALeb2f9mzR8G3AEMJNx3PNfMlsR55wPXxEWvM7O74vTDgclAN2Aq8DUz895pGmlrVTV/iL22LIy9tuzZqwsXfsh7bXHOtR1pf8kmAScAUyQtB5YCgwmN4V8EvpdmI5JKgBuBjxBylTMkTTGztxOL/RS428zuknQSoaPt8yT1IwzHdASh67VX4rprCAH4QkKPNVOBU4A/pzw2F23aWsUDLy/itn/M4731WwEY1r+UL48dxRmHea8tzrm2JW1D+EpJY4HzCMGrP1BOqABzr5lVpdzfkUC5mc0DkPQgcBqQDIAHAN+Ir6cBj8XXHwP+amar47p/BU6RNB3oZWYvxul3A6fjATC1tZXbuOv5hdz5/HzWxl5b9hvUk0tPHM2pBw7yXlucc21SQxrCbycUTd6xG/sbDCxOvM/cW0x6DTiTUEx6BtBTUv886w6OjyU5pu9C0kXARcC+FRV+2/L99Vu4/bn53Pfijl5bDhvah8tiry3eeN0515albQg/ACg1s0WJaRcDBwJPmtnjTZimbwK/ljQBeJZQ3FrdFBs2s1uBWyVNLysrG9sU2yxGi1aFXlt+98qOXls+NGYAl504mqNG9PPA55xrF9LmAO8g5KwuBZD0HcJ9vzXApZI+Z2YPpdjOUmDvxPshcVotM1tGyAEiqQfwaTNbK2kpO3e6PYTQQffS+DrvNl3w79hry5REry2nfGAQl544ig8O8V5bnHPtS9oAeAQ7D3p7CfBDM7tG0i8J9+zSBMAZwBhJIwhBajzwueQCMbe52sxqgKvYUeT6JPBDSX3j+48CV5nZaknrY2fdLwGfB36V8rjahVmL13LTtHKeejvRa8thg/ny2FGM8V5bnIFCRHgAACAASURBVHPtVNoA2A94D0DSgYTan5mA+Bgh6NTLzKokXU4IZiXAHWb2lqRJwEwzm0LI5V0vyQhFoJfFdVdL+j4hiAJMylSIIeRMJxOaQfwZrwCDmfHC3FXcOL2cf5avArzXFuecS0obAFexo5jxJGCZmc2J7zsBqasJmtlUQlOF5LRrE68fIbQ7zLVuzko4ZjaTcD+y3aupMf727vvcOK2cWYtDry09unTk3KOH8cXjh3uvLc45F6UNgE8DE2Px5BXsaJoAsB+wsKkT5hqmqrqGx18PvbbMfi/02tI39tryee+1xTnndpE2AP4/4F5Co/QZ7Nzw/RzguSZOl0tpy/Zqfv+vJdzyzDwWrQ69tgzq1ZWLThjJ+CP3prSz99rinHO5pG0I/x6hAXwuJwNbmixFLpVNW6u4/6XQa8v7G0KvLcP7l/LlcaM4/VDvtcU55+qz29kDM1vfFAlx6ayt3Mbk5xdw5z8XsG5z6LVl/7JeXDpuFKceVEZJB2/D55xzaXj5WJF4b/0Wbv/HPO57aRGVsdeWw4f15fITRzNu34HeeN055xrIA2CBW7Sqkt88O5dHZi5hW3XoteWEfQZy2bhRHOm9tjjnXKN5ACxQs5dv4Obp5Ux5bRk1BhL854GDuHTcaA4a0ru1k+ecc0XPA2CBeXXRGm6cNpen3wm9tnSMvbZcMnYUo/fo0cqpc865tsMDYAEwM56fu4obp5Xz/NzQa0uXjh0Y/x97c+EJIxnS13ttcc65ppY6AEoaCXwWGApkdydiZnZBUyasPaipMf76znvcNH0uryV6bTnvmGF88bgRDOzZpZVT6JxzbVfa4ZBOBx4mdHn2PrA1axFr4nS1aVXVNfzp9WXcNG0uc97fCEC/7p254PgRnHv0MHp3815bnHOuuaXNAX6fMPTQOWa2ovmS07Zt2V7NI68s4ZZn57J49WYAynrHXlv+YyjdOnvjdeecaylpA+BI4AoPfo2zeVs197y4gNv+MZ8VsdeWkQO6c8nY0GtL546p+xJ3zjnXRNIGwHeB/s2ZkLbs6kff4A+vhjF6DyjrxWUnjuaUAwd5ry3OOdeKGtIZ9s8lvWRm85ozQW3RO8vD6Aw/PuuDfObwId543TnnCkDaADiRkAN8R9IcYHXWfDOzsU2ZsLakYl2433fivnt48HPOuQKRNgBWA7ObMyFt1eZt1ayt3E6nEtG/e+fWTo5zzrko7XBI45o5HW3W8vVhpKg9e3Wlg9/zc865guHVD5tZpvizrHd23wHOOedaU+oAKKlM0k8lzZA0Nz7/WNKg5kxgsVu+LuQAy3p3a+WUOOecS0oVACXtA8wCvgpsBF6Oz18DZkkak3aHkk6RNFtSuaQrc8wfKmmapFclvS7p1Dj9HEmzEo8aSYfEedPjNjPz9kibnuZWURsAPQfonHOFJG0lmB8B64GjzGxBZqKkYcBTcf6Z9W1EUglwI/ARYAkwQ9IUM3s7sdg1wMNmdrOkA4CpwHAzuw+4L27nIOAxM5uVWO8cM5uZ8nhaTKYIdJAHQOecKyhpi0BPBL6TDH4AZraQ0ETixJTbORIoN7N5ZrYNeBA4LWsZA3rF172BZTm2c3Zct+At9xygc84VpLQBsDOwIc+8DXF+GoOBxYn3S+K0pInAuZKWEHJ/X8mxnf8CHsiadmcs/vyO8jS2k3SRpJnA4RUVFSmTvHsyRaCD/B6gc84VlLQBcBbwFUk7LR8DzaVxflM5G5hsZkOAU4F7kvuVdBRQaWZvJtY5x8wOAj4UH+fl2rCZ3WpmRwCvlJWVNWGS88vkAPfyHKBzzhWUtPcAJwGPE3qCeQioAAYBnwHGAB9PuZ2lwN6J90PitKQLgFMAzOwFSV2BAYRhmADGk5X7M7Ol8XmDpPsJRa13p0xTs9myvZpVm7bRsYPo38PH9nPOuUKSKgdoZn8BPkEo7ryaUJHlGkJN0E+Y2VMp9zcDGCNphKTOhGA2JWuZRcCHASTtTxh8d0V834EwKG/t/T9JHSUNiK87xXS+SQF4L9EI3ju+ds65wpJ6RPgYBP8iqRToC6wxs8qG7MzMqiRdDjwJlAB3mNlbkiYBM81sCnAFcJukrxMqxEwws8yAuycAi7M65O4CPBmDXwnwNHBbQ9LVXHbc//PiT+ecKzSpA2BGDHoNCnxZ608lVG5JTrs28fpt4Lg8604Hjs6atgk4vLHpaU7LPQA651zByhsAJV0L3G5my+LrupiZfb9pk1b8KrwCjHPOFay6coATgb8Q2uFNrGc7BngAzLKjEbw3gXDOuUKTNwCaWYdcr1163g2ac84VrrR9gQ6NlUxyzesoaWjTJqtt8HuAzjlXuNLm7OYDh+aZd3Cc77J4DtA55wpX2gBYVyO2TkBNE6SlTdlWVcPKjVsp6SD26OkB0DnnCk1dtUD7AP0SkwZLGpm1WDfgfGB5M6StqGUawe/Rs4s3gnfOuQJUVy3QrwHfJdTwNOCRPMspLucSvBG8c84VtroC4GPAAkKAuwO4DpibtcxW4G0ze71ZUlfEMk0g/P6fc84VprqaQbwGvAYgyYDHzWxVSyWs2NXWAO3lbQCdc64QpeoKzczuau6EtDW1vcD08Rygc84VotR9gUr6APAlYF/CCA1JZmYfbsqEFTtvA+icc4UtVQCMg9A+Q7gnOAZ4nTAixFDCqO7lzZS+ouX3AJ1zrrClbQf4Q+APwAcIlWIuMLPhwMmEIYiua5bUFbEdtUD9HqBzzhWitAHwg8C9hOYQEIIeZvZ3QvC7vumTVry2V9ewYuNWpNAO0DnnXOFJGwA7A5vMrAZYDZQl5s0GDmzqhBWz9zdsxSwEv04l3o+4c84VorS/zuXA4Pj6deCLkjpI6gB8Ae8JZifLfRgk55wreGlrgf4JGAfcT7gf+ASwHqgGegBfbY7EFatla2Mn2L28AoxzzhWqtO0AJyZePy3paODTQCnwFzN7qnmSV5y8CYRzzhW+1O0Ak8zsVeDVJk5Lm+HDIDnnXOFLOyDu0ZI+m2feZ2I7wVQknSJptqRySVfmmD9U0jRJr0p6XdKpcfpwSZslzYqP3yTWOVzSG3Gbv5TUqsMvLF8f2wD28XuAzjlXqNJWgrme0AYwl/1J2QxCUglwI/CfwAHA2ZIOyFrsGuBhMzsUGA/clJg318wOiY9LEtNvBi4kNNIfA5ySJj3NxXOAzjlX+NIGwIOBF/PMe5nQTjCNI4FyM5tnZtuAB4HTspYxoFd83RtYVtcGJZUBvczsRTMz4G7g9JTpaRY7OsL2AOicc4UqbQDsWseyJUD3lNsZDCxOvF/CjuYVGROBcyUtAaYCX0nMGxGLRp+R9KHENpfUs80WU1VdUzsY7p4eAJ1zrmClDYDvAJ/KM+9ThMbwTeVsYLKZDQFOBe6J7Q0rgKGxaPQbwP2SetWxnV1IukjSTODwioqKJkzyDis2bqXGYECPLnTu6I3gnXOuUKX9hf4NcKGkn0jaR1KppDGSfgJcwM736eqyFNg78X5InJZ0AfAwgJm9QMh9DjCzrZnxCM3sFcLgvPvE9YfUs03ierea2RHAK2VlZbkW2W0+DJJzzhWHVAHQzG4DfgZ8nZAb3AC8G9/fYGa3ptzfDGCMpBGSOhMquUzJWmYR8GEASfsTAuAKSQNjJRokjSRUdplnZhXA+lhTVcDngT+mTE+T8/t/zjlXHFK3AzSzb0q6mTACRH9gJfC0mc1rwDaqJF0OPEm4d3iHmb0laRIw08ymAFcAt0n6OqFCzAQzM0knAJMkbQdqgEvMbHXc9KXAZKAb8Of4aBVeA9Q554pDgxrCm9lcQtFjo5nZVELlluS0axOv3waOy7He74Hf59nmTAqkQ+6Ktd4PqHPOFYO8AVDSUKDCzLbH13Uys0VNmrIiVbHec4DOOVcM6soBLgCOJrTzW8COsQDzKWmaJBW35V4E6pxzRaGuAPgFdhR3fpH6A6AjGQC9CNQ55wpZXQGwNztydX8nFoc2f5KKV3WN1TaC36OXjwTvnHOFrK5mEDcAw+Pr+cChzZ6aIrdy41aqaoz+3TvTtZOXCDvnXCGrKwCuBQbF18KLQOtV4eMAOudc0airCPSfwF2SXovvb5a0Ps+yZmYfbtqkFZ/l6+IwSH7/zznnCl5dOcALgQcIjc6NECw75Xl0bt5kFgdvBO+cc8Ujbw7QzN4j9LCCpBrgIjN7uaUSVoyWexGoc84VjbQ9wYwgjMbg6uA5QOecKx6pAqCZLWzuhLQFFesy3aB5AHTOuUKX9x6gpGpJR8bXNfF9vkdVyyW5cNUOheSVYJxzruDVlQOcxI6R1ifhzSDqVJNoBO85QOecK3x1VYL5XuL1xBZJTRFbtWkb26uNvqWdvBG8c84VgbQjwu9CUj9Jh0vyPr9I1gD14k/nnCsGqQKgpGskXZ94fwJhhIiXgTmSxjRP8orHstpG8F786ZxzxSBtDvBcIDny+4+A14DTgfeA7zdxuoqOD4PknHPFJW07wMHAHABJA4EjgQ+b2XRJnYFfNlP6ioa3AXTOueKSNgdYzY7uzk4AthD6CgVYAfRr4nQVneW1bQD9HqBzzhWDtAHwLeBcST0Ig+M+kxgbcG/g/eZIXDHxHKBzzhWXtAFwEvBZYB3wYcI9wIxTgX+l3aGkUyTNllQu6coc84dKmibpVUmvSzo1Tv+IpFckvRGfT0qsMz1uc1Z87JE2PU3Fh0JyzrnikrYrtCcl7Q8cBswys7mJ2c8SKsTUS1IJcCPwEUIj+xmSppjZ24nFrgEeNrObJR0ATCUMzLsS+KSZLZN0IPAk4d5kxjlmNjNNOpqamXklGOecKzJpK8FgZvMJI8NnT7+lAfs7Eig3s3kAkh4ETgOSAdCAXvF1b2BZ3M+riWXeArpJ6mJmWxuw/2axetM2tlXX0LtbJ0o7p/5InXPOtaK07QBPk/SFxPthkl6QtEHSI/HeYBqDgcWJ90vYORcHMJFwv3EJIff3lRzb+TTwr6zgd2cs/vyOJKVMT5Pw+3/OOVd80t4DvAYYmHj/M2AIcCuhVujEJkzT2cBkMxtCuL94j6TadEr6AOEe5MWJdc4xs4OAD8XHebk2LOkiSTOBwysqmm50Jx8H0Dnnik/aADgKeB1AUjdCYPqGmV0BfBs4I+V2lhJqjWYMidOSLgAeBjCzF4CuwIC47yHAo8Dnk/chzWxpfN4A3E8oat2Fmd1qZkcAr5SVlaVMcv0q1nsO0Dnnik3aANgV2BxfH0u4d/hUfD8b2CvldmYAYySNiA3oxwNTspZZRKhpSqx40xVYIakP8ARwpZll2iAiqaOkTIDsBHwCeDNleppExdpMN2jeBtA554pF2gC4ADg+vj4NeMXM1sX3exCaR9TLzKqAywk1ON8h1PZ8S9IkSZ+Ki10BXCjpNeABYIKZWVxvNHBtVnOHLsCTkl4HZhFylLelPK4m4UWgzjlXfNJWWbwF+KmkM4BDgC8n5h3DzrU462RmUwmVW5LTrk28fhs4Lsd61wHX5dns4Wn33xy8EoxzzhWftO0AfyFpJXA08EszuzsxuydwZ3Mkrlgs93uAzjlXdBrSDvA+4L4c0y/OsXi7YWZUeD+gzjlXdBo9IK4L1lZuZ8v2Gnp27UiPLt4I3jnnikXqABjb0L0qqVJSdfajORNZyPz+n3POFae0PcF8HvgVoRlDV8I9v3uB9cBcQmfZ7dLy9V786ZxzxShtDvC/gevZUfvzJjM7HxhJaB+4qhnSVhRqc4C9PAfonHPFJG0AHEMY9aEmPjoDmNka4AfA15oldUXA2wA651xxShsANwMdYoP05YScX8ZG0vcE0+YsWxsC4F59PAA651wxSVtt8Q1CLyxPA/8Avi1pPlBF6Aj73WZJXRHwe4DOOVec0gbAW9mR6/sOIRA+F99vAE5v4nQVDa8F6pxzxSltTzAPJV6XxyGJjgFKgefNbGUzpa+gJUeC93uAzjlXXBrVctvMNhFyge3a+i1VVG6rpnvnEnp6I3jnnCsqeX+1JQ1tyIbMbNHuJ6e4ZLpAK+vTjRYehN4559xuqivbsgCwBmyrZPeSUnz8/p9zzhWvugLgF2lYAGx3au//eSN455wrOnkDoJlNbsF0FCXPATrnXPHK2xBewSclHVjHMgdJ+mTzJK3wLfdhkJxzrmjV1RPMecADwKY6ltkAPCDp7CZNVZGozQF6LzDOOVd06gqA5wJ3mtn8fAuY2QLgt8D5TZyuouBFoM45V7zqCoCHAU+l2MbTwBFNk5zisrx2JAgvAnXOuWJTVwDsCaxJsY01cdl2ZcOW7WzcWkW3TiX06uaN4J1zrtjUFQBXAsNSbGNoXDYVSadImi2pXNKVOeYPlTQtjj7/uqRTE/OuiuvNlvSxtNtsDssTxZ/eCN4554pPXQHwOdLd25vAjo6x6ySpBLgR+E/gAOBsSQdkLXYN8LCZHQqMB26K6x4Q338AOAW4SVJJym02Oa8A45xzxa2uAPhz4MOSbpDUOXumpE6Sfg6cBNyQcn9HAuVmNs/MtgEPAqdlLWNAr/i6N7Asvj4NeNDMtsaKOeVxe2m22eQy3aAN8vt/zjlXlOpqCP+CpCuA/wPOkfQUsDDOHgZ8BOgPXGFmL6bc32BgceL9EuCorGUmAk9J+grQHTg5sW5yP0viNFJsEwBJFwEXAftWVFSkTHJuXgPUOeeKW50jwpvZz4ETgRnAGcBV8XEGMBM40cx+0cRpOhuYbGZDgFOBeySlHbm+TmZ2q5kdAbxSVla2W9vyYZCcc6641Vt90cyeBZ6NQWhAnLzKzKobsb+lwN6J90PitKQLCPf4MrnQrnG/da1b3zabXCYHuJffA3TOuaKUOmdlZjVm9n58NCb4QchJjpE0It5XHA9MyVpmEfBhAEn7A12BFXG58ZK6SBoBjAFeTrnNJrejI2y/B+icc8WoRRuwmVmVpMuBJwnDJ91hZm9JmgTMNLMpwBXAbZK+TqgQM8HMDHhL0sPA20AVcFkmEOfaZnMfy7LMWIBeBOqcc0WpxVtwm9lUYGrWtGsTr98Gjsuz7g+AH6TZZnPauLWKDVuq6NKxA31KO7XUbp1zzjWhJqlc0t54I3jnnCt+HgAbYUcA9Pt/zjlXrDwANkKF3/9zzrmi5wGwEbwNoHPOFT8PgI2wzHuBcc65oucBsBGWZ/oB9XuAzjlXtDwANoL3A+qcc8XPA2AjLF/vAdA554qdB8AG2rytmrWV2+lc0oF+3XcZJco551yR8ADYQLXjAHojeOecK2oeABvIm0A451zb4AGwgWqHQfIA6JxzRc0DYANlKsB4EwjnnCtuHgAbyLtBc865tsEDYANVrPV7gM451xa0+HiAxc4bwTvn2pLKykqWLl1KVVXVTtM7duzI4MGDKS0tbaWUNT/PATbQjkbwfg/QOVf8VqxYsUvwA6iqqmLFihWtkKKW4znABtheXUPPrh3ZVlVDf28E75xrAzZvDvUaxowZQ8eOISRUVVUxZ86c2nltlQfABuhU0oFn/udEzMwbwTvn2gQzA6gNfgAlJSU7zWurvAi0ETz4OefasvbyG9fiAVDSKZJmSyqXdGWO+TdImhUf/5a0Nk4/MTF9lqQtkk6P8yZLmp+Yd0hLH5dzzrUVbT3nl9GiRaCSSoAbgY8AS4AZkqaY2duZZczs64nlvwIcGqdPAw6J0/sB5cBTic3/j5k90uwH4ZxzbYgkzIyqqqraYtDq6uraeW1ZS98DPBIoN7N5AJIeBE4D3s6z/NnAd3NMPwv4s5lVNksqnXOunejWrRuVlZXMmTMn57y2rKWLQAcDixPvl8Rpu5A0DBgB/D3H7PHAA1nTfiDp9ViE2iXPNi+SNBM4vKKiosGJd865tmbgwIE7VYDJ6NixIwMHDmyFFLWcQq4FOh54xMyqkxMllQEHAU8mJl8FLAc6A7cC3wImZW/QzG4FbpU0vaysbGxzJdw554pFaWkpY8aMae1ktIqWzgEuBfZOvB8Sp+WSK5cH8FngUTPbnplgZhUWbAXuJBS1Ouecc3m1dA5wBjBG0ghC4BsPfC57IUn7AX2BF3Js42xCji+5fJmZVSjcsT0deLOedIyeNWsW48aNa/gROOdcOzVr1iyA0a2djqbSogHQzKokXU4oviwB7jCztyRNAmaa2ZS46HjgQcuqiytpOCEH+UzWpu+TNBAQMAu4pJ6krFi3bh3PPPNMj/h+VmOPqQBlmoD4MRU2P6bC19aOB3b/mEYDbaZ/NLWX9h65SJoOYGbjWjclTcePqTj4MRW+tnY80DaPaXd4TzDOOefaJQ+Azjnn2qV2XQTqnHOu/fIcoHPOuXbJA6Bzzrl2yQOgc865dqndBsD6hmUqVJL2ljRN0tuS3pL0tTi9n6S/SpoTn/vG6ZL0y3icr0s6rHWPID9JJZJelfR4fD9C0ksx7Q9J6hynd4nvy+P84a2Z7nwk9ZH0iKR3Jb0j6ZhiP0+Svh6vuzclPSCpa7GdJ0l3SHpf0puJaQ0+L5LOj8vPkXR+axxLIi25jukn8dp7XdKjkvok5l0Vj2m2pI8lphfl72KjmVm7exAa4c8FRhL6D30NOKC105Uy7WXAYfF1T+DfwAHAj4Er4/QrgR/F16cCfyZ0EnA08FJrH0Mdx/YN4H7g8fj+YWB8fP0b4Mvx9aXAb+Lr8cBDrZ32PMdzF/Cl+Loz0KeYzxOh4/r5QLfE+ZlQbOcJOAE4DHgzMa1B5wXoB8yLz33j674FdkwfBTrG1z9KHNMB8TevC2HAgbnxN7Fofxcb/bm1dgJa6WI5Bngy8f4q4KrWTlcjj+WPhPEVZwNlcVoZMDu+vgU4O7F87XKF9CD0C/s34CTg8fiDszLxBa49Z4SehI6JrzvG5dTax5B1PL1jsFDW9KI9T+wYzaVf/NwfBz5WjOcJGJ4VLBp0XghdMt6SmL7TcoVwTFnzzgDui693+r3LnKe29LuY9tFei0BTD8tUyGKR0qHAS8CeZpYZ42k5sGd8XSzH+nPg/wE18X1/YK2ZVcX3yXTXHlOcvy4uX0hGELqMujMW694uqTtFfJ7MbCnwU2ARUEH43F+huM9TRkPPS8GfryxfJORkoe0c025rrwGw6EnqAfwe+G8zW5+cZ+HvW9E08JT0CeB9M3ultdPShDoSiqRuNrNDgU2EorVaRXie+hIGsB4B7AV0B05p1UQ1g2I7L/WRdDVQBdzX2mkpNO01ADZkWKaCI6kTIfjdZ2Z/iJPfUxgrMTNm4vtxejEc63HApyQtAB4kFIP+AugjKdNhezLdtccU5/cGVrVkglNYAiwxs5fi+0cIAbGYz9PJwHwzW2FhOLI/EM5dMZ+njIael2I4X0iaAHwCOCcGdijyY2pK7TUA1g7LFGusjQem1LNOQZAk4LfAO2b2s8SsKUCmJtr5hHuDmemfj7XZjgbWJYp6CoKZXWVmQ8xsOOFc/N3MzgGmAWfFxbKPKXOsZ8XlC+ofu5ktBxZL2jdO+jDwNkV8nghFn0dLKo3XYeaYivY8JTT0vDwJfFRS35gz/ig7D9Ld6iSdQrit8Ckzq0zMmgKMj7V0RwBjgJcp4t/FRmvtm5Ct9SDU7vo3odbT1a2dngak+3hC8czrhCFNZsVj6U+oRDIHeBroF5cXcGM8zjeAI1r7GOo5vnHsqAU6kvDFLAd+B3SJ07vG9+Vx/sjWTneeYzkEmBnP1WOE2oJFfZ6A7wHvEsbcvIdQk7CozhNhoO0KYDshp35BY84L4b5aeXx8oQCPqZxwTy/zO/GbxPJXx2OaDfxnYnpR/i429uF9gTrnnGuX2msRqHPOuXbOA6Bzzrl2yQOgc865dskDoHPOuXbJA6Bzzrl2yQOgaxBJp0t6NvY8v1nSQkmPxTZHmWUmSLJC6f0/Q1IHST+XVCGpRtJjeZbrI2lirhEZJE2X9Fwzp/N0Sd9opm1Pjh0ONGbdBZImN22KWk68Lr/Y2ulwhcObQbjUJH2V0EPLHYR2bZuAUcDHgX+b2f+Lyw2M0181s62tlNxdSPos8BBwBfACsMrM/p1jueGEjqwvNLPbs+ZNJ3T8fHwzpnMycLKZDWmGbY8CepnZq41Y91BgvZnNbep0tYSWOHeuuHSsfxHnan0TeMzMLkhM+ztwm6Ta0gQzW0HoCLrQ7B+ff25mNXUuWSQkdWnIn4zdCV6NCZrOFTIvAnUN0Y/QU/4ukgEluwg0FrtZnse4xHoHS5oiaU0sXv2npA+lSVgcyPOFuN66WCy7b2L+AmBifFsd9z0hx3aGE3J/EAK75VpW0smS/iWpUmFw2DNybKvBxxNzf+cDgxP7XhDnjYvvz5R0m6QVwHtx3mhJ90iaH/c1T9LNsZuunbafLAKVNDxu82JJk2Lx8FpJf5I0JGvdnYpAE+f5aEn3SVovaZnCALJds9YdKWlq/Lzel/R/ki5KU1Qu6WOSno/ndaPCgK3XNuSzjrm/scBxic91el37dW2fB0DXEC8D50v6H0n7NGC977NjvLHM459AJaF/SeL9tucJQfZC4NOEjpOflnR4XRtXuP/4BLAR+C/gy8CBwHOSMsO5nAFMjq8zaXgix+YqgDPj6+vzLDuKUBT8s7hsBfA7SaMTaWrs8XwfmErIQWf2nR1cf0Xoous8woC0EEZnWAz8N2GMvkmEvjqn1rGvpKuA0YTuvb4W93tvynXvIXSddSZwM3BZ3B4ACv1K/hX4IOHcTCCMKHF1fRuWNJLQH+V8wrn9FOFz755YJs1nfSnwKqFbusznemnK43NtVWv3xeaP4nkA+xB+QDLDxawk9EH40azlJsT5w/Ns55tANXB6YtrfgHeAzolpJXHaY/WkayahD8eOiWkjCP0i/iwx7TriaDf1bG94TP+XcsybHrc7JjFtj3g8326i45lMGEkie/q4mK5HUxxDR3b0G3to1rYX5DjW6TnOkQF7JaYtACbnOM/fy1r3ccI94cz7i+JyRyamiTDieN7rJC53VlymVx3LgBjkYAAABCtJREFUpPqs47l7rrW/R/4onIfnAF1qFiqMHEooSvoBoYPdM4AnJV2TZhuSPgn8CPiWmT0Wp3WL2/wdUCOpo8LwOSJ0THxCHdvrThhm6CHbMSgrZjafkMsc29DjTGGOmc1J7Ot9wvA5Q2OaGn08KT2aPUFSZ0nflvSupM2EIP2POHvf7OVzyM4pvhGfh6ZYNzsn/UbWekcDi8zs5cwEMzPCkF71mUU4lgclnSVpj+TMFvisXRvmAdA1iJlVm9mzZnaNmZ1MGAngDeC72febskk6GLgf+K2Z/TQxqx/hH/t3CD92ycflQF8lKtlk6Uv4scs1dNDyuO2mtjrHtK2E0Q9g944njVzHej3hHue9hFq5R7KjKLdrjuWzZR9TpmJNY9ftknifHF8v6b36Nmxm5YQi3Q6Eotblkl6UlPlj09yftWvDvBao2y1mtkzS7YR7YplxxXYhaRDwJ+BFdr33shaoIQw7c3ee/eSrtbmGUEQ2KMe8QeQOVs1td44njVxtl8YDd5vZdZkJknrsxj6aUgVwQI7pe6ZZ2cymAdMkdSEMwDsJeCJWnmnuz9q1YR4AXWqSyiz3IK37xeecNURjjcDHCJVUzkoWVQKY2SZJ/wAOBv7VkB+suO4rwGckTTSz6rjPYcCxhAojDZXJ/XRrxLq7dTyJ/Td036WEXE/SFxq4jebyIvAFSUdmikEliVBZJTULzT3+HgP7H4ERZjajAZ/1VqBno47AtUkeAF1DvCnpacL9ovlAL8IAmpcAD5vZojzr/Zxwn24CsH/47av1tpmtB74BPEu4n/hbQq5hQFyvxMyurCNd3yHch3pc0k1AD8LAreuA/2vEcb5HqEU4XtLrhAb/881sVQO2sTvH8zbQT9KXCRV8tpjZG3UsD/AXQg3dNwgDoZ5J+ANQCCYD3wL+IOlqQg3XLxGKryHk4HKSdAnhPt5UQi3XAYQapssIg/JC+s/6beBSSf9FqLW6wcxmN9ExuiLkZeOuIa4m5EwmAU8RelU5BriSUCU/n/34/+3cr0oFQRSA8W8R4wXfwWSwCDaTwWY0+AgGgw+hYlQUrnANBkGsGgy+gSAYfQG7IIiWYzgrWJS7ynWR8/1g2+7s/IE5O7MzA9PAGXkCy+drASAi7oBFMvActOnvA/Nk5/aliLgm/3vNABfAkFwBuBQRj10L2Y4iPjroG+AWWO2Yxo/LA4yAc2CHnFK+HOOVm+R2gW2yXQbAepc8T0pEvAEr5AriIXBKBrOj9panbx6/J7c87JJ1eEh+fC1HxEub/rh1vUeuGB2RbXr8+9LpP/MoNEm9aJrmCpiLiNm+86KanAKVNHFNHu79TO7XHABr5Kh9o898qTYDoKS/8ApskfsDp4AH8qCBk15zpdKcApUkleQiGElSSQZASVJJBkBJUkkGQElSSQZASVJJ7w+FQ2fG1CEsAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(nbSamplesConfig,scores)\n",
    "plt.xlabel('Size of the training set')\n",
    "plt.ylabel('Classification score')\n",
    "plt.legend(loc=\"lower right\", fontsize='small')\n",
    "plt.title('Model accuracy, against number of training samples')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Differenciating the multi-class logisitc loss\n",
    "We recall that logisitc loss L, reads:\n",
    "\\begin{align*}\n",
    "  L(W,b) = c(s(Wx+b))\n",
    "\\end{align*}\n",
    "\n",
    "Where:\n",
    "  * $W$ is the weight, a $m\\times n$ matrix, where each element is indexed: $W_{i,j}$\n",
    "  * $b$ is the bias vector of size n, equivalent to the number of class\n",
    "  * $c$ is the cross entropy function\n",
    "  * $s$ is the softmax function"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Jacobian\n",
    "\n",
    "We first recall that, for a function $f:\\underset{\\mathbb{R}^n}{x} \\mapsto \\underset{\\mathbb{R}^n}{f(x)}$, we can derive the matrix of partial derivative, called the jacobian matrix:\n",
    "\n",
    "\\begin{align*}\n",
    "  J=\\begin{pmatrix}\n",
    "    \\frac{\\partial f_0}{\\partial x_0} & \\frac{\\partial f_0}{\\partial x_1} & \\dots & \\frac{\\partial f_0}{\\partial x_{n-1}}\\\\\n",
    "    \\frac{\\partial f_1}{\\partial x_0} & \\frac{\\partial f_1}{\\partial x_1} & \\dots & \\frac{\\partial f_1}{\\partial x_{n-1}}\\\\\n",
    "  \\vdots & \\vdots & \\vdots & \\vdots \\\\\n",
    "  \\frac{\\partial f_{n-1}}{\\partial x_0} & \\frac{\\partial f_{n-1}}{\\partial x_1} & \\dots & \\frac{\\partial f_{n-1}}{\\partial x_{n-1}}\n",
    "\\end{pmatrix}\n",
    "\\end{align*}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Jacobian of cross entropy $c$ for a single sample\n",
    "\n",
    "We recall that, given a vector $v$ of size $n$, and a set of probability (for each label) p of size $n$ (the number of classes), we have:\n",
    "\\begin{align*}\n",
    "  c(v) = -\\sum_{k=0}^{n-1} p_k log(v_k)\n",
    "\\end{align*}\n",
    "\n",
    "Then the jacobian becomes a simple row vector (1 output, multiple input variables) which is very simple to compute:\n",
    "\\begin{align*}\n",
    "  J_{cross} & \\begin{pmatrix}\n",
    "    \\frac{\\partial c}{\\partial v_0} , \\dots, \\frac{\\partial c}{\\partial v_{n-1}} \\end{pmatrix} \\\\\n",
    "  = &\\begin{pmatrix} -\\frac{p_0}{v_0}, \\dots, -\\frac{p_{n-1}}{v_{n-1}} \\end{pmatrix}\n",
    "\\end{align*}\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Jacobian for the softmax\n",
    "For the case of softmax, given a vector $v$ of size $n$ (the number of classes), we have also an output of size n:\n",
    "\\begin{align*}\n",
    "  s(v)_i = s_i = \\frac{e^{v_i}}{\\sum_{k=0}^{n-1}e^{v_k}}\n",
    "\\end{align*}\n",
    "\n",
    "Let's now try to find $\\frac{\\partial s_i}{\\partial v_j}$:\n",
    "\n",
    "In the case $i\\neq j$, we have to differentiate a function $f$ of the form $f(x)=\\frac{a}{b+e^x}$, where $a=e^{v_i}$ and $b=\\sum_{k\\neq j}^{n-1}e^{v_k}$ which can be computed thanks to the composition rule ($(g\\circ f)' = g'\\circ f \\times f'$) between $\\frac{a}{x}$ and $b+e^x$: \n",
    "\\begin{align*}\n",
    "  f'(x) &= -\\frac{a}{\\left(b+e^x\\right)^2} \\times e^x\\\\\n",
    "  &= -\\frac{ae^x}{\\left(b+e^x\\right)^2}\n",
    "\\end{align*}\n",
    "From there we can express:\n",
    "\\begin{align*}\n",
    "  \\frac{\\partial s_i}{\\partial v_j} &= -\\frac{e^{v_i+v_j}}{\\left(\\sum_{k=0}^{n-1}e^{v_k}\\right)^2} \\\\\n",
    "  &= -s_i s_j\n",
    "\\end{align*}\n",
    "\n",
    "Otherwise we have to differentiate a function $f$ of the form $f(x)=\\frac{a e^x}{b+e^x}$, where $a=1$ and $b=\\sum_{k\\neq i}^{n-1}e^{v_k}$ which can be computed thanks to the multiplication rule ((uv)'=u'v+uv') between $\\frac{a}{b+e^x}$ and $e^x$:\n",
    "\\begin{align*}\n",
    "  f'(x) &= -\\frac{ae^x}{\\left(b+e^x\\right)^2} \\times e^x + \\frac{a}{b+e^x} \\times e^x\\\\\n",
    "  &= -\\frac{ae^{2x}}{\\left(b+e^x\\right)^2} + \\frac{ae^x}{b+e^x}\n",
    "\\end{align*}\n",
    "From there we can express:\n",
    "\\begin{align*}\n",
    "  \\frac{\\partial s_i}{v_i} &= -\\frac{e^{2v_i}}{\\left(\\sum_{k=0}^{n-1}e^{v_k}\\right)^2} + \\frac{e^{v_i}}{\\sum_{k=0}^{n-1}e^{v_k}} \\\\\n",
    "  &= \\frac{e^{v_i}}{\\sum_{k=0}^{n-1}e^{v_k}} \\left( 1 - \\frac{e^{v_i}}{\\sum_{k=0}^{n-1}e^{v_k}}\\right)\\\\\n",
    "  &= s_i (1-s_i)\n",
    "\\end{align*}\n",
    "\n",
    "Jacobian for the softmax is then a $n \\times n$ matrix of the form:\n",
    "\\begin{align*}\n",
    "  J_{soft}=\\begin{pmatrix}\n",
    "    s_0 (1-s_0) & -s_0 s_1 & \\dots & -s_0 s_{n-1}\\\\\n",
    "    -s_1 s_0 & s_1 (1-s_1) & -s_1 s_{n-1}\\\\\n",
    "  \\vdots & \\vdots & \\vdots & \\vdots \\\\\n",
    "  -s_{n-1} s_0 & -s_{n-1} s_1 & \\dots & s_{n-1} (1-s_{n-1})\n",
    "\\end{pmatrix}\n",
    "\\end{align*}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Jacobian of the linear part\n",
    "We recall, that our linear model, for a single sample (feature vector) x, reads $f(x) = Wx+b$, and we are interested in differentiating this expression with respect to $W_{ij}$ and $b_i$.\n",
    "\n",
    "\\begin{align*}\n",
    "  \\frac{\\partial f_i}{\\partial W_{kj}} &= \\begin{cases} x_k \\quad &\\text{ if } k=i \\\\ 0 & \\text{otherwise} \\end{cases}\n",
    "\\end{align*}\n",
    "\n",
    "\\begin{align*}\n",
    "  \\frac{\\partial f_i}{\\partial b_{j}} &= \\begin{cases} 1 \\quad &\\text{ if } i=j \\\\ 0 & \\text{otherwise} \\end{cases}\n",
    "\\end{align*}\n",
    "\n",
    "Jacobian for the linear part is a pretty big matrix, that has n rows and M columns where M is the total amount of parameters in both $W$ and $b$, and n is the number of classes we account for in the model:\n",
    "\n",
    "\\begin{align*}\n",
    "  J_{lin}=\\begin{pmatrix}\n",
    "    \\frac{\\partial f_0}{\\partial W_{00}} & \\frac{\\partial f_0}{\\partial W_{01}} & \\dots & \\frac{\\partial f_0}{\\partial b_{0}} & \\dots \\\\\n",
    "    \\frac{\\partial f_1}{\\partial W_{00}} & \\frac{\\partial f_1}{\\partial W_{01}} & \\dots & \\frac{\\partial f_1}{\\partial b_{0}} & \\dots \\\\\n",
    "  \\vdots & \\vdots & \\vdots & \\vdots \\\\\n",
    "  \\frac{\\partial f_{n-1}}{\\partial W_{00}} & \\frac{\\partial f_{n-1}}{\\partial W_{01}} & \\dots & \\frac{\\partial f_{n-1}}{\\partial b_{0}} & \\dots \\\\\n",
    "\\end{pmatrix}\n",
    "\\end{align*}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Composition Rule\n",
    "Let's recall the composition rule with 3 fnctions: $c,s$ and $l$ in $x$:\n",
    "\\begin{align*}\n",
    "  (c \\circ s \\circ l)'(x) = (c' \\circ s \\circ l)(x) \\times (s' \\circ l)(x) \\times l'(x)\n",
    "\\end{align*}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Composition rule with partial derivatives \n",
    "\n",
    "We recall the chain rule in higher dimensions, using the jacobian matrix (matrix of partial derivative):\n",
    "\n",
    "\\begin{align*}\n",
    "  J_{f\\circ g}(a) = J_{f}(g(a)) J_g(a)\n",
    "\\end{align*}\n",
    "\n",
    "Which, gives, in the case of 3 mappings:\n",
    "\n",
    "\\begin{align*}\n",
    "  J_{c\\circ s \\circ l}(a) &= J_{c}(s\\circ l(a)) J_s(c(a)) J_l(a) \\\\\n",
    "  &= J_{cross}(s\\circ l(a)) J_{soft}(c(a)) J_{lin}(a)\n",
    "\\end{align*}\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gradient descent and convexity\n",
    "We won't prove here that optimizing this multiclass logisitc regression in convex, but let's assume it is, and that we can use gradient descent on it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now doing this by hand\n",
    "def crossEntropy(label,y):\n",
    "  '''label and y are a nSample*nClass matrices'''\n",
    "  #print('Cross input is '+str(label.shape)+'  '+str(y.shape))\n",
    "  return -np.dot(label.flatten(),np.log(y).flatten())\n",
    "\n",
    "def softmax(x):\n",
    "  '''x is a nSample*nClass matrices'''\n",
    "  #print('softmax input is '+str(x.shape))\n",
    "  exp = np.exp(x)\n",
    "  ret = (exp.T/exp.T.sum(axis=0)).T\n",
    "  #print('softmax output should be same as input '+str(ret.shape))\n",
    "  return ret\n",
    "\n",
    "def linearModel(weight, bias, features):\n",
    "  '''weight a nClass*nFeatures matrix\n",
    "     bias a nClass vector\n",
    "     features is a nSample*nFeatures matrix\n",
    "     output should be a nSample*nClass matrix\n",
    "    '''\n",
    "  #print('linear input is '+str(weight.shape)+'  '+str(bias.shape)+'  '+str(features.shape))\n",
    "  ret=np.dot(weight,features.T).T+bias\n",
    "  #print('output should be compliant with labels '+str(ret.shape))\n",
    "  return ret\n",
    "\n",
    "def loss(weight,bias,features,labels):\n",
    "  return crossEntropy(labels, softmax(linearModel(weight,bias,features)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gradientWb(weight, bias, features, labels):\n",
    "  \" return gradient for W and b as a tuple\"\n",
    "  nSample=labels.shape[0]\n",
    "  nClass=labels.shape[1]\n",
    "  nFeatures=weight.shape[1]\n",
    "  gradW=np.zeros(weight.size)\n",
    "  gradB=np.zeros(bias.size)\n",
    "  for sample in range(nSample):\n",
    "    # computing local values for reuse\n",
    "    linLocal=linearModel(weight, bias, features[sample,:])\n",
    "    softLocal=softmax(linLocal)\n",
    "    # computing cross entropy jacobian for single sample\n",
    "    Jcross_sol = np.divide(-labels[sample,:], softLocal)\n",
    "    Jcross_sol = Jcross_sol.reshape(1,-1)\n",
    "    # computing softmax jacobian for single sample\n",
    "    Jsoft_l = -np.outer(softLocal, softLocal)\n",
    "    np.fill_diagonal(Jsoft_l, softLocal*(1-softLocal))\n",
    "    # computing linear jacobian separately for W and b\n",
    "    JlW = np.zeros((nClass,weight.size))\n",
    "    Jlb = np.identity(bias.size)\n",
    "    for i in range(nClass):\n",
    "      JlW[i,i*nFeatures:(i+1)*nFeatures]=features[sample,:]\n",
    "    Jprod=np.dot(Jcross_sol,Jsoft_l).reshape(1,-1)\n",
    "    #print('Jprod is '+str(Jprod.shape)+' jlW is '+str(JlW.shape))\n",
    "    gradW+=np.dot(Jprod,JlW).flatten()\n",
    "    gradB+=np.dot(Jprod,Jlb).flatten()\n",
    "  return gradW.reshape(weight.shape), gradB.reshape(bias.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "def makeMultiClassLabel(label,nbclass):\n",
    "  '''Make 2d 1-hot encoded matrix from 1d integer list'''\n",
    "  #nbclass=np.unique(label).size\n",
    "  ret=np.zeros((label.size,nbclass))\n",
    "  ret[(range(label.size),label)]=1 \n",
    "  return ret\n",
    "\n",
    "def selectData(nSamples,data,label):\n",
    "  # We want data to be: dimension 0: list of samples, dimension 1: features\n",
    "  x = data[:nSamples,:]\n",
    "  y = makeMultiClassLabel(label[:nSamples],10)\n",
    "  return x,y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "stepLength = 1e-4\n",
    "nbIter=500\n",
    "\n",
    "#Select input data\n",
    "nbSample=800\n",
    "features, labels = selectData(nbSample,train_dataset,train_labels)\n",
    "\n",
    "#Sizes\n",
    "nbClass=labels.shape[1]\n",
    "nbFeatures=features.shape[1]\n",
    "\n",
    "#Perform gradient descent\n",
    "w=np.random.rand(nbClass,nbFeatures)\n",
    "b=np.random.rand(nbClass)\n",
    "losslist = []\n",
    "for i in range(nbIter):\n",
    "  gW,gB=gradientWb(w, b, features, labels)\n",
    "  w-=stepLength*gW\n",
    "  b-=stepLength*gB\n",
    "  l=loss(w, b, features, labels)\n",
    "  #print('loss is '+str(l))\n",
    "  losslist.append(l)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZYAAAEiCAYAAADao/T1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO2deZgU1dX/P9+ZYRFxQ1xG0YCK+rpEXOKSGEXjFqNxjcEYFTWicYvG9xc1vkY0mmjiEo1LYiIR424SE6O44IJLIooL4oqgooIDKCgIyDIz5/fHvT1T03T31Mx0z9DT5/M89XTXubeqTlVX16l77rnnysxwHMdxnGJR1dUKOI7jON0LNyyO4zhOUXHD4jiO4xQVNyyO4zhOUXHD4jiO4xQVNyyO4zhOUXHD4jiO4xQVNyyO4zhOUalpS2VJ6wLrASsBnwLvm9nSUijmOI7jlCdqbeS9pB2AHwH7AhtmFS8FJgB3AreZ2RelUNJxHMcpH/IalmhQrgB2A14DngZeAT4BvgT6AYOAnYBvxc1+A1xpZotLq7bjOI6zolKoj+UpYCKwpZltY2anm9koM/u3mT1mZveY2eVmdiiwDjACOBT4f52gd5ciabgkk7RJV+tSbCQNjec2tKt1WVGRNE7SuMR6ya9ZvOeO7+A+RkrqFskBJa0ez2e7lPVvkTQtsT4wbr9RyZRsXaehUYeqLPnAeD8N7yLVOkwhw7KxmZ1pZm+1thMzW2xmd5vZ9sBNxVPPccqCl4Fd4mepGA50yLB0M1YHLgRSGRbgl8AhifWBcfsuMyzA0KhD9nO4jnA/PdjZChWLvJ33ZjazPTs0s1ntV8dxSoekHkC9FTmlt5nNB8YXc59OcTGzd0t9DEnVhO6F+o7sx8yWUOb3U6pwY0n9JW2YJTtJ0u8lHVAa1cofST+U9KqkxZI+lfRXSbVZdX4g6RVJCyTNl/SapJMS5V+TNFbSHElfSnpP0g3t1OciSS/H43wq6QlJO6fYTpLOkjRZ0lJJdZKuk7RqVj2TdImkMyS9L+kLSU9J2jKrXnWsVydpUdRj87j9yBT67BWv2WJJUyX9KI+rwySdIuk3kj4GlgCrS1pL0h8lvROP/5GkOyStn+NYwyS9LWmJpDckHZKjTk5XmKRDJY2Px/hc0r05/kfTJN0Wj/OWpIWSXpS0a6LOOGB34BvxOJZ0xeW5RttKeiZeoxmSLgCUo16NpPMS5/ixpCsl9c6q80tJ7ybu5WeTOsZ6J8b760tJn8Xf/uuJ8j6SLo/3xtL4eb4SrqDEtfxuvMc+jcttklaPdQYC78dN/pS4JsMLXI+m+yP+Tk/GorGJ7Ycm6o9Qy//uzZL6Ze3TJF0q6VxJ7xOCmbaW1FvS1ZJeV/hfz5T0b0mbJ7YdSWitACzL6JA5v1zno3TPk1bvp1ivaM+VnJhZqwtwP3BDYv0CoBGYAzQA30+zn+6yENwSBmxSoM6IWOcuYH9CZN1s4B2gb6yza7yOvwP2AvYBzgDOieV9gbnAw8CBhKbzcOCmdur9Z+BoYA/ggKjbUmDrRJ2hUe+hCdmvouw6QnTgWcAC4BmgKlHPgGnAI8B3gcMJD4CpQE2i3iXxvC8H9gbOjdfFgJGtnMMWBAPxDHAwcAQhuORDYFqi3sC4vxnAP+P5HkQIld8MuAY4jBCcMowQ3TgN6J3Yx15Rz38D34nX/kOCq2JcK9fs5CgbFX//7wNvxeuxSqLeNOCDePzDo56vAJ8DqyfO+WXgVWDnuGxR4Br1Bz6Lx/t+vE7/AT4CLKvuXcBC4BfxfE+Px/57os758ff+CcHAHQhcBHw3UeeKeL5/juXfIbifhsXymvibzQHOJAT8nA8sJgT8ZF/L94HfE/4TpxMChkbHOr0Ibi0j3JuZa7JWgWtyS+b+AFYFTonbn57YftVYfhmwDLgyHv84wn30PFCddb/PiOd1GLAfob95tXgdhsXrdQgwNv4m68ZtB8Q6Bnwjo0PWvTu8Lc+TNtxPRX2u5LzeKR9IHwOHJdZnAJfE79cCzxdLoXJYaMWwANXALODJLPmucbsz4vr/AnMLHGeHWP+rJTiHasKffTJwTUI+lMRDkhD9twS4JWv7H8Z6yYeLAVOAHgnZ4VH+9bi+BuEhdUPW/n5KOsNyByEysU9CVkt4QE1LyDJ/zpeJ0Y+tXIsNYv1DEvL/AG/S0njuHOuNK3DN+gLzgFFZxxlEMORnJmTTCA+cNXL87j9IyMYBz6b8bS+Nx9kgIVuZMPbMErJvxuMck7X9UVE+JK4/APyjwPE2IbxgXlWgztFxn7tlyc+Puq6ddS1HZ9W7Lv7GmUjWzO/7o5TX5Jas+yNznL2y6g2M5/KLLPk3Yv2Ds+73j4GVUtxffYAvgLMS8pFxHzU5dGgyLKR8nqS9nyjhcyWzpB153y+eGJK2AtYFRseyfxLeAJ1mNgPWBm5PCs3sWcLbxO5RNAFYIzZdD8g09RNMIbxp/DE2gzfoiFIKLqQnJc0B6glvZZtS+PfbGegJ3JYlvyvuY/cs+VgzW5ZYfy1+ZlxAWxMecvdmbfe3VCcR9BljZosyAjOrA/6bp/4/Lf6bkkj6cXQrLCCcx4exaLNYXg18DfibmTUmjjWe8OctxC6Et+LboxupRlINocXwNqGVlOQ5M/sssZ59zdrKLsB4M/soofdCQssryX6Eh/rfsvR8NJZn9JwA7B/dPrtK6pm1n70IbvVCgTv7Ee79/+Y4Vg/C75oku+P6NUJLZZ0CxygGexPOJfu3e55gGLJ/u4fN7MvsnUg6QtLzkj4n3F8LCS8c7XlWpn2eZGjtfirqcyUXaQ3LHELTDWBP4GMzmxLXe7RhP5VCxhdbl6NsZqbczJ4Cvkd4W74P+ETSY5K+GsvnEdxWHwM3AB9Gv+1hbVVIISxzDKG1cALhj/w1gnuld4FNc56LhQ7KOYnyDHOz1pfEz8wxMj7h2Vn10gZ91ObYttD2y/0Gkk4nXM/HCCHyO9L8YMvo2Z9wb+fab2u6rh0/HyMY7+SyNbBmVv0W18xC521Sl7ZSm0fHbNnahJeGhVk6Zq5vRs9fEfoDvkt0Z0n6i6T+WfWmF9BpbeArLH89XsjaR4bW7qNSkfntprK8rquwvJ657q8DgbsJrsgfEMb6fY3Q0m6P/qmeJwkK3k/FfK7kI21Kl8eAkfFGOpvQSsmwOcFqOs1kfth1c5StC7yUWTGzvxHeGPsSmueXAw9LGmBmjWY2ETgsvjXtAJwH3CNpGzN7vQ06HUZ4czo02aKQtAbh7SXNubyR2K6G8CfLfgC0RubPsXZyf6R/E62j+c+fJN/2y7VWCL7vx83s7IxA0qCsOp8SHia59rsOhe/5OfFzOC3PMUOpM1TUkV/vJHMI7qVv5tnPxwDxfrkcuFwhrdMBwFUE9873CdcKYH2CazUXcwj9JkfkKZ+WR97ZZH67fQgupXzlGfLdX1PNbHhGoBCRmG0A0pL6eZKWIj5XcpK2pfEzQjP+18C7hI67DEcBz3ZUkW7GZMLb4bCkMEbIfIXgL2+BmS0wsweAPxLeONfMKq+PbpgLCL/b/7RRpz4E33HTH0HSnrTubhlPcJcMy5J/n/BiMq6NerxGeEP+XpY8e72QPvtL6pMRxMiYb7RBhz4Eo5HkuOSKmTUQO0CzopZ2IvjAC/FfgvHYxMxezLHke/gWYgkh8CANzwE7J10cklYmdNQmeZjwFrtaHj0/zt6xmc00sz8TXja3iuLHCEEOIwro9DChZb4gz7E+LbBtLjJv4WmvSdrtxxLOZcM8er6fYt99CC9xSY4m9JWk0SGbNj9P0lKE50pOUrVYLIxN2TtP8V6EiI1KZD9J2eN95pnZWEm/IPgwbyP0T6xP6FSdQogUQtLFhLfIJwlvhwMIUWETzewThVDuEYQW4vuEvokzCA+t5zIHjGGKo5NvSDl4mBCNc4ukvxD6Vi4gBGLkxczmSroSOE/SQoI77X8IkV3P0sZBXGb2maTfAT+X9AXhobQdwT0H4U9diEsIAQGPSLqC4He/gPDHa23bDA8D50j6OcEVs2fcZzYXEvoA/inpj8BahJeqgmO8zGy+pP8HXC9pLeAhQmf++gR/+DgzuyOlrhneBE6R9H3Cy90XBQzU1YSop0djWOsSQkaMFv9TMxsn6U5Ci/kqwrVoJBjO/QnRie9I+hfBZfoy4S1+W0KfyR/jft6VdDXwU0mrEKJIGwguxrfN7G5C/8BxwOPxfnqV4IbbmOBiOzjZb5aCWYTWwzBJkwgvK++bWXaLIh/vEB7+x0uaG6/R5HgulwPXSdqMkIFkMcEo7g382cyebGXfDwMHx2vyAKFFkIm2S/Jm/Dxb0kNAg5m9mL0zM2tI8zxJS9rnSodI08MPPAFsnqdsU+CJUkUXrIgLzVFhuZbXE/V+SPgDLSH8Cf4K1CbKv0MIza2LdT4CbgbWi+WbEXy17xNu7k8ID/adEvtYOR73shR6nx739SXhbXwvwtvOuESdoSwfOitCiPFkQuulDrieGJ6ZqGfEaMGEbCDLh05WE/4UM6Mu44Cvx3o/SXEeexPSDS0B3gNOIvRRvZLjuMtFDRHeEG+M1/MLwp9/EDmi0oAj43kvIbi1DklzzaJ8f8JLw3xgEc0PgS0SdaYRErhm69hCF4LLY0zUt0VUWp5rtB2hP2Qx4eXhAoJRtKx6VYQw4ldj3Xnx+28ILRkI7u/xhHv4y3g9RpKI/ov1TgYmxWs1N16nXRLlveN2byfqTIiymqxrmR2tNTzKByZkBxMezsuy77Ec1+MWElFhUXZSvH/qs38/QgtjPMFgLSD0l1wHDCh0vyeu6SWEl8VFBOO0bfytb8n6H1xP6NNqzPw25PjPpHmepL2fSPFc6ejSanZjAEmNhBjrF3KUbQ+8YGbZzTynE5C0DyHaZ2MzK9R5ukIj6XBCpNhuZvZMG7ftS+hsfdDMTmitvuM4paUt87Hks0AbEyy60zXsTnCDlY1Rif0U3yGEcC4GticMkhxPiv46Sb8n9GN8TJgf6CeE8THXlEhlx3HaQF7DIuk4mjs0Dbgp+sSTrETowHu8NOo5rWFm53e1Du1gAWE8wKmE8R6zgXuA8yxNEzq4VC4n9E8tJfQN7GVmk0qjruM4baFQi6WR0AEHwceeXM8wh+CrvjzNwSSNIoQqzjazraLsbpoHDa0OfG5mQxTyAb1Fc/jieDM7OW6zPcFnuhLBN/gTMzOFXD53E3yU04AjrOVAIWcFwMzeIPjS27v9icXTxnGcYpO2j+VJ4Mdm9naHDibtRnhbvTVjWLLKryREVV0cDcsDeeq9QIhieJ5gWK41s4ck/YaQIuUySecS0hqc04pOrxCifaZ25Nwcx3EqjE2AT8xs2+yCtOHGexRDCzN7OhqM5ZAkwuCpPQvtI45ZWNVC7DWSbiVEhzxESDI4NFYdTYhKKWhYgLVWW2219YcMGbJcZlvHcRwnNxMnTmTevHk5ywr1sRxDiLKZE78XxMxubb+KQBj9O8uaU8UADIotivnA/8VoofVpmTpiepQBrGMhbxSEUNY0o7mnDhkyZP1x48Z1SHnHcZxKYujQoTz11FM5PT2FWiy3EPInzYnfC2FARw3LkcCdifU6wujXObFP5Z/KmtejoEKhzyWvn0/SCMIgoc3q6nKl4HEcx3HaQyHDMojmvE7ZeZSKSsxXcygh7BRoSpy2JH5/SdK7hMGYM2hOiEn8nhk9PktSrZnVRZdZrmSFmf3fRIh0G1dbW5udHdRxHMdpJ4WmJv4g1/cSsRch9UOTiyumwphrIZ3BRsBg4D0LKUbmK8x8+DxwDGFCIAipJI4lTNRzLPCvEuvtOI7jZNGWAZJI2oSQ/2d94oxq1oa5pGNeoqFAf0nTgQvN7GZCcrU7s6rvBlwsaRkh1PlkM8tk+TyF5nDjh+ICwaDcI+kEQvbZfJlUHcdxnBKRyrAozH99A8tn6GyQNBo41Zpz/ufFzI7MIx+eQ/Z34O956r9Ic2bVpHwOYcpTx3Ecp4tImzb/CkJ6/AsJscurxM+RBGPz21Io5ziO45QfaV1hw4CLzOxXCdl7wKVh+AlnEQYsOo7jOBVO2hZLL5qnEM3mecK8ChWJmXHhv17n7gkftl7ZcRynAmjL1MT7xM9s9iHM11KRvDZjHqOfC0Fz3/9aa5MxOo7jdH8KjbzfKLF6FfDXOL3pvYTZ29YhRF3tT5iApiJZUp920kLHcZzKoFCLZSot52AR8GPCLHFJGYQZ0ipyoq+aKrVeyXEcp4IoZFiOK1DmRKrUbFgaGo1qNzSO41Q4hUbej+5MRcqV+sZmV9iXyxro26tNY04dx3G6HWmjwpw8LK1v9hYuWlrfhZo4juOsGLhh6SDLGhItlqXZE2w6juNUHm5YOkjSsCxyw+I4juOGpaMsrW/Zx+I4jlPpuGHpIEvdFeY4jtMCNywdZFlDsvPeDYvjOE7q2Ng4Ev8IYEOgd1axmdkJxVSsXGjZx+JRYY7jOGnnYzkYuIfQwplNnDI4Qd655bs7HhXmOI7TkrQtll8C44CjzOyT0qlTfnjnveM4TkvSGpaNgLPdqCyP97E4juO0JG3n/dvAmqVUpFxJusLqGyrWI+g4jtNEWsPyM+DnWan0HVoalgZzw+I4jpPWFTaS0GJ5S9IUYG5WuZnZ7sVUrFxIjmNpbHTD4jiOk7bF0gBMBv4LfBLXk0uq2a4kjZI0W9LrCdlISTMkTYzL/omy8yRNlTRZ0r4J+X5RNlXSuQn5IEnPR/ndkko+ZfKyRBJKb7E4juOkbLGY2dAiHe8W4Drg1iz51WZ2RVIgaQtgGLAlsB7wmKRNY/H1wN7AdGCCpPvN7E3g8rivuyT9ATgBuLFIuudkmbdYHMdxWtCpI+/N7GmWd6Pl4yDgLjNbYmbvE2a03DEuU83sPTNbCtwFHCRJwJ7A3+L2o4GDi3oCOWjRee+GxXEcp+Cc97sBL5vZgvi9INFotJfTJB0DvEgIa/4MWB8Yn6gzPcoAPsqS70ToA/rczOpz1F8OSSOAEcBmdXV17VY82cfS4IbFcRynoCtsHLAz8EL8nu+pqVjW3jnvbyQMwLT4eSVwfDv3lRozuwm4SdK42tradgceJMexNHofi+M4TkHDsgfwZuJ7STCzWZnvkv4EPBBXZwAbJKoOiDLyyOcAq0uqia2WZP2SsazeWyyO4zhJCs15/1Su78VGUq2ZZXxRhwCZiLH7gTskXUXovB9MaD0JGCxpEMFwDAN+YGYm6UngcEK/y7HAv0qld4YWnffeYnEcx0mf3bgYSLoTGAr0lzQduBAYKmkIwRU2DTgJwMzekHQPodVUD5xqZg1xP6cBjxDcb6PM7I14iHOAuyRdArwC3Fzqc/I+FsdxnJYU6ry/H7jQzF5JsyNJvYFTgEVm9odcdczsyBzivA9/M7sUuDSHfAwwJof8PULUWKexpIUrrDOP7DiOs2JSKNx4GjA+Djg8Q9J2kloYIknrSTpY0s1AHWHcyMulU3fFI5kq311hjuM4hftYzpB0DXAmIaXLaoBJmk+Yj2V1oCehz+OFWO+2jLuqUkhO7uWuMMdxnFb6WMzsXeB0SWcDuxDGi6xHmEFyDiHr8dNm9kGpFV1RSbZYPKWL4zhO+pQuS4Gn4uIkWJh0hXmLxXEcp3NTunRHWrRY3LA4juO4YWkvDY3GWXdPbJk2311hjuM4bljaS3WVuO+VlgP7vcXiOI7jhqVD9O7R8vL5zMSO4zhuWDpEn54tYx8aGn2EpOM4jhuWDrBSj5YJnd0V5jiO0wbDImlbSf+Q9KmkeknbRfmvJO1XOhVXXFbu1dKweIPFcRwnpWGRtCvwHLA5cEfWdo3AycVXbcVnpWxXmEeFOY7jpG6xXEbIJrwl8NOsspeB7YqpVLmwck93hTmO42ST1rBsB9xoZsbyM0l+CqxVVK3KhOFfHwjAwDX7AD6OxXEcB9IblsVAnzxltcC84qhTXuyz5bqM+9+h/G7YtoC3WBzHcSC9YXkWOFNS0veTeYqeADxRVK3KiIH9V6ZHtQA3LI7jOJB+BskLgP8ArwJ/IxiVY+O0wdsDXyuNeuVBdVUwLO4KcxzHSdliMbNXgd2AWcD5hDlYTovFu5vZ5NKoVx5Uy1ssjuM4GVLPeW9mLwPfilMQ9wM+N7NFJdOsjKhqarF0sSKO4zgrAKkNSwYzWwx8XAJdyhZvsTiO4zSTyrBI+kUrVczMfpliP6OAA4DZZrZVlP0WOBBYCrwLHGdmn0saCLwFZNxs483s5LjN9sAtwErAGOAnZmaS+gF3AwOBacARZvZZmnPsCJk+FjcsjuM46VssIwuUZZ6mrRoWgjG4Drg1IRsLnGdm9ZIuB84Dzoll75rZkBz7uRE4EXieYFj2Ax4CzgUeN7PLJJ0b18/JsX1RqfLOe8dxnCbSdt5XZS9Af2A48DqwScr9PA3MzZI9amb1cXU8MKDQPiTVAqua2fg4YPNW4OBYfBAwOn4fnZCXFHeFOY7jNNPu7MZmNtfMbiW0Qq4vkj7HE1oeGQZJekXSU5K+GWXrA9MTdaZHGcA6ZlYXv88E1imSXgWpilfRDYvjOE5x0uZnQpE7hKTzgXrg9iiqAzY0s20J+cnukLRq2v3lST+TPN4ISS8C29fV1eWrloqaaFk8CaXjOE5xDMsBwCcd2YGk4XE/R0WDgJktMbM58ftLhI79TYEZtHSXDYgygFnRVZZxmc3Od0wzu8nMdgBeqq2t7Yj67gpzHMdJkDYqbFQOcU9gK2Br4ML2KhDncvkZYaDlooR8LWCumTVI2ggYDLxnZnMlzZe0M6Hz/hjg93Gz+4FjCdmYjwX+1V692kLGFdbohsVxHCd1VNieLO9WWgx8APyO5g7zgki6ExgK9Jc0nWCQzgN6AWMV3vwzYcW7ARdLWkac88XMMh3/p9AcbvwQzf0ylwH3SDoh6nZEyvPrEE3hxu4KcxzHSWdYzGxgMQ5mZkfmEN+cp+7fgb/nKXuR0FrKls8BvtURHdtDVXSF+QySjuM4Pud9UfAWi+M4TjN5WyyS2hTpFceoVCTeee84jtNMIVfYOAqE6yZQrFfdWsXuSmbkPYQO/OS64zhOpVHIsOzRaVp0A6qrREOj0WBGFW5YHMepXPIaFjN7qjMVKXeqJRowGhqNHhXbdnMcx/HO+6LRNJbFO/Adx6lwUs/HImlL4EfAZkDvrGIzs04P812R8A58x3GcQNqR9zsBTxHmOBkMTALWADYkJIGcWiL9yoam1Pk+lsVxnAonrSvsV8A/gC0JUWAnxEGTexGiwS4piXZlhI9lcRzHCaQ1LF8FbqM5/LgawMyeIBiVXxdftfIi4wqr9yaL4zgVTlrD0hNYaGaNhIm6kumAJ5MjvUqlUe2uMMdxHCC9YZlK82Rak4DjJVVJqgKOI0yqVdG4K8xxHCeQNirs34SsxHcQ+lseBOYDDUBf4IxSKFdONCeidMPiOE5lkza78cjE98fiXCiHAX2Ah83s0dKoVz40tVjcsDiOU+GkHseSxMxeAV4psi5ljbvCHMdxAqn6WCTdJ+lgST1KrVC5ksk76a4wx3EqnbSd95sRxrHMlHRDdIU5CbzF4jiOE0hlWMxsC+BrwF+BQ4H/SJoi6RdxPvqKp8pTujiO4wBtSEJpZi+Z2ZmEsOMDgQnAOcAUSc+USL+ywcexOI7jBNqc3djMGsxsjJn9ADgE+Bj4etE1KzPcFeY4jhNoc1RYdH0dDRwFbAzUAVcWWa+yw11hjuM4gbRRYWtIOlnSf4ApwP8C44F9gQ3M7GdpDyhplKTZkl5PyPpJGhv7bcZKWiPKJelaSVMlTZK0XWKbY2P9KZKOTci3l/Ra3OZaSZ0ynWOTK8xbLI7jVDhpXWEzgeuAhcCxwDpmdoyZPWbW5ifpLcB+WbJzgcfNbDDweFwH+DYhTf9gYARwIwRDBFwI7ATsCFyYMUaxzomJ7bKPVRJ8PhbHcZxAWsNyPrChme1jZreZ2aL2HtDMniYkskxyEDA6fh8NHJyQ32qB8cDqkmoJLaWxZjbXzD4DxgL7xbJVzWx8NHi3JvZVUppmkHTD4jhOhZM23PgKM/u4hHqsY2Z18ftMYJ34fX3go0S96VFWSD49h3w5JI2Q9CKwfV1dXa4qbSLjCqt3w+I4ToWzws15H1saJX86m9lNZrYD8FJtbW2r9VujOjZZPCrMcZxKZ0UxLLOiG4v4OTvKZwAbJOoNiLJC8gE55CWn2lO6OI7jACuOYbmfEBRA/PxXQn5MjA7bGZgXXWaPAPvEaLU1gH2AR2LZfEk7x2iwYxL7Kime3dhxHCfQruzGHUHSnYS5XfpLmk6I7roMuEfSCcAHwBGx+hhgf8JEY4sIk4phZnMl/ZIw+h/gYjPLBAScQog8Wwl4KC4lp2k+FneFOY5T4XS6YTGzI/MUfStHXQNOzbOfUcCoHPIX6YKpkptbLJ19ZMdxnBWLtAMkD5J0XGL9K5Kek/SFpL9J6ls6FcuDKk/p4jiOA6TvY/k/YK3E+lWEjvGbgN2AkcVVq/yo9qmJHcdxgPSGZWNgEoCklQj9Hj81s7OBnxOSUVY03nnvOI4TSGtYegNfxu9fJ/TNZOa5nwysV2S9yo6mJJTuCnMcp8JJa1imAbvG7wcBL5nZvLi+NjAv10aVRLWndHEcxwHSR4X9EbhC0iHAEODHibJdgDeLrVi54fOxOI7jBFIZFjO7RtKnwM7AtWZ2a6J4FeAvpVCunKjyznvHcRygDeNYzOx24PYc8pOKqlGZ4p33juM4gbTjWDaVtGNifSVJv5b0b0mnlU698qG5876LFXEcx+li0nbeXwccnli/FDibEA12taSco+MriZqmFosPvXccp7JJa1i2Af4DIKmKkNzxHDPbHriEMLtjReMpXRzHcQJpDctqwJz4fVtgDeBvcX0csFFx1So/MildltQ3eD+L4zgVTVrDMgvYJH7fB3jXzDIzOPYF6outWLmRSenyu8emcPD1/+libRzHcbqOtFFh9wO/lrQVMJwwriXD1sB7Rdar7Mi0WABem1Hx40Udx6lg0hqWcwlpXfYlGJlfJcq+S3N6l4ol02JxHMepdNIOkFwInJin7OtF1ahMqV5R5uJ0HMfpYto00ZekfoQULv2AuTCf87gAAB6FSURBVMBziZkbK5qkK8xxHKeSSW1YJF1CGLvSKyFeIukKM7ug6JqVGe4KcxzHCaQdeX8mYd6V24A9gP+Jn7cBP5d0Rsk0LBOqvcXiOI4DpG+xnAxcY2ZnJWSTgackLQBOAa4ttnLlRFVWi8XMkLdiHMepQNJ2OQ8EHsxT9mAsbzeSNpM0MbHMl3SmpJGSZiTk+ye2OU/SVEmTJe2bkO8XZVMlndsRvdpCdovFB0k6jlOppG2xzAG2Ah7LUbYlzaPy24WZTSbM84KkamAGcB9wHHC1mV2RrC9pC2BYPPZ6wGOSNo3F1wN7A9OBCZLuN7OSzxeT3XnfYNa2yAjHcZxuQtoWy33ALyUdLakGQFKNpCOBi4G/F1GnbxFG9n9QoM5BwF1mtsTM3gemAjvGZaqZvWdmS4G7Yt2Sk91577koHcepVNIalvOAicBo4EtJs4AvCfOzvEro2C8Ww4A7E+unSZokaZSkNaJsfeCjRJ3pUZZPvhySRkh6Edi+rq6uw0pnj2PxmSQdx6lUUhkWM/sC2I0wyv4qwuj7q4ADgN3NbEExlJHUMx7j3ii6EdiY4CarA64sxnEAzOwmM9sBeKm2trbD+6uuankpvY/FcZxKpS0zSBrwQFxKxbeBl81sVjzmrEyBpD8ljj0D2CCx3YAoo4C8pCzXYnHD4jhOhbKiJSI5koQbTFKyKXEI8Hr8fj8wTFIvSYOAwcALwARgsKRBsfUzLNYtOdnhxm5YHMepVPK2WCQ1AmmfjmZmHQqCkrQyIZrrpIT4N5KGRD2mZcrM7A1J9wBvElL2n2pmDXE/pwGPANXAKDN7oyN6pSU73LjR+1gcx6lQChmDi0lvWDpMTHS5Zpbs6AL1LyVMkZwtHwOMKbqCreAtFsdxnEBew2JmIztRj7Knbt7iFutuWBzHqVRWtD6WsqV/354t1t0V5jhOpeKGpUjsv3Utvzx4q6Z1b7E4jlOpuGEpEj2qqzh656+w8VorA25YHMepXNywFJlMdJiPvHccp1Jxw1JkMtFh3mJxHKdSSTvR16at13KgucXiSSgdx6lU0rZY3pb0uKTvZbIbO7lxV5jjOJVOWsNyPLAScDcwXdKvYioVJwt3hTmOU+mkzW58i5l9nZBl+O+EqYinSHpY0kGSvK8m0uQK8xaL4zgVSpsMgplNMrNTCbM2ngSsA/wD+DBOI7xOCXQsKzKGpb7BDYvjOJVJe1saA4Gvxs+lhKzDPwWmSjqkKJqVKZmZJL3F4jhOpZLasEjqKekoSU8DrwEHApcBG5jZfsBXgIcJE4BVLE2d997H4jhOhZIqwkvSlcAxwBqElPTfBcbEyb8AMLPPJF0DPF0KRcuFKo8KcxynwkkbOnw0MAr4g5m9X6De28BxHdaqjKmO2fMbvcXiOE6FktawDDCzpa1VMrNPgdEdU6m8cVeY4ziVTirDkjEqkjYBdgTWJ8wl/7yZvVs69cqPKu+8dxynwknbx9IbuIHgEqtOFDVIGk2YGnhJCfQrO2qiL6zeWyyO41QoaaPCrgCOAi4ENgFWiZ8jCcbmt6VQrhzxkfeO41Q6aftYhgEXmdmvErL3gEsVHqRnAWcUWbeyxEfeO45T6aRtsfQCXshT9jzQM09ZxVHd1GLpYkUcx3G6iLSG5TFgnzxl+wBPFEMZSdMkvSZpoqQXo6yfpLGSpsTPNaJckq6VNFXSJEnbJfZzbKw/RdKxxdAtLVVNafO9xeI4TmWS1hV2FfBXSSsD9wKzCHnCjgD2B34oaaNMZTN7rwM67RHDljOcCzxuZpdJOjeunwN8Gxgcl52AG4GdJPUj9AXtABjwkqT7zeyzDuiUmqYWi7vCHMepUNIalqfi54+BkxNyZZVnqKZ4HAQMjd9HA+MIhuUg4NY4+n+8pNUl1ca6Y81sLoCkscB+wJ1F1CkvVT6OxXGcCietYems0fQGPCrJgD+a2U3AOmZWF8tnElpKEMbSfJTYdnqU5ZO3QNIIYASwWV1dXXZxu6lxw+I4ToWTdoBkZ42m39XMZkhaGxgr6e0sPSwanQ4TjdZNksbV1tbuXox9go+8dxzHaVPa/NhhvqWkb8ZPtb5VesxsRvycDdxHGOU/K7q4iJ+zY/UZwAaJzQdEWT55p+Aj7x3HqXTakjb/R0AdMInQzzEJ+FjSCcVQRNLKklbJfCdEm70O3A9kIruOBf4Vv98PHBON3c7AvOgyewTYR9IaMYJsnyjrFKrjFfUWi+M4lUralC5HATcBjwO3Efo61iWMxr9J0iIz62jn+DrAfbERVAPcYWYPS5oA3BMN2AeESDSAMYSItKnAImI/kJnNlfRLYEKsd3GmI78z8LT5juNUOmk7738G3G5mR2fJR0v6KyFKq0OGJYYob5NDPgf4Vg65Aafm2dcoQpr/TqdpBklvsTiOU6GkdYVtRmip5OK2WO6Q7LzvYkUcx3G6iLSG5QtCJ3guBsRyh6RhccviOE5lktawPAT8StI3k0JJuwCXxHIHH3nvOI7Tlj6WnYFxkmYQosPWJbRWpsZyh+TI+y5WxHEcp4tIO0BypqQhwPHAN4F+wDRCKpdbzGxRyTQsM3pkJvpyy+I4ToXSqmGR1IMQ1jvJzK4Driu5VmVMTVXwLvoMko7jVCqt9rGY2TLgHmBgybXpBmRaLMu8xeI4ToWStvP+PWDtUirSXaiJQ+/rG7zF4jhOZZLWsPwGOF/SWqVUpjuQyW68zMONHcepUNJGhe1J6LB/X9J4QlRY8pXczKxTZ2pcUenhLRbHcSqctIblm8Ay4BNg47gk8adopCYTFeYtFsdxKpS04cYDS6xHtyETFbbMWyyO41QoqfpYJO0mqW+espUl7VZctcoXH8fiOE6lk7bz/klgizxlm8dyh0RUmI9jcRynQklrWArNFNkLaCiCLt2CHlU+jsVxnMombx+LpIHARgnRDjncYSsR0rx8WHTNypTkOJbXps/jR7dO4MIDt2T/rWu7WDPHcZzOoVDn/bHAhYSILwN+T8uWi8X1evJMuFWJZKLCljUaZ90zkVnzl3DK7S8z7bLvdLFmjuM4nUMhw3ILYW57AU8QjMebWXWWAO905tS/Kzo9MrnCGhrdHeY4TkWS17CY2QeEOeaRtAfwspn5hF6t0KMmExVmNPqcLI7jVCBpx7E8VWpFugtN41gaG/Exko7jVCJpx7H0lHShpLclLZLUkLXUd0QJSRtIelLSm5LekPSTKB8paYakiXHZP7HNeZKmSposad+EfL8omyrp3I7o1R6ax7EY5i0Wx3EqkLQpXX5L6GN5CPgHoW+lmNQDZ5vZy5JWAV6SNDaWXW1mVyQrS9oCGAZsCawHPCZp01h8PbA3MB2YIOl+M8vuGyoZzVFhjT49seM4FUlaw3I4cKGZXVoKJcysjpDYEjP7QtJbwPoFNjkIuMvMlhASY04FdoxlU83sPQBJd8W6nWZYmsax+ABJx3EqlLQDJPsCz5VSkQxx/My2wPNRdJqkSZJGSVojytYHPkpsNj3K8slzHWeEpBeB7evq6oqmf7LF0ujGxXGcCiStYfk3UPJ8YHEA5t+BM81sPnAjIZPyEEKL5spiHcvMbjKzHYCXamuLN3ixJtHH4q4wx3EqkbSusN8Dt0pqBMYAy41bybif2oukHgSjcruZ/SPuc1ai/E/AA3F1BrBBYvMBUUYBeaeQGceytKGRXkprtx3HcboPaQ1Lxg02kjAaPxfV7VVCkoCbgbfM7KqEvDb2vwAcArwev98P3CHpKkLn/WDgBcJgzsGSBhEMyjDgB+3Vqz00z8di9PAGi+M4FUhaw3I8pZ3M6xvA0cBrkiZG2c+BIyUNiceeBpwEYGZvSLqH0ClfD5xqZg0Akk4DHiEYulFm9kYJ9V6OzNTEDY1Gg/exOI5TgaQdIHlLKZUws2fJnUF5TIFtLgWWi1IzszGFtis1kqipEvWNxlJP6eI4TgXSpk4ASVWStpK0u6SVS6VUuZNxh3mLxXGcSiS1YZF0KjATmERISrlZlP9T0hmlUa88yXTgO47jVCJpU7qcCFwD/BM4gpZuq2eAw4qvWvmSabE4juNUImlfrX8KXGlmI4D7ssreJrZenEBmkKTjOE4lkvYJOIgQaZWLhcDqxVGne5BJ65LEE1I6jlMppDUsnwID85RtRicPQlzR6Vmz/GWtz9OR39BofLm0odQqOY7jdBppDcsDwC8kbZSQmaT+wFmEvhcnUrvaSsvJsmeTbGw0Trz1RTb++Ri2ufhRFizp0MwDjuM4KwxpDcv/EVLlvw48RhiweC3wFtAAXFwS7cqUgf37LCdbVt+yxfL2zC8Y+2bIWLO0vpE3P57fKbo5juOUmlSGxcw+BXYAfg30AN4lDK68DtjFzOaVTMMyZMN+yw/xyR4subje3V+O43RP0qZ0Ic53/8u4OAVYvU+Ppu9rrdKLT75YspwrbN6iZS3Wv1jcct1xHKdc8bjYErD/1rUM2WB1LjhgC1bqEXJzZhuWuQuXtlif96UbFsdxugepWyxOelZbqQf/PPUbANzx/AfA8obls0VuWBzH6Z54i6XE9IiDJRcvK9ximfjR552mk+M4Tilxw1Ji1ls9hB5/OHdRk2zW/MXcMO7dFvX+NfFjnnt3Tqfq5jiOUwrcsJSYwev0BWDKrAVNsvsnfpyz7p0vfNgpOjmO45SSdhsWSVtIOkzSesVUqLsxeO1VAHhn9hcAjJs8m0vHvNVUvtOgfk3fa3KkgnEcxyk30mY3vk7SHxLrhwKvAvcCb0r6Won0K3u2XG9VAMa9PZurx77D8L9MaCp78IxdueJ72zSt/+OVGTw4qW65fTiO45QTaVss3wb+m1i/iJDmZRvCXPMXFlmvbsPm667Cd75ay8KlDVzz+JQm+TcH92fzdVdlg359ePzs3ZvkNz3zXleo6TiOUzTShhvXEuacR9IAYEvgBDN7TdK1wM2lUa/8kcSV39uGQWuuzNszv2CL9Vbl21uty//UrtpUZ+O1+nL49gP420vTmfbpQhobjSp3izmOU6akNSyLgL7x++7AfODFuL4AWKXIenUreveo5n/3LTxlzW8P/yrPTvmUmfMXc+xfXmD3Tddi8DqrMHjtvtSu1hvJDY3jOOVBWsPyMnCqpA+BU4GxZpYZmDEI8I6BDiKJiw7akv+991WemfIpz0z5tKmsb68aBvbvw5or92LNvj1ZZ9XerNGnByv1qKbfyr1YpXcNK/eqpk/PGvr2CkufXtX0qKrylo/jOJ1OWsNyPvAwocP+c+DkRNnBhH6WFQZJ+xGmUq4G/mxml3WxSqnYd8t12f4razDmtTremfUFU2Yt4N1PFvDpgqW8PqPt2Y97VIteNdWs1LO6yeCs1LOaPk1LkK3Su4Y+PWvo07OalXu1NE41VaJPzxpW7V1Dn141VEv0qJbPkuk4Tl5SGRYzmyBpQ2BzYIqZJZ9yNwFTcm/Z+UiqBq4H9gamAxMk3W9mb3atZuno37cXx+wysIVs7sKlfDh3EZ8tXMonXyxh5vzFzP9yGQuXNjB34RLmf1nPgiX1LFxaz8Il9SxYXM/CpQ0sazCWNYSyT75YUlQ9e9ZU0au6iupqUVMlelRXBVlNFT2qq6iprqJntaipqqIm1qmSkELrrEpQpWZZVUJG1nrL+pntM9ulq5NZb66fXbdAnaqoF806V2WfR1VYV746Vc3r2W3IpJtTLeSJ74mSfF7RNPVT1cmzz+zS/NuU4Hzy7LM10nqQl/9VOr7PsN9iV0yva1o9+/aqoXfMaVgs2pLdeCHwUlImaU0ze7CoGnWcHYGpZvYegKS7gIOAsjAsuei3ck/6rdyzTds0NBqLltbz5bIGvlzawIJocBYta2DRkga+XNbAoqX1zfKloV4oD/JFSxswjEVLGpi/uJ5FS+tpNGNJfSNL4+I4TnlzycFb8cOdv1LUfaYyLJJOBFY3s9/G9a2Bh4BaSa8AB5jZzKJq1n7WBz5KrE8HdsquJGkEMALYrK6u+3URVVeJVXr3YJXePVqv3EYsY1waGqlvMOobGlnWaCytb2RJfQPL6o1ljaFsWUMjyxoaaTSjsREazbC4j0YL642WWU/UMTBy1GlslmXXMSOrPG6Xok4LvRpz6NWGOsBy59IY6zRkz1BtzQLLLcYSJS3kLeq0/H1y/24p9llgP5ZnJd82HTmfFofKs8/WyLe/5eq1Ya9p9xn2W/x9pt1rW/bZK8dU6h0lbYvldILLK8NVhL6Wy4EzCDNIjiiuaqXFzG4CbpI0rra2dvdWN3CakETvHtVFbz47jtM9SGtYvgK8DSBpNULI8cFmNkbSHMLMkisKM4ANEusDosxxHMfpBNK2gaqAjEN9V0J7bFxc/whYu7hqdYgJwGBJgyT1BIYB93exTo7jOBVDWsMyBfhO/D4M+K+ZZfLArwfMLbZi7cXM6oHTgEeAt4B7zOyNrtXKcRynckjrCrsC+KukY4E1gO8lyvYAJhVbsY5gZmOAMV2th+M4TiWSdhzLHXHU/U7ABDN7OlE8C3c1OY7jOJG2jGN5Fng2h9wzGzuO4zhNpDYskvoAxxMiwvoR+lWeBP5iZl+WRj3HcRyn3FC+wVQtKknrEqLANgU+AGYC6xLCkCcDQ81sVunULB2Spq+22mrrDxkypKtVcRzHKRsmTpzIvHnzZpjZgOyytIblVmBf4FAz+09C/nXg78AjZja8eCp3HjFzwFrA1HZsnrFGE4un0QqPn3Nl4OdcGXTknDcBPjGzbbML0hqWT4BzzGxUjrITgMvMbK12KFbWSBoHYGZDu1aTzsPPuTLwc64MSnXOacex9AU+zlM2neZJwBzHcZwKJ61hmQwcnafsh8R0L47jOI6T1hX2Q+BW4AngDsKMkesSRuHvBRxtZneUUE/HcRynTEhlWKApzfzFtMwLNgv4hZn9qQS6OY7jOGVIasMCIKkK2IzmcSyTzcxne3Icx3GaaNWwxAzB44FzzezRTtHKcRzHKVta7bw3s6XAIKC+9OqUF5L2kzRZ0lRJ53a1PsVC0ihJsyW9npD1kzRW0pT4uUaUS9K18RpMkrRd12nePiRtIOlJSW9KekPST6K8O59zb0kvSHo1nvNFUT5I0vPx3O6OL5ZI6hXXp8bygV2pf0eQVC3pFUkPxPVufc6Spkl6TdJESS9GWUnv7bRRYWOBfdpzgO6KpGrgeuDbwBbAkZK26FqtisYtwH5ZsnOBx81sMPB4XIdw/oPjMgK4sZN0LCb1wNlmtgWwM3Bq/C278zkvAfY0s20Ig+T2k7QzYVbYq81sE+Az4IRY/wTgsyi/OtYrV35CmFIjQyWc8x5mNsTMdojrpb23Lc77XWgBvklI5XIFYaKvjYGNkkua/XSnBdiFkHEgs34ecF5X61XE8xsIvJ5YnwzUxu+1hP41gD8CR+aqV64L8C9g70o5Z6AP8DIhe/mnQE2UN93jhPmNdonfa2I9dbXu7TjXAfFBuifwAKAKOOdpQP8sWUnv7bQtlqcI0/3+NH5/hzD5V3KpNNYnzJ6ZYXqUdVfWMbO6+H0msE783q2uQ3R3bAs8Tzc/5+gSmgjMJngl3gU+tzBZHrQ8r6ZzjuXzgDU7V+Oi8DvgZzTPiLsm3f+cDXhU0ksxuhdKfG+nzW58XFt37HRfzMwkpQ8nLBMk9SXkvjvTzOZLairrjudsZg3AEEmrA/cBm3exSiVF0gHAbDN7SdLQrtanE9nVzGZIWhsYK6nFgPZS3NtpJ/oaXcyDdhNmEFpxGQZEWXdllqRaM6uTVEt4y4Vuch0k9SAYldvN7B9R3K3POYOZfS7pSYIbaHVJNfENPXlemXOeLqkGWA2Y0yUKt59vAN+VtD/QG1gVuIbufc6Y2Yz4OVvSfcCOlPjezusKk1Ql6UBJWxWos7WkA9t60G7CBGBwjCjpSchC0J1n0rwfODZ+P5bQD5GRHxOjSXYG5iWa2GWBQtPkZuAtM7sqUdSdz3mt2FJB0kqEPqW3CHMsHR6rZZ9z5locDjxh0QlfLpjZeWY2wMwGEv6vT5jZUXTjc5a0sqRVMt8JQVivU+p7u0CHzzHAAmBQgToDY50j89XpzguwP6G/6V3g/K7Wp4jndSchbc8ygo/1BIJv+XFCf9pjQL9YV4TouHeB14Adulr/dpzvrgQ/9CRC+vCJ8bftzuf8VeCVeM6vEzJoQAjGeYEwjcS9QK8o7x3Xp8bysg7YAYYCD3T3c47n9mpc3sg8p0p9b+cdICnpUUKkwOk5KzTXuwbYzMyyw1Mdx3GcCqRQVNh2QJqR9o8BO7Ray3Ecx6kIChmWVQiDhVrjs1jXcRzHcQoalk8Jc9q3xoaxruM4juMUNCzP0hw1UIjhsa7jOI7jFDQsvwO+JenqTFK2JJJ6SPodITXC1aVS0HEcxykv8hoWM3sOOBs4gzBA6DZJl8blNkIY6mmE5H3jO0ddp5KRNFySSfo8k401UVYTy0bmqL9JifXaJWa/XRiPNyRPvZHJEc6SVo+yLsuOLGlI1KFfjrIW19Nx0lJw5L2Z/U7Sy8A5wCHASrHoS2AccJmZPVNSDR1neVYj3JMrylQFNxP+EwcCiwhjm3LxZ+DhxPrqwIWEl7SXS6lgAYZEHW4jTN6XZBeCbo7TJlpN6WJmTwNPK8we2T+K51jIM+Q4XcGjwOmSrjazWV2piJpnVb3UzJ4oVNfMplPiB3XMItDDwjxKHcI9EU57SZvdGDNrNLPZcXGj4nQll8TP/0tZfz1J/5S0QNIcSdfHNCYFkbSqpOskfSxpicKkbmfFhzeShgMNhP/RBdF1NK3A/ppcYTGL8vux6E9xW4v7zNQ/VNJ4SYui++9eSRtm7XNadFMfH5MLLgW+E8sukvSypPmSPpX0REzTkdl2OPCXuDolocPAWL6cK0xhcrvnJH0paV68rptl1Rkn6VlJe8XjL5L0uqRDsuptKuk+hUnlFkv6MJ5j2uS4zgpKasPiOCsQdcB1wAhJaULibyOk5TiUEGhyIq1MYBRbIg8SMntfSXBzPQxcBVwaqz1ISAcDwR22C8FlnPYcDo3ffx233SXuE0knE5JivknIU3USsBXwVCb3U4I9CFNaXESYoG1SlK8fz/cgQvTmbIL3YeuE/hkj/b2EDjlzQ0naL26zAPg+8OOo07OSslOrb0xI8HhVPM864N6s/q4Ho44/BvYluDaX4M+l8qerc9n44kvahfBwNGAToB/wOTAqltXEspE56v8haz/nE1oamxY41gFx2+FZ8j8THn798x23wD5Hhr9c0/rAuO2Psur1Jcz9MSpLPojQIjkzIZtG6NdZt5VjV0ddJwPX5LqmObbJvp4vEnJL1WTptAy4KiEbF2WDE7K14zX/eVzvH/f/3a6+r3wp/uJvBk5ZYmZzCS2JY7JdMTm4J2v9LsJb8Y4FttmNMBnUHVny24CehDf7UrELIaX77THarSa6hz4C3o66JRlvZjOzdxJdUU9KmkOYfnkZsCmhT6hNxMy42wF3W/OkWJjZ+8B/gN2zNpliZlMS9WYTWkwZV94c4D3gMkknShrcVp2cFRc3LE45czUhkuniVupld/Bn1gvNjNcPmGvLd4LPTJSXirXj52MEY5Bctmb5WQyXc13FEOYxBLfVCcDOwNcIWW57t0OnNQiZb3O5yWay/PXIjjCD0NLrDbHZFlL1v0hwBb4j6T1JP26Hbs4KhneSOWWLmS2Q9GtCy+W3BaquQ0gZnlyHwhMYzQX6SeqZZVzWTZSXisxkUsNpqXeGL7LWc6UoP4zQSjnUzJZlhHH8z+ft0OmzeJx1c5StSzuuh5m9R5z7A9iGMC7uBknTzOyhdujorCB4i8Upd24gGIhLCtQ5Imt9GMHN9XyBbZ4i/D++lyU/itDP8Vzb1MzJkviZHaH2X4Lx2MTMXsyxTE6x7z6EPo3kgMw9aXZFtaZDC8xsIfAS8D1J1Yl9fgX4OqFfpV1YYCIhAAFCQIBTxniLxSlrzGyJpIuBmwpU21/SbwnjX3YkDAi8NdkHkIOHCDnw/iBpLULLYX/gR8CvzawYiVdnEVonwyRNAhYC75vZHEn/D7g+HvshQmf++oS+jHFmlt33k83DwJnALZL+QuhbuYDlW2lvxs9TJY0muNsm5XABErd/EHhA0g2EIIOLom5XtuG8kfRVQtTY3YSIvWpCC60eKDgeyFnx8RaL0x34CyFaKR8/JDxY7yOkKfoTcEqhHZpZI2E8yGjCKP8H4/pPCVFlHSYe40eE/ovHCNNdHxjL/gh8l9DR/ldCf8lIwsvgxBT7foSQjukbwAPA8YRZYadm1Xs17vdAgiGdAKyXZ58PE67B6oSAiD8QpjPe1cw+TnnaGWYCHxKu5/2EWUvXAw4ws5fauC9nBSPvDJKO4ziO0x68xeI4juMUFTcsjuM4TlFxw+I4juMUFTcsjuM4TlFxw+I4juMUFTcsjuM4TlFxw+I4juMUFTcsjuM4TlH5/0QqcUc2CSbCAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy is 0.91\n"
     ]
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(np.arange(nbIter),losslist)\n",
    "plt.xlabel('Nb of iterations')\n",
    "plt.ylabel('Cross entropy loss value (in bits)')\n",
    "plt.title('Loss, along gradient descent iterations')\n",
    "plt.show()\n",
    "\n",
    "nbSample=test_labels.size\n",
    "nbClass=10\n",
    "features, labels = selectData(nbSample,test_dataset,test_labels)\n",
    "\n",
    "prediction = softmax(linearModel(w,b,features)).reshape(nbSample,nbClass).T\n",
    "#print('prediction '+str(prediction))\n",
    "#print('max pred '+str(prediction.max(axis=0)))\n",
    "prediction=(prediction==prediction.max(axis=0)).T\n",
    "#print('1 hot pred '+str(prediction))\n",
    "#print('labels '+str(labels))\n",
    "prediction=(prediction * labels)\n",
    "#print('labels*pred '+str(prediction))\n",
    "\n",
    "accuracy = np.count_nonzero(prediction)/prediction.shape[0]\n",
    "print('Accuracy is '+str(accuracy))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Simplified version\n",
    "\n",
    "Following the chain rule, and given the fact that probability labels are one-hot encoded, we should be able to give a more simple version.\n",
    "\n",
    "In this case, the Cross-entropy function jacobian for a single sample whose label is k looks like that:\n",
    "\\begin{align*}\n",
    "  \\begin{pmatrix} 0 & \\dots & 0 & -\\frac{1}{s_k} & 0 & \\dots & 0\\end{pmatrix}\n",
    "\\end{align*}\n",
    "\n",
    "Then, there is only one useful row (the $k^{th}$ row) in the softmax jacobian:\n",
    "\n",
    "\\begin{align*}\n",
    "  J_{soft}=\\begin{pmatrix}\n",
    "   0 & 0 & \\dots & \\dots & \\dots & 0\\\\\n",
    "   \\vdots & \\vdots & \\vdots & \\vdots & \\vdots & \\vdots & \\\\\n",
    "   -s_k s_0 & -s_k s_1 & \\dots & s_k(1-s_k) & \\dots & -s_k s_{n-1}\\\\\n",
    "   \\vdots & \\vdots & \\vdots & \\vdots & \\vdots & \\vdots & \\\\\n",
    "   0 & 0 & \\dots & \\dots & \\dots & 0\\\\\n",
    "\\end{pmatrix}\n",
    "\\end{align*}\n",
    "\n",
    "Moreover, we recall that there was only one useful "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gradientWbv1(weight, bias, features, labels):\n",
    "  \" return gradient for W and b as a tuple\"\n",
    "  nSample=labels.shape[0]\n",
    "  nClass=labels.shape[1]\n",
    "  nFeatures=weight.shape[1]\n",
    "  gradW=np.zeros_like(weight)\n",
    "  gradB=np.zeros_like(bias)\n",
    "  for sample in range(nSample):\n",
    "    # computing local values for reuse\n",
    "    linLocal=linearModel(weight, bias, features[sample,:])\n",
    "    softLocal=softmax(linLocal)\n",
    "    # Right class idx\n",
    "    rc = labels[sample,:].argmax()\n",
    "    src = softLocal[rc]\n",
    "    # Now jcross is a single value\n",
    "    Jcross_sol = -1/src\n",
    "    # Now Jsoft is a single value as well\n",
    "    Jsoft_l = src*(1-src)\n",
    "    # When looking carefully, we have Jcross_sol*Jsoft_l=src-1\n",
    "    prod = src-1\n",
    "    # Now we can directly update gradient\n",
    "    gradW[rc,:]+=features[sample,:]*prod\n",
    "    gradB[rc]+=prod\n",
    "  return gradW, gradB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gradientWbv2(weight, bias, features, labels):\n",
    "  \" Yeah, this one is a bit trickier to understand without pen and paper\"\n",
    "  softLocal = softmax(linearModel(weight, bias, features))\n",
    "  rc = labels.argmax(axis=1)\n",
    "  prod = softLocal[np.arange(labels.shape[0]),rc]-1\n",
    "  gradW = np.dot(labels.T,(features.T*prod).T)\n",
    "  gradB = np.dot(labels.T,prod)\n",
    "  return gradW, gradB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "outputs": [],
   "source": [
    "stepLength = 5e-5\n",
    "nbIter=80\n",
    "\n",
    "#Select input data\n",
    "nbSample=800\n",
    "features, labels = selectData(nbSample,train_dataset,train_labels)\n",
    "\n",
    "#Sizes\n",
    "nbClass=labels.shape[1]\n",
    "nbFeatures=features.shape[1]\n",
    "\n",
    "#Perform gradient descent\n",
    "w=np.random.rand(nbClass,nbFeatures)\n",
    "b=np.random.rand(nbClass)\n",
    "losslist = []\n",
    "for i in range(nbIter):\n",
    "  gW,gB=gradientWbv1(w, b, features, labels)\n",
    "  w-=stepLength*gW\n",
    "  b-=stepLength*gB\n",
    "  l=loss(w, b, features, labels)\n",
    "  #print('loss is '+str(l))\n",
    "  losslist.append(l)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZYAAAEiCAYAAADao/T1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nOzdeZhcVbX38e+v58zz0BkggQSQMQwyyRAVEBEE1IugIghX4Iog6r1XkCsggoICTqivqAiKCKigiAwyJQgS5gBhDkmAJJ157qTn9f6xd3WfVLo7p9NVXT2sz/PUU3X2GWpVdXWtOnvvs7fMDOeccy5XigodgHPOud7FE4tzzrmc8sTinHMupzyxOOecyylPLM4553LKE4tzzrmc8sTinHMupzyxOOecy6mSjmwsaSwwDugHrADmm1ldPgJzzjnXM2lrV95L2g/4T+AjwHZZq+uAZ4A/AreY2fp8BOmcc67naDOxxIRyDXAY8DLwGPACsBzYBAwHJgMHAB+Ou30fuNbMavIbtnPOue6qvTaWmcBsYDcz28vMzjOzG83s72b2kJndYWZXm9kngDHAWcAngP/pgrgLStLpkkzSlELHkmuSpsfXNr3QsXRXkmZImpFYzvt7Fj9zZ3TyGJdJ6hWDA0oaGl/PPim3v0nSgsTypLj/DnkLcusxTY8xFGWVT4qfp9MLFFqntZdYdjSzC8zsta0dxMxqzOx2M9sXuCF34TnXIzwPHBTv8+V0oFOJpZcZClwKpEoswHeAExPLk+L+BUsswPQYQ/b3cBXh8/SPrg4oV9psvDezJdtyQDNbuu3hOJc/kkqBBsvxkN5mtg6Ylctjutwys7fz/RySignNCw2dOY6Z1dLDP0+puhtLGilpu6yysyX9VNKx+Qmt55P0OUkvSqqRtELS7yVVZm3zGUkvSNogaZ2klyWdnVj/fkkPSlopaZOkeZJ+vo3xfFvS8/F5Vkh6RNKBKfaTpK9KekNSnaQqSddLGpy1nUm6QtL5kuZLWi9ppqTdsrYrjttVSdoY49gl7n9ZiniOiO9ZjaS5kv6zjaoOk/QlSd+XtBioBYZKGiXpl5LejM//nqRbJY1v5blOlvS6pFpJr0g6sZVtWq0Kk/QJSbPic6yR9KdW/o8WSLolPs9rkqolPSvpkMQ2M4DDgQ/E57FkVVwb79Hekv4V36NFkr4FqJXtSiRdlHiNiyVdK6kia5vvSHo78Vl+PBlj3O6L8fO1SdLq+Lc/OLG+v6Sr42ejLt5frERVUOK9/Hj8jK2It1skDY3bTALmx11+lXhPTm/n/Wj+fMS/06Nx1YOJ/acntj9Lm//v/kbS8KxjmqQrJV0oaT6hM9Mekiok/VDSHIX/6yWS/i5pl8S+lxHOVgDqMzFkXl9rr0fpvk+2+nmK2+Xse6VVZrbVG3A38PPE8reAJmAl0Ah8Os1xesuNUC1hwJR2tjkrbnMbcAyhZ90y4E1gYNzmkPg+/gg4AjgKOB/4Rlw/EFgF3A8cRzh1Ph24YRvj/jVwKvBB4NgYWx2wR2Kb6THu6Ymy78ay6wm9A78KbAD+BRQltjNgAfAA8HHgU4QvgLlASWK7K+Lrvho4Ergwvi8GXLaV17ArIUH8CzgBOInQueRdYEFiu0nxeIuAv8bXezyhq/zOwI+BTxI6p5xM6N24AKhIHOOIGOffgY/F9/5dQlXFjK28Z+fEshvj3//TwGvx/RiU2G4B8E58/k/FOF8A1gBDE6/5eeBF4MB427Wd92gksDo+36fj+/QE8B5gWdveBlQDl8TXe1587r8ktrk4/r2/QkhwxwHfBj6e2Oaa+Hp/Hdd/jFD9dHJcXxL/ZiuBCwgdfi4GaggdfrLfy/nATwn/E+cROgzdHLcpJ1RrGeGzmXlPRrXzntyU+XwAg4Evxf3PS+w/OK6/CqgHro3P/wXC5+gpoDjr874ovq5PAkcT2puHxPfh5Ph+nQg8GP8mY+O+E+I2BnwgE0PWZ/f0jnyfdODzlNPvlVbf75RfSIuBTyaWFwFXxMc/AZ7KVUA94cZWEgtQDCwFHs0qPyTud35c/m9gVTvPs1/cfs88vIZiwj/7G8CPE+XTSXxJEnr/1QI3Ze3/ubhd8svFgLeA0kTZp2L5wXF5GOFL6udZx/sa6RLLrYSeif0TZZWEL6gFibLMP+fzxN6PW3kvJsbtT0yUPwG8yubJ88C43Yx23rOBwFrgxqznmUxI5BckyhYQvnCGtfJ3/0yibAbweMq/7ZXxeSYmygYQrj2zRNmh8Xk+n7X/Z2P5tLh8D3BnO883hfAD87p2tjk1HvOwrPKLY6yjs97Lm7O2uz7+jTM9WTN/3/9M+Z7clPX5yDzPEVnbTYqv5ZKs8g/E7U/I+rwvBvql+Hz1B9YDX02UXxaPUdJKDM2JhZTfJ2k/T+TxeyVzS3vl/fD4wpC0OzAWuDmu+yvhF6BrsTMwGvhDstDMHif8mjg8Fj0DDIunrsdmTvUT3iL80vhlPA2e2JmgFKqQHpW0Emgg/Crbifb/fgcCZcAtWeW3xWMcnlX+oJnVJ5ZfjveZKqA9CF9yf8ra78+pXkSI514z25gpMLMq4N9tbP9Xi/9NSZL+K1YrbCC8jnfjqp3j+mLg/cCfzawp8VyzCP+87TmI8Kv4D7EaqURSCeGM4XXCWVLSk2a2OrGc/Z511EHALDN7LxF3NeHMK+lowpf6n7Pi/Gdcn4nzGeCYWO1ziKSyrOMcQahWb6/jztGEz/6/W3muUsLfNSm74fplwpnKmHaeIxeOJLyW7L/dU4TEkP23u9/MNmUfRNJJkp6StIbw+aom/ODYlu/KtN8nGVv7POX0e6U1aRPLSsKpG8CHgMVm9lZcLu3AcfqKTF1sVSvrlmTWm9lM4D8Iv5bvApZLekjSnnH9WkK11WLg58C7sd72kx0NSKFb5r2Es4UzCf/I7ydUr1S0s2urr8VCA+XKxPqMVVnLtfE+8xyZOuFlWdul7fRR2cq+7e2/xd9A0nmE9/MhQhf5/Wn5YsvEOZLw2W7tuFuLdXS8f4iQvJO3PYARWdtv9p5ZaLxNxtJRlW3EmF02mvCjoTorxsz7m4nzu4T2gI8Tq7Mk/VbSyKztFrYT02hge7Z8P57OOkbG1j5H+ZL5281ly1gHsWWcrX2+jgNuJ1RFfoZwrd/7CWfa2xJ/qu+ThHY/T7n8XmlL2iFdHgIuix+krxPOUjJ2IWRN1yLzhx3byrqxwHOZBTP7M+EX40DC6fnVwP2SJphZk5nNBj4ZfzXtB1wE3CFpLzOb04GYPkn45fSJ5BmFpGGEXy9pXssrif1KCP9k2V8AW5P55xidPB7pf4lW0fLPn9TW/lucrRDqvh82s69nCiRNztpmBeHLpLXjjqH9z/zKeH86m7/GjHyPUFFF23EnrSRULx3axnEWA8TPy9XA1QrDOh0LXEeo3vk04b0CGE+oWm3NSkK7yUltrF/QRnlXy/ztjiJUKbW1PqOtz9dcMzs9U6DQIzE7AaSV+vskrRx+r7Qq7ZnG/xJO478HvE1ouMv4LPB4ZwPpZd4g/Do8OVkYe8hsT6gv34yZbTCze4BfEn5xjsha3xCrYb5F+Lu9r4Mx9SfUHTf/I0j6EFuvbplFqC45Oav804QfJjM6GMfLhF/I/5FVnr3cXjzHSOqfKYg9Yz7QgRj6E5JG0heSC2bWSGwAzeq1dAChDrw9/yYkjylm9mwrt7a+fNtTS+h4kMaTwIHJKg5JAwgNtUn3E37FDmkjzsXZBzazJWb2a8KPzd1j8UOETg5ntRPT/YQz8w1tPNeKdvZtTeZXeNr3JO3+DxJey3ZtxDk/xbH7E37EJZ1KaCtJE0O2Dn+fpJWD75VWpTpjsXBtypFtrD6C0GOjLzpaUvb1PmvN7EFJlxDqMG8htE+MJzSqvkXoKYSkywm/Ih8l/DqcQOgVNtvMlit05T6LcIY4n9A2cT7hS+vJzBPGboo3J38hteJ+Qm+cmyT9ltC28i1CR4w2mdkqSdcCF0mqJlSnvY/Qs+txOngRl5mtlvQj4JuS1hO+lPYhVM9B+KduzxWEDgEPSLqGUO/+LcI/3tb2zbgf+IakbxKqYj4Uj5ntUkIbwF8l/RIYRfhR1e41Xma2TtL/AD+TNAq4j9CYP55QHz7DzG5NGWvGq8CXJH2a8ONufTsJ6oeEXk//jN1aawkjYmz2f2pmMyT9kXDGfB3hvWgiJM5jCL0T35T0N0KV6fOEX/F7E9pMfhmP87akHwJfkzSI0Iu0kVDF+LqZ3U5oH/gC8HD8PL1IqIbbkVDFdkKy3SyFpYSzh5MlvUT4sTLfzLLPKNryJuHL/wxJq+J79EZ8LVcD10vamTACSQ0hKR4J/NrMHt3Kse8HTojvyT2EM4JMb7ukV+P91yXdBzSa2bPZBzOzxjTfJ2ml/V7plDQt/MAjwC5trNsJeCRfvQu6442WXmGt3eYktvsc4R+olvBP8HugMrH+Y4SuuVVxm/eA3wDj4vqdCXW18wkf7uWEL/YDEscYEJ/3qhRxnxePtYnwa/wIwq+dGYltprNl11kRuhi/QTh7qQJ+RuyemdjOiL0FE2WT2LLrZDHhn2JJjGUGcHDc7ispXseRhOGGaoF5wNmENqoXWnneLXoNEX4h/iK+n+sJ//yTaaVXGnBKfN21hGqtE9O8Z7H8GMKPhnXARlq+BHZNbLOAMIBrdoybxUKo8rg3xrtZr7Q23qN9CO0hNYQfD98iJEXL2q6I0I34xbjt2vj4+4QzGQjV37MIn+FN8f24jETvv7jdOcBL8b1aFd+ngxLrK+J+rye2eSaWlWS9l9m9tU6P5ZMSZScQvpzrsz9jrbwfN5HoFRbLzo6fn4bsvx/hDGMWIWFtILSXXA9MaO/znnhPryD8WNxISE57x7/1TVn/Bz8jtGk1Zf42tPI/k+b7JO3niRTfK529bXV0YwBJTYQ+1k+3sm5f4Gkzyz7Nc11A0lGE3j47mll7jafdmqRPEXqKHWZm/+rgvgMJja3/MLMzt7a9cy6/OjIfS1sZaEdCRneFcTihGqzHJJXYTvExQhfOGmBfwkWSs0jRXifpp4R2jMWE+YG+Qrg+5sd5Ctk51wFtJhZJX6ClQdOAG2KdeFI/QgPew/kJz22NmV1c6Bi2wQbC9QDnEq73WAbcAVxkaU6hQ5XK1YT2qTpC28ARZvZSfsJ1znVEe2csTYQGOAh17MnljJWEuuqrcx+a663M7BVCXfq27v/F3EXjnMu1tG0sjwL/ZWav5z+kriXpBUJvn7mFjsU553qQKcByM9s7e0WqxNKbSVo4ZMiQ8dOmTSt0KM4512PMnj2btWvXLjKzCdnr2mtj+Tyhl83K+LhdZva7TsZZKHOnTZs2fsaMGYWOwznneozp06czc+bMVmt62mtjuYkwftLK+Lg9BvTUxOKccy6H2kssk2kZ1yl7HCXnnHOuVe1NTfxOa4+dc8659nTkAkkkTSGM/zOeOKOadcFc0s4553qOVIlFYf7rn7PlCJ2Nkm4GzrWWMf/bO86NhCG3l5nZ7rHsdlomvxkKrDGzaQrzWr9GyzDcs8zsnLjPvoR2n36EMW6+YmamMCf17YSxdhYAJ9nmE94455zLs7TD5l9DGB7/UkLf5UHx/jJCsvlByuPcRBgVtZmZfdrMppnZNOAvwJ2J1W9n1mWSSvQL4IvA1HjLHPNCwjwbUwmjAVyYMi7nnHM5kjaxnAx828y+a2bzzKw63l8JXE6YJW2rzOwx2pgYSpIIkwD9sb1jxLk3BpvZrDj8x+8Io5wCHE/LlMk3J8qdc851kbSJpZyWKUSzPUWYV6GzDgWWWsuUxwCTJb0gaaakzCx349l8CtSFsQxgjIX5zyEMyd7mrISSzpL0LLBvVVVrM362r7HJOO3Gpzn++sdpaurbF5k651xS2sTyEGGqztYcRZivpbNOYfOzlSrCLG57A18DbpU0OO3B4tlMm9/4ZnaDme0HPFdZWdnWZm0qLhLPvbOaFxeuZV1N9mSEzjnXd7V35f0OicXrgN/H6U3/RJi9bQyh6uoYwgQ02yzOu/wJwvDpAMTOALXx8XOS3iZMKraIMNNixgRaZkFcKqnSzKpildmyzsS1NSMGlrGhtoGV1XUM7Z+LkzbnnOv52usVNpfNf/EL+C/CLHHJMggzpHVmoq8jCFOYNldxxSldV1mYlnMHQiP9PAtT5a6TdCChGu7zwE/jbncDpwFXxfu/dSKmrRoxoIx3Vm5k5YY6dhyVz2dyzrmeo73E8oV21m2TOL/2dGCkpIXApWb2G0LngOxG+8OAyyXVE4bsP8fMMg3/X6Klu/F98QYhodwh6UzgHcIZVd4MH1AOwKrqrfa0ds65PqO9K+9vbmvdtjKzU9ooP72Vsr8Quh+3tv2zhAnGsstXAh/uXJTpjRwYqr9WbKjrqqd0zrluL23jvWvF8AEhsayq9sTinHMZnlg6YcTAUBW2coNXhTnnXIYnlk7IVIWt9DMW55xr5omlEzJVYSu9jcU555p5YumEEbFX2ErvFeacc808sXTCiIHeeO+cc9lSz8cSL1I8CdgOqMhabWZ2Zi4D6wmG9W9JLE1NRlGRtrKHc871fmnnYzkBuINwhrOMONRKQp8chbGspIjBFSWsq2lgzab65jYX55zry9KesXwHmAF81syW5y+cnmfkwHLW1TSwckOtJxbnnCN9G8sOwDWeVLY0wrscO+fcZtImlteBEfkMpKfyLsfOObe5tInlf4FvZg2l72i5+t4HonTOuSBtG8tlhDOW1yS9xZbTC5uZHZ7LwHqKEQN8IErnnEtKm1gagTfyGUhPNcIHonTOuc2kSixmNj3PcfRYzQNRelWYc84BfuV9p3lVmHPOba69Oe8PA543sw3xcbvM7LGcRtZDtDTee2JxzjlovypsBnAg8HR83NbV9YrrOjPnfY/V0t3Yq8Kccw7aTywfBF5NPHatGNa/FAnWbKqnobGJkmKvXXTO9W3tzXk/s7XHbnMlxUUM7VfK6o31rN5Yz6hB5YUOyTnnCqpLf15LulHSMklzEmWXSVokaXa8HZNYd5GkuZLekPSRRPnRsWyupAsT5ZMlPRXLb5fUJYN3eTuLc861aDOxSLpb0t5pDySpQtLXJJ3TzmY3AUe3Uv5DM5sWb/fG4+0KnAzsFvf5uaRiScXAz4CPArsCp8RtAa6Ox5oCrAa6ZCj/Ed7O4pxzzdo7Y1kAzIpnAOdL2kfSZlVnksZJOkHSb4Aqwhf5820dMPYcy75qvy3HA7eZWa2ZzQfmAvvH21wzm2dmdcBtwPGSBHwI+HPc/2bghJTP1Sk+EKVzzrVoM7GY2fmEM4KnCUO6PAPUSFolqUrSJuA94E7CWcUFwJ5m9vQ2xPFlSS/FqrJhsWx8PH7GwljWVvkIYI2ZNWSV513zFMV+xuKcc+1feW9mbwPnSfo6cBBwADCOMIPkSsKox4+Z2TudiOEXhPleLN5fC5zRieOlIuks4Cxg56qqqk4dq7nLsZ+xOOdc6iFd6oCZ8ZZTZrY081jSr4B74uIiYGJi0wmxjDbKVwJDJZXEs5bk9q097w3ADZJmVFZWdmoAzZFeFeacc80KftGFpMrE4olApsfY3cDJksolTQamEqrlngGmxh5gZYQG/rvNzIBHgU/F/U8D/tYVr6F5vDCvCnPOudSjG+eEpD8C04GRkhYClwLTJU0jVIUtAM4GMLNXJN1BuEizATjXzBrjcb4MPEC42v9GM3slPsU3gNskXQG8APymK17XcB/h2DnnmnVpYjGzU1opbvPL38yuBK5spfxe4N5WyucReo11qeaqMB+I0jnnCl8V1hsMz/QK8zMW55zzxJILQ/uVUiRYu6me+samQofjnHMF5YklB4qK1NzOstrPWpxzfVzqxCJpb0l3SlohqUHSPrH8u5JaG6alT8lcJOkTfjnn+rpUiUXSIcCTwC7ArVn7NQHtjQ/WJ7QM6+Jdjp1zfVvaM5arCN17dwO+lrXueWCfXAbVE3mXY+ecC9J2N94H+ISZmaTsmSRXAKNyG1bPM3KgV4U55xykP2OpAfq3sa4SWJubcHquljMWrwpzzvVtaRPL48AFcS6UjMyZy5nAIzmNqgca4RdJOucckL4q7FvAE8CLhPlODDhN0nXAvsD78xNezzHCRzh2zjkg5RmLmb0IHAYsBS4GBHw5rj7czN7IT3g9hw9E6ZxzQeqxwszseeDDkiqA4YRJtTbmLbIexs9YnHMu6PAglGZWAyzOQyw92ujBFQAsXVeDmRFmSnbOub4nVWKRdMlWNjEz+04O4umxBpaXMKRfKWs31bOquq65asw55/qatGcsl7WzLtM7rE8nFoBxQ/uxdlM9i9Zs8sTinOuz0jbeF2XfgJHA6YQZH6fkMcYeY/zQfgAsWr2pwJE451zhbPNEX2a2CvidpBHAz4BjchZVDzVhWEwsazyxOOf6rlwMm5/pitznjRsaGvA9sTjn+rJcJJZjgeU5OE6PN35oGPXGq8Kcc31Z2l5hN7ZSXAbsDuwBXJrLoHqq8V4V5pxzqdtYPkRL76+MGuAd4EfAzbkMqqfKVIUt9sTinOvDUiUWM5uUiyeLZz7HAsvMbPdY9gPgOKAOeBv4gpmtkTQJeA3IDBczy8zOifvsC9wE9APuBb4Sh/QfDtwOTAIWACeZ2epcxJ7GyAHllJUUsXpjPRvrGuhfts19I5xzrsfq6jnvbwKypzF+ENjdzPYE3gQuSqx728ymxVtylspfAF8EpsZb5pgXAg+b2VTg4bjcZYqK5F2OnXN9Xps/qSV1qKeXmT2WZpt4JpIs+2dicRbwqfaOIakSGGxms+Ly74ATgPuA44HpcdObgRnAN9LEnyvjhlYwf0U1i9ZsYuqYQV351M451y20V1czgy3bVVqjuF3x1jZM4QxCVVbGZEkvAOuA/zOzfwHjgYWJbRbGMoAxZlYVHy8BxrQZtHQWcBawc1VVVVubdVjzGYu3szjn+qj2EssHuywKQNLFQAPwh1hUBWxnZitjm8pfJe2W9nhtTKOcXH8DcIOkGZWVlYd3JvYk73LsnOvr2kwsZjazq4KQdDqhUf/DZmbx+WuB2vj4OUlvAzsBi4AJid0nxDKApZIqzawqVpkt66KX0My7HDvn+rqubrzfgqSjgf8FPp6c30XSqMxUyJJ2IDTSz4tVXeskHagwNv3ngb/F3e4GTouPT0uUdxnvcuyc6+tS94eN1VD/CewMVGStNjP7cIpj/JHQuD5S0kLChZUXAeXAg3EOk0y34sOAyyXVA03AOXF8MoAv0dLd+L54A7gKuEPSmYRrbE5K+/pyZYJXhTnn+ri0V94fAMwkXBsyFXgJGAZsR2g8n5vmOGZ2SivFv2lj278Af2lj3bOEq/6zy1cCW01w+TR2SAUSLFlXQ31jE6XFBT8pdM65LpX2W++7wJ3AboReYGfGiyaPIPQGuyIv0fVAZSVFjB5UTpOF2SSdc66vSZtY9gRuoaX7cTGAmT1CSCrfy31oPZdfJOmc68vSJpYyoNrMmoBVQGVi3Ru0Ui3Vl40fFttZvAHfOdcHpU0sc2m5CPEl4AxJRZKKgC8QLkZ0kfcMc871ZWl7hf2d0JvrVkJ7yz8IV8M3AgOB8/MRXE81wa++d871YWlHN74s8fghSQcCnwT6A/dnjffV52UuklzobSzOuT5om8Z1N7MXgBdyHEuvMS6esXhVmHOuL0rVxiLpLkknSCrNd0C9QXIgyjhCjXPO9RlpG+93JlzHskTSz2NVmGvDoIpSBleUUFPfxKrqukKH45xzXSpVYjGzXYH3A78HPgE8IektSZfEcbxclnHegO+c66NSjzdiZs+Z2QWEbsfHAc8QJtF6S9K/8hRfjzVhmLezOOf6pg4PZGVmjWZ2r5l9BjgRWAwcnPPIerhMO4v3DHPO9TUd7hUWq75OBT4L7EiYkOvaHMfV4/m8LM65virt6MbDgE8TEsqBwEbgLsLw9Q+bd33awjgfL8w510elPWNZQhh48hHCBFp3JiflclvabngYL+zdVf42Oef6lrSJ5WLgVjNbnM9gepMdRg0EYN6KahqbjOIiFTgi55zrGmm7G1/jSaVjBpaXMHZwBXUNTSxc7Wctzrm+w6c3zKMpo8NZy9xlGwociXPOdR1PLHm046gBALy93BOLc67v8MSSR5kzlreXVRc4Euec6zqeWPJox9iAP9fPWJxzfUiXJxZJN0paJmlOomy4pAfj+GMPxutmUPATSXMlvSRpn8Q+p8Xt35J0WqJ8X0kvx31+Iqlg3bGSbSx+qY9zrq9IO2z+8ZK+kFjeXtKTktZL+rOkgR14zpuAo7PKLiRcaDkVeDguA3wUmBpvZwG/iM8/HLgUOADYH7g0k4ziNl9M7Jf9XF1m1KByBpWXsHZTPSt9lGPnXB+R9ozl/4BRieXrgAnADcBhwGVpn9DMHgNWZRUfD9wcH98MnJAo/50Fs4ChkiqBjwAPmtkqM1sNPAgcHdcNNrNZcTSA3yWO1eUksWNzO4tXhznn+oa0iWVH4CUASf2AY4CvmdnXgW8SBqPsjDFmVhUfLwHGxMfjgfcS2y2MZe2VL2ylfAuSzpL0LLBvVVVVa5vkhLezOOf6mrSJpQLIDHp1MOGK/cw8928A43IVUDzTyHuDhJndYGb7Ac9VVlbm7Xl2HB27HHvPMOdcH5E2sSwADomPjweeM7O1cXk0sLa1nTpgaazGIt4vi+WLgImJ7SbEsvbKJ7RSXjBT/IzFOdfHpE0svwQui1VHXwJ+k1h3EPBqJ+O4mzC4JfH+b4nyz8feYQcCa2OV2QPAUZKGxUb7o4AH4rp1kg6MvcE+nzhWQXgbi3Our0k1CKWZ/VjSCsKQ+T8xs98lVg8Cfpv2CSX9EZgOjJS0kNC76yrgDklnAu8AJ8XN7yW058wlDNX/hRjPKknfIcxiCXC5mWU6BHyJ0POsH3BfvBXMdsP7U1osFq3ZxMa6BvqXdXgKHOec61FSf8uZ2R+AP7RSfnZHntDMTmlj1Ydb2daAc9s4zo3Aja2UPwvs3pGY8qm0uIjtRwxg7rINzFteze7jhxQ6JOecy6u017HsJGn/xHI/Sd+T9HdJX85feL1Dpp3FxwxzzvUFadtYrgc+laPoL94AACAASURBVFi+Evg6oTfYDyW1elbhgpaeYZ5YnHO9X9rEshfwBICkIkKj+DfMbF/gCsJV8a4NzUO7+BmLc64PSJtYhgAr4+O9gWHAn+PyDGCH3IbVu2QukvRrWZxzfUHaxLIUmBIfHwW8bWaZK98HAg25Dqw3ySSW+SuqaWhsKnA0zjmXX2kTy93A9yRdQ2hb+VNi3R7AvFwH1psMKC+hckgFdY1NLFy9aes7OOdcD5Y2sVwI3EMY/PFu4LuJdR+nZXgX14bmMcO8Ad8518ulvUCymjAUfWvrDs5pRL3UlNEDeXzuCt5evoEjmsfYdM653qdDl4HHeVAOAoYThr5/MnHFu2tHZmiXN5asL3AkzjmXX6kTi6QrCO0r5YniWknXmNm3ch5ZL7NnvOJ+9sI1BY7EOefyK+2V9xcQ5l25Bfgg8L54fwvwTUnn5y3CXuJ9lYMpKyli3vJq1m6sL3Q4zjmXN2kb788BfmxmXzSzmWb2Rrz/IvATwsCPrh1lJUXsNm4wAC/6WYtzrhdLm1gmAf9oY90/4nq3FdMmDgVg9nueWJxzvVfaxLKStkcM3o2Wq/JdOzyxOOf6grSJ5S7gO5JOlVQCIKlE0inA5cBf8hVgb7L3xGFASCxhRgDnnOt90iaWi4DZwM3AJklLgU2E+VleJDTsu62YOLwfwweUsaq6jvdW+RX4zrneKe0FkuslHQZ8DDiUlutYZgL3mf/8TkUSe00YwqNvLGf2wjVsN6J/oUNyzrmc68gMkkYY1uWe/IXT+02bOCwklnfX8PG9xhU6HOecy7m0VWEuR6Ztl2nAX13gSJxzLj/aPGOR1ASkreIyM+vQ8DB91bQJIbHMWbyOuoYmyko8tzvnepf2ksHlpE8snSJpZ+D2RNEOwCXAUMLgl8tj+TfN7N64z0XAmUAjcL6ZPRDLjwZ+DBQDvzazq7riNaQ1pH8pO4wcwLwV1by+ZB17xkTjnHO9RZuJxcwu66ogzOwNYBqApGJgEaGL8xeAH5rZNcntJe0KnEy4hmYc8JCkneLqnwFHAguBZyTdbWavdskLSWnaxKHMW1HN7PfWeGJxzvU63bEe5sOEGSrfaWeb44HbzKzWzOYDc4H9422umc0zszrgtrhtt9LczvKuXyjpnOt9umNiORn4Y2L5y5JeknSjpGGxbDzwXmKbhbGsrfItSDpL0rPAvlVVVTkLPg2/At8515t1q8QiqYwwI2Vm6uNfADsSqsmqgGtz9VxmdoOZ7Qc8V1lZmavDprLL2DjS8Qof6dg51/t0q8QCfBR43syWApjZUjNrNLMm4FeEqi4IbTATE/tNiGVtlXcrZSVF7B5HOvb5WZxzvU13SyynkKgGk5Q8lTgRmBMf3w2cLKlc0mRgKvA08AwwVdLkePZzcty229l7u1Cr9+wCn4DTOde7pJ3oa6etb9U5kgYQenPdmSj+vqSXJb1EmFjsqwBm9gpwB/AqcD9wbjyzaQC+DDwAvAbcEbftdg7ecQQAT8xdUeBInHMut9Je1Pi6pEeB/wfcFb/Ac8rMqoERWWWntrP9lcCVrZTfC9yb6/hy7YAdRlBcJF5cuJZ1NfUMrigtdEjOOZcTaavCzgD6ES5iXCjpu7EKym2jgeUl7D1xKI1NxlPzvDrMOdd7pEosZnaTmR1M6J31F8JUxG9Jul/S8ZK6W1tNj/CBKSMBrw5zzvUuHUoIZvaSmZ1LuNr9bGAMoU3kXUmXSRqThxh7rUOmhsTyuCcW51wvsq1nGpOAPeN9HaG31teAuZJOzElkfcC0iUMZUFbM3GUbWLK2ptDhOOdcTqROLJLKJH1W0mPAy8BxwFXARDM7Gtie0EPrurxE2guVFhdxwA7eO8w517uk7W58LeFCw5uB9YSr43c0s6vNbAWAma0mjCq8fZ5i7ZW8ncU519uk7W58KnAj8P/ioI9teZ0wIrFL6ZApLe0sZoakAkfknHOdkzaxTIijBbcrnr3c3LmQ+padxgxk5MBylq2vZe6yDUwdM6jQITnnXKek7W5cByBpiqTPSPqfeL9jfsPr/SRxyJTQzuK9w5xzvUHaNpYKSTcShkm5Bbg63r8u6deSyvMYY6/n7SzOud4kba+wa4DPApcCU4BB8f4yQvvLD/IRXF+RSSyz5q2ivrGpwNE451znpE0sJwPfNrPvxtkZq+P9lcDlwGfyF2LvN25oP3YYNYANtQ0++ZdzrsdLm1jKCcPSt+YpoCw34fRd03caDcB9Ly8pcCTOOdc5aRPLQ8BRbaw7CngkN+H0XR/bM0w9c+/LVTQ1WYGjcc65bZe2u/F1wO/jnCl/ApYSxgk7CTgG+JykHTIbm9m8XAfa2+2z3VDGD+3HojWbeO7d1bx/0vBCh+Scc9skbWKZGe//CzgnUa6s9RnFnQmqL5LEx/as5IbH5nHPi4s9sTjneqy0icWvpu8Cx8bEcu+cJVxy3G4UF/lV+M65nidVYjEzv5q+C+wxfgjbDe/Pu6s28tT8lRy848hCh+Sccx3WoWHzFewm6dB47z+pc0gSx8ZG/HteqipwNM45t206Mmz+fwJVwEvAjHi/WNKZ+Qmtbzp2z3EA3D9nCQ1+saRzrgdKO6TLZ4EbCPOwnEHoCXZGXL5B0il5i7CPeV/lIHYYNYBV1XX8++2VhQ7HOec6LO0Zy/8CfzCzI83sZjN7IN4fBdwKfCMXwUhaIOllSbMlPRvLhkt6UNJb8X5YLJekn0iaK+klSfskjnNa3P4tSaflIrauEqrDwlnLPS8tLnA0zjnXcWkTy86EQSdbc0tcnysfNLNpZrZfXL4QeNjMpgIPx2WAjwJT4+0s4BcQEhFhTLMDgP2BSzPJqKc4Lraz3D9nCXUNXh3mnOtZ0iaW9cCENtZNiOvz5Xha5ni5GTghUf47C2YBQyVVAh8BHjSzVXFWyweBo/MYX85NHTOIXcYOYl1NA/981Yd4cc71LGkTy33AdyUdmiyUdBBwRVyfCwb8U9Jzks6KZWPMLNNFagnhin+A8cB7iX0XxrK2yjcj6axY3bZvVVX364H12QO2A+D3T75T4Eicc65jOtLGshaYIeldSU9Jegd4HFgX1+fCIWa2D6Ga61xJhyVXmpkRkk+nmdkNsbrtucrKylwcMqdO2Hs8A8qKeWr+Kt5cms8TQuecy620M0guAaYBXwGeJCSTWcB5wN5mtjQXwZjZoni/DLiL0EayNFZxEe+Xxc0XARMTu0+IZW2V9yiDKko5cZ9wouVnLc65nmSriUVSqaTjCVVS15vZp2PvsE+b2c/NbGMuApE0QNKgzGPCqMlzgLuBTM+u04C/xcd3A5+PvcMOBNbGKrMHgKMkDYuN9kfFsh7n1AMnAXDn8wvZUNtQ2GCccy6lrSYWM6sH7gAm5TmWMcDjkl4kzP3yDzO7H7gKOFLSW8ARcRngXmAeMBf4FfClGO8q4DvAM/F2eSzrcXYeO4j9Jw+nuq6Ru17ocSddzrk+Ku0glPOA0fkMJA61v1cr5SuBD7dSbsC5bRzrRuDGXMdYCKceuD1Pz1/F759cwOcO2A4fRcc5192lbbz/PnCxpFH5DMZt6SO7jWXUoHLeXLqBp+f3yBMv51wfk/aM5UPAcGC+pFmEMcOSvbPMzHrUFe49RVlJEae8fyI/eWQuv5/1DgfsMKLQITnnXLvSJpZDgXpgObBjvCX5XLp5dMoB2/GzGW9z/5wlLF6ziXFD+xU6JOeca1Pa7saTzGxyO7cdtn4Ut60qh/TjmD0qaWgyfvvE/EKH45xz7Uo7uvFhkga2sW5A9oWMLvfOPizk7lufepe1m+oLHI1zzrUtbeP9o8CubazbJa53ebT7+CF8YMoIqusa+cNTfsGkc677SptY2uvjWg405iAWtxVnHxaatn77xAJq6v0td851T2023kuaBCTbTvZrpTqsH2HCr3dzHpnbwqFTR/K+ysG8VrWOv76wiJP3367QITnn3Bba6xV2GmFek8zAjz9l8zMXi8sNtHGhosstSZxz+A585bbZ3PCveZy030SKivyCSedc99JeYrmJMLe9gEcIyePVrG1qgTd76pApPdExe1Ty/fvfYN7yah56bSlH7Ta20CE559xm2kwsZvYO8A6ApA8Cz5uZj99eYKXFRZx5yGQuv+dVfvboXI543xg/a3HOdStpr2OZ6Uml+/j0+ycyalA5Ly5cy1+eX1jocJxzbjNpr2Mpk3SppNclbZTUmHXzMd270IDyEr55zC4AXHXf66zd6Ne1OOe6j7RDuvyA0MZyH3AnoW3FFdAJ08bzx6fe4+kFq7j2wTe4/PjdCx2Sc84B6RPLp4BLzezKfAbj0pPE5Sfsxsd+8ji3zHqHk/abyO7jhxQ6LOecS32B5EDClMSuG9ll7GBOO2gSTQaX/G0OTU0+FqhzrvDSJpa/Az4eWDd0wZFTGTWonOffXcNtz7xX6HCccy51YvkpcIqkSyTtJ2mH7Fs+g3RtG1xRysXHvA+Ab/1tDnd6LzHnXIGlbWPJVINdRrgavzXFnY7GbZPjp41j7rINXP/oXL52x4tsqG3g8wdNKnRYzrk+Km1iOQOfzKvbksR/f2RnBlWU8L37XueSv73C+poGvjR9RyS/eNI517VSJRYzuymfQUiaCPwOGENIYDeY2Y8lXQZ8kTBzJcA3zezeuM9FwJmEkZXPN7MHYvnRwI8JZ1C/NrOr8hl7d3L24TsyuF8p37zrZX7wwBu8vWwDlxy3K0P7lxU6NOdcH5L2jAUASUWEeVlGAM+aWXWO4mgAvm5mz0saBDwn6cG47odmdk1WHLsCJwO7AeOAhyTtFFf/DDgSWAg8I+luM8se46zXOmX/7RhYXsJ//+lF7nxhEY+9tYIrTtiNo3evLHRozrk+Im3jPZLOBZYALxEGpdw5lv9V0vmdCcLMqszs+fh4PfAaML6dXY4HbjOzWjObD8wF9o+3uWY2z8zqgNvitn3KcXuN476vHMr+k4azYkMt59zyPF/6w3MsWrOp0KE55/qAtEO6fJFQvfRX4CQ2Hz7/X8AncxVQnAdmb+CpWPRlSS9JulHSsFg2Hkj2rV0Yy9oqb+15zpL0LLBvVVVVrsLvNnYYNZDbzjqQy4/fjf5lxdz78hI++IMZfPvvr7B8vQ+c4JzLn7RnLF8DrjWzs4C7sta9Tjx76aw4kdhfgAvMbB3wC2BHYBpQBVybi+cBMLMbzGw/4LnKyt5ZTVRUJD5/0CQeuOAwPr7XOOoam/jtEws47PuP8v37X2fZ+ppCh+ic64XSJpbJwANtrKsGhnY2EEmlhKTyBzO7E8DMlppZo5k1Ab8iVHUBLAImJnafEMvaKu/TJg7vz09O2Zv7vnIoR7xvDJvqG/n5jLf5wFWP8LXbZzNn0dpCh+ic60XSJpYVwKQ21u1MJ7+8FfrE/gZ4zcyuS5QnTyVOBObEx3cDJ0sqlzQZmAo8DTwDTJU0WVIZoYH/7s7E1pu8r3Iwvz5tP+780sF8ZLcxNDQZd76wiGN/+jif/MW/ue3pd1lX4yMlO+c6J22vsHuASyTNIE7+BZikkcBXCW0vnfEB4FTgZUmzY9k3CVf7TyN0QV4AnA1gZq9IuoMwo2UDcK6ZNQJI+jLh7KoYuNHMXulkbL3OPtsN45en7sd7qzZy878XcPsz7/HcO6t57p3VXHr3Kxy56xhO3Hs8H5gykopSv+7VOdcxMtv6dY8xgTxBqGZ6ijBu2L+BXYBlwMFm1iPrUyTNOPzwww+fMWNGoUMpmOraBu6bs4S/PLeQJ+etbC7vX1bM9J1H8ZHdxjJ9p9EM6V9awCidc93J9OnTmTlz5kwzm569Lu0Fkisk7QdcAHwEeDvuez3hOpN1OYzXdbEB5SV8at8JfGrfCSxcvZG/vrCI++Ys4ZXF67j35SXc+/ISigR7ThjKYVNHcsjUUUybOJSyktS91Z1zfUiqM5bezM9Y2vbeqo08+OpS/vnqEp57ZzX1jS2flYrSIqZNHMr+k0ew/6Th7DVxCIMq/IzGub6i02csrm+aOLw/ZxwymTMOmUx1bQNPzV/JY2+u4Im5K3hr2QZmzVvFrHmrAJBg6uiBTJs4lGkTh7HH+CHsNHYg5SXeRuNcX+OJxaUyoLyED+0yhg/tMgaAlRtqeWbBap5ZsIpnFqzi1cXreHPpBt5cuoE7ng1D95cUiZ3GDGK3cYPZpXIwu4wdxC5jBzFiYHkhX4pzLs88sbhtMmJgOUfvPpajdx8LQE19I68sXsfs99bw4ntreGXxWuatqObVqnW8WrV5E9zIgWVMGT2QKaMHMnX0IHYcNZAdRg1g7OAKiop8NGbnejpPLC4nKkqL2Xf7Yey7/bDmsuraBl6LieX1Jet5vWodbyxZz4oNdazY0FKN1nKMIiaNGMDkkQPYbkR/th8+gO1H9Ge74f0ZO6SC0mLvLOBcT+CJxeXNgPIS9ps0nP0mDW8ua2oyqtbV8NbS9cxdtoG3lm5g/opq5q2oZsWG2pCAlqzf4lhFgsoh/ZgwrB/jh/Vj/NB+jIu3yiEVjBlcweCKEp9/xrluYJsTSxy6/n3Ak2a2OHchud6sqEiMHxoSw/SdR2+2bl1NPfOXV7NgZTXvrtzIO6s28s7Kahau3sSSdTUsWrMpjNA8v/Vj9y8rZuzgCkYNKmfM4ApGDypn9OByRg0qZ9TAUD5yYBlD+5dR7FVuzuVNqsQi6XqgxMzOicufAG4nXN2+TtKRZvZM/sJ0fcHgilL2mjiUvSZuOfRcbUMjVWtqeG/1RqrWhCSzeM0mFq/dxJK1NVStrWFjXSPz4tlPe4oEwweEJDN8QLiNGFDG8AHlDB9QyrABZQzvHxLQsAGlDOtf5iMQONcBac9YPgp8O7H8beIwL4QRhy8Fjs1taM61KC8pZtLIAUwaOaDV9WbG+toGlq6tYdn6Wpauq2HpulqWr69l+YZaVmTuN9SyZmM9K+Lj9M9fxLD+ZQztX8qQfqXN98nb4MytopQh/UoYVFHKoIoS+pUWexWd61PSJpZKwlhdSJpAmLnxTDN7WdJPCANIOlcwkhhcEb7Up44Z1O629Y1NrK6uY/mGWlZV17Gquo6VG8L96o3htqq6jtXV9azeWMeajfXUNjSxZF0NS9Z1fKqBkiIxqKKEgRUlDCoPyWZQRQkDy0sYUB7KB5bFx7FsQHkxA8tL6F8WyvqXFzOgrISK0iJPUq7bS5tYNgID4+PDgXXAs3F5A9D+f7Jz3UhpcRGjB1cwenBFqu3NjI11jazeWMfaTfWs3VjPmk31rNlYz7qa+lC2qZ51m+pZV9MQ7uPj9TUhKa3eWM/qjfVA52bxlKB/aTH9ykroX1acuJXQr6yYfqVhOfm4orRluV9pMRVlxVSUtJRVlBZRUVocb0WUFXvycp2TNrE8D5wr6V3gXODBOEcKhLlaet8UjM5FkuJZRAkThm19+2y1DY2sr2lgQ00D62saWF9bz/qaBqprG9iQuTUvN1Jd20B1XVjeWNcYHzeysa6Bmvomqusaqa5rzP0LjSRiwimmoiQknfJMAippSUTlJS0Jqby0iPK4rrykZV15SdHmj0s33yZTVl5SREmRPKH1EmkTy8XA/cCLwBrgnMS6EwhzoTjnWlFeUkz5wGJG5mDEgcYmo7qugU11jWysC8lmY11j83JNfbjfVN/IproGNtVnypualzfVN1FTH7bdFLetqW+itr6RmoZG6hstHjt/yas1RQrvVVlzwgkJqKw487iIsuRycRFlJfGWfByXy0uKKI3lmfuy4mSZNisvKVbz+pK4rrS4yHsQboO0oxs/I2k7wjD5b2WNZnwD8FY+gnPOba64qKUtKV8aGpuoaWhJPjUxEdU2tDxuvm9opLa+aYv72oaWRFXX0BQSV9y/riE8ro3PEZabaGiymPi6NqFtTZGgpLiI0iJRWlJESVFISiXForQoJKGS5nuFbWNZabEojmUlRXG7IjVvWxz3Ky7KLKs5mRUrs29yfVHzdiVFoiixX3FzeRHFRTRvW6RwjCJtvm1m/4HlJTnv9Zj6OhYzqwaeS5ZJGmFm/8hpRM65giopLmJgcREDy7v2+umGxibqGptaElNDS9LJXk7e1zU0UteYedxEbeJxfWPm3qiNy81lTUZ9ZpvGJhoajbrE+sxykxGOB9DFZ3Fd4TvH78apB03K6THTXsfyRWComf0gLu8B3AdUSnoBONbMluQ0Mudcn1JSXERJcRH9ywodSQszo7HJaGiymJSMhsaQlBoyy00hCdU3hrOuhkxZ5nEsb2xKbNNkNCYfx20bm8Kxm5LlTU00NkFj4piNZjQ2ZtYZjRbXNxpNFvZtajLq43JjU8vryDxujOX9ynL/AyLtEc8jVHllXEdoa7kaOB+4HDgrt6E551xhKVYjlRTjF8l2QNrEsj3wOoCkIYQuxyeY2b2SVgLfy1N8zjnnepi0w8UWAZnuxYcABsyIy+8Bo1vZxznnXB+UNrG8BXwsPj4Z+LeZbYzL44BVre7lnHOuz0mbWK4BLpC0AvgM8NPEug8CL+U6sM6QdLSkNyTNlXRhoeNxzrm+JO11LLfGq+4PAJ4xs8cSq5cCd+cjuG0hqRj4GXAksBB4RtLdZvZqYSNzzrm+oSPXsTwOPN5K+aU5jajz9gfmmtk8AEm3AccDnlicc64LpE4skvoDZxB6hA0ntKs8CvzWzDo3sl5ujSd0KMhYSDjT2oykswhdpHeuqvKhzpxzLlfSXiA5ltALbCfgHWAJsAPwSeA8SdPNbGm+gswHM7sBuEHSwqVLlw6cPn16oUNyzrkeY/bs2QBTWluX9ozl+8Aw4FAzeyJTKOlg4C+ECyVP71SUubMImJhYnhDL2rJ87dq1zJw5c+42PNe0eD97G/bNt+4aW3eNC7pvbN01Lui+sXXXuKD7xtbRuKYAy1tbITPb6t6SlgPfMLMbW1l3JnCVmY1KGUxeSSoB3gQ+TEgozwCfMbNX8vBcMwDMbHquj91Z3TW27hoXdN/Yumtc0H1j665xQfeNLZdxpT1jGQgsbmPdQlomASs4M2uQ9GXgAaAYuDEfScU551zr0iaWN4BTCXOyZPsccbiX7sLM7gXuLXQczjnXF6WtCvsc8DvgEeBWwoyRYwlX4R8BnGpmt+YxTueccz1EqsQCzd1zL2fzccGWApeY2a/yEJtzzrkeKHViAZBUBOxMy3Usb5hZU/t7Oeec60u2mlgklQGzgAvN7J9dEpVzzrkea6uDUJpZHTAZaMh/OD1LdxnsUtKNkpZJmpMoGy7pQUlvxfthBYptoqRHJb0q6RVJX+kO8UmqkPS0pBdjXN+O5ZMlPRX/prfHH1ZdTlKxpBck3dPN4log6WVJsyU9G8u6y2dtqKQ/S3pd0muSDip0bJJ2ju9V5rZO0gWFjisR31fj53+OpD/G/4tOf9bSjm78IHBURw/emyUGu/wosCtwiqRdCxTOTcDRWWUXAg+b2VTg4bhcCA3A181sV+BA4Nz4PhU6vlrgQ2a2F+HCsKMlHUi42PeHZjYFWA2c2cVxZXwFeC2x3F3iAvigmU0zs/3icqH/lhk/Bu43s12AvQjvX0FjM7M34ns1DdgX2AjcVei4ACSNJ8wAvJ+Z7U64PONkcvFZM7Ot3oBDCUO5XEOY6GtHwpAuzbc0x+lNN+Ag4IHE8kXARQWMZxIwJ7H8BlAZH1cS2sO6w/v2N8LI090mPqA/8DxhTLkVQElrf+MujGcC4cvmQ8A9gLpDXPG5FwAjs8oK/rcEhgDzidX73Sm2RCxHAU90l7hoGVdxOOHSk3uAj+Tis5b2jGUmYZiUr8XHbxIm/0re+prWBrscX6BYWjPGzDKjay4BxhQyGABJk4C9gafoBvHF6qbZwDLCWfnbwBozy1T7Fupv+iPgf2mZtXVEN4kLwuyx/5T0XOwpCt3gb0morl8O/DZWIf5a0oBuElvGycAf4+OCx2VmiwgnC+8SLiFZCzxHDj5raS+Q/EJHD+y6DzMzSem7/+WBpIGEceUuMLN1kprXFSo+M2sEpkkaSqie2KWrY8gm6VhgmZk9J2l6oeNpxSFmtkjSaOBBSZtdHF3Az1oJsA9wnpk9JenHZFUvFfL/ILZTfJxQs7GZQsUV23WOJyTlNcCf2LJKfZuknejr5lw8WS/T0cEuu9pSSZVmViWpkvCrvCAklRKSyh/M7M7uFp+ZrZH0KOG0f6ikkviLrRB/0w8AH5d0DFABDCa0HRQ6LqD5Vy5mtkzSXYT5j7rD33IhsNDMnorLfyYklu4QG4S22OetZRT47hDXEcB8M1sOIOlOwuev05+1NqvCJBVJOk7S7u1ss4ek4zr6pL3EM8DU2IOijHCa221m0iTEclp8fBqhbaPLKZya/AZ4zcyuS6wqaHySRsUzFST1I7T7vEaYY+hThYrLzC4yswlmNonwmXrEzD5b6LgAJA2QNCjzmNBmMIdu8FkzsyXAe5J2jkUfJkzuV/DYolNoqQaD7hHXu8CBkvrH/9PMe9b5z1o7DTufBzYAk9vZZlLc5pSubnjqDjfgGEJ709vAxQWM44+EOtJ6wi+3Mwn18g8T2r8eAoYXKLZDCPXyLxGG454d37eCxgfsCbwQ45pDGEECQmeUp4G5hKqB8gL+XacD93SXuGIML8bbK5nPfKH/lon4pgHPxr/pXwlTfRQ8NmAAsBIYkigreFwxjm8TxnqcA/weKM/FZ63NCyQl/ZPQU+G8Vjdo2e7HwM5mlpO6Oeeccz1be73C9gHSXGn/ELDfVrdyzjnXJ7SXWAYRLo7ZmtVxW+ecc67dxLIC2D7FMbaL2zrnnHPtJpbHaem10J7T47bOOedcu4nlR8CHJf2wtUHIJJVK+hFh2Ikf5itA55xzPUubicXMngS+ThikbKGkWyRdGW+3ELq1fpkwwOCsrgnX9WWSTpdkktZkjwYrqSSuu6yV7afkOa6D4miw1fH5prWx3WXJK6zjaLyXSdonn/G1R9K0GMPwpJuoVwAABghJREFUVtZt9n46l1a7V96b2Y8kPQ98AzgR6BdXbQJmAFeZ2b/yGqFzWxpC+EwWbKqCLL8h/E8cRxi99s02tvs1cH9ieShwKeFH2vP5DLAd02IMtxAm70s6iBCbcx2y1SFdzOwx4DGF2SNHxuKVFsZZcq4Q/gmcJ+mH1jJERkGoZVbVK83skfa2NbOF5PmLOl5BXWphHqVO8ZoIt63Sjm6MmTWZ2bJ486TiCumKeP9/KbcfJ+mvkjZIWinpZ3EYl3ZJGizpekmLJdUqTOr21fjljaTTgUbC/9G3YtXRgnaO11wVFkd6nh9X/Srua/GYme0/IWmWpI2x+u9PkrbLOuaCWE19RhwQsg74WFz3bUnPK0wutULSIwpzzmT2PR34bVx8KxHDpLh+i6owhcntnpS0SdLa+L7unLXNDEmPSzoiPv9GhYmkTszabidJdylMUlcj6d34GtMOjuu6qdSJxblupAq4HjhLUpou8bcQhqf4BKGjyReBX7S3QzwT+QdhZO9rCdVc9wPXAVfGzf5BGLIGQnXYQYQq47Sv4RPx8ffivgfFYyLpHMLAna8Sxm06G9gdmJkZryvhg4QpLb5NGJ32pVg+Pr7e4wm9N5cRah/2SMSfSdL/kYghM5z7ZiQdHff5/+2dT4hXVRTHP4cZRMJFidlgi4rUXQVBohUGLkVdFJZQhNS0mNqIEIQiTBIYSEYLhylDs//Woo3TGAgZ9JcmMBeBGaMQ6BjMJKQLGem0+N5fvrnz5s38fr3B3y/OBy6Pe9999573+P3efeeee8+5DDwB9CWZvjYFjSpyN3KeuS/d5wXg08zeNZRk7ENxQF5CAdjivdTp3Aj/NJEitZLQy9GB5Sg40SXgYDrXnc71l9QfzNrZiTSNlRV9bUjXbs3K30YvvyUz9VvRZr/+cv/m70zX9mb1FqHYGAez8ruQRrKtUHYO2XV6Zum7K8l6Gnij7JmWXJM/zxHk26o7k2kS2FcoO5HKVhTKlqZnviPll6T2N93o31Wk+lN8GQQdibtPIE3i6XwqpoRPsvzH6Kt4VcU1a1GgrQ+z8veBBejLfr5Yg9zlf5BWu3Wn6aHfkcPAtVn9713efaeQpqK+NLNxFCJ6EliJbEJNkbwZ3w8c8etBoHD3s8A3wCPZJWfc/Uyh3h9IY2pM5Y0Do8CrZvacma1oVqagfYmBJehkXkcrmXbPUi838DfyVZHxFgMTPt0IPlY4P18sTcfjaDAopnuQZ9wi06au0hLmz9G01bPAauAB5Jl4YQsy3YJCJJdNk40x/XnkK8xAmt5CSGqbQhWMoKnAX81s1Mz6WpAtaDPCSBZ0LO5+2cz2IM1lb0XV25Cb92IeqgMYTQCLzWxBNrj0FM7PF+PpuJWpcjf4K8uXuSh/DGkpj7r7ZKMw7f+51IJMf6Z+ekrO9dDC83D3UaRxGnAf2hc3YGbn3H24BRmDNiE0lqDTGUADxCsVdR7P8lvQNNcPJXUbfIX+H5uz8ieRneO75sQs5Wo65ivUvkWDx3J3HylJp+fQ9k3IplHckLmO61NRs8kwBXe/guKhbzazrkKbdwAPIrtKS7g4iRYggBYEBB1MaCxBR+PuV81sN/BWRbX1ZrYX7X9ZhTYEvlu0AZQwjHzgDZrZrUhzWA/0AnvcvQ7HqxeRdrLFzE4BV1Co2HEzexHYn/oeRsb825Et44S757afnGPANuAdMzuEbCu7mK6l/ZKOL5jZYTTddqpkCpB0/RBw1MwG0CKDl5NsrzVx35jZvWjV2BG0Yq8LaWjXgMr9QEH7ExpL8H/gEFqtNBNPoRfrZ8hN0QHg+aoG3f1vtB/kMNrlP5Ty29Gqsv9M6qMX2S+Oo3DXG9O5N4FNyND+HrKX9KOPwZNzaPsL5I7pIeAo8AyKCvtbVu/n1O5GNJD+CCyboc1j6BncjBZEDKJwzg+7+/k53naDMRQadzsK0/tR6neDu//UZFtBmzFjBMkgCIIgaIXQWIIgCIJaiYElCIIgqJUYWIIgCIJaiYElCIIgqJUYWIIgCIJaiYElCIIgqJUYWIIgCIJaiYElCIIgqJV/ABLvY9tyE6yFAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy is 0.89\n"
     ]
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(np.arange(nbIter),losslist)\n",
    "plt.xlabel('Nb of iterations')\n",
    "plt.ylabel('Cross entropy loss value (in bits)')\n",
    "plt.title('Loss, along gradient descent iterations')\n",
    "plt.show()\n",
    "nbSample=test_labels.size\n",
    "nbClass=10\n",
    "features, labels = selectData(nbSample,test_dataset,test_labels)\n",
    "prediction = softmax(linearModel(w,b,features)).reshape(nbSample,nbClass).T\n",
    "prediction = prediction.argmax(axis=0)\n",
    "labels = labels.argmax(axis=1)\n",
    "prediction = np.logical_and(prediction, labels)\n",
    "accuracy = np.count_nonzero(prediction)/prediction.size\n",
    "print('Accuracy is '+str(accuracy))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "default_view": {},
   "name": "1_notmnist.ipynb",
   "provenance": [],
   "version": "0.3.2",
   "views": {}
  },
  "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.12.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}