Exercise_04_1_solution.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "rF2trPuyzm9C"
   },
   "source": [
    "# Exercise 4.1 - Solution\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "id": "ipcsUFDUzm9C"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "MCJe_ITJzm9G"
   },
   "source": [
    "**Disclaimer**\n",
    "\n",
    "The book mistakently refers to the page for Exercise 4.2 when introducting Exercise 4.1, etc. Of course, these numbers should match: Book Exercise 4.1 is discussed under Exercise 4.1 \n",
    "\n",
    "**Simple Network**\n",
    "\n",
    "We continue with the dataset first encountered in Exercise 3.2. Please refer to the discussion there for an introduction to the data and the learning objective.\n",
    "\n",
    "Here, we manually implement a simple network architecture"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "NopU99AT9G7s",
    "outputId": "d7e8848e-b9c0-4eb4-8f18-5acda9d8c343"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/bin/sh: wget: command not found\r\n"
     ]
    }
   ],
   "source": [
    "# The code snippet below is responsible for downloading the dataset\n",
    "# - for example when running via Google Colab.\n",
    "#\n",
    "# You can also directly download the file using the link if you work\n",
    "# with a local setup (in that case, ignore the !wget)\n",
    "\n",
    "!wget https://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-white.csv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "5ONqeI5Uzm9H",
    "outputId": "d31ba8d4-cf0a-4f25-8a93-9091c0dd041a"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "data: (4898, 12)\n",
      "First example:\n",
      "Features: [6.300e+00 2.700e-01 1.800e-01 7.700e+00 4.800e-02 4.500e+01 1.860e+02\n",
      " 9.962e-01 3.230e+00 4.700e-01 9.000e+00]\n",
      "Quality: 5.0\n"
     ]
    }
   ],
   "source": [
    "# Before working with the data, \n",
    "# we download and prepare all features\n",
    "\n",
    "# load all examples from the file\n",
    "data = np.genfromtxt('winequality-white.csv',delimiter=\";\",skip_header=1)\n",
    "\n",
    "print(\"data:\", data.shape)\n",
    "\n",
    "# Prepare for proper training\n",
    "np.random.shuffle(data) # randomly sort examples\n",
    "\n",
    "# take the first 3000 examples for training\n",
    "# (remember array slicing from last week)\n",
    "X_train = data[:3000,:11] # all features except last column\n",
    "y_train = data[:3000,11]  # quality column\n",
    "\n",
    "# and the remaining examples for testing\n",
    "X_test = data[3000:,:11] # all features except last column\n",
    "y_test = data[3000:,11] # quality column\n",
    "\n",
    "print(\"First example:\")\n",
    "print(\"Features:\", X_train[0])\n",
    "print(\"Quality:\", y_train[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "jiwnyNHpzm9L"
   },
   "source": [
    "# Solutions\n",
    "\n",
    "The goal is to implement the training of a neural network with one input layer, one hidden layer, and one output layer using gradient descent. We first (below) define the matrices and initialise with random values. We need W, b, W' and b'. The shapes will be:\n",
    "  * W: (number of hidden nodes, number of inputs) named `W`\n",
    "  * b: (number of hidden nodes) named `b`\n",
    "  * W': (number of hidden nodes) named `Wp`\n",
    "  * b': (one) named `bp`\n",
    "\n",
    "Your tasks are:     \n",
    "   * Implement a forward pass of the network as `dnn` (see below)\n",
    "   * Implement a function that uses one data point to update the weights using gradient descent. You can follow the `update_weights` skeleton below\n",
    "   * Now you can use the code below (training loop and evaluation) to train the network for multiple data points and even over several epochs. Try to find a set of hyperparameters (number of nodes in the hidden layer, learning rate, number of training epochs) that gives stable results. What is the best result (as measured by the loss on the training sample) you can get?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(50, 11)\n"
     ]
    }
   ],
   "source": [
    "# Initialise weights with suitable random distributions\n",
    "hidden_nodes = 50 # number of nodes in the hidden layer\n",
    "n_inputs = 11 # input features in the dataset\n",
    "\n",
    "# See section 4.3 of the book for more information on\n",
    "# how to initialise network parameters\n",
    "W = np.random.randn(hidden_nodes,11)*np.sqrt(2./n_inputs)\n",
    "b = np.random.randn(hidden_nodes)*np.sqrt(2./n_inputs)\n",
    "Wp = np.random.randn(hidden_nodes)*np.sqrt(2./hidden_nodes)\n",
    "bp = np.random.randn((1))\n",
    "\n",
    "print(W.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# You can use this implementation of the ReLu activation function\n",
    "def relu(x):\n",
    "    return np.maximum(x, 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def dnn(x,W,b,Wp,bp):\n",
    "    # SOLUTION\n",
    "    # sum_i W'_ki*Relu(sum_j W_ij*x_j + b_i) + b'_k \n",
    "    return np.dot(Wp, relu(np.dot(W,x) + b)) + bp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update_weights(x,y, W, b, Wp, bp):\n",
    "    \n",
    "    lr = 0.00005\n",
    "\n",
    "    # SOLUTION\n",
    "\n",
    "    # Calculate the network output\n",
    "    phi = dnn(x,W,b,Wp,bp)\n",
    "\n",
    "    # Use the formulas derived to calculate the gradient for each of W,b,Wp,bp\n",
    "    delta_bp = 2 * (phi - y)\n",
    "    delta_Wp = 2 * (phi - y) * relu(np.dot(W,x) + b)\n",
    "    delta_b  = 2 * (phi - y) * Wp * np.heaviside(np.dot(W,x) + b, 0.5)\n",
    "    delta_W  = 2 * (phi - y) * np.outer(Wp * np.heaviside(np.dot(W,x) + b, 0.5), x)\n",
    "                \n",
    "    # Update the weights/bias following the rule:  X_new = X_old - learning_rate * gradient    \n",
    "    bp -= lr * delta_bp\n",
    "    Wp -= lr * delta_Wp\n",
    "    b  -= lr * delta_b\n",
    "    W  -= lr * delta_W\n",
    " \n",
    "    \n",
    "    return -1 # no return value needed, you can modify the weights in-place"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Training loop and evaluation below"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoch: 0 Train Loss: 5.645403582072304 Test Loss: 5.789360719845873\n",
      "Epoch: 1 Train Loss: 1.1533433753334428 Test Loss: 1.1637637243130456\n",
      "Epoch: 2 Train Loss: 0.8064459925389651 Test Loss: 0.8256775797887149\n",
      "Epoch: 3 Train Loss: 0.7194986944212736 Test Loss: 0.7463420451468886\n",
      "Epoch: 4 Train Loss: 0.6823412559002283 Test Loss: 0.6977324105881756\n",
      "Epoch: 5 Train Loss: 0.6664846138051684 Test Loss: 0.6876245778967164\n",
      "Epoch: 6 Train Loss: 0.6587998893746104 Test Loss: 0.6813185867359125\n",
      "Epoch: 7 Train Loss: 0.6542220249941378 Test Loss: 0.6763715413742601\n",
      "Epoch: 8 Train Loss: 0.6515229736431492 Test Loss: 0.6726972058370031\n",
      "Epoch: 9 Train Loss: 0.6500603789067503 Test Loss: 0.6696083273661432\n",
      "Epoch: 10 Train Loss: 0.6493454499784376 Test Loss: 0.6681520553512201\n",
      "Epoch: 11 Train Loss: 0.6555439871036638 Test Loss: 0.6725116282689468\n",
      "Epoch: 12 Train Loss: 0.6550890293925202 Test Loss: 0.6699690448726485\n",
      "Epoch: 13 Train Loss: 0.6526835138767458 Test Loss: 0.6670276424175429\n",
      "Epoch: 14 Train Loss: 0.6524323996393876 Test Loss: 0.6657859446027024\n",
      "Epoch: 15 Train Loss: 0.6508383266412704 Test Loss: 0.6634289610832423\n",
      "Epoch: 16 Train Loss: 0.6475738523938723 Test Loss: 0.6588708333747422\n",
      "Epoch: 17 Train Loss: 0.6448067370447748 Test Loss: 0.656277910633608\n",
      "Epoch: 18 Train Loss: 0.6429866949447117 Test Loss: 0.6541361924738887\n",
      "Epoch: 19 Train Loss: 0.6452836868194535 Test Loss: 0.6603206884328453\n",
      "Epoch: 20 Train Loss: 0.6390966116784169 Test Loss: 0.6484080376897273\n",
      "Epoch: 21 Train Loss: 0.643608872593549 Test Loss: 0.6568905242756689\n",
      "Epoch: 22 Train Loss: 0.6365683454156014 Test Loss: 0.6459070379047909\n",
      "Epoch: 23 Train Loss: 0.6344200247151789 Test Loss: 0.6405648728904523\n",
      "Epoch: 24 Train Loss: 0.6307892457380325 Test Loss: 0.6502294953157811\n",
      "Epoch: 25 Train Loss: 0.6271675389745466 Test Loss: 0.6394346181103574\n",
      "Epoch: 26 Train Loss: 0.625209796186921 Test Loss: 0.6367168482951644\n",
      "Epoch: 27 Train Loss: 0.6254154882759166 Test Loss: 0.6365726656263813\n",
      "Epoch: 28 Train Loss: 0.6233292217982007 Test Loss: 0.6342570747926862\n",
      "Epoch: 29 Train Loss: 0.6194417930925457 Test Loss: 0.6365037775275856\n",
      "Epoch: 30 Train Loss: 0.6153231814680605 Test Loss: 0.6317002972628465\n",
      "Epoch: 31 Train Loss: 0.6126120664039123 Test Loss: 0.6287678460188807\n",
      "Epoch: 32 Train Loss: 0.6129166430079448 Test Loss: 0.6284321403254162\n",
      "Epoch: 33 Train Loss: 0.6104762499877989 Test Loss: 0.6257073921750971\n",
      "Epoch: 34 Train Loss: 0.6082635861912785 Test Loss: 0.6230077493420016\n",
      "Epoch: 35 Train Loss: 0.6080842723192383 Test Loss: 0.6227811633866328\n",
      "Epoch: 36 Train Loss: 0.6035597556499765 Test Loss: 0.6178078657920026\n",
      "Epoch: 37 Train Loss: 0.6016812959797005 Test Loss: 0.6154486388653907\n",
      "Epoch: 38 Train Loss: 0.6000701547553112 Test Loss: 0.6136773857670896\n",
      "Epoch: 39 Train Loss: 0.5996881694114135 Test Loss: 0.6129554259443293\n",
      "Epoch: 40 Train Loss: 0.59857725747447 Test Loss: 0.6122699599476488\n",
      "Epoch: 41 Train Loss: 0.5970114764108893 Test Loss: 0.6105162649916903\n",
      "Epoch: 42 Train Loss: 0.5952629839015641 Test Loss: 0.6094251196258232\n",
      "Epoch: 43 Train Loss: 0.5940349734884387 Test Loss: 0.608027638521175\n",
      "Epoch: 44 Train Loss: 0.5934948491666505 Test Loss: 0.6075563960754147\n",
      "Epoch: 45 Train Loss: 0.5953789625146336 Test Loss: 0.6078046455326416\n",
      "Epoch: 46 Train Loss: 0.5948276458373926 Test Loss: 0.607165242559313\n",
      "Epoch: 47 Train Loss: 0.5934935232783387 Test Loss: 0.6065483328684542\n",
      "Epoch: 48 Train Loss: 0.59473349370651 Test Loss: 0.6076145117172599\n",
      "Epoch: 49 Train Loss: 0.5978180426888854 Test Loss: 0.6139122064546286\n",
      "Best loss: 0.6065483328684542 Final loss: 0.6139122064546286\n",
      "Correlation coefficient: 0.47560172145940866\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAEGCAYAAABvtY4XAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAZYUlEQVR4nO3dcXCkd33f8fd39x7hR+LonkF2I9nHYceVCz7uZG/tC9d4MAYLCrjyYTdc4+mETHItzWRIPVEmV5iCM9BLZgND0k6Z2GYSOgRjMGdNSFuEk+AkpPbByrI5O1ilJtiHRHwitlzbt3Dr1bd/7EpIut3V7mqfXem3n9eM5nS/5/f8fr/V89uPnn32t4/M3RERkfCkuj0AERFJhgJeRCRQCngRkUAp4EVEAqWAFxEJ1I5uD2C117zmNb5nz55uD0NEZNuYnp7+obsPVtu2pQJ+z5495PP5bg9DRGTbMLOnam3TJRoRkUAp4EVEAqWAFxEJlAJeRCRQCngRkUAp4EVEApXoMkkz+wDwy4ABd7r7J5PsT0TCMDkzR25qlvnFAkOZmImxEYA1ZdddPsjXnlhYU2d8dLilthvZbzuypG4XbGZXAJ8HrgbOAl8B3u/u36m1Tzabda2DF+ltkzNzHD1+kkKxtFIWpQwMiqXaeRVHaY4d2ls3rKu13ch+W5mZTbt7ttq2JC/R/FPgIXc/4+4vA38J3JRgfyISgNzU7JoABigued1wBygUS+SmZptuu5H9tqskA/4x4Foze7WZ9QP/Arh4fSUzO2JmeTPLLywsJDgcEdkO5hcLie1ba/tm+tzKEgt4d/828DvA/ZQvzzwKvFyl3h3unnX37OBg1dspiEgPGcrEie1ba/tm+tzKEl1F4+6fdvcr3f1a4Fmg5vV3ERGAibER4ii9pixKGVHa6u4XR+mVN2ObabuR/barpFfRXODup81sN3AI+Jkk+xOR7W/5zc4kVtHUanu7vsG6kcRW0QCY2V8DrwaKwG3u/uf16msVjYhIc+qtokn0DN7dfzbJ9kVEpDZ9klVEJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCtSPJxs3sPwC/BDhwEnifu/8oyT5FtprJmTlyU7PMLxYYysRMjI0wPjrcUt23feIBvnP6pZbHEqXg5SUYysRcd/kgX3tiYaWvPa+O+d/ffRb3n9Q34IKdfTzzwtm67e7qj3CHxUKRtBkl9zVltRjw8wd289Hxvfz8nQ/yN08+u7ItbVDymrvWddkFA/zf0y/R7O5xlOJHxSUy/REv/qhIcal6vUwc8ZEb37BybJo5xp1k7i3+BDdq2GwY+DrwencvmNkXgP/p7n9Ua59sNuv5fD6R8Yh0w+TMHEePn6RQLK2UxVGaY4f2nhMAG9XdbLhvZRc28Etkq4lSRu6WfQANH+MkmNm0u2erbUv6Es0OIDazHUA/MJ9wfyJbSm5qds0TH6BQLJGbmm26bqjhDmy7cAcoLjm5qdmmjnGnJRbw7j4H/C7wNPAD4Hl3/+r6emZ2xMzyZpZfWFhIajgiXTG/WGi4vJm6sjXMLxa29HFLLODNbBfwL4HXAUPAgJndur6eu9/h7ll3zw4ODiY1HJGuGMrEDZc3U1e2hqFMvKWPW5KXaN4K/J27L7h7ETgOvCnB/kS2nImxEeIovaYsjtJMjI00XfeyCwaSG2iXXbizr9tDaFqUMibGRpo6xp2WZMA/DRwws34zM+B64NsJ9iey5YyPDnPs0F6GMzEGDGfimm++bVT3/tvevOmQj1KstH3rgd1r+jp46fmYra1vNBa+u/ojMnEEQLrSyOqyWgy49cBuTnzwbRy89Pw129JWfZ9GXHbBAK3sHkcpjPLYozrpmIkjcrfsY3x0uKlj3GmJraIBMLPbgZ8DXgZmgF9y9x/Xqq9VNCIizam3iibRdfDu/mHgw0n2ISIi1emTrCIigVLAi4gESgEvIhIoBbyISKAU8CIigVLAi4gESgEvIhIoBbyISKAU8CIigVLAi4gESgEvIhIoBbyISKAU8CIigVLAi4gESgEvIhIoBbyISKAU8CIigVLAi4gESgEvIhIoBbyISKAU8CIigVLAi4gESgEvIhIoBbyISKB2JNWwmY0A96wqugT4T+7+yaT6FGnU5MwcualZ5hcLDGViJsZGGB8dbmsf13zsfp554eyasuFVfU3OzHH7lx/nuTPFtvYrybjsggEWXjjLYqH+8RroS/Oxm/auzKd6cy3peWju3rbGanZilgbmgGvc/ala9bLZrOfz+cTHI71tcmaOo8dPUiiWVsriKM2xQ3vb9uSqFu6r+3rPVcPc881TFEvJP/+k89Ip4+O37AOoOdfqbWtmHprZtLtnq21L7Ax+neuBJ+uFu0in5KZm1zypAArFErmp2bYFfK1wX+7r7hOnKHXg5Eq6o7Tk5KZmAWrOtXrb2jUPOxXw7wXurrbBzI4ARwB2797doeFIL5tfLDRVngSFe/jqzadWtzUr8TdZzawPuBH4YrXt7n6Hu2fdPTs4OJj0cEQYysRNlSchbdaxvqQ7hjJx3bnWiXnYiVU07wAedvdnOtCXyIYmxkaIo/SasjhKMzE20rY+LtzZV3NbHKU5fM3FRGmFfKjSKWNibKTuXOvEPOzEJZrD1Lg8I9INy9c3k1y9cOKDb9twFU32tedrFc020uoqGqg/17btKhoz6wdOAZe4+/Mb1dcqGhGR5nRtFY27nwFenWQfIiJSnT7JKiISKAW8iEigFPAiIoFSwIuIBEoBLyISKAW8iEigFPAiIoFSwIuIBEoBLyISKAW8iEigFPAiIoFSwIuIBEoBLyISKAW8iEigFPAiIoFSwIuIBEoBLyISqIYD3sxekeRARESkvTYMeDO72sxOAt+p/H+fmf2XxEcmIiKb0sgZ/O8D7wL+AcDdHwWuS3JQIiKyeY0EfMrdn1pXVkpiMCIi0j47GqhzysyuBtzM0sCvAv8n2WGJiMhmNXIG/37gNmA38AxwoFImIiJb2IZn8O5+GnhvB8YiIiJttGHAm9mdgK8vd/cjDeybAe4Crqi08Yvu/mAL45TATM7MkZuaZX6xwFAmZmJshPHR4a718aHJk/zxiafxykxPGyx5lYkP7EgZLy9V2yIhe8WOFCmDQnEJgF39ER9+9xsanredmPPrNXIN/s9WfX8ecBNwqsH2fw/4irvfbGZ9QH+T45MATc7McfT4SQrF8nv1c4sFjh4/CdC2Cd9MHx+aPMlnH3p6TVmpTn4r3HvTj19eWvP/584Umbj3UWDjeduJOV/Nhtfg3f2eVV+fAQ4Br99oPzN7FXAt8OlKO2fdfXGzA5btLzc1uzLRlxWKJXJTs13p4+4TjZ6viKxVLHlD87YTc76aVm5V8DrgtQ3UuwRYAP7QzGbM7C4zG1hfycyOmFnezPILCwstDEe2m/nFQlPlSfdRcp2RS+sambedmPPVNPJJ1ufM7NnK1yJwP/AfG2h7B3Al8Cl3HwVeAn5zfSV3v8Pds+6eHRwcbHL4sh0NZeKmypPuI23Wtn6l9zQybzsx56upG/BmZsA+YLDytcvdL3H3LzTQ9veB77v7icr/76Uc+NLjJsZGiKP0mrI4SjMxNtKVPg5fc3Hb+pXeEqWtoXnbiTlfTd03Wd3dzew+d7+q2Ybd/e/N7JSZjbj7LHA98LetDlTCsfymUpIrCprp46PjewG0ikbq2swqmk7M+WrMN7j+aGafAu5094ebbtxsP+Vlkn3Ad4H3uftztepns1nP5/PNdiMi0rPMbNrds9W21TyDN7Md7v4y8M+BXzazJylfRzfKJ/cbXm5x90eAqh2LiEiy6l2i+Qbla+bjHRqLiIi0Ub2ANwB3f7JDYxERkTaqF/CDZnZbrY3u/okExiMiIm1SL+DTwCupnMmLiMj2Ui/gf+Duv9WxkYiISFvV+6CTztxFRLaxegF/fcdGISIibVcz4N392U4ORERE2quVu0mKiMg2oIAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQNX7m6ybZmbfA14ASsDL7p5Nsj/pvMmZOXJTs8wvFhjKxEyMjTA+OtxyvUb3W19+3eWD/OmjP2CxUARgV3/EO9/4U9zzjacpLv2kvQt39nH6xbO4t+1HIJtw8NLzeeTU87x0tlR1ewpYdfjoj1L850NvPGfuTM7McfuXH+e5M8WVsoG+NO7OmcoEiKMU50VpFs8Um5qD25l5gjO9EvBZd/9hI/Wz2azn8/nExiPtNTkzx9HjJykUf/LkjKM0xw7tXfPEabReo+2/56phvjQ9t6ZcekfK4BP/av/K3JmcmWPi3kcplprLskbm4HZgZtO1Tp51iUZalpuaPSdkC8USuanZluo12v7dJ04p3HvYkrNm7uSmZpsOd2hsDm53SQe8A181s2kzO1KtgpkdMbO8meUXFhYSHo600/xioaHyRus1ur2k6ys9b/Xc2GgeNdpOiJIO+IPufiXwDuBXzOza9RXc/Q53z7p7dnBwMOHhSDsNZeKGyhut1+j2tFkDo5OQrZ4bG82jRtsJUaIB7+7zlX9PA/cBVyfZn3TWxNgIcZReUxZHaSbGRlqq12j7h6+5+Jxy6R0pY83cmRgbIUo3/0u/kTm43SUW8GY2YGY7l78HbgAeS6o/6bzx0WGOHdrLcCbGgOFMXPVNq0brNdr+R8f3nlN+64HdZOJoZd9d/RG3HthNtG6GX7izD70A2DoOXno+A321f1mvD6j+KLXmDVYoz5PczfvY1R+tqTvQl6Z/1QSIoxS7+qOm5uB2l9gqGjO7hPJZO5SXY37O3T9Wbx+tohERaU69VTSJrYN39+8C+5JqX0RE6tMySRGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFA7ku7AzNJAHphz93cl3Z/0lsmZOXJTs8wtFkibUXJf2WYG8Y4UheISQ5mY6y4f5GtPLDC/WOAfxRFmsHimyFAmZmJsBIDc1OzK9mJpiZfOlpoajwG+Ya2ytEFpg8oDfWleOlvCDLzRhuuIoxTHDr2R8dFhJmfmuP3Lj/PcmWLVuga86dLzeeTU82t+DgN9aT52017GR4eBtccgZbDkP9nfgeF1P/vln/fy/u22PJ5O9LXVmbdj1tTrwOw2IAu8aqOAz2azns/nEx2PhGNyZo6jx09SKDYXwtVEaQOH4lKyz4etIAX86wO7ueebpyhu9BumhnTK+Pgt+wBaOgZxlObYob1tD95qcyKpvrYKM5t292y1bYleojGzi4B3Ancl2Y/0ptzUbFvCHaBY8p4Id4Al4O4TrYc7QGnJyU3NtnwMCsUSuanZlvuvpdp4kuprO0j6Es0ngd8AdtaqYGZHgCMAu3fvTng4EpL5xUK3h7Btldrwyn2zP/8kjl+tNnt1riR2Bm9m7wJOu/t0vXrufoe7Z909Ozg4mNRwJEBDmbjbQ9i20mabbmMoE2/qGCRx/Gq12atzJclLNAeBG83se8DngbeY2WcT7E96zMTYCHGUbktbUdqIUpsPve0gBRy+5uLy+w4tSqeMibGRlo9BHKVX3thup2rjSaqv7SCxgHf3o+5+kbvvAd4L/IW735pUf9J7xkeHOXZoL8OVs7P1Z6Vm0B+lMMorOW49sJvhTIwBmThiV3+0si138z5yt+xbs32gr/ngaiYyG8nX5TG04YQbKK+i+cTP7eej43vJ3byPXf1RzboGHLz0/HN+DgN9aT5+yz7GR4fPOQarf0cuf7v+Zz+ciRN703P1eJLuaztIfBUNgJm9Gfh1raIREWmveqtoEl8HD+DuDwAPdKIvEREp0ydZRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFA7kmrYzM4D/gp4RaWfe939w+3uZ3JmjtzULPOLBYYyMRNjI4yPDre7G9lCmjnm6+ted/kgX3ti4Zx9N9NmtbrNzsvN1M/0R7jD84Ximu8383yYnJnj9i8/znNnigBk4oiP3PiGltvSc7Q7zN2TadjMgAF3f9HMIuDrwAfc/aFa+2SzWc/n8w33MTkzx9HjJykUSytlcZTm2KG9mkCBauaYV6u7Xhylec9Vw3xpeq7lNtfXbXZetqN+vcfX7PNhcmaOiXsfpVhamw1Rysjdsq/ptvQcTZaZTbt7ttq2xC7ReNmLlf9Gla+2/jbJTc2eM8kLxRK5qdl2diNbSDPHvFrd9QrFEnefOLWpNtfXbXZetqN+La08H3JTs+eEO0BxyVtqS8/R7kn0GryZpc3sEeA0cL+7n6hS54iZ5c0sv7Cw0FT784uFpspl+2vmmDc6D0o1XsU20+bq8mbnZbvKa2ln/Xa1pedoZyQa8O5ecvf9wEXA1WZ2RZU6d7h71t2zg4ODTbU/lImbKpftr5lj3ug8SJttus3V5c3Oy3aV19LO+u1qS8/RzujIKhp3XwQeAN7eznYnxkaIo/SasjhKMzE20s5uZAtp5phXq7teHKU5fM3Fm2pzfd1m52U76tfSyvNhYmyEKH3uL70oZS21pedo9yS5imYQKLr7opnFwFuB32lnH8tv0ugd+t7RzDGvVrfWKprsa89vuc31dZudl5ut3+5VNMv127GKRs/R7kpyFc0bgc8AacqvFL7g7r9Vb59mV9GIiPS6eqtoEjuDd/dvAaNJtS8iIvXpk6wiIoFSwIuIBEoBLyISKAW8iEigFPAiIoFKbJlkK8xsAXiq2+Oo4jXAD7s9iC7RY+9NvfrYt+Pjfq27V70NwJYK+K3KzPK11pmGTo9dj72XhPa4dYlGRCRQCngRkUAp4BtzR7cH0EV67L2pVx97UI9b1+BFRAKlM3gRkUAp4EVEAqWAr8PMzjOzb5jZo2b2uJnd3u0xdVLlTy7OmNmfdnssnWRm3zOzk2b2iJn11P2rzSxjZvea2RNm9m0z+5luj6kTzGykcryXv/6fmf1at8e1WYndLjgQPwbe4u4vmlkEfN3M/pe7P9TtgXXIB4BvA6/q9kC64Dp3324feGmH3wO+4u43m1kf0N/tAXWCu88C+6F8YgPMAfd1dVBtoDP4Orzsxcp/o8pXT7wrbWYXAe8E7ur2WKQzzOxVwLXApwHc/Wzlz232muuBJ919K36qvikK+A1ULlM8ApwG7nf3E90eU4d8EvgNYKnbA+kCB75qZtNmdqTbg+mgS4AF4A8rl+buMrOBbg+qC94L3N3tQbSDAn4D7l5y9/3ARcDVZnZFt8eUNDN7F3Da3ae7PZYuOejuVwLvAH7FzK7t9oA6ZAdwJfApdx8FXgJ+s7tD6qzKZakbgS92eyztoIBvUOWl6gPA27s8lE44CNxoZt8DPg+8xcw+290hdY67z1f+PU35OuzV3R1Rx3wf+P6qV6n3Ug78XvIO4GF3f6bbA2kHBXwdZjZoZpnK9zHwVuCJ7o4qee5+1N0vcvc9lF+u/oW739rlYXWEmQ2Y2c7l74EbgMe6O6rOcPe/B06Z2Uil6Hrgb7s4pG44TCCXZ0CraDbyU8BnKu+qp4AvuHtPLRnsQRcC95kZlJ8fn3P3r3R3SB31q8AfVy5VfBd4X5fH0zFm1g+8Dfi33R5Lu+hWBSIigdIlGhGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngJRhmVqrcCfAxM/tiZdlbq229efkummZ2o5nV/ERn5Q6M/76FPj5iZr/e6hhFNqKAl5AU3H2/u18BnAX+3eqNVtb0nHf3P3H3365TJQM0HfAiSVPAS6j+GvhpM9tTua/5fwMeBi42sxvM7EEze7hypv9KADN7e+U+6F8HDi03ZGa/YGb/tfL9hWZ2X+VvBDxqZm8Cfhu4tPLqIVepN2Fm3zSzb63+OwJm9kEzmzWzPwNGEEmQAl6CY2Y7KN9T5GSlaAT476tuoPUh4K2VG4rlgdvM7DzgTuDdwM8C/7hG878P/KW776N8n5bHKd+Q68nKq4cJM7sBuIzyPWz2A1eZ2bVmdhXlWz+MUv4F8s/a/NBF1tCtCiQkceXWzlA+g/80MAQ8teqPtBwAXg/8TeV2BH3Ag8DlwN+5+3cAKjdXq3ar4LcA/wbKdxoFnjezXevq3FD5mqn8/5WUA38ncJ+7n6n08SeberQiG1DAS0gKlVs7r6iE+Euriyjf1//wunr7ad8fczHgmLv/wbo+fq2NfYhsSJdopNc8BBw0s5+G8g2mzOyfUL5L6OvM7NJKvcM19v9z4P2VfdOVv4L0AuWz82VTwC+uurY/bGYXAH8F3GRmceWOle9u82MTWUMBLz3F3ReAXwDuNrNvUQ78y939R5QvyfyPypustf5c2weA68zsJDANvMHd/4HyJZ/HzCzn7l8FPgc8WKl3L7DT3R8G7gEeAb5E+TKSSGJ0N0kRkUDpDF5EJFAKeBGRQCngRUQCpYAXEQmUAl5EJFAKeBGRQCngRUQC9f8BeH9u/DsJ36gAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAEICAYAAABVv+9nAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAeb0lEQVR4nO3dfbAddZ3n8fe3+/Q55948QnLBSNTEh8UnEDBSOLhTovIQYB0Vl1GHKXWYjeVayOwOjrBVOuXUWuvMVrmsOz4sKuqWDyMDsj6hE3DIqKOCCUaMgBtAlBggl0jg5uE8dX/3j+5zcwk3yeUmfTv59edVdar79Hn4fbty8+lf/06f3zF3R0REwhNVXYCIiJRDAS8iEigFvIhIoBTwIiKBUsCLiARKAS8iEqhSA97MFpvZ9WZ2j5ndbWavLLM9ERHZq1Hy+/9P4Lvu/mYzawKjB3ry0qVLfcWKFSWXJCISjg0bNjzq7mPTPVZawJvZQuAPgXcAuHsP6B3oNStWrGD9+vVllSQiEhwz+83+HitziOa5wDjwOTP7mZl9xszmTVPcGjNbb2brx8fHSyxHRKReygz4BnAa8El3PxXYBVy575Pc/Rp3X+Xuq8bGpj3LEBGRWSgz4LcAW9z9tuL+9eSBLyIic6C0MXh3f9jMHjSzE939V8BrgbvKak9E6qnf77NlyxY6nU7VpZSq3W6zfPlykiSZ8WvKvormMuBLxRU09wPvLLk9EamZLVu2sGDBAlasWIGZVV1OKdyd7du3s2XLFlauXDnj15Ua8O6+EVhVZhsiUm+dTifocAcwM5YsWcLTvRBF32QVkaNeyOE+NJt9DCPg/+Xv4N5bqq5CROSIEkbA//BquO/WqqsQkRrasWMHn/jEJ572684//3x27NhRQkV7hRHwSRv6e6quQkRqaH8Bn6bpAV930003sXjx4rLKAsq/imZuNEZgEPYlUiJyZLryyiu57777OOWUU0iShPnz57Ns2TI2btzIXXfdxRve8AYefPBBOp0Ol19+OWvWrAH2Ts2yc+dOVq9ezate9Sp+9KMfccIJJ/D1r3+dkZGRQ64tjIBXD15EgA9985fctfWJw/qeL37mQv76371kv49/5CMfYdOmTWzcuJF169ZxwQUXsGnTpsnLGa+99lqOPfZY9uzZwyte8QouuugilixZ8qT32Lx5M1/5ylf49Kc/zcUXX8wNN9zAJZdccsi1BxHw3miT9fcQV12IiNTe6aef/qRr1T/2sY9x4403AvDggw+yefPmpwT8ypUrOeWUUwB4+ctfzgMPPHBYagki4Dc+3GXxnu3M/PJ/EQnRgXrac2XevL1zKq5bt45bbrmFH//4x4yOjvLqV7962m/ctlqtyfU4jtmz5/CMSATxIWvfmsSpxuBFZO4tWLCAiYmJaR97/PHHOeaYYxgdHeWee+7hJz/5yZzWFkQPvmctonRn1WWISA0tWbKEM888k5e+9KWMjIxw/PHHTz523nnn8alPfYqTTz6ZE088kTPOOGNOawsi4AdRi0a2veoyRKSmvvzlL0+7vdVq8Z3vfGfax4bj7EuXLmXTpk2T26+44orDVlcQQzR5wHerLkNE5IgSRMCnCngRkacII+DjNk1XwIuITBVEwGeNNokf8Pe8RURqJ4iA97hNwgCyA8/9ICJSJ2EEfKOdr2i6AhGRSUEEPMOA14RjIjLHZjtdMMDVV1/N7t27D3NFe4UR8Ekx65p68CIyx47kgA/ii04k6sGLSDWmThd89tlnc9xxx3HdddfR7XZ54xvfyIc+9CF27drFxRdfzJYtW0jTlA984AM88sgjbN26lbPOOoulS5dy662H/0eLggh4S0YBSHu7NaOkSJ1950p4+BeH9z2fcRKs/sh+H546XfDatWu5/vrruf3223F3Xv/61/P973+f8fFxnvnMZ/Ltb38byOeoWbRoER/96Ee59dZbWbp06eGtuRDEEE3UzIdo+p3yTnVERA5m7dq1rF27llNPPZXTTjuNe+65h82bN3PSSSdxyy238P73v58f/OAHLFq0aE7qCaIHHxVj8N3OLtoV1yIiFTpAT3suuDtXXXUV73rXu57y2IYNG7jpppu46qqrOOecc/jgBz9Yej1B9OAb7XyIpt9VD15E5tbU6YLPPfdcrr32WnbuzGe3/d3vfse2bdvYunUro6OjXHLJJVxxxRXccccdT3ltGYLowcfNPOAHGqIRkTk2dbrg1atX87a3vY1XvvKVAMyfP58vfvGL3Hvvvbzvfe8jiiKSJOGTn/wkAGvWrGH16tUsW7ZMH7LuT6NVBLx68CJSgX2nC7788sufdP95z3se55577lNed9lll3HZZZeVVlcQQzRJe+9VNCIikgsj4FvDgNcXnUREhoII+GY7/5HbTAEvUkvuXnUJpZvNPpY6Bm9mDwATQAoM3H1VGe20iiEaBbxI/bTbbbZv386SJUsws6rLKYW7s337dtrtp3ch+Fx8yHqWuz9aZgMjzQYdT3DNRSNSO8uXL2fLli2Mj49XXUqp2u02y5cvf1qvCeIqmnYS0aGJ9zUXjUjdJEnCypUrqy7jiFT2GLwDa81sg5mtme4JZrbGzNab2frZHoFbSUyHJjZQD15EZKjsgD/T3U8DVgPvMbM/3PcJ7n6Nu69y91VjY2OzamQkiel4E9MQjYjIpFID3t23FsttwI3A6WW0k8RGlyam6YJFRCaVFvBmNs/MFgzXgXOATSW1RdeaWKqAFxEZKvND1uOBG4vLlhrAl939u2U11rcWcdot6+1FRI46pQW8u98PvKys999XP2oRpxqDFxEZCuKbrACDqEmcqQcvIjIUUMC3aSjgRUQmBRPwadQiUcCLiEwKJ+DjFokr4EVEhoIJ+Cxu0/Re1WWIiBwxwgn4RpsWPajBtKEiIjMRTMB7YyRf0bdZRUSAgAKeRjFPsuajEREBAgp4S4qAVw9eRAQIKOBJ8l910o9+iIjkggn4YQ++391dcSUiIkeGYAI+auY9+N4eBbyICAQV8PlVNP3OroorERE5MgQT8HHRg9cQjYhILpiAb7SKHrwCXkQECCjgk1begx8o4EVEgIACvtGeB0CqgBcRAQIK+GY7H6JJe7oOXkQEggr4ogevgBcRAQIMeH2TVUQkF0zAt1tNeh7j6sGLiAABBfxIEtOhCX19yCoiAgEFfDuJ6NLUbJIiIoVwAr4R0/EmpoAXEQECCvgoMrqmgBcRGQom4AF61iRKFfAiIhBcwLeI027VZYiIHBGCCvi+tdSDFxEpBBXwg6hFI1MPXkQE5iDgzSw2s5+Z2bfKbitVwIuITJqLHvzlwN1z0A6DuEWigBcRAUoOeDNbDlwAfKbMdoayuE3iCngRESi/B3818FdAtr8nmNkaM1tvZuvHx8cPqbEsbtH03iG9h4hIKEoLeDO7ENjm7hsO9Dx3v8bdV7n7qrGxsUNqM2uM0KQH7of0PiIiISizB38m8HozewD4B+A1ZvbFEtuDRpuYDNJ+qc2IiBwNSgt4d7/K3Ze7+wrgLcA/u/slZbUH4I12vjLQlMEiIkFdB+9J/rN99PVlJxGRxlw04u7rgHVlt2ON4e+y7iYuuzERkSNcUD34qJkHfK+zq+JKRESqF1bAJ8OA1xi8iEhYAa8evIjIpKACvtEeBaDf1e+yiogEFfBxMw/4gXrwIiJhBXyjlQ/RDLoagxcRCSrgm615QH6ZpIhI3QUV8EkxBp/21IMXEQkq4JvtvAefqQcvIhJWwLdGhgGvHryISFAB3261SN3wvgJeRCSogG81G3RoggJeRCSsgG8nUR7wA80mKSISVMA344guTUwBLyISVsCbmQJeRKQQVMAD9KxFlCrgRUSCC/i+NYnTbtVliIhUbkYBb2bPM7NWsf5qM3uvmS0ut7TZ6UctYvXgRURm3IO/AUjN7PnAZ4GVwJdLq+oQ9K1FnKkHLyIy04DP3H0AvBG42t3/E7CsvLJmL41bJAp4EZEZB3zfzN4KvB34VrEtKaekQ5PGbRoKeBGRGQf8O4FXAh9291+b2Urgi+WVNXtp1CJxBbyISGMmT3L3u4D3ApjZMcACd/9ImYXNVtZo0/Re1WWIiFRuplfRrDOzhWZ2LPBz4HNm9tFyS5sdj9u0UMCLiMx0iGaRuz8BvAn4nLu/HHhdeWXNnjfaJAwgS6suRUSkUjMN+IaZLQMuZu+HrEckb+S/y6oZJUWk7mYa8H8D/BNwn7v/1MyeC2wur6zZsyQPeM0JLyJ1N9MPWf8R+Mcp9+8HLiqrqEOStAHod3fTnF9xLSIiFZrph6zLzexGM9tmZo+Y2Q1mtrzs4mYjauY/vN3ds6viSkREqjXTIZrPAd8AngmcAHyz2LZfZtY2s9vN7Odm9ksz+9ChlTozUTMfoul39MPbIlJvMw34MXf/nLsPitvngbGDvKYLvMbdXwacApxnZmccQq0zMgz4ngJeRGpupgH/qJldYmZxcbsE2H6gF3huZ3E3KW5+CLXOSDzswXcV8CJSbzMN+D8jv0TyYeAh4M3k0xccUHEw2AhsA25299umec4aM1tvZuvHx8dnXvl+JK18DL7f0Ri8iNTbjALe3X/r7q939zF3P87d30D+paeDvS5191OA5cDpZvbSaZ5zjbuvcvdVY2MHG/U5uEZrHgBpTz14Eam3Q/lFp/880ye6+w5gHXDeIbQ3I0k7H6IZaIhGRGruUALeDvig2djwV5/MbIR8aoN7DqG9GRn24LOevugkIvU2oy867cfBPjBdBnzBzGLyA8l17l76NAetkXwMXgEvInV3wIA3swmmD3IDRg70Wne/Ezh19qXNTnNEPXgREThIwLv7grkq5HBpF1fRaC4aEam7QxmDPyKNtBp0PIGBAl5E6i24gG83Yjo0YdCpuhQRkUoFF/BRZHRpYgp4Eam54AIeoGtNIgW8iNRckAHfo0WUKuBFpN6CDPh+pIAXEQkz4K1JnHarLkNEpFJBBvwgatHIFPAiUm+BBnxbAS8itRdkwKdxi0QBLyI1F27AuwJeROotyIDP4jZN71VdhohIpYIMeG+0aaKAF5F6CzLgaYzQpgde+m98i4gcsYIM+CwppqrXdAUiUmNBBrw12gCk+tEPEamxMAO+mffge52dFVciIlKdMAO+kQd8d8/uiisREalOkAEfFT34bkcBLyL1FWTAx838d1n7nV0VVyIiUp0wA76V9+AV8CJSZ0EGfKOV9+AHXQ3RiEh9BR7wukxSROoryIBvFgGfqgcvIjUWZMAn7XkAZPqik4jUWJAB3xwpevB9BbyI1FeQAd8uevDe0xCNiNRXkAHfHCkCvq/JxkSkvkoLeDN7lpndamZ3m9kvzezystraV7uV0PMYBhqiEZH6apT43gPgL939DjNbAGwws5vd/a4S2wSgGUdM0AT14EWkxkrrwbv7Q+5+R7E+AdwNnFBWe1OZGT2amOaDF5Eam5MxeDNbAZwK3DbNY2vMbL2ZrR8fHz9sbXathaUKeBGpr9ID3szmAzcAf+HuT+z7uLtf4+6r3H3V2NjYYWu3Z00i9eBFpMZKDXgzS8jD/Uvu/rUy29pX31rEmQJeROqrzKtoDPgscLe7f7Ssdvanby3itDvXzYqIHDHK7MGfCfwp8Boz21jczi+xvScZRE0amQJeROqrtMsk3f2HgJX1/gcziNs0Bk8Z8hcRqY0gv8kKkEYtEvXgRaTGgg34LG6TuAJeROor3IBvtGl6r+oyREQqE27Ax21aKOBFpL6CDXhPRmh5F9yrLkVEpBLBBjyNNrE5nqoXLyL1FGzAW9IGoN/Rj36ISD0FHPAjAHT27Kq4EhGRagQb8CT577L2Ogp4EamnYAM+buY9+P4eDdGISD0FH/Bd9eBFpKYCDvh8iGbQVQ9eROop2IBvtBXwIlJvwQZ80sqHaAY9BbyI1FOwAd9oDXvweyquRESkGsEGfNKeB4CrBy8iNRVswLdG8oDPeurBi0g9hRvwxYesCngRqatwA77owfugU3ElIiLVCDbg2602qRv01YMXkXoKNuDjOKJDEwYKeBGpp2ADHqBrTWyg32UVkXoKOuB7tIjUgxeRmgo74K2JperBi0g9BR3wfWsRp7qKRkTqKeyAj1rEmXrwIlJPQQf8Y8nxPKf7//DOE1WXIiIy54IO+B0n/zkLfCe/XfvxqksREZlzQQf8Wa+7gNvsZBZt/N/6wpOI1E5pAW9m15rZNjPbVFYbB9NOYrae/B4WZ4+x9dZrqipDRKQSZfbgPw+cV+L7z8hrz3sTd/iJtG//exj0qi5HRGTOlBbw7v594Pdlvf9MLRxpsvmF7+bYwTbGf/SFqssREZkzlY/Bm9kaM1tvZuvHx8dLaeOsC97CL/y52A+vhnRQShsiIkeaygPe3a9x91XuvmpsbKyUNo5bOMKdKy9laW8Lj6//ailtiIgcaSoP+Lnyby98B7/KltNb998hy6ouR0SkdLUJ+Gcvnc9tJ7yDsT2/Zted36i6HBGR0pV5meRXgB8DJ5rZFjO7tKy2ZuoVF17Kr7Pj2XXLfwP3qssRESlVmVfRvNXdl7l74u7L3f2zZbU1Uy864VjWjV3CcTvvoXvP2qrLEREpVW2GaIZOOn8Nv/Ml7LzpA7D9vqrLEREpTe0CftXznsFXF69h3sSv8f+1Cr/hP8C2e6ouS0TksKtdwAOc/9b/yJ8t+izXDFbT+cU38E+cgX/1T+GhO6suTUTksDE/gj5sXLVqla9fv35O2soy55t3buUz/7Sesye+xp8naxn13fD8s+GlF8ELzoF5S+akFhGR2TKzDe6+atrH6hrwQ71Bxld/+ls+e8vPubDzTS5tf49j0t/jGCw/HTvxXPg3q+G4F4HZnNYmInIwCvgZ2N0b8Ll/fYDP//B+nrH7V7w2voOz4428xO4HoDO6DI5ZQWPeYhojx0B7EbQX5svmPEhGIRkpbsV6YwQaTWi0i1sL4hbEjUr2UUTCo4B/Gtyd3/5+Nxt+8xgbfvMYv/71vTx7+7/yquhOxuxxFrKbRbabhbab+eyeXRsW4VET4gTiKcsoycM/TrAogagBcQOLE7AYLIKoWA5vWUo26JINung/XzLo4p5B1MCtgUdxsZ4vrXhfogSbbK8x/ftbhEXR5OssjrEof741EqJGqzhwNZ+8jBrFLZ6y3sj3o3i/fD3Ol2k3n7O/v7tYFrekDSPHwsgxMHpsvp6M5GdTWQq9ndDdWSwn8m3tRTCyGNqL89eLBOxAAa+u5D7MjOcsmcdzlszjTactB05ionMhv9jyOA8+3mHbRIdtT3QZn+jy6BO72Dmxg7SzC+/tIk47tOkxYt18SY8mfVrWp8mAFj1axf0GKU0GNEhJGNC0Acnk/ZQGA2I6JKQkNiDCicjypTkRTkzGwI0eCT0a9LxRrLfJMGKcmJQGKTE9GqQ0LL+f37KinYyEAQaYTWlnss2MBlnxXsV9Syv7N+pbk8xiWtnBf8Slb026jYX0kgUM4jZZ1CKNW/my0cajJkm2h9ZgJ8lgJ83+BI3BBI3+TjxqkDYXkbUWkbUX4a3F+YGjOS8/4MX5AdKiBGsk+YGQ6YfxrDiYRnFcHCiLA5xn4Gl+YMpSyAb5fYvzA+XwrG/qehTnnYIoP0hPHkQnhxCL5fD+1ANs1CheG095fXEg1hBkcBTwM7CgnfAHz1960OcN0ozd/ZTd3ZRdvQF7eim9NKPbz+ilGb1BcUtTBqmTZs4gc3qZsyvNSDMnc8jcyYr1tFh3HHdwKJb5/SQ2RpKYdhIz0oxpN/JlEuX/WTOgO3mSVrxn5k+6DTInLebnGT51eGLnvrcmH9bjTpZmZOkAH3TxQRcGHRh0yQY9bNDFPcXTAWSDyeVkeGUpeIoV6+YZadSkH7UZRC3SuD0ZxpZ2aXYfo9l/nPbgcdqDJ5ifPkGUDehE8+jGo/TiUbrxPPrxPDKLSfoTtAZP0B5MMJo+wWhnggWd3bQnD7BT1umxixYPMcqEjzLBciZ8hAlGaZCyqLOLRbaLRTbBIh5moe1ilO7kwS5mQLPCg93hlFmMF7epJv8mLCaLW3jcIoua+NSDz9QzNYuLM8QYK87SrDgztOKMLT8Ykp8lYpMHF4uKs0YMs3y7mUEUYRYTRU9+v/zM1op1m/J+U85Co6nPKW55a0W79uSD2+Soxr6jG/s+f5rXD+/jxX9Uzw/iFMu0N82Z6m48bmGrP3KI/4JPpYA/jBpxxMI4YmE7qboUmcLd6aUZWQaDbO8y9fwAZ6mTpBkLU6efZsXNGaT5gbmfOr9PMx6Zsn144B0eILOiBz61Dzz8P+/ukKV4NiAbDHBPIR2QpQNSdzppRDczupnRS41OCmmaFgfO/GZpsZ71iLIBpAPMB1g2IPIBeHEGNmx72L47ESmWDTBPiT1fJqRTDlL5mV1+Jpdi7J2Mb+r+5Ged/cmz0Fax3mQXsWXEZMXZ4t73HZ4FxmTElmHFmWeEY3l3pag7vz/cnre797GpZ5P5ex05Q8uzlbmxhyZ7aPGoLeGFCniRp8/MaDWGvdL4gM+tAy8ObMMD3L5nc/v7WC51pzfI6A7SvWejg4xumk0e7DJ30mzvmWdWtPGks7/MJ89Eh/UMZZNnqcVzi07wvme1nmVknuFpsczSJ93Psiw/UKUZ7vnNsjQ/s3QAz7flFRSddSeyiMjAIiM2I44jzAx3z9/PwTMHMijqw4dtMNmWW0w0eeYRFWciEcRN0rgNySgWN2k0YhqRsaDd4IUl/Fsr4EVqxsxoxKb//DVQy2+yiojUgQJeRCRQCngRkUAp4EVEAqWAFxEJlAJeRCRQCngRkUAp4EVEAnVEzSZpZuPAb2b58qXAo4exnKOF9rtetN/1MpP9fo67j033wBEV8IfCzNbvb8rMkGm/60X7XS+Hut8aohERCZQCXkQkUCEF/DVVF1AR7Xe9aL/r5ZD2O5gxeBERebKQevAiIjKFAl5EJFBHfcCb2Xlm9iszu9fMrqy6njKZ2bVmts3MNk3ZdqyZ3Wxmm4vlMVXWeLiZ2bPM7FYzu9vMfmlmlxfbg95vADNrm9ntZvbzYt8/VGxfaWa3Ffv+VTNrVl3r4WZmsZn9zMy+VdwPfp8BzOwBM/uFmW00s/XFtln/rR/VAW9mMfBxYDXwYuCtZvbiaqsq1eeB8/bZdiXwPXd/AfC94n5IBsBfuvuLgDOA9xT/xqHvN0AXeI27vww4BTjPzM4A/hb4H8W+PwZcWmGNZbkcuHvK/Trs89BZ7n7KlOvfZ/23flQHPHA6cK+73+/uPeAfgD+quKbSuPv3gd/vs/mPgC8U618A3jCnRZXM3R9y9zuK9Qny//QnEPh+A3huZ3E3KW4OvAa4vtge3L6b2XLgAuAzxX0j8H0+iFn/rR/tAX8C8OCU+1uKbXVyvLs/BHkYAsdVXE9pzGwFcCpwGzXZ72KoYiOwDbgZuA/Y4e6D4ikh/s1fDfwVkBX3lxD+Pg85sNbMNpjZmmLbrP/Wj/bf3bVptum6zwCZ2XzgBuAv3P2JvFMXPndPgVPMbDFwI/Ci6Z42t1WVx8wuBLa5+wYze/Vw8zRPDWaf93Gmu281s+OAm83snkN5s6O9B78FeNaU+8uBrRXVUpVHzGwZQLHcVnE9h52ZJeTh/iV3/1qxOfj9nsrddwDryD+HWGxmw85ZaH/zZwKvN7MHyIdcX0Peow95nye5+9ZiuY38gH46h/C3frQH/E+BFxSfsDeBtwDfqLimufYN4O3F+tuBr1dYy2FXjL9+Frjb3T865aGg9xvAzMaKnjtmNgK8jvwziFuBNxdPC2rf3f0qd1/u7ivI/z//s7v/CQHv85CZzTOzBcN14BxgE4fwt37Uf5PVzM4nP8LHwLXu/uGKSyqNmX0FeDX5FKKPAH8N/F/gOuDZwG+Bf+/u+34Qe9Qys1cBPwB+wd4x2f9CPg4f7H4DmNnJ5B+qxeSdsevc/W/M7LnkvdtjgZ8Bl7h7t7pKy1EM0Vzh7hfWYZ+LfbyxuNsAvuzuHzazJczyb/2oD3gREZne0T5EIyIi+6GAFxEJlAJeRCRQCngRkUAp4EVEAqWAl1oxs7SYqW94O2yTlJnZiqkzfYpU7WifqkDk6drj7qdUXYTIXFAPXoTJebj/tph//XYze36x/Tlm9j0zu7NYPrvYfryZ3VjM1f5zM/uD4q1iM/t0MX/72uIbqCKVUMBL3YzsM0Tzx1Mee8LdTwf+nvzb0RTr/8fdTwa+BHys2P4x4F+KudpPA35ZbH8B8HF3fwmwA7io5P0R2S99k1Vqxcx2uvv8abY/QP7jGvcXk5s97O5LzOxRYJm794vtD7n7UjMbB5ZP/bp8MZ3xzcUPM2Bm7wcSd/+v5e+ZyFOpBy+yl+9nfX/Pmc7U+VFS9DmXVEgBL7LXH09Z/rhY/xH5rIYAfwL8sFj/HvBumPxRjoVzVaTITKl3IXUzUvxC0tB33X14qWTLzG4j7/i8tdj2XuBaM3sfMA68s9h+OXCNmV1K3lN/N/BQ6dWLPA0agxdhcgx+lbs/WnUtIoeLhmhERAKlHryISKDUgxcRCZQCXkQkUAp4EZFAKeBFRAKlgBcRCdT/B+hvQNY3rKw3AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# The code below implements the training.\n",
    "# If you correctly implement  dnn and update_weights above, \n",
    "# you should not need to change anything below. \n",
    "# (apart from increasing the number of epochs)\n",
    "\n",
    "train_losses = []\n",
    "test_losses = []\n",
    "\n",
    "# How many epochs to train\n",
    "n_epochs = 50\n",
    "\n",
    "# Loop over the epochs\n",
    "for ep in range(n_epochs):\n",
    "        \n",
    "    # Each epoch is a complete over the training data\n",
    "    for i in range(X_train.shape[0]):\n",
    "        \n",
    "        # pick one example\n",
    "        x = X_train[i]\n",
    "        y = y_train[i]\n",
    "\n",
    "        # use it to update the weights\n",
    "        update_weights(x,y,W,b,Wp,bp)\n",
    "    \n",
    "    # Calculate predictions for the full training and testing sample\n",
    "    y_pred_train = [dnn(x,W,b,Wp,bp)[0] for x in X_train]\n",
    "    y_pred = [dnn(x,W,b,Wp,bp)[0] for x in X_test]\n",
    "\n",
    "    # Calculate aver loss / example over the epoch\n",
    "    train_loss = sum((y_pred_train-y_train)**2) / y_train.shape[0]\n",
    "    test_loss = sum((y_pred-y_test)**2) / y_test.shape[0] \n",
    "    \n",
    "    # print some information\n",
    "    print(\"Epoch:\",ep, \"Train Loss:\", train_loss, \"Test Loss:\", test_loss)\n",
    "    \n",
    "    # and store the losses for later use\n",
    "    train_losses.append(train_loss)\n",
    "    test_losses.append(test_loss)\n",
    "    \n",
    "    \n",
    "# After the training:\n",
    "    \n",
    "# Prepare scatter plot\n",
    "y_pred = [dnn(x,W,b,Wp,bp)[0] for x in X_test]\n",
    "\n",
    "print(\"Best loss:\", min(test_losses), \"Final loss:\", test_losses[-1])\n",
    "\n",
    "print(\"Correlation coefficient:\", np.corrcoef(y_pred,y_test)[0,1])\n",
    "plt.scatter(y_pred_train,y_train)\n",
    "plt.xlabel(\"Predicted\")\n",
    "plt.ylabel(\"True\")\n",
    "plt.show()\n",
    "\n",
    "# Prepare and loss over time\n",
    "plt.plot(train_losses,label=\"train\")\n",
    "plt.plot(test_losses,label=\"test\")\n",
    "plt.legend()\n",
    "plt.xlabel(\"Epoch\")\n",
    "plt.ylabel(\"Loss\")\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Hint 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We want a network with one hidden layer. As activiation in the hidden layer $\\sigma$ we apply element-wise ReLu, while no activation is used for the output layer. The forward pass of the network then reads:\n",
    "$$\\hat{y}=\\mathbf{W}^{\\prime} \\sigma(\\mathbf{W} \\vec{x}+\\vec{b})+b^{\\prime}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Hint 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the regression problem the objective function is the mean squared error between the prediction and the true label $y$:\n",
    "$$\n",
    "L=(\\hat{y}-y)^{2}\n",
    "$$\n",
    "\n",
    "Taking the partial derivatives - and diligently the applying chain rule - with respect to the different objects yields:\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial L}{\\partial b^{\\prime}} &=2(\\hat{y}-y) \\\\\n",
    "\\frac{\\partial L}{\\partial b_{k}} &=2(\\hat{y}-y) \\mathbf{W}_{k}^{\\prime} \\theta\\left(\\sum_{i} \\mathbf{W}_{i k} x_{i}+b_{k}\\right) \\\\\n",
    "\\frac{\\partial L}{\\partial \\mathbf{W}_{k}^{\\prime}} &=2(\\hat{y}-y) \\sigma\\left(\\sum_{i} \\mathbf{W}_{i k} x_{i}+b_{k}\\right) \\\\\n",
    "\\frac{\\partial L}{\\partial \\mathbf{W}_{k m}} &=2(\\hat{y}-y) \\mathbf{W}_{m}^{\\prime} \\theta\\left(\\sum_{i} \\mathbf{W}_{i k} x_{i}+b_{m}\\right) x_{k}\n",
    "\\end{aligned}\n",
    "\n",
    "Here, $\\Theta$ denotes the Heaviside step function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "collapsed_sections": [],
   "name": "Exercise 4",
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}