Neural_networks_with_keras_solved.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import division, print_function, generators\n",
    "\n",
    "import tensorflow as tf\n",
    "import keras\n",
    "\n",
    "import numpy as np\n",
    "import keras.datasets.mnist as mnist\n",
    "import keras.datasets.cifar10 as cifar10\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Schedule\n",
    "\n",
    "- What is MNIST? What is CIFAR? (10 minutes)\n",
    "\n",
    "- Loading images and converting them to numpy (5 minutes)\n",
    "\n",
    "- What is Keras? Creating a basic model, Callbacks, saving/loading a model, changing optimizer, loss  (30 minutes)\n",
    "\n",
    "- Mini-Project: Train a feed forward NN on MNIST/CIFAR10 with Keras (No scaling) (15 minutes)\n",
    "\n",
    "- Preprocessing images: Scaling (5 minutes)\n",
    "\n",
    "- Mini-Project: Train a feed forward NN on MNIST/CIFAR10 with Keras (10 minutes)\n",
    "\n",
    "- Preprocessing images: Data augmentation and generators with Keras (10 minutes)\n",
    "\n",
    "- Mini-Project: Train a feed forward NN on MNIST/CIFAR10 using data augmentation with Keras (15 minutes)\n",
    "\n",
    "- Break (20 minutes)\n",
    "\n",
    "- Convolutional NN, Pooling, Flattening, etc  (20 minutes)\n",
    "\n",
    "- How to build a CNN in Keras (10 minutes)\n",
    "\n",
    "- Mini-Project: Train CNN on MNIST/CIFAR10 (30 minutes)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## What is CIFAR10?\n",
    "\n",
    "More information at\n",
    "https://www.cs.toronto.edu/~kriz/cifar.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# load dataset\n",
    "(x_train, y_train), (x_test, y_test) = cifar10.load_data()\n",
    "\n",
    "# This is needed to know what each class means\n",
    "label_id_to_class_name = {0:'airplane', 1:'automovile', 2:'bird', 3:'cat', 4:'deer', 5:'dog', \n",
    "                          6:'frog', 7:'horse', 8:'ship', 9:'truck'}\n",
    "\n",
    "print(x_train.shape[0], 'Number of train samples')\n",
    "print(x_test.shape[0], 'Number of test samples')\n",
    "print('x_train shape:', x_train.shape)\n",
    "\n",
    "# plot images\n",
    "for image_id in range(0, 5):\n",
    "    plt.imshow(x_train[image_id])\n",
    "    plt.title(\"The true label is %s\" % label_id_to_class_name[int(y_train[image_id])])\n",
    "    plt.show()\n",
    "    \n",
    "# This makes sure the image has the correct order in the axis for Tensorflow, it would be different for Theano backend\n",
    "x_train = x_train.reshape(x_train.shape[0], 32, 32, 3)\n",
    "x_test = x_test.reshape(x_test.shape[0], 32, 32, 3)\n",
    "\n",
    "# Convert values to floats, originally they are integers\n",
    "x_train = x_train.astype('float32')\n",
    "x_test = x_test.astype('float32')\n",
    "\n",
    "# Convert values of labels from 0 to 9 to categorical (one_hot encoding)\n",
    "y_train = keras.utils.to_categorical(y_train, 10)\n",
    "y_test = keras.utils.to_categorical(y_test, 10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## What is MNIST?\n",
    "\n",
    "More information at http://yann.lecun.com/exdb/mnist/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "60000 Number of train samples\n",
      "10000 Number of test samples\n",
      "x_train shape: (60000, 28, 28)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP8AAAEICAYAAACQ6CLfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEg5JREFUeJzt3X2QXXV9x/H3xyC1hECIaIgxEIM0DjAYpyE4TqbI0PBM\nIco4RplioURbYnHGUpk41mAbhyrQIQOjCQVMrEaYIiUwjiCP0WFMWSFACCKIYUhcEmlYSMKTyX77\nxz2hN8ve3717n87N/j6vmZ2993zPwzc3+9lz7jn37E8RgZnl5x1lN2Bm5XD4zTLl8JtlyuE3y5TD\nb5Yph98sUw5/F0haJOk/y+6jWZJC0gcbmG9qMe8+TWwjuaykJyR9fKTrtdoc/jaQtL3qa1DSa1XP\nP9vmbX1P0r+2c517g4g4KiLuH8kyVb9Qqv9/vtahFvc6I/4NbW8XEfvvfixpA/C3EXF31bRF3epF\n0j4RsbNb29tLjPdr8nbe83fPvpJWSNpWHMLO3F2Q9D5Jt0j6g6TfSfqH4VYgaT7wWeCfir3Y7cX0\nDZK+IukxYIekfYYeqg89YpB0hqS1kgYkPSjpmEb+EZJOl/SIpFckPV/jF9v5kn4vqV/SP1Yt+w5J\nl0r6raT/lXSzpAkNbneDpL8sHs+S1Ff0sFnSVY2sw/bk8HfPXwE/AsYDq4BroBII4HbgUWAycCLw\nJUknD11BRCwDfgB8KyL2j4gzq8rzgNNpYC8n6SPADcDngXcDS4FVkv6kgX/HDuCvi3/H6cDfSTp7\nyDwnAEcAJwFf2R1a4IvA2cDxwPuAl4BrG9jmUFcDV0fEAcDhwM115n9O0kZJN0o6uIntjUoOf/f8\nIiJ+EhG7gO8DHy6mHwu8JyK+ERFvRsSzwHXAp0e4/iUR8XxEvNbAvPOBpRGxJiJ2RcRy4A3go/UW\njIj7I+LxiBiMiMeAlVTCXO2yiNgREY8DN1L5xQTwBeCrEbExIt4AFgHnNHGC8I/AByUdHBHbI+KX\nNeZ7kcrrexjw58A4Kr88DYe/m16oevwq8K7ih/4w4H3F4feApAFgITBxhOt/fgTzHgZ8ecg2p1DZ\nGydJOk7SfcVblJepBHro3rS6l+eq1nsYcGvVNp8EdjHyf+sFwJ8Bv5b0kKQzhpup+MXQFxE7I2Iz\nsAA4SdK4EW5vVHL4y/c88LuIGF/1NS4iTqsxf63bMIdOfxXYr+r5IUO2uXjINveLiJUN9PtDKm9b\npkTEgcB3AQ2ZZ0rV40OB31dt99Qh231XRGxqYLtviYinI2Ie8F7g34D/kjS2kUWL7/65xy9CL/gf\nYFtxwu5PJY2RdLSkY2vMvxmY1sB61wKfKdZ3Cnseml8HfKHYi0vS2OJEXiN7xHHA1oh4XdIs4DPD\nzPM1SftJOgr4G+CmYvp3gcWSDgOQ9B5JZzWwzT1IOlfSeyJiEBgoJg8OM99xkqYXJxrfDSwB7o+I\nl0e6zdHI4S9ZcQ7gDGAG8Dsq71P/AziwxiLXA0cWh87/nVj1xcCZVMLxWeCteSOiD7iQyknHl4Bn\ngM812PLfA9+QtA34Z4Y/2fZAsc57gCsi4q5i+tVUjhruKpb/JXBcg9utdgrwhKTtxTo/XeNcxzTg\np8A2YB2V8xrzhpkvS/If8zDLk/f8Zply+M0y5fCbZcrhN8tUV2/skeSzi2YdFhFDP3cxrJb2/JJO\nkfSUpGckXdrKusysu5q+1CdpDPAbYA6wEXgImBcR6xPLeM9v1mHd2PPPAp6JiGcj4k0qd6yN+NNa\nZlaOVsI/mT1v4NhYTNuDpPnFvdd9LWzLzNqs4yf8invQl4EP+816SSt7/k3seffW+4tpZrYXaCX8\nDwFHSPqApH2p/PGJVe1py8w6renD/ojYKWkBcCcwBrghIp5oW2dm1lFdvavP7/nNOq8rH/Ixs72X\nw2+WKYffLFMOv1mmHH6zTDn8Zply+M0y5fCbZcrhN8uUw2+WKYffLFMOv1mmHH6zTDn8Zply+M0y\n5fCbZcrhN8uUw2+WKYffLFMOv1mmHH6zTDn8Zply+M0y5fCbZcrhN8uUw2+WKYffLFMOv1mmHH6z\nTDU9RLftHcaMGZOsH3jggR3d/oIFC2rW9ttvv+Sy06dPT9YvuuiiZP2KK66oWZs3b15y2ddffz1Z\nv/zyy5P1yy67LFnvBS2FX9IGYBuwC9gZETPb0ZSZdV479vwnRMSLbViPmXWR3/ObZarV8Adwt6Rf\nSZo/3AyS5kvqk9TX4rbMrI1aPeyfHRGbJL0X+JmkX0fE6uoZImIZsAxAUrS4PTNrk5b2/BGxqfi+\nBbgVmNWOpsys85oOv6SxksbtfgycBKxrV2Nm1lmtHPZPBG6VtHs9P4yIn7alq1Hm0EMPTdb33Xff\nZP1jH/tYsj579uyatfHjxyeX/eQnP5msl2njxo3J+pIlS5L1uXPn1qxt27Ytueyjjz6arD/wwAPJ\n+t6g6fBHxLPAh9vYi5l1kS/1mWXK4TfLlMNvlimH3yxTDr9ZphTRvQ/djdZP+M2YMSNZv/fee5P1\nTt9W26sGBweT9fPPPz9Z3759e9Pb7u/vT9ZfeumlZP2pp55qetudFhFqZD7v+c0y5fCbZcrhN8uU\nw2+WKYffLFMOv1mmHH6zTPk6fxtMmDAhWV+zZk2yPm3atHa201b1eh8YGEjWTzjhhJq1N998M7ls\nrp9/aJWv85tZksNvlimH3yxTDr9Zphx+s0w5/GaZcvjNMuUhuttg69atyfoll1ySrJ9xxhnJ+iOP\nPJKs1/sT1ilr165N1ufMmZOs79ixI1k/6qijatYuvvji5LLWWd7zm2XK4TfLlMNvlimH3yxTDr9Z\nphx+s0w5/GaZ8v38PeCAAw5I1usNJ7106dKatQsuuCC57Lnnnpusr1y5Mlm33tO2+/kl3SBpi6R1\nVdMmSPqZpKeL7we10qyZdV8jh/3fA04ZMu1S4J6IOAK4p3huZnuRuuGPiNXA0M+vngUsLx4vB85u\nc19m1mHNfrZ/YkTsHuzsBWBirRklzQfmN7kdM+uQlm/siYhInciLiGXAMvAJP7Ne0uylvs2SJgEU\n37e0ryUz64Zmw78KOK94fB5wW3vaMbNuqXvYL2kl8HHgYEkbga8DlwM3S7oAeA74VCebHO1eeeWV\nlpZ/+eWXm172wgsvTNZvuummZH1wcLDpbVu56oY/IubVKJ3Y5l7MrIv88V6zTDn8Zply+M0y5fCb\nZcrhN8uUb+kdBcaOHVuzdvvttyeXPf7445P1U089NVm/6667knXrPg/RbWZJDr9Zphx+s0w5/GaZ\ncvjNMuXwm2XK4TfLlK/zj3KHH354sv7www8n6wMDA8n6fffdl6z39fXVrF177bXJZbv5szma+Dq/\nmSU5/GaZcvjNMuXwm2XK4TfLlMNvlimH3yxTvs6fublz5ybrN954Y7I+bty4pre9cOHCZH3FihXJ\nen9/f7KeK1/nN7Mkh98sUw6/WaYcfrNMOfxmmXL4zTLl8Jtlytf5Lenoo49O1q+66qpk/cQTmx/M\neenSpcn64sWLk/VNmzY1ve29Wduu80u6QdIWSeuqpi2StEnS2uLrtFaaNbPua+Sw/3vAKcNM//eI\nmFF8/aS9bZlZp9UNf0SsBrZ2oRcz66JWTvh9UdJjxduCg2rNJGm+pD5Jtf+Ym5l1XbPh/w4wDZgB\n9ANX1poxIpZFxMyImNnktsysA5oKf0RsjohdETEIXAfMam9bZtZpTYVf0qSqp3OBdbXmNbPeVPc6\nv6SVwMeBg4HNwNeL5zOAADYAn4+IujdX+zr/6DN+/Phk/cwzz6xZq/e3AqT05ep77703WZ8zZ06y\nPlo1ep1/nwZWNG+YydePuCMz6yn+eK9Zphx+s0w5/GaZcvjNMuXwm2XKt/Raad54441kfZ990hej\ndu7cmayffPLJNWv3339/ctm9mf90t5klOfxmmXL4zTLl8JtlyuE3y5TDb5Yph98sU3Xv6rO8HXPM\nMcn6Oeeck6wfe+yxNWv1ruPXs379+mR99erVLa1/tPOe3yxTDr9Zphx+s0w5/GaZcvjNMuXwm2XK\n4TfLlK/zj3LTp09P1hcsWJCsf+ITn0jWDznkkBH31Khdu3Yl6/396b8WPzg42M52Rh3v+c0y5fCb\nZcrhN8uUw2+WKYffLFMOv1mmHH6zTNW9zi9pCrACmEhlSO5lEXG1pAnATcBUKsN0fyoiXupcq/mq\ndy193rzhBlKuqHcdf+rUqc201BZ9fX3J+uLFi5P1VatWtbOd7DSy598JfDkijgQ+Clwk6UjgUuCe\niDgCuKd4bmZ7ibrhj4j+iHi4eLwNeBKYDJwFLC9mWw6c3akmzaz9RvSeX9JU4CPAGmBiROz+fOUL\nVN4WmNleouHP9kvaH7gF+FJEvCL9/3BgERG1xuGTNB+Y32qjZtZeDe35Jb2TSvB/EBE/LiZvljSp\nqE8Ctgy3bEQsi4iZETGzHQ2bWXvUDb8qu/jrgScj4qqq0irgvOLxecBt7W/PzDql7hDdkmYDPwce\nB3bfI7mQyvv+m4FDgeeoXOrbWmddWQ7RPXFi+nTIkUcemaxfc801yfqHPvShEffULmvWrEnWv/3t\nb9es3XZben/hW3Kb0+gQ3XXf80fEL4BaKztxJE2ZWe/wJ/zMMuXwm2XK4TfLlMNvlimH3yxTDr9Z\npvynuxs0YcKEmrWlS5cml50xY0ayPm3atKZ6aocHH3wwWb/yyiuT9TvvvDNZf+2110bck3WH9/xm\nmXL4zTLl8JtlyuE3y5TDb5Yph98sUw6/Waayuc5/3HHHJeuXXHJJsj5r1qyatcmTJzfVU7u8+uqr\nNWtLlixJLvvNb34zWd+xY0dTPVnv857fLFMOv1mmHH6zTDn8Zply+M0y5fCbZcrhN8tUNtf5586d\n21K9FevXr0/W77jjjmR9586dyXrqnvuBgYHkspYv7/nNMuXwm2XK4TfLlMNvlimH3yxTDr9Zphx+\ns0wpItIzSFOAFcBEIIBlEXG1pEXAhcAfilkXRsRP6qwrvTEza1lEqJH5Ggn/JGBSRDwsaRzwK+Bs\n4FPA9oi4otGmHH6zzms0/HU/4RcR/UB/8XibpCeBcv90jZm1bETv+SVNBT4CrCkmfVHSY5JukHRQ\njWXmS+qT1NdSp2bWVnUP+9+aUdofeABYHBE/ljQReJHKeYB/ofLW4Pw66/Bhv1mHte09P4CkdwJ3\nAHdGxFXD1KcCd0TE0XXW4/CbdVij4a972C9JwPXAk9XBL04E7jYXWDfSJs2sPI2c7Z8N/Bx4HBgs\nJi8E5gEzqBz2bwA+X5wcTK3Le36zDmvrYX+7OPxmnde2w34zG50cfrNMOfxmmXL4zTLl8JtlyuE3\ny5TDb5Yph98sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y1S3h+h+EXiu6vnBxbRe1Ku99Wpf4N6a\n1c7eDmt0xq7ez/+2jUt9ETGztAYSerW3Xu0L3FuzyurNh/1mmXL4zTJVdviXlbz9lF7trVf7AvfW\nrFJ6K/U9v5mVp+w9v5mVxOE3y1Qp4Zd0iqSnJD0j6dIyeqhF0gZJj0taW/b4gsUYiFskrauaNkHS\nzyQ9XXwfdozEknpbJGlT8dqtlXRaSb1NkXSfpPWSnpB0cTG91Ncu0Vcpr1vX3/NLGgP8BpgDbAQe\nAuZFxPquNlKDpA3AzIgo/QMhkv4C2A6s2D0UmqRvAVsj4vLiF+dBEfGVHultESMctr1DvdUaVv5z\nlPjatXO4+3YoY88/C3gmIp6NiDeBHwFnldBHz4uI1cDWIZPPApYXj5dT+eHpuhq99YSI6I+Ih4vH\n24Ddw8qX+tol+ipFGeGfDDxf9XwjJb4Awwjgbkm/kjS/7GaGMbFqWLQXgIllNjOMusO2d9OQYeV7\n5rVrZrj7dvMJv7ebHREzgFOBi4rD254UlfdsvXSt9jvANCpjOPYDV5bZTDGs/C3AlyLilepama/d\nMH2V8rqVEf5NwJSq5+8vpvWEiNhUfN8C3ErlbUov2bx7hOTi+5aS+3lLRGyOiF0RMQhcR4mvXTGs\n/C3ADyLix8Xk0l+74foq63UrI/wPAUdI+oCkfYFPA6tK6ONtJI0tTsQgaSxwEr039Pgq4Lzi8XnA\nbSX2sodeGba91rDylPza9dxw9xHR9S/gNCpn/H8LfLWMHmr0NQ14tPh6ouzegJVUDgP/SOXcyAXA\nu4F7gKeBu4EJPdTb96kM5f4YlaBNKqm32VQO6R8D1hZfp5X92iX6KuV188d7zTLlE35mmXL4zTLl\n8JtlyuE3y5TDb5Yph98sUw6/Wab+D/SyD4qiDGDaAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x125c9eb38>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP8AAAEICAYAAACQ6CLfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEk1JREFUeJzt3XuwXXV5xvHvA0oEQiIXDZFwEQ3TUQcOEmL+YEqUSjFg\nE0tBIjVYnQZbsehYK2qA1MJYGaATbEUCpAaJSYCABEeLShB0qpYjRC6CEjFAwiEh3HICCpK8/WOv\n6M7x7N/eZ9/Wzvk9n5k9Z5/1rst7NnlY970UEZhZfnYpuwEzK4fDb5Yph98sUw6/WaYcfrNMOfxm\nmXL4u0DSfEnXlt1HsySFpDc3MN4hxbivamIZyWklPSBp+kjna7U5/G0gaUvVa5uk31b9fnqbl/U1\nSRe0c547g4h4a0T8YKTTSTpO0kOSXpR0u6SDO9DeTsnhb4OIGLv9BTwGvLdq2JJu9tLMWne0krQf\ncCNwLrAP0A8sL7WpHuLwd89ukq6RNFhswk7ZXpD0BkkrJD0l6TeS/mm4GUiaC5wO/EuxVXFLMXyt\npM9Iuhd4QdKrhm6qD91ikHSSpNWSnpP0v5IOb+SPkHSipHskbZb0uKT5w4z2YUlPSBqQ9M9V0+4i\n6RxJv5b0tKTrJO3T4HLXSvqL4v1USf1FDxskXVpjsr8GHoiI6yPid8B84AhJf9bIMkc7h797/gpY\nBrwWWAn8J1QCAdwC/Bw4ADgO+ISkvxw6g4hYCCwBLiq2Kt5bVZ4NnAi8NiJeSTUi6UhgEXAmsC9w\nBbBS0pgG/o4XgDnF33Ei8A+SZg0Z553AZOB44DPbQwt8HJgFHAu8AXgW+K8GljnUAmBBRIwD3gRc\nV2O8t1L5XAGIiBeANcXw7Dn83fOjiPh2RGwFvg4cUQw/GnhdRHwhIl6OiEeAK4HTRjj/yyLi8Yj4\nbQPjzgWuiIifRsTWiFgMvARMqzdhRPwgIu6LiG0RcS+wlEqYq/1rRLwQEfcB/03lf0wAHwU+HxHr\nIuIlKmviv2liV+X3wJsl7RcRWyLiJzXGGws8P2TYZmCvES5vVHL4u+fJqvcvAq8p/tEfDLyh2Px+\nTtJzwOeACSOc/+MjGPdg4FNDlnkglbVxkqR3FAfOnpL0PJVA75fo5dGq+R4M3FS1zAeBrYz8b/0I\ncBjwkKS7JJ1UY7wtwLghw8YDgyNc3qjk8JfvceA3EfHaqtdeETGjxvi1bsMcOvxFYI+q3/cfsswL\nhyxzj4hY2kC/36Cy23JgRIwHvgpoyDgHVr0/CHiiarnvGbLc10TE+gaW+wcR8XBEzAZeD3wJuEHS\nnsOM+gB/3MKiGOdNxfDsOfzl+z9gsDhgt7ukXSW9TdLRNcbfABzawHxXAx8o5ncCO26aXwl8tFiL\nS9KexYG8RjaH9wKeiYjfSZoKfGCYcc6VtIektwJ/xx+PsH8VuHD76TZJr5M0s4Fl7kDS30p6XURs\nA54rBm8bZtSbgLdJOlnSa4DzgZ9HxEMjXeZo5PCXrDgGcBLQB/wG2ARcRWXzdDhXA28pNp2/mZj1\n2cB7qYTjdOAP40ZEP/D3VA46PkvlINiHGmz5H4EvSBoEzmP4g213FPO8Dbg4Ir5bDF9AZavhu8X0\nPwHe0eByq50APCBpSzHP04Y71hERTwEnAxdS+TunMvJjKaOW/GUeZnnymt8sUw6/WaYcfrNMOfxm\nmerqTSCSfHTRrMMiYuh1F8Nqac0v6QRJv5S0RtI5rczLzLqr6VN9knYFfgW8G1gH3AXMjohfJKbx\nmt+sw7qx5p8KrImIRyLiZSp3rI34ai0zK0cr4T+AHW/gWFcM24GkucW91/0tLMvM2qzjB/yKe9AX\ngjf7zXpJK2v+9ex499akYpiZ7QRaCf9dwGRJb5S0G5UbJla2py0z67SmN/sj4hVJZwG3ArsCiyLC\n90mb7SS6elef9/nNOq8rF/mY2c7L4TfLlMNvlimH3yxTDr9Zphx+s0w5/GaZcvjNMuXwm2XK4TfL\nlMNvlimH3yxTDr9Zphx+s0w5/GaZcvjNMuXwm2XK4TfLlMNvlimH3yxTDr9Zprr6iG4bfY466qhk\n/ayzzqpZmzNnTnLaa665Jln/8pe/nKzffffdyXruvOY3y5TDb5Yph98sUw6/WaYcfrNMOfxmmXL4\nzTLlp/RaUl9fX7K+atWqZH3cuHHtbGcHzz//fLK+7777dmzZvazRp/S2dJGPpLXAILAVeCUiprQy\nPzPrnnZc4ffOiNjUhvmYWRd5n98sU62GP4DvS/qZpLnDjSBprqR+Sf0tLsvM2qjVzf5jImK9pNcD\n35P0UETcWT1CRCwEFoIP+Jn1kpbW/BGxvvi5EbgJmNqOpsys85oOv6Q9Je21/T1wPHB/uxozs85q\nZbN/AnCTpO3z+UZE/E9burKumTo1vbG2YsWKZH38+PHJeuo6ksHBweS0L7/8crJe7zz+tGnTatbq\n3etfb9mjQdPhj4hHgCPa2IuZdZFP9ZllyuE3y5TDb5Yph98sUw6/WaZ8S+8osMcee9Ssvf3tb09O\ne+211ybrkyZNStaLU701pf591TvddtFFFyXry5YtS9ZTvc2bNy857Re/+MVkvZc1ekuv1/xmmXL4\nzTLl8JtlyuE3y5TDb5Yph98sUw6/Wab8iO5R4IorrqhZmz17dhc7GZl61yCMHTs2Wb/jjjuS9enT\np9esHX744clpc+A1v1mmHH6zTDn8Zply+M0y5fCbZcrhN8uUw2+WKZ/n3wkcddRRyfqJJ55Ys1bv\nfvt66p1Lv+WWW5L1iy++uGbtiSeeSE57zz33JOvPPvtssv6ud72rZq3Vz2U08JrfLFMOv1mmHH6z\nTDn8Zply+M0y5fCbZcrhN8uUv7e/B/T19SXrq1atStbHjRvX9LK/853vJOv1vg/g2GOPTdZT981f\nddVVyWmfeuqpZL2erVu31qy9+OKLyWnr/V31njlQprZ9b7+kRZI2Srq/atg+kr4n6eHi596tNGtm\n3dfIZv/XgBOGDDsHuC0iJgO3Fb+b2U6kbvgj4k7gmSGDZwKLi/eLgVlt7svMOqzZa/snRMRA8f5J\nYEKtESXNBeY2uRwz65CWb+yJiEgdyIuIhcBC8AE/s17S7Km+DZImAhQ/N7avJTPrhmbDvxI4o3h/\nBnBze9oxs26pe55f0lJgOrAfsAE4H/gmcB1wEPAocGpEDD0oONy8stzsP+yww5L1888/P1k/7bTT\nkvVNmzbVrA0MDNSsAVxwwQXJ+g033JCs97LUef56/+6XL1+erJ9++ulN9dQNjZ7nr7vPHxG1rvI4\nbkQdmVlP8eW9Zply+M0y5fCbZcrhN8uUw2+WKX91dxuMGTMmWU99fTXAjBkzkvXBwcFkfc6cOTVr\n/f39yWl33333ZD1XBx10UNktdJzX/GaZcvjNMuXwm2XK4TfLlMNvlimH3yxTDr9Zpnyevw2OPPLI\nZL3eefx6Zs6cmazXe4y22XC85jfLlMNvlimH3yxTDr9Zphx+s0w5/GaZcvjNMuXz/G1w6aWXJutS\n+puU652n93n85uyyS+1127Zt27rYSW/ymt8sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y5TP8zfo\npJNOqlnr6+tLTlvvcdArV65sqidLS53Lr/ffZPXq1e1up+fUXfNLWiRpo6T7q4bNl7Re0uri1dq3\nVZhZ1zWy2f814IRhhv9HRPQVr2+3ty0z67S64Y+IO4FnutCLmXVRKwf8Pi7p3mK3YO9aI0maK6lf\nUvqhcWbWVc2G/3LgUKAPGAAuqTViRCyMiCkRMaXJZZlZBzQV/ojYEBFbI2IbcCUwtb1tmVmnNRV+\nSROrfn0fcH+tcc2sN9U9zy9pKTAd2E/SOuB8YLqkPiCAtcCZHeyxJ6SeY7/bbrslp924cWOyvnz5\n8qZ6Gu3GjBmTrM+fP7/pea9atSpZ/+xnP9v0vHcWdcMfEbOHGXx1B3oxsy7y5b1mmXL4zTLl8Jtl\nyuE3y5TDb5Yp39LbBS+99FKyPjAw0KVOeku9U3nz5s1L1j/96U8n6+vWratZu+SSmhelArBly5Zk\nfTTwmt8sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y5TP83dBzl/Nnfpa83rn6d///vcn6zfffHOy\nfvLJJyfrufOa3yxTDr9Zphx+s0w5/GaZcvjNMuXwm2XK4TfLlM/zN0hSUzWAWbNmJetnn312Uz31\ngk9+8pPJ+rnnnluzNn78+OS0S5YsSdbnzJmTrFua1/xmmXL4zTLl8JtlyuE3y5TDb5Yph98sUw6/\nWaYaeUT3gcA1wAQqj+ReGBELJO0DLAcOofKY7lMj4tnOtVquiGiqBrD//vsn65dddlmyvmjRomT9\n6aefrlmbNm1actoPfvCDyfoRRxyRrE+aNClZf+yxx2rWbr311uS0X/nKV5J1a00ja/5XgE9FxFuA\nacDHJL0FOAe4LSImA7cVv5vZTqJu+CNiICLuLt4PAg8CBwAzgcXFaIuB9GVsZtZTRrTPL+kQ4Ejg\np8CEiNj+nKknqewWmNlOouFr+yWNBVYAn4iIzdXXs0dESBp2x1fSXGBuq42aWXs1tOaX9GoqwV8S\nETcWgzdImljUJwIbh5s2IhZGxJSImNKOhs2sPeqGX5VV/NXAgxFxaVVpJXBG8f4MIP1VqmbWU1Tv\nNJWkY4AfAvcB24rBn6Oy338dcBDwKJVTfc/UmVd6YT3slFNOqVlbunRpR5e9YcOGZH3z5s01a5Mn\nT253Ozv48Y9/nKzffvvtNWvnnXdeu9sxICLS95gX6u7zR8SPgFozO24kTZlZ7/AVfmaZcvjNMuXw\nm2XK4TfLlMNvlimH3yxTdc/zt3VhO/F5/tStq9dff31y2qOPPrqlZdf7avBW/humbgcGWLZsWbK+\nM3/t+GjV6Hl+r/nNMuXwm2XK4TfLlMNvlimH3yxTDr9Zphx+s0z5PH8bTJw4MVk/88wzk/V58+Yl\n662c51+wYEFy2ssvvzxZX7NmTbJuvcfn+c0syeE3y5TDb5Yph98sUw6/WaYcfrNMOfxmmfJ5frNR\nxuf5zSzJ4TfLlMNvlimH3yxTDr9Zphx+s0w5/GaZqht+SQdKul3SLyQ9IOnsYvh8SeslrS5eMzrf\nrpm1S92LfCRNBCZGxN2S9gJ+BswCTgW2RMTFDS/MF/mYdVyjF/m8qoEZDQADxftBSQ8CB7TWnpmV\nbUT7/JIOAY4EfloM+rikeyUtkrR3jWnmSuqX1N9Sp2bWVg1f2y9pLHAHcGFE3ChpArAJCODfqOwa\nfLjOPLzZb9ZhjW72NxR+Sa8GvgXcGhGXDlM/BPhWRLytznwcfrMOa9uNPap8dezVwIPVwS8OBG73\nPuD+kTZpZuVp5Gj/McAPgfuAbcXgzwGzgT4qm/1rgTOLg4OpeXnNb9Zhbd3sbxeH36zzfD+/mSU5\n/GaZcvjNMuXwm2XK4TfLlMNvlimH3yxTDr9Zphx+s0w5/GaZcvjNMuXwm2XK4TfLlMNvlqm6X+DZ\nZpuAR6t+368Y1ot6tbde7QvcW7Pa2dvBjY7Y1fv5/2ThUn9ETCmtgYRe7a1X+wL31qyyevNmv1mm\nHH6zTJUd/oUlLz+lV3vr1b7AvTWrlN5K3ec3s/KUveY3s5I4/GaZKiX8kk6Q9EtJaySdU0YPtUha\nK+m+4rHjpT5fsHgG4kZJ91cN20fS9yQ9XPwc9hmJJfXWE49tTzxWvtTPrtced9/1fX5JuwK/At4N\nrAPuAmZHxC+62kgNktYCUyKi9AtCJP05sAW4Zvuj0CRdBDwTEf9e/I9z74j4TI/0Np8RPra9Q73V\neqz8hyjxs2vn4+7boYw1/1RgTUQ8EhEvA8uAmSX00fMi4k7gmSGDZwKLi/eLqfzj6boavfWEiBiI\niLuL94PA9sfKl/rZJfoqRRnhPwB4vOr3dZT4AQwjgO9L+pmkuWU3M4wJVY9FexKYUGYzw6j72PZu\nGvJY+Z757Jp53H27+YDfnzomIvqA9wAfKzZve1JU9tl66Vzt5cChVJ7hOABcUmYzxWPlVwCfiIjN\n1bUyP7th+irlcysj/OuBA6t+n1QM6wkRsb74uRG4icpuSi/ZsP0JycXPjSX38wcRsSEitkbENuBK\nSvzsisfKrwCWRMSNxeDSP7vh+irrcysj/HcBkyW9UdJuwGnAyhL6+BOS9iwOxCBpT+B4eu/R4yuB\nM4r3ZwA3l9jLDnrlse21HitPyZ9dzz3uPiK6/gJmUDni/2vg82X0UKOvQ4GfF68Hyu4NWEplM/D3\nVI6NfATYF7gNeBj4PrBPD/X2dSqPcr+XStAmltTbMVQ26e8FVhevGWV/dom+SvncfHmvWaZ8wM8s\nUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y9T/A/k0FxVjH2hNAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x125e7fd68>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP8AAAEICAYAAACQ6CLfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEQBJREFUeJzt3X2wVPV9x/H3B42ND6ighKJSxASdMQ7iDJJoHSWDsQS1\nmpnqSLCSxvGamhqdSac6ZprQWDuaMenYaU2CDxWVmDhFlGS0CTKJNlOlgCXIg4pBUJCHGKWCWgzw\n7R97sMv17tm9u2f37L2/z2tm5+6e33n43jP3c3/nYXd/igjMLD1Dyi7AzMrh8JslyuE3S5TDb5Yo\nh98sUQ6/WaIc/g6QNEvSg2XX0SxJIekTDcx3fDbvgU1sI3dZSaskTe7veq02h78AknZWPfZKeq/q\n9YyCt3WfpL8vcp0DQUR8MiJ+2ezykr6R/XM5t8CyBjSHvwARcdi+B/AqcGHVtLmdrKWZXnewk/Rx\n4BJgc9m1dBOHv3MOknS/pB3ZIezEfQ2SjpE0T9JvJb0i6at9rUBSDzAD+JvsqOIn2fT1km6QtAJ4\nR9KBvQ/Vex8xSLpA0nJJ2yX9p6TxjfwSks6X9N+S3pb0mqRZfcz2JUmvS9os6a+rlh0i6UZJv5H0\nO0kPSxre4HbX7+u1JU2StDSrYauk79ZZ/F+AG4D3G9lWKhz+zvlT4EfAkcAC4J+hEgjgJ8CvgWOB\nKcD1kv6k9woiYjYwF/h2dlRxYVXzdOB84MiI2J1XiKTTgHuBq4GjgB8ACyT9QQO/xzvAFdnvcT7w\nl5Iu7jXPZ4BxwHnADVWH2tcCFwPnAMcAb1EJZn/dAdwREYcDHwcerjWjpEuAXRHxeBPbGdQc/s75\nVUQ8HhF7gAeAU7PppwMjIuJbEfF+RKwD7gIu6+f6/ykiXouI9xqYtwf4QUQsjog9ETEH2AV8ut6C\nEfHLiHg+IvZGxArgISphrvZ3EfFORDwP/CuVf0wAXwa+HhEbI2IXMAv4syZOVX4PfELS0RGxMyKe\n7WsmSUOBfwCu6+f6k+Dwd86WqufvAh/N/ujHAMdkh9/bJW0HbgJG9nP9r/Vj3jHA13ptczSV3jiX\npE9J+kV2ivI/VAJ9dE4tG6rWOwaYX7XNNcAe+v+7XgmcCLwgaYmkC2rMNwt4ICLW93P9SXD4y/ca\n8EpEHFn1GBoR02rMX+tjmL2nvwscUvX6D3tt85Ze2zwkIh5qoN4fUjltGR0RRwDfB9RrntFVz/8I\neL1qu5/rtd2PRsSmBrb7gYhYGxHTgY8BtwH/JunQPmadAnxV0hZJW7K6HpZ0Q3+2N1g5/OX7L2BH\ndsHuYEkHSDpF0uk15t8KnNDAepcDX8jWN5X9D83vAr6c9eKSdGh2IW9oA+sdCrwZEf8raRLwhT7m\n+VtJh0j6JPAXwI+z6d8HbpE0BkDSCEkXNbDN/Ui6XNKIiNgLbM8m7+1j1inAKcCE7PE6lesczVxn\nGHQc/pJl1wAuoPLH+QrwBnA3cESNRe4BTs4OnR/NWfV1wIVUwjED+GDeiFgKXEXlouNbwMvAFxss\n+RrgW5J2AN+g74ttT2XrXATcHhE/z6bfQeWo4efZ8s8Cn2pwu9WmAqsk7czWeVlf1zoi4ncRsWXf\ng8opxlsRsbOJbQ468pd5mKXJPb9Zohx+s0Q5/GaJcvjNEtXRD4FI8tVFszaLiN7vu+hTSz2/pKmS\nXpT0sqQbW1mXmXVW07f6JB0AvAR8FtgILAGmR8TqnGXc85u1WSd6/knAyxGxLiLep/KJtX6/W8vM\nytFK+I9l/w9wbMym7UdST/bZ66UtbMvMCtb2C37ZZ9Bngw/7zbpJKz3/Jvb/9NZx2TQzGwBaCf8S\nYJyksZIOovLlEwuKKcvM2q3pw/6I2C3pr4CfAQcA90bEqsIqM7O26uin+nzOb9Z+HXmTj5kNXA6/\nWaIcfrNEOfxmiXL4zRLl8JslyuE3S5TDb5Yoh98sUQ6/WaIcfrNEOfxmiXL4zRLl8JslyuE3S5TD\nb5Yoh98sUQ6/WaIcfrNEOfxmiXL4zRLl8JslyuE3S5TDb5Yoh98sUQ6/WaIcfrNEOfxmiXL4zRLV\n9BDdZt1uypQpNdvmzp2bu+w555yT2/7iiy82VVM3aSn8ktYDO4A9wO6ImFhEUWbWfkX0/J+JiDcK\nWI+ZdZDP+c0S1Wr4A3hS0jJJPX3NIKlH0lJJS1vclpkVqNXD/rMiYpOkjwELJb0QEU9XzxARs4HZ\nAJKixe2ZWUFa6vkjYlP2cxswH5hURFFm1n5Nh1/SoZKG7nsOnAesLKowM2uvVg77RwLzJe1bzw8j\n4t8LqaoNzj777Nz2o446Krd9/vz5RZZjHXD66afXbFuyZEkHK+lOTYc/ItYBpxZYi5l1kG/1mSXK\n4TdLlMNvliiH3yxRDr9ZopL5SO/kyZNz28eNG5fb7lt93WfIkPy+a+zYsTXbxowZk7tsdgt7UHPP\nb5Yoh98sUQ6/WaIcfrNEOfxmiXL4zRLl8JslKpn7/FdccUVu+zPPPNOhSqwoo0aNym2/6qqrarY9\n+OCDucu+8MILTdU0kLjnN0uUw2+WKIffLFEOv1miHH6zRDn8Zoly+M0Slcx9/nqf/baB5+677256\n2bVr1xZYycDkRJglyuE3S5TDb5Yoh98sUQ6/WaIcfrNEOfxmiRo09/nHjx+f2z5y5MgOVWKdcsQR\nRzS97MKFCwusZGCq2/NLulfSNkkrq6YNl7RQ0trs57D2lmlmRWvksP8+YGqvaTcCiyJiHLAoe21m\nA0jd8EfE08CbvSZfBMzJns8BLi64LjNrs2bP+UdGxObs+Rag5gm1pB6gp8ntmFmbtHzBLyJCUuS0\nzwZmA+TNZ2ad1eytvq2SRgFkP7cVV5KZdUKz4V8AzMyezwQeK6YcM+uUuof9kh4CJgNHS9oIfBO4\nFXhY0pXABuDSdhbZiGnTpuW2H3zwwR2qxIpS770ZY8eObXrdmzZtanrZwaJu+CNieo2mKQXXYmYd\n5Lf3miXK4TdLlMNvliiH3yxRDr9ZogbNR3pPOumklpZftWpVQZVYUW6//fbc9nq3Al966aWabTt2\n7GiqpsHEPb9Zohx+s0Q5/GaJcvjNEuXwmyXK4TdLlMNvlqhBc5+/VUuWLCm7hAHp8MMPz22fOrX3\nd7/+v8svvzx32fPOO6+pmva5+eaba7Zt3769pXUPBu75zRLl8JslyuE3S5TDb5Yoh98sUQ6/WaIc\nfrNE+T5/Zvjw4aVt+9RTT81tl5Tbfu6559ZsO+6443KXPeigg3LbZ8yYkds+ZEh+//Hee+/VbFu8\neHHusrt27cptP/DA/D/fZcuW5banzj2/WaIcfrNEOfxmiXL4zRLl8JslyuE3S5TDb5YoRUTnNia1\nbWN33nlnbvvVV1+d217v892vvvpqv2tq1Pjx43Pb693n3717d822d999N3fZ1atX57bXuxe/dOnS\n3PannnqqZtvWrVtzl924cWNu+7Bhw3Lb672HYbCKiPw/mEzdnl/SvZK2SVpZNW2WpE2SlmePaa0U\na2ad18hh/31AX1/H8o8RMSF7PF5sWWbWbnXDHxFPA292oBYz66BWLvhdK2lFdlpQ8+RLUo+kpZLy\nTw7NrKOaDf/3gBOACcBm4Du1ZoyI2RExMSImNrktM2uDpsIfEVsjYk9E7AXuAiYVW5aZtVtT4Zc0\nqurl54GVteY1s+5U9/P8kh4CJgNHS9oIfBOYLGkCEMB6IP8megdcc801ue0bNmzIbT/zzDOLLKdf\n6r2H4NFHH81tX7NmTc22Z599tqmaOqGnpye3fcSIEbnt69atK7Kc5NQNf0RM72PyPW2oxcw6yG/v\nNUuUw2+WKIffLFEOv1miHH6zRCXz1d233XZb2SVYL1OmTGlp+Xnz5hVUSZrc85slyuE3S5TDb5Yo\nh98sUQ6/WaIcfrNEOfxmiUrmPr8NPvPnzy+7hAHNPb9Zohx+s0Q5/GaJcvjNEuXwmyXK4TdLlMNv\nliiH3yxRDr9Zohx+s0Q5/GaJcvjNEuXwmyXK4TdLlMNvlqhGhugeDdwPjKQyJPfsiLhD0nDgx8Dx\nVIbpvjQi3mpfqZYaSbntJ554Ym57Nw9P3g0a6fl3A1+LiJOBTwNfkXQycCOwKCLGAYuy12Y2QNQN\nf0Rsjojnsuc7gDXAscBFwJxstjnAxe0q0syK169zfknHA6cBi4GREbE5a9pC5bTAzAaIhr/DT9Jh\nwDzg+oh4u/p8LCJCUtRYrgfoabVQMytWQz2/pI9QCf7ciHgkm7xV0qisfRSwra9lI2J2REyMiIlF\nFGxmxagbflW6+HuANRHx3aqmBcDM7PlM4LHiyzOzdmnksP+PgT8Hnpe0PJt2E3Ar8LCkK4ENwKXt\nKdFSFdHnmeQHhgzx21RaUTf8EfEroNYN19YGWDez0vhfp1miHH6zRDn8Zoly+M0S5fCbJcrhN0uU\nh+i2AeuMM87Ibb/vvvs6U8gA5Z7fLFEOv1miHH6zRDn8Zoly+M0S5fCbJcrhN0uU7/Nb16r31d3W\nGvf8Zoly+M0S5fCbJcrhN0uUw2+WKIffLFEOv1mifJ/fSvPEE0/ktl9yySUdqiRN7vnNEuXwmyXK\n4TdLlMNvliiH3yxRDr9Zohx+s0Sp3hjokkYD9wMjgQBmR8QdkmYBVwG/zWa9KSIer7Ou/I2ZWcsi\noqEvQmgk/KOAURHxnKShwDLgYuBSYGdE3N5oUQ6/Wfs1Gv667/CLiM3A5uz5DklrgGNbK8/Mytav\nc35JxwOnAYuzSddKWiHpXknDaizTI2mppKUtVWpmhap72P/BjNJhwFPALRHxiKSRwBtUrgPcTOXU\n4Et11uHDfrM2K+ycH0DSR4CfAj+LiO/20X488NOIOKXOehx+szZrNPx1D/tV+QrVe4A11cHPLgTu\n83lgZX+LNLPyNHK1/yzgP4Dngb3Z5JuA6cAEKof964Grs4uDeetyz2/WZoUe9hfF4Tdrv8IO+81s\ncHL4zRLl8JslyuE3S5TDb5Yoh98sUQ6/WaIcfrNEOfxmiXL4zRLl8JslyuE3S5TDb5Yoh98sUZ0e\novsNYEPV66Ozad2oW2vr1rrAtTWryNrGNDpjRz/P/6GNS0sjYmJpBeTo1tq6tS5wbc0qqzYf9psl\nyuE3S1TZ4Z9d8vbzdGtt3VoXuLZmlVJbqef8Zlaesnt+MyuJw2+WqFLCL2mqpBclvSzpxjJqqEXS\neknPS1pe9viC2RiI2yStrJo2XNJCSWuzn32OkVhSbbMkbcr23XJJ00qqbbSkX0haLWmVpOuy6aXu\nu5y6StlvHT/nl3QA8BLwWWAjsASYHhGrO1pIDZLWAxMjovQ3hEg6G9gJ3L9vKDRJ3wbejIhbs3+c\nwyLihi6pbRb9HLa9TbXVGlb+i5S474oc7r4IZfT8k4CXI2JdRLwP/Ai4qIQ6ul5EPA282WvyRcCc\n7PkcKn88HVejtq4QEZsj4rns+Q5g37Dype67nLpKUUb4jwVeq3q9kRJ3QB8CeFLSMkk9ZRfTh5FV\nw6JtAUaWWUwf6g7b3km9hpXvmn3XzHD3RfMFvw87KyImAJ8DvpId3nalqJyzddO92u8BJ1AZw3Ez\n8J0yi8mGlZ8HXB8Rb1e3lbnv+qirlP1WRvg3AaOrXh+XTesKEbEp+7kNmE/lNKWbbN03QnL2c1vJ\n9XwgIrZGxJ6I2AvcRYn7LhtWfh4wNyIeySaXvu/6qqus/VZG+JcA4ySNlXQQcBmwoIQ6PkTSodmF\nGCQdCpxH9w09vgCYmT2fCTxWYi376ZZh22sNK0/J+67rhruPiI4/gGlUrvj/Bvh6GTXUqOsE4NfZ\nY1XZtQEPUTkM/D2VayNXAkcBi4C1wJPA8C6q7QEqQ7mvoBK0USXVdhaVQ/oVwPLsMa3sfZdTVyn7\nzW/vNUuUL/iZJcrhN0uUw2+WKIffLFEOv1miHH6zRDn8Zon6P7b2eWhzcmVkAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x125f8dc18>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP8AAAEICAYAAACQ6CLfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEGZJREFUeJzt3X+wVOV9x/H3BxNDRSsgBsEQiI39w2QMzhCSmTItGQlF\nxUJmOozEtNg4uaZNaZyJUxxjG0zrjI2JLTWa5FrlR2pInKCBZJwmyvijDqMVLUFUiAahcoOXWKWC\nsSFyv/1jD3a53j27d3fPnr33+bxm7tzd85x9ni87fO75sXvOo4jAzNIzpuwCzKwcDr9Zohx+s0Q5\n/GaJcvjNEuXwmyXK4e8ASSsl/WvZdTRLUkh6fwPrzcjWfUcTY+S+VtLTkuYOt1+rzeFvA0mHq34G\nJL1R9fzSNo+1RtLft7PPkSAiPhARDw7nNZJOlPR9SXuyPyxzi6luZHL42yAiTj72A/wXcHHVsjs7\nWUszW91R7hHgU8BLZRfSbRz+zjlR0jpJh7Jd2FnHGiRNlbRB0i8lvSDpr4bqQFIPcCnw19lexQ+z\n5XskrZC0HXhd0jsG76oP3mOQtFDSNkkHJW2RdG4j/whJF0n6T0mvSXpR0sohVvu0pF9I2i/pqqrX\njpF0taSfS/pvSXdJmtjguHskzcsez5a0NauhX9JNQ70mIo5ExD9FxCPA0UbGSYnD3zl/BHwXGA9s\nAr4OlUAAPwR+CpwJnA9cKekPB3cQEb3AncBXsr2Ki6ualwIXAeMj4s28QiSdB9wBXAGcBnwL2CTp\nXQ38O14H/jT7d1wE/LmkxYPW+RhwNjAfWHEstMByYDHwB8BU4FXglgbGHGwVsCoifhv4HeCuJvpI\nnsPfOY9ExL0RcRT4NvChbPmHgdMj4svZlmo3cBtwyTD7/+eIeDEi3mhg3R7gWxHxWEQcjYi1wK+B\nj9Z7YUQ8GBFPRcRARGwH1lMJc7XrIuL1iHgKWE3lDxPAZ4EvRsS+iPg1sBL44yYOVX4DvF/SpIg4\nHBGPDvP1hsPfSdXHnL8Cxmb/6acDU7Pd74OSDgLXAJOH2f+Lw1h3OvCFQWNOo7I1ziXpI5IeyA5R\n/odKoCfl1LK3qt/pwD1VYz5LZXd8uP/Wy4HfBXZKelzSwmG+3nD4u8GLwAsRMb7q55SIuLDG+rUu\nwxy8/FfASVXPzxg05vWDxjwpItY3UO93qBy2TIuIU4FvAhq0zrSqx+8FflE17gWDxh0bEX0NjPuW\niHguIpYC7wb+Afi+pHHD6cMc/m7wH8Ch7ITdb0k6QdIHJX24xvr9wFkN9LsN+GTW3wKO3zW/Dfhs\nthWXpHHZibxTGuj3FOCViPhfSbOBTw6xzt9IOknSB4A/A76XLf8mcL2k6QCSTpe0qIExjyPpU5JO\nj4gB4GC2eKDGuu+SNDZ7eqKksZIG/7FKksNfsuwcwEJgJvAC8DLwL8CpNV5yO3BOtuv8g5yuPw9c\nTCUclwJvrRsRW4HPUDnp+CrwPHBZgyX/BfBlSYeAv2Xok20PZX1uBr4aET/Jlq+istfwk+z1jwIf\naXDcaguApyUdzvq8JOdcxy7gDSonU3+cPZ7exJijjnwzD7M0ectvliiH3yxRDr9Zohx+s0R19CIQ\nST67aFawiGjoo8yWtvySFkjaJel5SVe30peZdVbTH/VJOgH4GfBxYB/wOLA0Ip7JeY23/GYF68SW\nfzbwfETsjogjVK5YG/a3tcysHK2E/0yOv4BjX7bsOJJ6smuvt7Ywlpm1WeEn/LJr0HvBu/1m3aSV\nLX8fx1+99Z5smZmNAK2E/3HgbEnvk3QilZtPbGpPWWZWtKZ3+yPiTUl/SeVKqROAOyLi6bZVZmaF\n6uhVfT7mNyteR77kY2Yjl8NvliiH3yxRDr9Zohx+s0Q5/GaJcvjNEuXwmyXK4TdLlMNvliiH3yxR\nDr9Zohx+s0Q5/GaJcvjNEuXwmyXK4TdLlMNvliiH3yxRDr9Zohx+s0Q5/GaJcvjNEuXwmyXK4TdL\nlMNvliiH3yxRDr9Zohx+s0Q1PUW3WdGuvfba3Pbrrrsut33MmNrbtrlz5+a+9qGHHsptHw1aCr+k\nPcAh4CjwZkTMakdRZla8dmz5PxYRL7ehHzPrIB/zmyWq1fAHcL+kJyT1DLWCpB5JWyVtbXEsM2uj\nVnf750REn6R3A/dJ2hkRD1evEBG9QC+ApGhxPDNrk5a2/BHRl/0+ANwDzG5HUWZWvKbDL2mcpFOO\nPQbmAzvaVZiZFauV3f7JwD2SjvXznYj4t7ZUZUm47LLLcttXrFiR2z4wMND02BE+Am06/BGxG/hQ\nG2sxsw7yR31miXL4zRLl8JslyuE3S5TDb5YoX9JrpZk+fXpu+9ixYztUSZq85TdLlMNvliiH3yxR\nDr9Zohx+s0Q5/GaJcvjNEuXP+a1Q8+bNq9m2fPnylvreuXNnbvvChQtrtvX397c09mjgLb9Zohx+\ns0Q5/GaJcvjNEuXwmyXK4TdLlMNvlih/zm8tmTNnTm776tWra7adeuqpLY1944035rbv3bu3pf5H\nO2/5zRLl8JslyuE3S5TDb5Yoh98sUQ6/WaIcfrNE+XN+a8myZcty26dOndp03w8++GBu+7p165ru\n2xrY8ku6Q9IBSTuqlk2UdJ+k57LfE4ot08zarZHd/jXAgkHLrgY2R8TZwObsuZmNIHXDHxEPA68M\nWrwIWJs9XgssbnNdZlawZo/5J0fE/uzxS8DkWitK6gF6mhzHzArS8gm/iAhJkdPeC/QC5K1nZp3V\n7Ed9/ZKmAGS/D7SvJDPrhGbDvwk49hnPMmBje8oxs05RRP6euKT1wFxgEtAPfAn4AXAX8F5gL7Ak\nIgafFByqL+/2jzCTJk3Kba93//uBgYGabQcPHsx97ZIlS3LbH3jggdz2VEWEGlmv7jF/RCyt0XT+\nsCoys67ir/eaJcrhN0uUw2+WKIffLFEOv1mifElv4mbMmJHbvmHDhsLGvvnmm3Pb/VFesbzlN0uU\nw2+WKIffLFEOv1miHH6zRDn8Zoly+M0S5c/5E7dgweB7sx7v3HPPban/zZs312xbtWpVS31ba7zl\nN0uUw2+WKIffLFEOv1miHH6zRDn8Zoly+M0SVffW3W0dzLfu7rjFi/OnUVyzZk1u+7hx43Lbt2zZ\nktued/vterf9tuY0eutub/nNEuXwmyXK4TdLlMNvliiH3yxRDr9Zohx+s0T5ev5RIO/e+0Xedx9g\n9+7due3+LL971d3yS7pD0gFJO6qWrZTUJ2lb9nNhsWWaWbs1stu/Bhjqdi//GBEzs59721uWmRWt\nbvgj4mHglQ7UYmYd1MoJv+WStmeHBRNqrSSpR9JWSVtbGMvM2qzZ8H8DOAuYCewHvlZrxYjojYhZ\nETGrybHMrABNhT8i+iPiaEQMALcBs9tblpkVranwS5pS9fQTwI5a65pZd6r7Ob+k9cBcYJKkfcCX\ngLmSZgIB7AGuKLBGq2PFihU12wYGBgod+4Ybbii0fytO3fBHxNIhFt9eQC1m1kH+eq9Zohx+s0Q5\n/GaJcvjNEuXwmyXKl/SOADNnzsxtnz9/fmFjb9y4Mbd9165dhY1txfKW3yxRDr9Zohx+s0Q5/GaJ\ncvjNEuXwmyXK4TdLlKfoHgEOHDiQ2z5hQs27qNX16KOP5rZfcMEFue2HDx9uemwrhqfoNrNcDr9Z\nohx+s0Q5/GaJcvjNEuXwmyXK4TdLlK/nHwFOO+203PZWbs9966235rb7c/zRy1t+s0Q5/GaJcvjN\nEuXwmyXK4TdLlMNvliiH3yxRjUzRPQ1YB0ymMiV3b0SskjQR+B4wg8o03Usi4tXiSh29Vq9ends+\nZkxxf6O3bNlSWN/W3Rr5X/Um8IWIOAf4KPA5SecAVwObI+JsYHP23MxGiLrhj4j9EfFk9vgQ8Cxw\nJrAIWJutthZYXFSRZtZ+w9qflDQDOA94DJgcEfuzppeoHBaY2QjR8Hf7JZ0MbACujIjXpP+/TVhE\nRK3780nqAXpaLdTM2quhLb+kd1IJ/p0RcXe2uF/SlKx9CjDkXSYjojciZkXErHYUbGbtUTf8qmzi\nbweejYibqpo2Acuyx8uA/OlczayrNLLb/3vAnwBPSdqWLbsGuAG4S9LlwF5gSTEljnz1ptieN29e\nbnu9S3aPHDlSs+2WW27JfW1/f39uu41edcMfEY8Ate4Dfn57yzGzTvE3/MwS5fCbJcrhN0uUw2+W\nKIffLFEOv1mifOvuDhg/fnxu+xlnnNFS/319fTXbrrrqqpb6ttHLW36zRDn8Zoly+M0S5fCbJcrh\nN0uUw2+WKIffLFEOv1miHH6zRDn8Zoly+M0S5fCbJcrhN0uUw2+WKIffLFG+nr8Ddu7cmdteb5rs\nOXPmtLMcM8BbfrNkOfxmiXL4zRLl8JslyuE3S5TDb5Yoh98sUYqI/BWkacA6YDIQQG9ErJK0EvgM\n8Mts1Wsi4t46feUPZmYtiwg1sl4j4Z8CTImIJyWdAjwBLAaWAIcj4quNFuXwmxWv0fDX/YZfROwH\n9mePD0l6FjiztfLMrGzDOuaXNAM4D3gsW7Rc0nZJd0iaUOM1PZK2StraUqVm1lZ1d/vfWlE6GXgI\nuD4i7pY0GXiZynmAv6NyaPDpOn14t9+sYG075geQ9E7gR8CPI+KmIdpnAD+KiA/W6cfhNytYo+Gv\nu9svScDtwLPVwc9OBB7zCWDHcIs0s/I0crZ/DvDvwFPAQLb4GmApMJPKbv8e4Irs5GBeX97ymxWs\nrbv97eLwmxWvbbv9ZjY6OfxmiXL4zRLl8JslyuE3S5TDb5Yoh98sUQ6/WaIcfrNEOfxmiXL4zRLl\n8JslyuE3S5TDb5aoTk/R/TKwt+r5pGxZN+rW2rq1LnBtzWpnbdMbXbGj1/O/bXBpa0TMKq2AHN1a\nW7fWBa6tWWXV5t1+s0Q5/GaJKjv8vSWPn6dba+vWusC1NauU2ko95jez8pS95Tezkjj8ZokqJfyS\nFkjaJel5SVeXUUMtkvZIekrStrLnF8zmQDwgaUfVsomS7pP0XPZ7yDkSS6ptpaS+7L3bJunCkmqb\nJukBSc9IelrS57Plpb53OXWV8r51/Jhf0gnAz4CPA/uAx4GlEfFMRwupQdIeYFZElP6FEEm/DxwG\n1h2bCk3SV4BXIuKG7A/nhIhY0SW1rWSY07YXVFutaeUvo8T3rp3T3bdDGVv+2cDzEbE7Io4A3wUW\nlVBH14uIh4FXBi1eBKzNHq+l8p+n42rU1hUiYn9EPJk9PgQcm1a+1Pcup65SlBH+M4EXq57vo8Q3\nYAgB3C/pCUk9ZRczhMlV06K9BEwus5gh1J22vZMGTSvfNe9dM9Pdt5tP+L3dnIiYCVwAfC7bve1K\nUTlm66bPar8BnEVlDsf9wNfKLCabVn4DcGVEvFbdVuZ7N0RdpbxvZYS/D5hW9fw92bKuEBF92e8D\nwD1UDlO6Sf+xGZKz3wdKructEdEfEUcjYgC4jRLfu2xa+Q3AnRFxd7a49PduqLrKet/KCP/jwNmS\n3ifpROASYFMJdbyNpHHZiRgkjQPm031Tj28ClmWPlwEbS6zlON0ybXutaeUp+b3ruunuI6LjP8CF\nVM74/xz4Yhk11KjrLOCn2c/TZdcGrKeyG/gbKudGLgdOAzYDzwH3AxO7qLZvU5nKfTuVoE0pqbY5\nVHbptwPbsp8Ly37vcuoq5X3z13vNEuUTfmaJcvjNEuXwmyXK4TdLlMNvliiH3yxRDr9Zov4PJI5t\n+G9VwgkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x129cf7940>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP8AAAEICAYAAACQ6CLfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEddJREFUeJzt3XuwVfV5xvHvo0TxGiUgRYMalP6hjtWRYGbqiBmNGq9E\nW5Voq2iDsSmNMzrVIcbQqJ00k5gh05gE0IKRSDBoRMfRoBNv04aKliiKCV6wchG0SEWiIcLbP/bC\nHI5n//Y++7Y25/d8Zs6cc9a7Li8bHtbaa+21fooIzCw/O5XdgJmVw+E3y5TDb5Yph98sUw6/WaYc\nfrNMOfwdIGmqpDvK7qNRkkLSoXXMd3Ax76AGtpFcVtLzkk7o73qtOoe/BSS92+Nrq6T3evx+YYu3\nNUvSja1c544gIg6PiEf7u5ykv5P0UvF38aCk/dvQ3g7J4W+BiNhz2xfwP8CZPabN6WQvjex1B6ri\nSOFfgLOBIcCrwJ1l9tRNHP7O2UXS7ZI2FoewY7YVJO0vab6kNyW9Kukf+1qBpEnAhcA/FXuy+4rp\nKyRdI+lZYJOkQb0P1XsfMUg6Q9ISSRsk/YekI+v5Q0g6XdJ/S3pH0uuSpvYx26WSVktaI+nqHsvu\nJOlaSS9L+l9J8yQNqXO7KySdVPw8VtLiooe1km6ustgZwM8j4vmI2AzcABwv6ZB6tjnQOfydcxYw\nF9gHWAD8G1QCAdwH/AY4ADgRuFLSKb1XEBHTgTnAt4ujijN7lCcApwP7RMQHqUYkHQ3cBlwOfAL4\nMbBA0q51/Dk2AX9b/DlOB66QNL7XPJ8FRgMnA9dsCy0wGRgPjAP2B94GflDHNnubBkyLiL2BQ4B5\ndS6n4vsRDWxzwHH4O+fJiHggIrYAPwH+opj+aWBYRHwzIjZHxCvADOCCfq7/+xHxekS8V8e8k4Af\nR8SiiNgSEbOBPwCfqbVgRDwaEc9FxNaIeJbKYfS4XrP9c0RsiojngH+n8h8TwJeBr0XEyoj4AzAV\n+KsG3qr8EThU0tCIeDcifl1lvgeBv5Z0pKTdgOuBAHbv5/YGJIe/c97o8fPvgcHFP/qDgP2Lw+8N\nkjYAU4Dh/Vz/6/2Y9yDgql7bHEllb5wk6VhJvyreovwflUAPTfTyWo/1HgTc02Oby4At9P/Pehnw\n58CLkp6SdEZfM0XEw1T+g5kPrCi+NgIr+7m9AcnhL9/rwKsRsU+Pr70i4rQq81e7DbP39N+z/R7u\nz3pt86Ze29w9Iuo5GfZTKm9bRkbEx4Ef8afD6W1G9vj5QGB1j+1+vtd2B0fEqjq2+6GIWB4RE4D9\ngH8Ffi5pjyrz/iAiRkfEcCr/CQwClvZnewOVw1++/wI2FifsdpO0s6QjJH26yvxrgVF1rHcJ8MVi\nfaey/aH5DODLxV5ckvYoTuTtVcd69wLWR8T7ksYCX+xjnq9L2l3S4cBE4GfF9B8BN0k6CEDSMEln\n17HN7Ui6SNKwiNgKbCgmb+1jvsHFaylJBwLTqZwreLu/2xyIHP6SFecAzgCOonIp6i1gJvDxKovc\nChxWHDr/IrHqrwJnUgnHhcCH80bEYuBLVE46vg28BFxSZ8t/D3xT0kYq76H7Otn2WLHOR4DvRMQv\ni+nTqBw1/LJY/tfAsXVut6dTgeclvVus84Iq5zoGUzlSeZfKf7L/CXy9ge0NSPLDPMzy5D2/WaYc\nfrNMOfxmmXL4zTLV0ZtAJPnsolmbRUTvz130qak9v6RTJf22uGXy2mbWZWad1fClPkk7A78DPkfl\n45JPARMi4oXEMt7zm7VZJ/b8Y4GXIuKV4nbJuVTumzazHUAz4T+A7W/gWFlM246kScW914ub2JaZ\ntVjbT/gV96BPBx/2m3WTZvb8q9j+7q1PFtPMbAfQTPifAkZL+pSkXag8fGJBa9oys3Zr+LA/Ij6Q\n9A/AQ8DOwG0R8XzLOjOzturoXX1+z2/Wfh35kI+Z7bgcfrNMOfxmmXL4zTLl8JtlyuE3y5TDb5Yp\nh98sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y5TDb5Yph98sUw6/WaYcfrNMOfxmmXL4zTLl8Jtl\nyuE3y5TDb5Yph98sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y9SgZhaWtALYCGwBPoiIMa1oysza\nr6nwFz4bEW+1YD1m1kE+7DfLVLPhD+BhSU9LmtTXDJImSVosaXGT2zKzFlJENL6wdEBErJK0H7AQ\nmBwRjyfmb3xjZlaXiFA98zW154+IVcX3dcA9wNhm1mdmndNw+CXtIWmvbT8DJwNLW9WYmbVXM2f7\nhwP3SNq2np9GxIMt6crM2q6p9/z93pjf85u1XUfe85vZjsvhN8uUw2+WKYffLFMOv1mmWnFjj3Wx\nY489Nlm/6KKLkvVx48Yl64cffni/e9rm6quvTtZXr16drB933HHJ+h133FG1tmjRouSyOfCe3yxT\nDr9Zphx+s0w5/GaZcvjNMuXwm2XK4TfLlO/qGwDOP//8qrVp06Yllx06dGiyXtyyXdWjjz6arA8b\nNqxq7bDDDksuW0ut3u66666qtQsuuKCpbXcz39VnZkkOv1mmHH6zTDn8Zply+M0y5fCbZcrhN8uU\n7+fvAoMGpf8axoxJD348Y8aMqrXdd989uezjj1cdYAmAG264IVl/8sknk/Vdd921am3evHnJZU8+\n+eRkvZbFiz1CXIr3/GaZcvjNMuXwm2XK4TfLlMNvlimH3yxTDr9ZpnydvwvUenb+zJkzG173woUL\nk/XUswAA3nnnnYa3XWv9zV7HX7lyZbI+e/bsptY/0NXc80u6TdI6SUt7TBsiaaGk5cX3fdvbppm1\nWj2H/bOAU3tNuxZ4JCJGA48Uv5vZDqRm+CPicWB9r8lnA9uOqWYD41vcl5m1WaPv+YdHxJri5zeA\n4dVmlDQJmNTgdsysTZo+4RcRkXowZ0RMB6aDH+Bp1k0avdS3VtIIgOL7uta1ZGad0Gj4FwAXFz9f\nDNzbmnbMrFNqPrdf0p3ACcBQYC3wDeAXwDzgQOA14LyI6H1SsK91ZXnYX+ue+ClTpiTrtf6Obrnl\nlqq16667Lrlss9fxa1m2bFnV2ujRo5ta97nnnpus33tvnvukep/bX/M9f0RMqFI6sV8dmVlX8cd7\nzTLl8JtlyuE3y5TDb5Yph98sU76ltwWuv/76ZL3WpbzNmzcn6w899FCyfs0111Stvffee8llaxk8\neHCyXuu23AMPPLBqrdYQ2zfeeGOynuulvFbxnt8sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y1TN\nW3pburEd+JbeffbZp2rtxRdfTC47dOjQZP3+++9P1sePb98jEg899NBkfc6cOcn6Mccc0/C258+f\nn6xfeumlyfqmTZsa3vZAVu8tvd7zm2XK4TfLlMNvlimH3yxTDr9Zphx+s0w5/GaZ8nX+Ou23335V\na6tXr25q3aNGjUrW33///WR94sSJVWtnnXVWctkjjjgiWd9zzz2T9Vr/flL1c845J7nsfffdl6xb\n33yd38ySHH6zTDn8Zply+M0y5fCbZcrhN8uUw2+WKV/nr1Pqfv7UMNQAw4YNS9ZrPb++nX9HtT6j\nUKu3ESNGJOtvvvlmw8taY1p2nV/SbZLWSVraY9pUSaskLSm+TmumWTPrvHoO+2cBp/Yx/XsRcVTx\n9UBr2zKzdqsZ/oh4HFjfgV7MrIOaOeE3WdKzxduCfavNJGmSpMWSFjexLTNrsUbD/0NgFHAUsAb4\nbrUZI2J6RIyJiDENbsvM2qCh8EfE2ojYEhFbgRnA2Na2ZWbt1lD4JfW8RvMFYGm1ec2sOw2qNYOk\nO4ETgKGSVgLfAE6QdBQQwArg8jb22BU2bNhQtVbrufq1nss/ZMiQZP3ll19O1lPj1M+aNSu57Pr1\n6XO5c+fOTdZrXauvtbyVp2b4I2JCH5NvbUMvZtZB/nivWaYcfrNMOfxmmXL4zTLl8JtlqubZfqtt\n0aJFyXqtW3rLdPzxxyfr48aNS9a3bt2arL/yyiv97sk6w3t+s0w5/GaZcvjNMuXwm2XK4TfLlMNv\nlimH3yxTvs6fud122y1Zr3Udv9ZjxX1Lb/fynt8sUw6/WaYcfrNMOfxmmXL4zTLl8JtlyuE3y5SH\n6LakLVu2JOu1/v2kHu2dGr7bGteyIbrNbGBy+M0y5fCbZcrhN8uUw2+WKYffLFMOv1mm6hmieyRw\nOzCcypDc0yNimqQhwM+Ag6kM031eRLzdvlatHU455ZSyW7CS1LPn/wC4KiIOAz4DfEXSYcC1wCMR\nMRp4pPjdzHYQNcMfEWsi4pni543AMuAA4GxgdjHbbGB8u5o0s9br13t+SQcDRwOLgOERsaYovUHl\nbYGZ7SDqfoafpD2B+cCVEfGO9KePD0dEVPvcvqRJwKRmGzWz1qprzy/pY1SCPyci7i4mr5U0oqiP\nANb1tWxETI+IMRExphUNm1lr1Ay/Krv4W4FlEXFzj9IC4OLi54uBe1vfnpm1Sz2H/X8J/A3wnKQl\nxbQpwLeAeZIuA14DzmtPi9ZOo0aNKrsFK0nN8EfEk0C1+4NPbG07ZtYp/oSfWaYcfrNMOfxmmXL4\nzTLl8JtlyuE3y5SH6M7cE088kazvtFN6/1BrCG/rXt7zm2XK4TfLlMNvlimH3yxTDr9Zphx+s0w5\n/GaZ8nX+zC1dujRZX758ebJe63kAhxxySNWah+gul/f8Zply+M0y5fCbZcrhN8uUw2+WKYffLFMO\nv1mmFNHnKFvt2ViVIb2se11yySXJ+syZM5P1xx57rGpt8uTJyWVfeOGFZN36FhHVHrW/He/5zTLl\n8JtlyuE3y5TDb5Yph98sUw6/WaYcfrNM1bzOL2kkcDswHAhgekRMkzQV+BKw7absKRHxQI11+Tr/\nDmbvvfdO1ufNm5esn3TSSVVrd999d3LZiRMnJuubNm1K1nNV73X+eh7m8QFwVUQ8I2kv4GlJC4va\n9yLiO402aWblqRn+iFgDrCl+3ihpGXBAuxszs/bq13t+SQcDRwOLikmTJT0r6TZJ+1ZZZpKkxZIW\nN9WpmbVU3eGXtCcwH7gyIt4BfgiMAo6icmTw3b6Wi4jpETEmIsa0oF8za5G6wi/pY1SCPyci7gaI\niLURsSUitgIzgLHta9PMWq1m+CUJuBVYFhE395g+osdsXwDSj4E1s65Sz6W+44AngOeAbeMxTwEm\nUDnkD2AFcHlxcjC1Ll/qG2BqXQq86aabqtauuOKK5LJHHnlksu5bfvvWskt9EfEk0NfKktf0zay7\n+RN+Zply+M0y5fCbZcrhN8uUw2+WKYffLFN+dLfZAONHd5tZksNvlimH3yxTDr9Zphx+s0w5/GaZ\ncvjNMlXP03tb6S3gtR6/Dy2mdaNu7a1b+wL31qhW9nZQvTN29EM+H9m4tLhbn+3Xrb11a1/g3hpV\nVm8+7DfLlMNvlqmywz+95O2ndGtv3doXuLdGldJbqe/5zaw8Ze/5zawkDr9ZpkoJv6RTJf1W0kuS\nri2jh2okrZD0nKQlZY8vWIyBuE7S0h7ThkhaKGl58b3PMRJL6m2qpFXFa7dE0mkl9TZS0q8kvSDp\neUlfLaaX+tol+irldev4e35JOwO/Az4HrASeAiZERFeMwCBpBTAmIkr/QIik44F3gdsj4ohi2reB\n9RHxreI/zn0j4pou6W0q8G7Zw7YXo0mN6DmsPDAeuIQSX7tEX+dRwutWxp5/LPBSRLwSEZuBucDZ\nJfTR9SLicWB9r8lnA7OLn2dT+cfTcVV66woRsSYinil+3ghsG1a+1Ncu0Vcpygj/AcDrPX5fSYkv\nQB8CeFjS05Imld1MH4b3GBbtDWB4mc30oeaw7Z3Ua1j5rnntGhnuvtV8wu+jjouIo4DPA18pDm+7\nUlTes3XTtdq6hm3vlD6Glf9Qma9do8Pdt1oZ4V8FjOzx+yeLaV0hIlYV39cB99B9Q4+v3TZCcvF9\nXcn9fKibhm3va1h5uuC166bh7ssI/1PAaEmfkrQLcAGwoIQ+PkLSHsWJGCTtAZxM9w09vgC4uPj5\nYuDeEnvZTrcM215tWHlKfu26brj7iOj4F3AalTP+LwNfK6OHKn2NAn5TfD1fdm/AnVQOA/9I5dzI\nZcAngEeA5cDDwJAu6u0nVIZyf5ZK0EaU1NtxVA7pnwWWFF+nlf3aJfoq5XXzx3vNMuUTfmaZcvjN\nMuXwm2XK4TfLlMNvlimH3yxTDr9Zpv4fb4v8973Rs5MAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x129da0400>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# load dataset\n",
    "(x_train, y_train), (x_test, y_test) = mnist.load_data()\n",
    "\n",
    "print(x_train.shape[0], 'Number of train samples')\n",
    "print(x_test.shape[0], 'Number of test samples')\n",
    "print('x_train shape:', x_train.shape)\n",
    "\n",
    "# plot images as gray scale\n",
    "for image_id in range(0, 5):\n",
    "    plt.imshow(x_train[image_id], cmap=plt.get_cmap('gray'))\n",
    "    plt.title(\"The true label is %s\" % str(y_train[image_id]))\n",
    "    plt.show()\n",
    "    \n",
    "# This makes sure the image has the correct order in the axis for Tensorflow, it would be different for Theano backend\n",
    "x_train = x_train.reshape(x_train.shape[0], 28, 28, 1)\n",
    "x_test = x_test.reshape(x_test.shape[0], 28, 28, 1)\n",
    "\n",
    "# Convert values to floats, originally they are integers\n",
    "x_train = x_train.astype('float32')\n",
    "x_test = x_test.astype('float32')\n",
    "\n",
    "# Convert values of labels from 0 to 9 to categorical (one_hot encoding)\n",
    "y_train = keras.utils.to_categorical(y_train, 10)\n",
    "y_test = keras.utils.to_categorical(y_test, 10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "print(x_train.shape, y_train.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Why use Keras?\n",
    "\n",
    "\"Keras is a high-level neural networks API, written in Python and capable of running on top of TensorFlow, CNTK, or Theano. \n",
    "It was developed with a focus on enabling fast experimentation. Being able to go from idea to result with the least possible delay is key to doing good research.\n",
    "\n",
    "Use Keras if you need a deep learning library that:\n",
    "\n",
    "- Allows for easy and fast prototyping (through user friendliness, modularity, and extensibility).\n",
    "\n",
    "- Supports both convolutional networks and recurrent networks, as well as combinations of the two.\n",
    "\n",
    "- Runs seamlessly on CPU and GPU.\" (Chollet, F.)\n",
    "\n",
    "More information at https://keras.io Descriptions about the functions and documentation are also taken from this website \n",
    "\n",
    "In summary, it makes your life way easier if you don't require to go to the level of granularity of Tensorflow"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## How does it work?\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To train a model we have 3 main steps:\n",
    "\n",
    "   - Define your architecture (number of layers, type of layers, activations, etc)\n",
    "    \n",
    "   - Compile your model (Define optimizer, callbacks, etc)\n",
    "    \n",
    "   - Train your model (Fit model to your data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define your architecture"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from keras.models import Sequential\n",
    "from keras.layers import Dense, Flatten, Dropout, Activation\n",
    "\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dropout(0.5)) # drop a unit with  50% probability.\n",
    "model.add(Dense(150, kernel_initializer='orthogonal'))\n",
    "model.add(Activation('sigmoid'))\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax')) # last layer, this has a softmax to do the classification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The sequential model is one layer after the other, no fancy connections. You can also use the Functional API if you want to have connections which skip layers and so on.\n",
    "\n",
    "To add a layer all we have to do is call the add method and then add a layer object.\n",
    "\n",
    "The first layer needs to have the input shape specified. After this the sizes are inferred automatically except for the output layer.\n",
    "\n",
    "We need to use the flatten layer to flatten the input since it is an image, which in this case is a 3 dimensional matrix. It has height, width and depth. In the case of a grayscale image like the ones from MNIST it has 1 dimensional depth. If it is a color image it has 3 layers of depth (Red, Green, Blue).\n",
    "\n",
    "A dense layer is a layer in which all units are connected to all units in the next layer. This is the most usual type of layer. You can specify things like the number of units and how the weights in the units are initialized (kernel_initializer). You can also specify an activation function, by default a linear function is used (No activation).\n",
    "\n",
    "If you want to use dropout in your model you can just add it as an extra layer in between layers or activations. Same thing for batch normalization.\n",
    "\n",
    "At the end we create a layer which will help us with the classification. For this we use a softmax layer in which the number of units will match the number of classes we have in our data.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Dense Layer\n",
    "\n",
    "keras.layers.Dense(units, activation=None, use_bias=True, \n",
    "                   kernel_initializer='glorot_uniform', bias_initializer='zeros', \n",
    "                   kernel_regularizer=None, bias_regularizer=None, \n",
    "                   activity_regularizer=None, \n",
    "                   kernel_constraint=None, bias_constraint=None)\n",
    "                   \n",
    "Arguments\n",
    "\n",
    "- units: Positive integer, dimensionality of the output space.\n",
    "- activation: Activation function to use (see activations). If you don't specify anything, no activation is applied (ie. \"linear\" activation: a(x) = x).\n",
    "- use_bias: Boolean, whether the layer uses a bias vector.\n",
    "- kernel_initializer: Initializer for the kernel weights matrix.\n",
    "- bias_initializer: Initializer for the bias vector.\n",
    "- kernel_regularizer: Regularizer function applied to the kernel weights matrix. (L1, L2)\n",
    "- bias_regularizer: Regularizer function applied to the bias vector.\n",
    "- activity_regularizer: Regularizer function applied to the output of the layer (its \"activation\").\n",
    "- kernel_constraint: Constraint function applied to the kernel weights matrix. (non-negative, etc)\n",
    "- bias_constraint: Constraint function applied to the bias vector."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compile the model\n",
    "\n",
    "Here we specify the loss, in our case categorical crossentropy. We can add an extra metric we want to measure, like accuracy\n",
    "\n",
    "We also specify the optimizer. Some examples are Adam and SGD"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sgd = keras.optimizers.SGD(lr=0.01, decay=1e-6, momentum=0.9)\n",
    "\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Compile function\n",
    "\n",
    "Arguments\n",
    "\n",
    "- optimizer: String (name of optimizer) or optimizer object. (SGD, RMSprop, Adam, Adagrad, etc)\n",
    "- loss: String (name of objective function) or objective function. You can also have multiple loss functions\n",
    "- metrics: List of metrics to be evaluated by the model during training and testing. \n",
    "- sample_weight_mode: If you need to do timestep-wise sample weighting (2D weights), set this to \"temporal\". None defaults to sample-wise weights (1D).\n",
    "- weighted_metrics: List of metrics to be evaluated and weighted by sample_weight or class_weight during training and testing."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Callbacks\n",
    "\n",
    "A callback is a set of functions to be applied at given stages of the training procedure. You can use callbacks to get a view on internal states and statistics of the model during training.\n",
    "\n",
    "Some examples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from keras.callbacks import EarlyStopping\n",
    "\n",
    "early_stopping = EarlyStopping(monitor='val_loss', min_delta=0, patience=5, verbose=0, mode='auto')\n",
    "my_callbacks = [early_stopping]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Train the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "history = model.fit(x_train, y_train, validation_split=0.2, epochs=10, batch_size=100, verbose=1, callbacks=my_callbacks)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Fit function\n",
    "\n",
    "Arguments\n",
    "\n",
    "- x: Numpy array of training data. \n",
    "- y: Numpy array of target (label) data.\n",
    "- batch_size: Integer or None. Number of samples per gradient update. If unspecified, it will default to 32.\n",
    "- epochs: Integer. Number of epochs to train the model. An epoch is an iteration over the entire x and y data provided.\n",
    "- verbose: 0, 1, or 2. Verbosity mode. 0 = silent, 1 = progress bar, 2 = one line per epoch.\n",
    "- callbacks: List of keras.callbacks.Callback instances. List of callbacks to apply during training.\n",
    "- validation_split: Float between 0 and 1. Fraction of the training data to be used as validation data.\n",
    "- validation_data: tuple (x_val, y_val) or tuple (x_val, y_val, val_sample_weights) on which to evaluate the loss and any model metrics at the end of each epoch.\n",
    "- shuffle: Boolean (whether to shuffle the training data before each epoch)\n",
    "- class_weight: Optional dictionary mapping class indices (integers) to a weight (float) value, used for weighting the loss function (during training only).\n",
    "- sample_weight: Optional Numpy array of weights for the training samples, used for weighting the loss function (during training only).\n",
    "- initial_epoch: Epoch at which to start training (useful for resuming a previous training run).\n",
    "- steps_per_epoch: Total number of steps (batches of samples) before declaring one epoch finished and starting the next epoch. When training with input tensors such as TensorFlow data tensors, the default None is equal to the number of samples in your dataset divided by the batch size, or 1 if that cannot be determined.\n",
    "- validation_steps: Only relevant if steps_per_epoch is specified. Total number of steps (batches of samples) to validate before stopping."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### History object\n",
    "\n",
    "Fit function returns a history object. This object stores information about your model while it was training. For example if you would like to plot the training loss and validation loss across the number of epochs you could do this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "print(\"Keys for history object\", history.history.keys())\n",
    "\n",
    "train_loss = history.history['loss']\n",
    "valid_loss = history.history['val_loss']\n",
    "epochs = list(range(1, len(train_loss)+1))\n",
    "\n",
    "plt.plot(epochs, train_loss, label=\"train loss\")\n",
    "plt.plot(epochs, valid_loss, label=\"validation loss\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Loss')\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "train_accuracy = history.history['acc']\n",
    "valid_accuracy = history.history['val_acc']\n",
    "\n",
    "epochs = list(range(1, len(train_accuracy)+1))\n",
    "\n",
    "plt.plot(epochs, train_accuracy, label=\"train accuracy\")\n",
    "plt.plot(epochs, valid_accuracy, label=\"validation accuracy\")\n",
    "\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Accuracy')\n",
    "plt.legend(loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Making predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "predictions = model.predict(x_test)\n",
    "\n",
    "print(predictions.shape)\n",
    "print(np.argmax(predictions, axis=1)[0:10])\n",
    "print(np.argmax(y_test, axis=1)[0:10])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Arguments\n",
    "\n",
    "- x: the input data, as a Numpy array.\n",
    "- batch_size: Integer. If unspecified, it will default to 32.\n",
    "- verbose: verbosity mode, 0 or 1.\n",
    "- steps: Total number of steps (batches of samples) before declaring the prediction round finished. Ignored with the default value of None.\n",
    "\n",
    "Returns\n",
    "\n",
    "A Numpy array of predictions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Loading and saving a model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from keras.models import load_model\n",
    "\n",
    "model_path = \"my_new_model.h5\"\n",
    "model.save(model_path)\n",
    "del model  # deletes the existing model\n",
    "model = load_model(model_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Mini-Project: Train you own fully connected neural network to classify handwritten digits\n",
    "\n",
    "- Train a logistic regression to create a benchmark\n",
    "- Train a neural network and compare. Start with a few layers, then experiment with more layers and different parameters.\n",
    "- Are you overfitting? Underfitting? How can you improve your model? Try other hyperparameters or adding regularization. Make some plots to understand the behaviour of your model.\n",
    "- Take a look at samples for which your model is predicting an incorrect label, what do you think is happening?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from keras.models import Sequential\n",
    "from keras.layers import Dense, Flatten, Dropout"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic regression (Benchmark)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])\n",
    "\n",
    "# Train the model\n",
    "log_history = model.fit(x_train, y_train, validation_split=0.2, epochs=10, batch_size=100, verbose=2)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Shallow Neural Network"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])\n",
    "\n",
    "# Train the model\n",
    "nn_history = model.fit(x_train, y_train, validation_split=0.2, epochs=10, batch_size=100, verbose=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "log_train_loss = log_history.history['loss']\n",
    "log_valid_loss = log_history.history['val_loss']\n",
    "nn_train_loss = nn_history.history['loss']\n",
    "nn_valid_loss = nn_history.history['val_loss']\n",
    "epochs = list(range(1, len(nn_train_loss)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_loss, label=\"log reg train loss\")\n",
    "plt.plot(epochs, log_valid_loss, label=\"log reg validation loss\")\n",
    "plt.plot(epochs, nn_train_loss, label=\"NN train loss\")\n",
    "plt.plot(epochs, nn_valid_loss, label=\"NN validation loss\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Loss')\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()\n",
    "\n",
    "log_train_accuracy = log_history.history['acc']\n",
    "log_valid_accuracy = log_history.history['val_acc']\n",
    "nn_train_accuracy = nn_history.history['acc']\n",
    "nn_valid_accuracy = nn_history.history['val_acc']\n",
    "epochs = list(range(1, len(log_train_accuracy)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_accuracy, label=\"log reg train accuracy\")\n",
    "plt.plot(epochs, log_valid_accuracy, label=\"log reg validation accuracy\")\n",
    "plt.plot(epochs, nn_train_accuracy, label=\"NN train accuracy\")\n",
    "plt.plot(epochs, nn_valid_accuracy, label=\"NN validation accuracy\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Accuracy')\n",
    "plt.legend(loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "predictions = model.predict(x_test)\n",
    "\n",
    "predicted_labels = np.argmax(predictions, axis=1)  # transform back from one_hot encoding\n",
    "true_labels = np.argmax(y_test, axis=1)\n",
    "misclassified_samples_positions = np.where(predicted_labels != true_labels)[0]\n",
    "\n",
    "for image_id in misclassified_samples_positions[0:5]:\n",
    "    print(image_id)\n",
    "    print(x_test[image_id].shape)\n",
    "    plt.imshow(x_test[image_id].reshape(28, 28), cmap=plt.get_cmap('gray'))\n",
    "    plt.title(\"The true label is %s but was classified as %s\" % (str(true_labels[image_id]), \n",
    "                                                                 str(predicted_labels[image_id])\n",
    "                                                                ))\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Scaling matters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Scale the inputs\n",
    "x_train /= 255.0  # The image is in grayscale and has values between 0 and 255\n",
    "x_test /= 255.0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Mini-Project: Train you own fully connected neural network to classify handwritten digits (Now with scaled data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Train a logistic regression to create a benchmark\n",
    "- Train a neural network and compare. Start with a few layers, then experiment with more layers and different parameters.\n",
    "- Are you overfitting? Underfitting? How can you improve your model? Try other hyperparameters or adding regularization. Make some plots to understand the behaviour of your model."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Extra credit\n",
    "\n",
    "- Try shuffling the order of the features before training (move the 5th pixel to a new random position) but do the same shuffling for all samples and retrain you model. How does it perform compared to the original model?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic regression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# Train the model\n",
    "log_history = model.fit(x_train, y_train, validation_split=0.2, epochs=10, batch_size=100, verbose=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Shallow neural network"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Train the model\n",
    "nn_history = model.fit(x_train, y_train, validation_split=0.2, epochs=10, batch_size=100, verbose=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "log_train_loss = log_history.history['loss']\n",
    "log_valid_loss = log_history.history['val_loss']\n",
    "nn_train_loss = nn_history.history['loss']\n",
    "nn_valid_loss = nn_history.history['val_loss']\n",
    "epochs = list(range(1, len(log_train_loss)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_loss, label=\"log reg train loss\")\n",
    "plt.plot(epochs, log_valid_loss, label=\"log reg validation loss\")\n",
    "plt.plot(epochs, nn_train_loss, label=\"NN train loss\")\n",
    "plt.plot(epochs, nn_valid_loss, label=\"NN validation loss\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Loss')\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()\n",
    "\n",
    "log_train_accuracy = log_history.history['acc']\n",
    "log_valid_accuracy = log_history.history['val_acc']\n",
    "nn_train_accuracy = nn_history.history['acc']\n",
    "nn_valid_accuracy = nn_history.history['val_acc']\n",
    "epochs = list(range(1, len(log_train_accuracy)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_accuracy, label=\"log reg train accuracy\")\n",
    "plt.plot(epochs, log_valid_accuracy, label=\"log reg validation accuracy\")\n",
    "plt.plot(epochs, nn_train_accuracy, label=\"NN train accuracy\")\n",
    "plt.plot(epochs, nn_valid_accuracy, label=\"NN validation accuracy\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Accuracy')\n",
    "plt.legend(loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "predictions = model.predict(x_test)\n",
    "\n",
    "predicted_labels = np.argmax(predictions, axis=1)  # transform back from one_hot encoding\n",
    "true_labels = np.argmax(y_test, axis=1)\n",
    "misclassified_samples_positions = np.where(predicted_labels != true_labels)[0]\n",
    "\n",
    "for image_id in misclassified_samples_positions[0:5]:\n",
    "    print(image_id)\n",
    "    print(x_test[image_id].shape)\n",
    "    plt.imshow(x_test[image_id].reshape(28, 28), cmap=plt.get_cmap('gray'))\n",
    "    plt.title(\"The true label is %s but was classified as %s\" % (str(true_labels[image_id]), \n",
    "                                                                 str(predicted_labels[image_id])\n",
    "                                                                ))\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Reload data before procedding"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# load dataset\n",
    "(x_train, y_train), (x_test, y_test) = mnist.load_data()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# This makes sure the image has the correct order in the axis for Tensorflow, it would be different for Theano backend\n",
    "x_train = x_train.reshape(x_train.shape[0], 28, 28, 1)\n",
    "x_test = x_test.reshape(x_test.shape[0], 28, 28, 1)\n",
    "\n",
    "# Convert values to floats, originally they are integers\n",
    "x_train = x_train.astype('float32')\n",
    "x_test = x_test.astype('float32')\n",
    "\n",
    "# Convert values of labels from 0 to 9 to categorical (one_hot encoding)\n",
    "y_train = keras.utils.to_categorical(y_train, 10)\n",
    "y_test = keras.utils.to_categorical(y_test, 10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data augmentation and generators\n",
    "\n",
    "Augmenting the data allows us to use more data for training by taking artificial modifications from the samples we already have. One example of this is rotating the image. This will allow our model to learn some rotational invariance to the data and prevent overfitting.\n",
    "\n",
    "To train a model when a dataset is too large to load in RAM memory we can use generators. Generators allow us to only pass a small portion of the data at a time so it doesn't use all the memory from the computer. In addition once it loops through all the data we can set the generator so it starts again from the beginning. This allows us to loop infinitely through the data.\n",
    "\n",
    "Keras gives us a nice tool for augmenting data which uses a generator. It uses a generator because the modifications are done on the fly and it would take too much memory to save the augmented dataset. Keras also provides us with a method to fit the model when we use a generator, fit_generator(). This method requires us to specify how many batches are equivalent to an epoch so the model knows when to do updates, otherwise since the generator loops through the data, it would go on training infinitely."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from keras.preprocessing.image import ImageDataGenerator\n",
    "\n",
    "data_generator = ImageDataGenerator(rescale=1./255,\n",
    "                                    #featurewise_center=True,\n",
    "                                    #featurewise_std_normalization=True,\n",
    "                                    rotation_range=10,\n",
    "                                    width_shift_range=0.1,\n",
    "                                    height_shift_range=0.1,\n",
    "                                    #horizontal_flip=True,\n",
    "                                    zoom_range=0.1\n",
    "                                    )\n",
    "\n",
    "# we git the model to the data. This needed for calculating mean and std.\n",
    "data_generator.fit(x_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### What does this look like?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "augmented_data_generator = data_generator.flow(x_train[0:40000], y_train[0:40000], batch_size=100)\n",
    "augmented_batch, labels = augmented_data_generator.next()\n",
    "\n",
    "print(augmented_batch.shape)\n",
    "print(augmented_batch[0].shape)\n",
    "print(labels[0])\n",
    "\n",
    "# plot images as gray scale\n",
    "for image_id in range(0, 5):\n",
    "    plt.imshow(augmented_batch[image_id].reshape(28, 28), cmap=plt.get_cmap('gray'))\n",
    "    plt.title(\"The true label is %s\" % str(np.argmax(labels[image_id])))\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Training a model using a generator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "history = model.fit_generator(data_generator.flow(x_train[0:40000], y_train[0:40000], batch_size=100), \n",
    "                              steps_per_epoch=40000/100,\n",
    "                              validation_data=(x_train[40000:], y_train[40000:]), \n",
    "                              epochs=10, verbose=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Arguments\n",
    "\n",
    "- generator: A generator. The output of the generator must be either\n",
    "a tuple (inputs, targets)\n",
    "a tuple (inputs, targets, sample_weights). All arrays should contain the same number of samples. The generator is expected to loop over its data indefinitely. An epoch finishes when steps_per_epoch batches have been seen by the model.\n",
    "- steps_per_epoch: Total number of steps (batches of samples) to yield from generator before declaring one epoch finished and starting the next epoch. It should typically be equal to the number of samples of your dataset divided by the batch size.\n",
    "- epochs: Integer, total number of iterations on the data. \n",
    "- verbose: Verbosity mode, 0, 1, or 2.\n",
    "- callbacks: List of callbacks to be called during training.\n",
    "- validation_data: This can be either\n",
    "A generator for the validation data\n",
    "A tuple (inputs, targets)\n",
    "A tuple (inputs, targets, sample_weights).\n",
    "- validation_steps: Only relevant if validation_data is a generator. Number of steps to yield from validation generator at the end of every epoch. It should typically be equal to the number of samples of your validation dataset divided by the batch size. Optional for Sequence: if unspecified, will use the len(validation_data) as a number of steps.\n",
    "- class_weight: Dictionary mapping class indices to a weight for the class.\n",
    "- max_queue_size: Maximum size for the generator queue\n",
    "- workers: Maximum number of processes to spin up\n",
    "- use_multiprocessing: if True, use process based threading.\n",
    "- initial_epoch: Epoch at which to start training (useful for resuming a previous training run).\n",
    "\n",
    "Returns\n",
    "\n",
    "A History object."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Mini-Project: Train you own fully connected neural network to classify handwritten digits (Now with augmented data)\n",
    "\n",
    "- Train a logistic regression to create a benchmark\n",
    "- Train a neural network and compare. Start with a few layers, then experiment with more layers and different parameters.\n",
    "- Are you overfitting? Underfitting? How can you improve your model? Try other hyperparameters or adding regularization. Make some plots to understand the behaviour of your model."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic regression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# Train the model\n",
    "log_history = model.fit_generator(data_generator.flow(x_train[0:40000], y_train[0:40000], batch_size=100), \n",
    "                              steps_per_epoch=40000/100,\n",
    "                              validation_data=(x_train[40000:], y_train[40000:]), \n",
    "                              epochs=10, verbose=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Shallow neural network"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Flatten(input_shape=(28, 28, 1)))  # Images are a 3D matrix, we have to flatten them to be 1D\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(200, kernel_initializer='normal', activation='relu'))\n",
    "model.add(Dense(10, kernel_initializer='normal', activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Train the model\n",
    "nn_history = model.fit_generator(data_generator.flow(x_train[0:40000], y_train[0:40000], batch_size=100), \n",
    "                              steps_per_epoch=40000/100,\n",
    "                              validation_data=(x_train[40000:], y_train[40000:]), \n",
    "                              epochs=10, verbose=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "log_train_loss = log_history.history['loss']\n",
    "log_valid_loss = log_history.history['val_loss']\n",
    "nn_train_loss = nn_history.history['loss']\n",
    "nn_valid_loss = nn_history.history['val_loss']\n",
    "epochs = list(range(1, len(log_train_loss)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_loss, label=\"log reg train loss\")\n",
    "plt.plot(epochs, log_valid_loss, label=\"log reg validation loss\")\n",
    "plt.plot(epochs, nn_train_loss, label=\"NN train loss\")\n",
    "plt.plot(epochs, nn_valid_loss, label=\"NN validation loss\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Loss')\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()\n",
    "\n",
    "log_train_accuracy = log_history.history['acc']\n",
    "log_valid_accuracy = log_history.history['val_acc']\n",
    "nn_train_accuracy = nn_history.history['acc']\n",
    "nn_valid_accuracy = nn_history.history['val_acc']\n",
    "epochs = list(range(1, len(log_train_accuracy)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_accuracy, label=\"log reg train accuracy\")\n",
    "plt.plot(epochs, log_valid_accuracy, label=\"log reg validation accuracy\")\n",
    "plt.plot(epochs, nn_train_accuracy, label=\"NN train accuracy\")\n",
    "plt.plot(epochs, nn_valid_accuracy, label=\"NN validation accuracy\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Accuracy')\n",
    "plt.legend(loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Convolutional Neural Networks"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Convolutional neural networks are similar to the previous networks we saw with a main difference, they start with the assumption that the input will be an image and optimize the architecture for that assumption. Now they are used in many other contexts besides images, but this is what they were created for.\n",
    "\n",
    "Convolutional layer define a kernel (weight matrix) which is then multiplied element by element with a section of the input of the same size. Sum all the resulting values. Move the kernel a number of pixels equal to a defined stride size and repeat until you go throught the whole image. Afterwards apply an activation function to each one of the values. Now repeat all of this for each filter in your layer.\n",
    "\n",
    "<img src=\"images/conv_layer.gif\" alt=\"convolutional layer\">\n",
    "Animation from Karpathy (http://cs231n.github.io/convolutional-networks/)\n",
    "\n",
    "A key aspect is the fact that the weights from the kernel are the same when applied in all sections of the input for the same filter. This means the convolutional layer has way less free parameters than a fully connected layer. This also gives the model location invariance since the filters will activate in the same way if the object is in the top left of the image or at the bottom right.\n",
    "\n",
    "Images are 3D objects which have height, length and depth. In the case of grayscale images the depth is 1, however in color images you have 3 channels, red, green and blue. So convolutional layers take as input a 3D tensor and outputs a 3D tensor.\n",
    "\n",
    "<img src=\"images/cnn.jpeg\" alt=\"convolutional layer\">\n",
    "Image from Karpathy (http://cs231n.github.io/convolutional-networks/)\n",
    "\n",
    "Notice the convolutional operation will reduce the size of the input. Depending on your architecture you might want to keep the size constant, you can do this by includding padding in the borders of the input. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Pooling layers\n",
    "\n",
    "The idea of the pooling layer is to reduce the size of the input and also help regularize the model. The most famous type of pooling is max pooling. In this case you select an area which is the size of the filter and from it you will only pass through the maximum values among the inputs. Then move the filter a specific stride size and repeat.\n",
    "\n",
    "<img src=\"images/maxpool.jpeg\" alt=\"convolutional layer\">\n",
    "Image from Karpathy (http://cs231n.github.io/convolutional-networks/)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Best Practices"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In practice it is best to use blocks of stacked convolutional layers with rectified linear units followed by a max pooling operation. Then at the end have a small number of fully connected layers and the output layer.\n",
    "\n",
    "Each convolutional layer usually has a small filter size (2, 2), (3, 3) and use max pooling with stride (2, 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Convolutional Neural Networks in Keras"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from keras.layers import Conv2D, MaxPooling2D"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Convolutional layer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Conv2D(filters, kernel_size, strides=(1, 1), padding='valid', \n",
    "                    activation=None, use_bias=True, \n",
    "                    kernel_initializer='glorot_uniform', bias_initializer='zeros', \n",
    "                    kernel_regularizer=None, bias_regularizer=None, activity_regularizer=None, \n",
    "                    kernel_constraint=None, bias_constraint=None)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Conv2D(32, (3, 3), padding='same')  # Keeps output of the same size as input\n",
    "Conv2D(32, (3, 3)) # Reduces the size of the output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Arguments\n",
    "\n",
    "- filters: Integer, the dimensionality of the output space (i.e. the number of output filters in the convolution).\n",
    "- kernel_size: An integer or tuple/list of 2 integers, specifying the width and height of the 2D convolution window. Can be a single integer to specify the same value for all spatial dimensions.\n",
    "- strides: An integer or tuple/list of 2 integers, specifying the strides of the convolution along the width and height. Can be a single integer to specify the same value for all spatial dimensions. Specifying any stride value != 1 is incompatible with specifying any dilation_rate value != 1.\n",
    "- padding: one of \"valid\" or \"same\" (case-insensitive).\n",
    "- activation: Activation function to use (see activations). If you don't specify anything, no activation is applied (ie. \"linear\" activation: a(x) = x).\n",
    "- use_bias: Boolean, whether the layer uses a bias vector.\n",
    "- kernel_initializer: Initializer for the kernel weights matrix (see initializers).\n",
    "- bias_initializer: Initializer for the bias vector (see initializers).\n",
    "- kernel_regularizer: Regularizer function applied to the kernel weights matrix (see regularizer).\n",
    "- bias_regularizer: Regularizer function applied to the bias vector (see regularizer).\n",
    "- activity_regularizer: Regularizer function applied to the output of the layer (its \"activation\"). (see regularizer).\n",
    "- kernel_constraint: Constraint function applied to the kernel matrix (see constraints).\n",
    "- bias_constraint: Constraint function applied to the bias vector (see constraints).\n",
    "\n",
    "Input shape\n",
    "\n",
    "4D tensor with shape: (samples, channels, rows, cols) if data_format='channels_first' or 4D tensor with shape: (samples, rows, cols, channels) if data_format='channels_last'.\n",
    "\n",
    "Output shape\n",
    "\n",
    "4D tensor with shape: (samples, filters, new_rows, new_cols) if data_format='channels_first' or 4D tensor with shape: (samples, new_rows, new_cols, filters) if data_format='channels_last'. rows and cols values might have changed due to padding."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "### Pooling layer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "MaxPooling2D(pool_size=(2, 2), strides=None, padding='valid')  # This will reduce the size of the output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Arguments\n",
    "\n",
    "- pool_size: integer or tuple of 2 integers, factors by which to downscale (vertical, horizontal). (2, 2) will halve the input in both spatial dimension. If only one integer is specified, the same window length will be used for both dimensions.\n",
    "- strides: Integer, tuple of 2 integers, or None. Strides values. If None, it will default to pool_size.\n",
    "- padding: One of \"valid\" or \"same\" (case-insensitive).\n",
    "\n",
    "Input shape\n",
    "\n",
    "4D tensor with shape: (batch_size, rows, cols, channels)\n",
    "\n",
    "Output shape\n",
    "\n",
    "4D tensor with shape: (batch_size, pooled_rows, pooled_cols, channels)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Mini-Project: Train you own convolutional neural network to classify handwritten digits (With or without augmented data)\n",
    "\n",
    "- Train a logistic regression to create a benchmark\n",
    "- Train a convolutional neural network and compare it with the logistic regression and fully connected neural network. Start with a few layers, then experiment with more layers and different parameters.\n",
    "- Are you overfitting? Underfitting? How can you improve your model? Try other hyperparameters or adding regularization. Make some plots to understand the behaviour of your model."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Convolutional neural network"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "from keras.layers import Conv2D, MaxPooling2D, Activation\n",
    "# Define model architecture\n",
    "model = Sequential()\n",
    "model.add(Conv2D(32, (3, 3), padding='same', input_shape=(28, 28, 1), activation='relu'))\n",
    "model.add(Conv2D(32, (3, 3), activation='relu'))\n",
    "model.add(MaxPooling2D(pool_size=(2, 2)))\n",
    "model.add(Dropout(0.25))\n",
    "\n",
    "model.add(Conv2D(64, (3, 3), padding='same', activation='relu'))\n",
    "model.add(Conv2D(64, (3, 3), activation='relu'))\n",
    "model.add(MaxPooling2D(pool_size=(2, 2)))\n",
    "model.add(Dropout(0.25))\n",
    "\n",
    "model.add(Flatten())\n",
    "model.add(Dense(512, activation='relu'))\n",
    "model.add(Dropout(0.5))\n",
    "model.add(Dense(10, activation='softmax'))\n",
    "\n",
    "# Compile model\n",
    "model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoch 1/10\n",
      " - 158s - loss: 0.4550 - acc: 0.8518 - val_loss: 0.3061 - val_acc: 0.9806\n",
      "Epoch 2/10\n",
      " - 156s - loss: 0.1381 - acc: 0.9584 - val_loss: 0.1904 - val_acc: 0.9881\n",
      "Epoch 3/10\n",
      " - 156s - loss: 0.0993 - acc: 0.9694 - val_loss: 0.1709 - val_acc: 0.9892\n",
      "Epoch 4/10\n",
      " - 156s - loss: 0.0813 - acc: 0.9752 - val_loss: 0.1744 - val_acc: 0.9891\n",
      "Epoch 5/10\n",
      " - 11744s - loss: 0.0710 - acc: 0.9783 - val_loss: 0.1543 - val_acc: 0.9903\n",
      "Epoch 6/10\n",
      " - 222s - loss: 0.0662 - acc: 0.9795 - val_loss: 0.1180 - val_acc: 0.9925\n",
      "Epoch 7/10\n",
      " - 173s - loss: 0.0574 - acc: 0.9821 - val_loss: 0.1387 - val_acc: 0.9914\n",
      "Epoch 8/10\n",
      " - 173s - loss: 0.0535 - acc: 0.9840 - val_loss: 0.1165 - val_acc: 0.9927\n",
      "Epoch 9/10\n",
      " - 170s - loss: 0.0525 - acc: 0.9837 - val_loss: 0.1327 - val_acc: 0.9917\n",
      "Epoch 10/10\n",
      " - 170s - loss: 0.0478 - acc: 0.9856 - val_loss: 0.1061 - val_acc: 0.9934\n"
     ]
    }
   ],
   "source": [
    "# Train the model\n",
    "cnn_history = model.fit_generator(data_generator.flow(x_train[0:40000], y_train[0:40000], batch_size=100), \n",
    "                              steps_per_epoch=40000/100,\n",
    "                              validation_data=(x_train[40000:], y_train[40000:]), \n",
    "                              epochs=10, verbose=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'log_history' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-12-3ca785f0c1bc>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mlog_train_loss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mlog_history\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mhistory\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'loss'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0mlog_valid_loss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mlog_history\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mhistory\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'val_loss'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0mnn_train_loss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnn_history\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mhistory\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'loss'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mnn_valid_loss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnn_history\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mhistory\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'val_loss'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mcnn_train_loss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mcnn_history\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mhistory\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'loss'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'log_history' is not defined"
     ]
    }
   ],
   "source": [
    "log_train_loss = log_history.history['loss']\n",
    "log_valid_loss = log_history.history['val_loss']\n",
    "nn_train_loss = nn_history.history['loss']\n",
    "nn_valid_loss = nn_history.history['val_loss']\n",
    "cnn_train_loss = cnn_history.history['loss']\n",
    "cnn_valid_loss = cnn_history.history['val_loss']\n",
    "epochs = list(range(1, len(log_train_loss)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_loss, label=\"log reg train loss\")\n",
    "plt.plot(epochs, log_valid_loss, label=\"log reg validation loss\")\n",
    "plt.plot(epochs, nn_train_loss, label=\"NN train loss\")\n",
    "plt.plot(epochs, nn_valid_loss, label=\"NN validation loss\")\n",
    "plt.plot(epochs, cnn_train_loss, label=\"CNN train loss\")\n",
    "plt.plot(epochs, cnn_valid_loss, label=\"CNN validation loss\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Loss')\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()\n",
    "\n",
    "log_train_accuracy = log_history.history['acc']\n",
    "log_valid_accuracy = log_history.history['val_acc']\n",
    "nn_train_accuracy = nn_history.history['acc']\n",
    "nn_valid_accuracy = nn_history.history['val_acc']\n",
    "cnn_train_accuracy = cnn_history.history['acc']\n",
    "cnn_valid_accuracy = cnn_history.history['val_acc']\n",
    "epochs = list(range(1, len(log_train_accuracy)+1))\n",
    "\n",
    "plt.plot(epochs, log_train_accuracy, label=\"log reg train accuracy\")\n",
    "plt.plot(epochs, log_valid_accuracy, label=\"log reg validation accuracy\")\n",
    "plt.plot(epochs, nn_train_accuracy, label=\"NN train accuracy\")\n",
    "plt.plot(epochs, nn_valid_accuracy, label=\"NN validation accuracy\")\n",
    "plt.plot(epochs, cnn_train_accuracy, label=\"CNN train accuracy\")\n",
    "plt.plot(epochs, cnn_valid_accuracy, label=\"CNN validation accuracy\")\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Accuracy')\n",
    "plt.legend(loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [py3Keras2_env]",
   "language": "python",
   "name": "Python [py3Keras2_env]"
  },
  "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.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}