unit06_svm/demo06_mnist_svm.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Demo 6: SVM for MNIST Digit Recognition\n",
    "\n",
    "In this demo, you will learn to:\n",
    "* Load and display images\n",
    "* Formulate image classification problems\n",
    "* Explain the limitations of linear classifiers for image classification\n",
    "* Build a simple SVM image classifier \n",
    "* Save and load results using `pickle`.\n",
    "\n",
    "For data, we will use the classic [MNIST](https://en.wikipedia.org/wiki/MNIST_database) data set used to recognize hand-written digits.  The dataset was originally produced in the 1980s and is now widely-used in machine learning classes as a simple image classification problem.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Loading the Data\n",
    "\n",
    "First, we load the standard packages."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the MNIST dataset is so widely-used, the `sklearn` module has the dataset built in.  It can be loaded as follows.  We see that the data is stored in an array with 1787 images, each 8 x 8."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1797, 8, 8)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn import datasets, linear_model, preprocessing\n",
    "digits = datasets.load_digits()\n",
    "images = digits.images\n",
    "labels = digits.target\n",
    "images.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The above is fairly low-resolution.  We can get a more high-resolution 28 x 28 image data from the following command.  But, this will take a few minutes to load for some reason.  In this data format, each image is stored as a 784 = 28x28 row vector."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.datasets import fetch_mldata\n",
    "mnist = fetch_mldata('MNIST original')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(70000, 784)"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mnist.data.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have 70000 image.  Each image is `28 x 28`, which is stored as a 784 pixel vector.  The following code plots a digit.  To do this, it reshapes the 784 dim vector to a `28 x 28` array and then uses the `plt.imshow` command.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWQAAABeCAYAAAAUjW5fAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAC9NJREFUeJzt3XmIldUbB/DvmFquaDaaWo5poVhm\nk4qhZqkoJhomlYqWmkuOpuCWK4oLJVIhFgpihAtlQ44aSuH4h6OM0jJOLklYSpO4RJaE+3r7Q57n\nPm/3vXO3ufeee+f7+ef3cO59557e3t/pvGd5Tk4gEAAREaVfjXRXgIiI7mGDTETkCDbIRESOYINM\nROQINshERI5gg0xE5Ag2yEREjmCDTETkCDbIRESOYINMROSImrF8OScnh/usIwgEAjnxXMd7G5UL\ngUAgN9aLeG+jEte9BXh/oxTV/WUPmTJJRborkMV4b5MrqvvLBpmIyBFskImIHBHTGDIRJW7dunUa\nP/DAAxq/+eab6agOOYQ9ZCIiR7BBJiJyBIcsiFJk0aJFAIDRo0dr2dWrV9NVHXIQe8hERI5Iew+5\nZcuWAIBZs2Zp2ZQpUzSuVasWAODAgQNadu7cuUr/5ooVKwAAN27c0LJff/1V4+vXrydQY6LoffDB\nBxrPmDEDAHD37l0t+/bbb1NeJ3IXe8hERI5gg0xE5IicQCD6bejJ2LO+bNkyAMCCBQuq+k97nDp1\nSmOZSJk5c6aWFRcXV8nvMJdFUpUFAoEusV6Ujntbt25dAMBvv/2mZQ8//DAA4IsvvtCykSNHprZi\n4cV1b4HMeXZv3rwJACgtLdWy3r17p+rno7q/7CETETmCDTIRkSPSvsqiqKgIAPDYY4/FfG3fvn01\nfuihhzS+7777Qr7bpk2bsL8NANOmTQMAfPbZZzHXg+j/OnfuDABo3ry5lp0/fx5A8Fmj5OvXr5/G\nNWrc63/afyeuYQ+ZiMgRae8hl5eXAwBGjRpVZX/z3XffBeDtdY8bN07jmjXv/WPXq1dPy2S9czaZ\nNGkSAGDEiBFaNmjQII0vXboU9lqZgAK8byL5+fkAgO7du2vZc889p/Hly5cBBNeXR/qdbLVx40YA\ngJ00l3v/999/p6VOmU4mSgFg8eLFAIA5c+akqzpJwR4yEZEj2CATETki7UMWybBy5UoAwPPPP69l\nb731VqXXyKt2pnvqqac0/uSTTwAEJzMAYMuWLRrLeux27dpp2eOPPx5yTe3ataP+/fr16wMA7r//\nfi2rLkMWdrKoSZMmAIA///xTy8rKylJep0xnhxX37NmjsUyaFhYWapnf/f33339Dytq2batxixYt\nND579mxila0C7CETETkia3rIdlJOTl5YvXq17+fi4sWLGn/33XdJrF3q2OV9d+7cAeDt7b700ktR\n/Z1bt25pbCehZOfZL7/8omU2nWR1ZieR5U1BJp8oPvbNtlu3bhrLW1ekt44nn3xSY/v/AyET/K5g\nD5mIyBFskImIHOFWf70SdideQUEBAO+kXevWrTWWialI7OuKTMKcPHkykWqm3ddff63x2LFjAQCv\nvPKKlj377LMaS7IVm+xG8k7/9ddfWnb48OGQ32natKnGHLK45/333w8pc2GiKBPVqVMHADB37lzf\nzysqKlJZnZRhD5mIyBFskImIHOHMkMVrr72mcW5ursYvv/wyAO+rth2+SESDBg00nj17dkg9Mt3n\nn3/u+V9KLvuMxuP1118HEHzmAe8Kl02bNgEAfvzxx4R+JxP0798fQPhEQAcPHozq7/z8888ay9FZ\nfqstXOFuzYiIqpm09JCbNWum8b59+wB4d89E+i+YJGyxh5ieOXNG4+3btwMAPvzwQy3r0aOHxrK7\nJycneLhHor2b6mbgwIG+5deuXQPgXcdcHcmhun6HmNo18Rs2bNB46NChALzpY+2BqG+//TYAbzKn\nn376qYpqnH6tWrXSeP369SGf2x2fkkAsku+//17jRHrIjRo10vj27dsaV/UOX/aQiYgcwQaZiMgR\nKR2y6NOnDwDgm2++0TLp8v/xxx9aduLECY1lcN8maXnvvfcAAB9//HHUv71161aNZcjDDlnEc2JJ\ndSR5kletWuX7ubyq+yV1yUb2Nduua5et+H73QQ72BYDhw4drvGvXLgDeoTY7vLFjxw4A3oM5s2HI\nQoZoPv30Uy2TfQFWSUmJxtE+X36nB1mSMxwA/vnnn5DPZYK1S5fg+aT24Nqnn346qnpEiz1kIiJH\npLSHfOTIEQDe3WRTp04F4F3q8+WXX2rctWtXAN7kP/Gkc7SJSWzPmGIzePBgAEDDhg19P1+wYEEq\nq5N28+bN09imHN27d2/Ya+zSSrsjcsqUKQC8b4uWJMPq1KlTXHV1VXFxMQDgxRdfDPnMTtxLuljA\nO1n3f7///rvG0n4A/omEwu0EFDI5vX//fi174403Kr0mEewhExE5gg0yEZEjUjpkceHCBQD+u+HW\nrVvne409JSBW9qSLJUuWaOw3ZBHtusbqyB5Yunz58pDP7Wu3TExVF/b5nDBhgsZPPPFE2GskVzIA\nHDp0SGO/oYr27dtrLBNddlI8U9mJyRdeeCHs9+wwkOxkjMSegGOH1mQyX9bKA0BeXp7G0j6lE3vI\nRESOYINMROQIZ5ILJcOiRYs0lvXM4fhtcaV7hgwZorEkfpLjoQBg+vTpKa+TK+z6drv655133gHg\nXYUhZGs/4M0l/cwzzwDwri22+blllUW0+b5dJikTAKC0tBQA0LNnz0qvsdvx5R5t3rxZy2S7tX02\njx07prHcNztM5MIwhcUeMhGRI7KyhyxHe9tdOJaclDF58mQts7tvKLgjDwDmz58f8rmdyGN6z3vs\nhOawYcMAABMnTtQymbi2CYXsRKD0Gn/44Qctsz08ScqVDc+q7cX27dsXgHcCTib9bPpN2xu2BxRX\nxu9tYvfu3bFVNoXYQyYicgQbZCIiR2TNkIVd2ymJbx588EHf78rrkt0Oef369STWLvOsXr1aY/va\nKNtXI03AVHeS1GbFihVaJkmz7LZqSbgFADNmzAAAdOjQQcvset1t27YB8KYWyAYyWWdPR/nqq6/S\nVZ20Yg+ZiMgRbJCJiByR8JBFvXr1NO7YsSMAoHHjxlp2/PhxjSsqKhL9OQDAypUrAQC9evXSMllZ\nAQCPPPJIyDV2jWhBQQGAYO5eCho7diwA4NVXX/X9XA6XPHXqVMrqlCnsundZOWCPK5OsZjbvrl2T\nLPmOBwwYoGU2D7hkg6PEdO7cOd1VCIs9ZCIiR+RIwo2ovpyTE/Jlu45STvCwWfptIo8rV64A8K7D\nrFu3LgBvrlKbU1eS2cgBkADQtGlTqU+l9T169KjG9tDEWE4aiVUgEIgr2bLfvU0V+0YjeadtcpzT\np09rLIlb7L/XFCoLBAJdIn/NKx33Vp5Xm7RK8hjbhDl+7Brd/Px8je2a5CSI694C6X12I7GHxAp7\nIpFN3pRkUd1f9pCJiBzBBpmIyBEJD1mUlZVpbF+vEmFfM2rUiO6/GXbSrrCwEEBwizSQugmRTByy\nWLp0qcYLFy4M+XzMmDEab9y4MRVVCidjhiz8yFDbnDlztMxOnsozbNeAh8sTngRZOWRh2wAZFuWQ\nBRERRZRwD7lBgwYaFxUVAQgu+UmW8vJyAN6j6Hfu3Kmx7LqzO/XOnDmT1DqJTOohy7JBewJFnTp1\nAHgn7eySwmiPX0+SjO4hOy4re8g2zafsLrUHp7Zu3Vpju8QwCdhDJiLKJGyQiYgckfBOPTuZNn78\neADAwIEDtcwmT0mEXbssr9h2vaafVA1TZBI7xLR27VoAwWEK66OPPtI4zcMURHGzE3gyZGHXgUda\nE55q7CETETmCDTIRkSMSXmVBXq6vsvBbN26fAVm5MnPmzFRUJ1ZcZZE8WbnK4tFHH9VYDje1Q3R5\neXka23zMScBVFkREmSRrTgyh8GTyDggmuLHWrFmjsaM9Y6K42MRYubm5aaxJdNhDJiJyBBtkIiJH\ncMiiGrD5qW2yppKSEgDA4sWLU14nIgrFHjIRkSO47K2Kub7sLcNx2VvyZOWyN4dw2RsRUSZhg0xE\n5IhYJ/UuAKhIRkWyRF7kr4TFextZvPeX9zYyPrvJFdX9jWkMmYiIkodDFkREjmCDTETkCDbIRESO\nYINMROQINshERI5gg0xE5Ag2yEREjmCDTETkCDbIRESO+A+M1cj+cCZTOwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10bf67a90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def plt_digit(x):\n",
    "    nrow = 28\n",
    "    ncol = 28\n",
    "    xsq = x.reshape((nrow,ncol))\n",
    "    plt.imshow(xsq,  cmap='Greys_r')\n",
    "    plt.xticks([])\n",
    "    plt.yticks([])\n",
    "\n",
    "# Convert data to a matrix\n",
    "X = mnist.data\n",
    "y = mnist.target\n",
    "\n",
    "# Select random digits\n",
    "nplt = 4\n",
    "nsamp = X.shape[0]\n",
    "Iperm = np.random.permutation(nsamp)\n",
    "\n",
    "# Plot the images using the subplot command\n",
    "for i in range(nplt):\n",
    "    ind = Iperm[i]\n",
    "    plt.subplot(1,nplt,i+1)\n",
    "    plt_digit(X[ind,:])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Trying a Logistic Regression Classifier\n",
    "\n",
    "To classify the digits, we will first use a logistic classifier.  We select a small number of samples for training. Generally, you would use more training samples, but the optimizer is very slow.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "y = mnist.target\n",
    "ntr = 5000\n",
    "Xtr = X[Iperm[:ntr],:]\n",
    "ytr = y[Iperm[:ntr]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we run the optimizer.  I have placed the `verbose=10` option so that you can see the progress.  It may not appear in the browser but in the command line where you launched jupyter notebook.  This can take several minutes and will likely say that it ran out of iterations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/yao/anaconda3/lib/python3.6/site-packages/sklearn/linear_model/logistic.py:718: UserWarning: lbfgs failed to converge. Increase the number of iterations.\n",
      "  warnings.warn(\"lbfgs failed to converge. Increase the number \"\n",
      "[Parallel(n_jobs=1)]: Done   1 out of   1 | elapsed:    6.5s remaining:    0.0s\n",
      "[Parallel(n_jobs=1)]: Done   1 out of   1 | elapsed:    6.5s finished\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,\n",
       "          intercept_scaling=1, max_iter=500, multi_class='multinomial',\n",
       "          n_jobs=1, penalty='l2', random_state=None, solver='lbfgs',\n",
       "          tol=0.0001, verbose=10, warm_start=False)"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "logreg = linear_model.LogisticRegression(verbose=10, multi_class='multinomial', \n",
    "                                         solver='lbfgs', max_iter=500)\n",
    "logreg.fit(Xtr,ytr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Before continuing, since it takes so long to run the optimizer, let's save the results in a file.  You can use `pickle` module for this purpose."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pickle\n",
    "\n",
    "with open( \"mnist_logreg.p\", \"wb\" ) as fp:\n",
    "    pickle.dump( [logreg, X, y, Iperm],  fp)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then, we can recover the objects via the `pickle.load` command"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "with open( \"mnist_logreg.p\", \"rb\" ) as fp:\n",
    "    logreg, X, y, Iperm = pickle.load(fp)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, we can measure the accuracy on the test data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuaracy = 0.863077\n"
     ]
    }
   ],
   "source": [
    "Xts = X[Iperm[ntr:],:]\n",
    "yts = y[Iperm[ntr:]]\n",
    "yhat = logreg.predict(Xts)\n",
    "acc = np.mean(yhat == yts)\n",
    "print('Accuaracy = {0:f}'.format(acc))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should get an accuracy of around 86%.  This may sounds OK, but it is actually very poor.  Had we increased the number of training samples, it would have improved to about 92%, but that is still very bad.  To illustrate, let's plot some of the errors.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAACbCAYAAABhy5/wAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAFYRJREFUeJzt3XuQVGV6x/HfCwooVxHUiMYLF+W6\nUKgolkhQ3Gy8rhEXIYIEL0gNSyEU8VIqqEERQ4QIUYNam8UFhRQagoIKARVNZUULLyAqBMQFRMRh\ncWS4OG/+OMdy3nN6prtn+nJ63u+nqgue0+97ztPd75x++vTb5xhrrQAAAHzVqNgJAAAAFBPFEAAA\n8BrFEAAA8BrFEAAA8BrFEAAA8BrFEAAA8BrFEAAA8FpRiyFjzFZjzKXFzCHKGDPSGGONMTfnYd3W\nGNOpDv2aGGMWh8+XNcYMzHVupSxJ48gYM8gY874x5s/GmC3GmFvzsI3V9Rmfxpjxxpj/M8ZUGGM2\nGmO65DK/UpSwMdTbGLPOGPND+G/vPGyjTvuisO+xxpi5xpg9xph9xpg3c51fqUrYOLrSGPOxMeZ7\nY8w7xphuedhGnR+vMeZpY8wmY0yVMeamHKeWtUQfGTLGHFXg7R0n6S5JnxRyuxl6W9LfSdpV7ERK\nTaHGkTHmaElLJD0lqbWk30iaaYz5RSG2n4mwiBot6XJJLSRdIWlPUZMqAQUcQ00kvSxpvqTjJP1O\n0svh8qR4WlJbSV3DfycUN53SUcBx1FnS85LGSGojaamk/yz0e2oa6yWNlfR+sRORJFlri3KT9HtJ\nVZIOSPpe0mRJp0uyCnbWX0p6U9JASV9F+m6VdGn4/0aS7pS0WdK3kl6U1LaOOT2p4MVZLenmWtq1\nlvSMpJ2S/iTpIUmNw/s6SVojaZ+CN5kXwuVvho+tIny8v6ljjl9JGlis1y1ptySNI0knhts9ttqy\nP0q6oYb2NW5TUjMFb4jfSioP13OipH+U9KOkyvDxPpFFfo0kbZd0SbFftyTdEjaGLgv3Kabasi8l\n/XUN7Qu6L5J0lqQ/S2pV7NctabeEjaMyScuqxY3CvFL+7UtqKumxMMevFbwXHhPe107Sf4X7ob2S\n3grXF3u8dXze3pZ0U7Ffv6IdGbLW3qjgib/SWtvCWvtotbsvVvCp45cZrOq3kq4J+5ws6TtJc366\n0xhTXsvtzmrtzpN0joJBkM7vJB1RsLPpo2AH9tPXFg9Kek3Bp7pTJP1L+HgHhPf/Iny8Lxhj/jJN\nfsMyyMVrSRpH1tqvJS2QNMoY09gYc4Gk0xT8sWe7zZEK3uhOlXS8gk94B6y19yjYGZWFj7cszO/D\nWvKbG67zlPDWwxizPfyqbKoxJtFHiPMtSWNIUndJH9rwXSL0Ybg8lULvi/pJ2iZpavg12UfGmL/N\n4Llp8BI2jkx4UyTuUcM2p0vqIqm3grHUQdJ94X0TFXwIb6/gA9ndwcNN/Xgzfc9NmiQdMqtuirW2\nQpKMMena3qbgjeGrsP0USV8aY2601h6x1rZJtwJjTGNJcyWNs9ZW1bZNY8yJkn4lqY219oCkCmPM\nP0u6VcHXI4cVvAGeHOZU0xuhrLVfKjiEifwo6DgKLZA0T9KsML7dWrs9220qGEfHS+pkrf1Q0rra\nNmqt7ZVBbqeE/14mqaeCsfeagh3dv2XQ30eFHkMtFBzJqW6fpJbRhkXaF52i4A31PxS8UV8gaZkx\nZoO1dmMG/X1V6HH0uqRHTDC/9B1J/yCpiaRjow1NkNAtknpZa/eGy6ZJ+oOCaSOHJf2FpNOstV8o\n+DBWoyz2lYmS1E+ENb15pHKapCU/VZ6SNir4GuHELNYxVsGnsXcz3N7RknZW2+ZTkk4I75+soAL/\nX2PMJ8aYv88iD+RWQceRMeZsSS9IGqFgx9Nd0mRjzOV12ObvJa2QtNAYs8MY82g4J6k+DoT/Pmqt\nLbfWblUwdv+mnuttyAq9L/peUqvIslaS9tewvULviw4oeHN8yFp7yFq7RtJ/KyiwUbOCjiNr7acK\nji4/oeAr1HaSNij44BPVXkGRtK7aNpeHyyVphqQvJL1mgh+FJPboTn0UuxiyGSyvULVqNjyK077a\n/dsl/cpa26barZm19k9h++9rud0druMSSb82xuwyxuyS1F/SPxljnkiR23ZJByW1q7a9Vtba7pJk\nrd1lrb3FWnuyggp/rqnhVxvhoena8hue5vlDICnjqIekTdbaFdbaKmvtJknLFHx6T6XGbVprD1tr\np1pruykYj1coKLJSPt7wza6m/H766neTpEO1PF8+S8oY+kRSL+MePuil1D/qKMa+6MManicEkjKO\nZK1dbK3tYa09XtL9CoqsP6bIbY+CIrd7te21tta2CNez31o70Vp7pqQrJd1hjLmkpsebaX6JU9uE\nonzfJP2PpFurxacreHKPqrastaQfFPz65WgFL+oR/TzZbIKCCc+nhXF7SVdnmUcbSSdVu70j6Q5J\nrWto/7KCr0FaKSgoO0q6OLxviKRTwv93VzDIzgjjXZIuq+Nz1VTBpNqvFHwKa6Zqkyx9viVoHHVU\n8Ml+kIJP5B0VfKK6pYb2NW5T0l8p+CqrsYJf7KxXOMlQ0kJJ0+r4XP27gsmQLRV85fGppNHFfg2L\nfUvQGGqiYE7O+PBvviyMm9TQvqD7ovBxfyHpXgXTLC5UcNTq7GK/hkm4JWUchf36hvuP9gqOWP+h\nlrazFEzUPiGMO0j6Zfj/KxTMIzIK5jDuVPgjnujjrcNYbyZprYKv6ZpJalS0167IA+dqBROwyiVN\nSjVwwnY3hS/A7rDdVrkz7+9Q8Kl3v4IZ+HV6o6i2vdVK/2uyf1VQmOyT9IGkoeF9jyr4Vcf3YS7V\n/zDGhI+jXNL1Wea0NXxuqt9OL+brl5RbksaRpOslfRyu4ysFExNT/oHXtk1JN4TLKxT8umP2T49H\nwTyNzxRMrJydZX6tFBRT+xV8Ar1PFNVJG0N9FMwRO6DgZ8d9amlbjH1Rd0nvhmNzg6RfF/v1S8ot\nYePo7bD/XgVfnzavpW0zSdMkbVHwa8GNkn4b3jchzK8iHGf31vR4s8xvteLvaQOL9dqZMCkAAAAv\nFXvOEAAAQFFRDAEAAK9RDAEAAK9RDAEAAK9RDAEAAK9ldTkOYww/PfOAtTbt+eLrijHkjT3W2vbp\nm9UN48gP7IuQAxntizgyBCAfthU7AQBQhvsiiiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1\niiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEA\nAOA1iiEAAOA1iiEAAOA1iiEAAOA1iiEAAOC1o4qdQENw0UUXOfHKlSudeNGiRbE+w4cPz2tOAEpH\n27ZtY8uWLVvmxG3atHFiY0ysz9SpU514xYoVTrx37966pgg0aBwZAgAAXqMYAgAAXqMYAgAAXmPO\nUJYaNYrXj8OGDXPio45yn9aLL744rzkBKC3ROUJr166NtenSpUut60g1Z2j+/PlOvGnTJieeN29e\nrM/MmTNr3Q78MmTIECd+/PHHnbhZs2axPh999JETDx482IkPHz6co+zyhyNDAADAaxRDAADAaxRD\nAADAa8Zam3ljYzJv3EDNmDEjtmzixIm19hk3blxs2Zw5c3KWU65Za+OTEXKEMZRax44dnXjAgAFO\nfPbZZ8f6fPrpp078/PPPO/GhQ4dylF2drLPWnpOvlZfaOLr22mud+IEHHnDirl27Zr3OVHOG0u3P\njxw5Elv25JNPOvH48eOzziVf2BcV3v79+504Ok+2oqIi1qdp06ZOfM899zjxE088kaPs6iSjfRFH\nhgAAgNcohgAAgNcohgAAgNcohgAAgNc46WIa/fv3d+KxY8dmvY4tW7bkKh2UgD59+jhxdEJqv379\nYn3OOOMMJ27SpEnW242u47777st6HciPxYsXO3E2P1zJpegJYSWprKzMiXfs2OHE06dPz2tOyJ8O\nHTo48auvvurE0YuKS/H908GDB5148+bNsT7HH3+8E5966qlZ5ZkEHBkCAABeoxgCAABeoxgCAABe\nY85QRM+ePZ14xYoVTnzMMcekXUdlZaUTb9u2rf6JIZFuueWW2LLoCcaOPvrotOvZs2ePE48cOTJt\nn1mzZjlxqgsoIv9uvfVWJ540aVKsTaoTJNZXPtYpSZMnT3bi5cuXx9qsX78+L9tGbr3yyitOvGTJ\nEieOnnBTkqqqqpx4zZo1Thy9kKsUv1Drt99+m1WeScCRIQAA4DWKIQAA4DWKIQAA4DWKIQAA4DWv\nJ1BPmTIltmzYsGFO3KJFCyfO5GRp0av+btiwIfvkkEgnn3yyE0+bNi3WJjphevfu3U4cHWOStGrV\nqqxziZ5ALV8GDx7sxKlOCLls2bKC5FJsEydOjC2LjoFUJzaM7jfyddLFXKw3+vp26dIl1oYJ1MUX\nPcnqggULYm2iP95ZvXq1E+/cuTPW591333Xidu3aOfHGjRuzSbNkcGQIAAB4jWIIAAB4jWIIAAB4\nrUHPGWrevLkTR+cITZgwIdanUSO3PiwvL3fiVCe2a9q0qRO/9NJL2aSJBLvhhhuceMaMGU7csmXL\nWJ8XXnjBiaMn5YvOKSukrl27OvEFF1wQaxP9uzjppJOcOPocSA13zlDbtm2d+Oabb461STVHqJRF\n95upTsy3aNGiQqWDUHScPf30006c6uKoF110kRNHT6j4xhtvxPqcf/75Tjxv3jwnPnLkSPpkSxBH\nhgAAgNcohgAAgNcohgAAgNcazJfd0e9GJem5555z4jPPPDPtej755BMnfuqpp5x49uzZadfBOThK\nU6rv3J999lknjs4PGzp0aKzPiy++mNvEMjR8+HAn7tWrV6zNbbfd5sSHDh2KtZkzZ44TP/zww2n7\nNBQdOnRw4uiFmlOdc6dYOnfuHFsWHcNz58514o4dO8b6pJvzlOpCw9Ftf/7557WuA/XXr18/Jx44\ncKATP/TQQ7E+o0ePduJRo0Y58bnnnhvr89577znx+PHjs0mzZHFkCAAAeI1iCAAAeI1iCAAAeI1i\nCAAAeK1kJ1APGjTIiZcvXx5rk25i4L59+2LLxo0b58Q9evTIOrfFixdn3QfFd+mll8aWRSdMR+Vr\nsnTv3r2dODrxWZJGjBjhxNETghpjYn2iF2G87LLLYm0qKioyzrOhiU6Ij56kspCmT5/uxN98840T\nb968OdYnuqxbt25OPGvWrFifsrKyWvOInoRRik+qTbcO5F7jxo2d+P7778/JeqNj5IcffsjJepOO\nI0MAAMBrFEMAAMBrFEMAAMBrJTtnqG/fvk5cl4sltm7dOrZs1apVdc4Jpe2zzz6LLYtelDAXF+WM\nzneTpBtvvNGJR44cmXY90Tlv06ZNc+L3338/1ufVV1/NJEWEUs27KtR67rzzTieOvnbRi0pLUvv2\n7Z148uTJWW83k1zHjh1bayzF87388suzzgU/q6ysdOLoBZ9T7ZuiJwkeM2aMEx977LGxPvPnz69r\niiWNI0MAAMBrFEMAAMBrFEMAAMBrxlqbeWNjMm9cYFdddVVs2VlnnVVrn+OOOy62bMiQIU7cqVMn\nJ/7xxx9jfe666y4nnjFjRq3bTTprbW4mSqSQ5DGUytKlS504Ou/hmmuuifUpLy934rvvvtuJU53b\nJ2rv3r1OvHbt2libCRMmOPGWLVvSrreA1llrz8nXyvM1jr7++msnbteuXU7WG52Hk81+N9N15mu9\nmawzk3PRtGzZMutc2BfVLHoB3bZt28baRMdz9FxVqS4qnup8ayUuo30RR4YAAIDXKIYAAIDXKIYA\nAIDXKIYAAIDXGswE6nxZt26dE6ealN2iRYtCpVMQTFr8WefOnZ041YkMo6IXTI36+OOPY8uiF1l8\n6623nPi7775Lu92EKckJ1FVVVU6ciwnJUsOfQH3ttdfGlkUn677zzjtZ58K+KLf27NnjxI899lis\nzSOPPFKodAqFCdQAAADpUAwBAACvUQwBAACvleyFWvPl+uuvd+I+ffo4cSYnF0Npil7EUJIefPBB\nJ27evHna9ezcudOJoxddfeONN+qQHQph9+7dThy98GkxRXNLNWeoWPmmmge3efPmImSC6m6//XYn\nbtWqlRN/8MEHhUwn0TgyBAAAvEYxBAAAvEYxBAAAvMacoYjevXvXen9lZWWBMkEupZrrs3z5cic+\n//zzY20aN27sxPv27XPiVOcUeu2115yYOUKlY9CgQU68cuVKJz7hhBMKmY5j//79adskaY4Tiq9v\n375OfOTIESfesGFDIdNJNI4MAQAAr1EMAQAAr1EMAQAAr1EMAQAArzGBOmLIkCG13r906dICZYJc\nSjWJuV+/fmn7LVmyxIlHjx7txAsXLqxfYkiU6ITSAQMGOPGmTZvqtN5UJ0jMVqdOndKuMx8Xal2/\nfn2szezZs52YEywm04gRI5x47dq1Trx9+/ZCppNoHBkCAABeoxgCAABeoxgCAABe83rO0NVXXx1b\n1rFjRyc+dOiQE0cv3InSkOqEitH5FV988UWszdChQ5348OHDTpzqJHy7du2qS4pIoOjrneqkq02b\nNs16vbmY25Ov9a5evdqJr7vuulibvXv31ns7yK0LL7wwtiw6XqNzHvEzjgwBAACvUQwBAACvUQwB\nAACveT1naOrUqWnbRM/DsGXLlnylgzxas2ZNbFn0HDLR87hI0o4dO5y4qqrKiVu3bh3r88wzz9Ql\nRSTQ1q1bnfi5556LtRkzZkyBsqm/8vLy2LIFCxY4cVlZWaHSQQ717NkztqxRI/d4x8GDBwuVTsnh\nyBAAAPAaxRAAAPAaxRAAAPAaxRAAAPCaVxOoBw8e7MTdu3dP2+f111/PVzoooFQnjrviiivS9ote\nzDU6oXbRokWxPkyyb7juvffe2LKZM2c68aRJk2JtBg4c6MT5uKCqJE2ZMsWJ33vvPSc+cOBArE/0\nRwIoDc2aNXPi6A9CpPgYiU6oxs94ZgAAgNcohgAAgNcohgAAgNdMNt9dG2Pyc3XBAnn77beduH//\n/rE227Ztc+Ju3bo5carv3Bsaa218MkKOlPoYQsbWWWvPydfKGUd+YF9Us+hFolNdIHrhwoVOPGzY\nsLzmlFAZ7Ys4MgQAALxGMQQAALxGMQQAALxGMQQAALzm1UkXmzRp4sSVlZWxNmPHjnViHyZMAwBK\ny6hRo9K2SXUCUKTGkSEAAOA1iiEAAOA1iiEAAOA1r+YMnXfeecVOAQCAemvTpk3aNuXl5QXIpGHg\nyBAAAPAaxRAAAPAaxRAAAPCaVxdqRWa4OCJygAu1ot7YFyEHuFArAABAOhRDAADAaxRDAADAaxRD\nAADAaxRDAADAaxRDAADAaxRDAADAaxRDAADAa9leqHWPpG35SASJcVqe188Y8gPjCPXFGEIuZDSO\nsjoDNQAAQEPD12QAAMBrFEMAAMBrFEMAAMBrFEMAAMBrFEMAAMBrFEMAAMBrFEMAAMBrFEMAAMBr\nFEMAAMBr/w/yQ4xVU6zQSgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10bf81400>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "nplt = 4\n",
    "Ierr = np.where(yts != yhat)[0]\n",
    "\n",
    "plt.figure(figsize=(10, 4))\n",
    "for i in range(nplt):        \n",
    "    plt.subplot(1,nplt,i+1)        \n",
    "    ind = Ierr[i]    \n",
    "    plt_digit(Xts[ind,:])        \n",
    "    title = 'true={0:d} est={1:d}'.format(yts[ind].astype(int), yhat[ind].astype(int))\n",
    "    plt.title(title)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Obviously, some of these digits are very easy to classify for a human.  We can get more fine-grained analysis of the digit errors by computing the confusion matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.928  0.     0.006  0.007  0.003  0.032  0.011  0.005  0.009  0.002]\n",
      " [ 0.     0.969  0.007  0.005  0.001  0.003  0.001  0.002  0.016  0.   ]\n",
      " [ 0.015  0.013  0.832  0.046  0.02   0.011  0.015  0.014  0.029  0.006]\n",
      " [ 0.006  0.006  0.033  0.843  0.003  0.051  0.004  0.011  0.038  0.013]\n",
      " [ 0.003  0.005  0.011  0.005  0.878  0.004  0.013  0.015  0.012  0.049]\n",
      " [ 0.016  0.006  0.016  0.042  0.02   0.767  0.026  0.006  0.064  0.014]\n",
      " [ 0.013  0.004  0.027  0.002  0.021  0.029  0.894  0.001  0.011  0.   ]\n",
      " [ 0.003  0.007  0.012  0.014  0.012  0.005  0.     0.893  0.007  0.052]\n",
      " [ 0.005  0.02   0.033  0.043  0.011  0.078  0.013  0.01   0.77   0.017]\n",
      " [ 0.005  0.008  0.005  0.017  0.049  0.011  0.001  0.053  0.018  0.834]]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x10c035e10>"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAScAAAD8CAYAAAA11GIZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAE1JJREFUeJzt3WuMXVd5xvH/47EdY+cCiUEttkMc\n6qRY6cUwMoRU3BKKAyiRKlolFFQQar40EC4tDVCFKlWpSimXDwF1moRWkJIWE1UWcjG3RGqrEjyx\no4BtIlyHxpNLcw9pQmLPnKcfzglMJjNz9sR7z1ln9vOTtnT2OdvveTX2vF5r7bXXkm0iIkqzbNAJ\nRETMJsUpIoqU4hQRRUpxiogipThFRJFSnCKiSClOEVGkFKeIKFKKU0QUaXkTQU8+eZk3rB+pPe4d\nPzih9pjRpWXN/D/lTqeRuG33JI9zxE/pWGK86fVr/OBDU5WuveW2p3bZ3nYs37dQjRSnDetH2Llz\nbe1x33Xqb9UeEwAd09/x4mrocaNlq9c0ErfzsycbiYsbKHpqqCPRqVYAFuJmf+eYYzz40BTf33Vq\npWtHfvnH9f9C99FIcYqI8hnoUG7LNsUpoqWMOer6W3V1SXGKaLG0nCKiOMZMFbxkUopTRIt1SHGK\niMIYmCq4OFW6dyppm6TbJR2UdHnTSUXE4ujgSscg9G05SRoBrgLeCEwAuyXtsL2/6eQiojkGjhY8\n5lSl5bQVOGj7kO0jwPXAhc2mFRFNM2aq4jEIVYrTOuDwtPOJ3nvPIOkSSeOSxh98qNzbkxHRY5iq\neAxCleI027Mdz0rX9pjtUdujp5yc54kjStedIV7tGIQqd+smgA3TztcDdzeTTkQsHjE1a9ujDFWK\n025gk6SNwF3ARcDbG80qIhrXHRAf4uJke1LSpcAuYAS41va+xjOLiEZ15zkNcXECsL0T2NlwLhGx\nyDrD3HKKiKVpSbScImLpMWKq4JW6U5wiWizduogojhFHXP9a/3VJcYpoqe4kzJZ16+74wQmNbEaw\n6+5ba48J8KZ1WxqJO0w6jz8+6BQGTsubaUW4gQ0O6pIB8Ygoji2m3LKWU0QMh05aThFRmu6AeLkl\noNzMIqJRrRwQj4jhMJV5ThFRmswQj4hidXK3LiJK033wN8UpIgpjxNE8vhIRpbHJJMyIKJEyCTMi\nymPKbjmVm1lENG6KZZWOfiRtk3S7pIOSLp/l81Ml3Shpr6TbJL25X8wUp4iWMqLjasd8JI0AVwHn\nA5uBiyVtnnHZnwH/YnsL3R2cPt8vv3TrIlqquzVULSVgK3DQ9iEASdcDFwL7Z3zdib3XJ1Fh78sU\np4jWWtCmmmsljU87H7M91nu9Djg87bMJ4JUz/vyfA9+U9F5gDXBevy9McYpoKbOgGeIP2B6d47PZ\nKpxnnF8M/IPtv5V0NvAlSWfZnnO38xSniBaraSXMCWDDtPP1PLvb9h5gG4Dt/5K0ClgL3DdX0AyI\nR7SULTpeVunoYzewSdJGSSvpDnjvmHHNncC5AJJeBqwC7p8vaFpOES3VHRA/9sdXbE9KuhTYBYwA\n19reJ+lKYNz2DuBDwN9L+kDvq99le2bX7xlSnCJaq741xG3vBHbOeO+Kaa/3A+csJGYzxUmg5fWH\n3nbqXONxx+b3D/yk9pjXn/1rtccE6DzxRCNxfXSykbgaaejB0rnHUZ8zrWjm18GTRxsIWk+IbKoZ\nEUXKkikRUZynZ4iXKsUposWywUFEFMeGo50Up4goTLdbl+IUEQWqaYZ4I/qWTUkbeuuwHJC0T9Jl\ni5FYRDTr6akEx7pkSlOqtJwmgQ/Z3iPpBOAWSd/qTaqKiKE15N062/cA9/RePybpAN0lElKcIobc\nkllDXNJpwBbg5iaSiYjF071btwS2hpJ0PPA14P22fzrL55cAlwCsYnVtCUZEM5bEJExJK+gWputs\n3zDbNb1V8cYATlx2cg1P/kRE04a6WydJwDXAAdufbj6liFgMpT/4W2Wo/hzgncAbJN3aO/pu6xIR\n5atpsblGVLlb9x/MvkZwRAwxW0wO81SCiFi6Su7WpThFtFTpY04pThEtluIUEcVZEvOcImJpGup5\nTs+JwVNTjYRuwldecWbtMd+9d0/tMQGuOWNjI3FHTjm5kbhTDz3cSFxU/12mZQ1sygENbfJQxwYH\nhsksNhcRJUq3LiKKkzGniCiWU5wiokTtGxCPiOLZGXOKiCKJqdyti4gSZcwpIoqTZ+siokzujjuV\nKsUposVyty4iiuMMiEdEqdKti4gilXy3rtw2XUQ0yu4WpypHP5K2Sbpd0kFJl89xze9J2i9pn6R/\n6hczLaeIFqtjKoGkEeAq4I3ABLBb0g7b+6ddswn4CHCO7Yclvahf3LScIlrMrnb0sRU4aPuQ7SPA\n9cCFM675Q+Aq2w93v9f39Qua4hTRUkZ0OssqHX2sAw5PO5/ovTfdGcAZkv5T0vckbesXNN26iBZb\nwM26tZLGp52P2R7rvZ6tbzgz9HJgE/A6YD3w75LOsv3IXF+Y4hTRVl7Q3boHbI/O8dkEsGHa+Xrg\n7lmu+Z7to8Adkm6nW6x2z/WF6dZFtJkrHvPbDWyStFHSSuAiYMeMa/4VeD2ApLV0u3mH5guallNE\ni9Uxz8n2pKRLgV3ACHCt7X2SrgTGbe/offbbkvYDU8Cf2H5wvrjNFacmpp6qmQljnZ89WXvMa848\nvfaYAB87tLeRuH/50i2NxG1k55GGdJ54opG4y046sfaYeuTYf64GOp16fqds7wR2znjvimmvDXyw\nd1SSllNEWxkoeIZ4ilNEi+XZuogoU4pTRJSn2nNzg5LiFNFmaTlFRHEMruluXRNSnCJardziVHmG\nuKQRSXslfb3JhCJiEdUzQ7wRC3l85TLgQFOJRMQADHtxkrQeeAtwdbPpRMSieXoSZpVjAKqOOX0W\n+DBwwlwXSLoEuARgFauPPbOIaFzJkzD7tpwkvRW4z/Yt811ne8z2qO3RFRxXW4IR0aCOqh0DUKXl\ndA5wgaQ3A6uAEyV92fY7mk0tIpqmYW452f6I7fW2T6O7Tst3U5giloCqg+EDKmCZ5xTRWoMb7K5i\nQcXJ9k3ATY1kEhGLr+BuXVpOEW3WGXQCc0tximirLDYXEaUq+W5dilNEmxVcnLI1VEQUqZmWk4RW\nrKw9rCeP1h4TmtkhZNma42uPCfCJl72ykbi/ueepRuLednYz/8Q6T9Wf78iLXlh7TIDOgw/VHtNT\nU7XESbcuIspjBvZoShUpThFtlpZTRJQo3bqIKFOKU0QUKcUpIkojp1sXEaXK3bqIKFFaThFRphSn\niChOxpwiolgpThFRIhW82FxWJYiIIqXlFNFmBXfr0nKKaCv/YiJmv6MfSdsk3S7poKTL57nubZIs\nabRfzBSniDarYd86SSPAVcD5wGbgYkmbZ7nuBOB9wM1VUktximizejbV3AoctH3I9hHgeuDCWa77\nC+CTwJNVUktximgp0b1bV+XoYx1weNr5RO+9X3yXtAXYYPvrVfPLgHhEWy1sEuZaSePTzsdsj/Ve\nz/aA3s8jS1oGfAZ410LSS3GKaLPqxekB23MNYk8AG6adrwfunnZ+AnAWcJMkgF8Cdki6wPb0gvcM\nKU4RbVbPVILdwCZJG4G7gIuAt//8K+xHgbVPn0u6Cfjj+QoTNFacDG5g6qmbmZShFfX/GKYeebT2\nmAA67rhG4t52zupG4v7VgRsbifunG+vfhWbq/gdrjwlAp56dUp6hpl+FOp6tsz0p6VJgFzACXGt7\nn6QrgXHbO55L3LScItqspiJneyewc8Z7V8xx7euqxExximgrl/1sXYpTRJsV/PhKilNEi2U9p4go\nU4pTRBSn2qMpA1Pp8RVJz5e0XdKPJB2QdHbTiUVEs0R9qxI0oWrL6XPAN2y/TdJKoJlJMRGxqIZ6\nzEnSicBr6D0X03vq+EizaUXEoii4OFXp1p0O3A98UdJeSVdLWtNwXhGxGOpZMqURVYrTcuDlwBds\nbwEeB5610p2kSySNSxo/6qdqTjMialfjSphNqFKcJoAJ20+vXredbrF6Bttjtkdtj65QM89/RUTN\nhrnlZPte4LCkM3tvnQvsbzSriFgUNS0214iqd+veC1zXu1N3CHh3cylFxGIZ6rt1ALZvBfrulhAR\nQ6TwSZiZIR7RZilOEVGap2eIlyrFKaLF1Cm3OqU4RbRVxpwiolTp1kVEmVpXnAyenKw/rmbbu+/Y\n+Uj9zzFreUN1f6qBnTyATgM/A2hmlxSAT9zx/dpjfvT0ZnIdOeXk2mPqkZF64rSuOEXEcEhxioji\nZPeViChR5jlFRLka2kW7DilOES2WllNElCeTMCOiVBkQj4gipThFRHlMBsQjokwZEI+IMqU4RURp\nMgkzIspkZ7G5iChUubUpxSmizUru1lXZ8TciliIDHVc7+pC0TdLtkg5KunyWzz8oab+k2yR9R9JL\n+sVMcYposxq2I5c0AlwFnA9sBi6WtHnGZXuBUdu/DmwHPtkvtRSniBaTqx19bAUO2j5k+whwPXDh\n9Ats32j7id7p94D1/YJmzCmixWq6W7cOODztfAKYb83j9wD/1i9oilNEWy1sVYK1ksannY/ZHuu9\nnm1x/1kjS3oHMAq8tt8XNlecGtiMQCtX1h4TmtmMQCua+dFO/fT/Gom7fN2LG4k7de//NhL3oy89\nu/aYr9hztPaYAHu21v935qljf2K3OwmzcnV6wPboHJ9NABumna8H7n7W90nnAR8DXmv7qX5fmDGn\niDbrVDzmtxvYJGmjpJXARcCO6RdI2gL8HXCB7fuqpJZuXUSLLaDlNCfbk5IuBXYBI8C1tvdJuhIY\nt70D+BvgeOCr6vaq7rR9wXxxU5wi2qrGlTBt7wR2znjvimmvz1tozBSniNbKs3URUaosNhcRxcmm\nmhFRrIJbTpWmEkj6gKR9kn4o6SuSVjWdWEQsghqerWtK3+IkaR3wProP7Z1F91bhRU0nFhHNU6dT\n6RiEqt265cDzJB0FVjPL7M+IGDKmygTLgenbcrJ9F/Ap4E7gHuBR29+ceZ2kSySNSxo/St+Z6REx\nYMLI1Y5BqNKtewHd5Q82Ai8G1vQe3nsG22O2R22PruC4+jONiPrZ1Y4BqDIgfh5wh+37bR8FbgBe\n3WxaEbEoCi5OVcac7gReJWk18DPgXGB8/j8SEcUrfMypb3GyfbOk7cAeYJLucptj8/+piBgGg7oT\nV0Wlu3W2Pw58vOFcImJRDa7LVkVmiEe0lUlxiohCldurS3GKaLNBzWGqIsUpos1SnCKiODbUsFFC\nU4aqOPnIkaGJq+Urao8JMPKCkxqJO3n3vY3EpTPVSNiRtafUHnPP1sdqjwnwzh/+d+0xf/w7NT0i\nlpZTRBQpxSkiimMga4hHRHkMzphTRJTGZEA8IgqVMaeIKFKKU0SUJw/+RkSJDAz7kikRsUSl5RQR\n5cnjKxFRIoMzzykiipQZ4hFRpIw5RURx7Nyti4hCpeUUEeUxnmpmva06pDhFtFXhS6ZU2Y48IpYq\nd6odfUjaJul2SQclXT7L58dJ+ufe5zdLOq1fzBSniJYy4I4rHfORNAJcBZwPbAYulrR5xmXvAR62\n/SvAZ4C/7pdfilNEW9l1tZy2AgdtH7J9BLgeuHDGNRcC/9h7vR04V5LmC5oxp4gWq2lAfB1weNr5\nBPDKua6xPSnpUeAU4IG5gjZSnB7j4Qe+3fnq/1S4dC3zJFeg6vk2s1HMQn5aS/dnC3B/c4lUsKBc\nv3VGIzm85FgDPMbDu77t7WsrXr5K0vi08zHbY73Xs7WAZvYFq1zzDI0UJ9svrHKdpHHbo03k0IRh\nyneYcoXhyneYcp2P7W01hZoANkw7Xw/cPcc1E5KWAycBD80XNGNOEXGsdgObJG2UtBK4CNgx45od\nwB/0Xr8N+K49/wzQjDlFxDHpjSFdCuwCRoBrbe+TdCUwbnsHcA3wJUkH6baYLuoXd9DFaaz/JUUZ\npnyHKVcYrnyHKddFYXsnsHPGe1dMe/0k8LsLiak+LauIiIHImFNEFGlgxanfdPdSSNog6UZJByTt\nk3TZoHOqQtKIpL2Svj7oXOYj6fmStkv6Ue9nfPagc5qPpA/0/h38UNJXJK0adE5L1UCKU8Xp7qWY\nBD5k+2XAq4A/KjjX6S4DDgw6iQo+B3zD9q8Cv0HBOUtaB7wPGLV9Ft3B374Du/HcDKrlVGW6exFs\n32N7T+/1Y3R/edYNNqv5SVoPvAW4etC5zEfSicBr6N7JwfYR248MNqu+lgPP683VWc2z5/NETQZV\nnGab7l70LzxA70nqLcDNg82kr88CHwbKXeaw63S6c72/2OuCXi1pzaCTmovtu4BPAXcC9wCP2v7m\nYLNaugZVnBY8lX3QJB0PfA14v+2fDjqfuUh6K3Cf7VsGnUsFy4GXA1+wvQV4HCh5/PEFdFv4G4EX\nA2skvWOwWS1dgypOVaa7F0PSCrqF6TrbNww6nz7OAS6Q9BO63eU3SPryYFOa0wQwYfvpluh2usWq\nVOcBd9i+3/ZR4Abg1QPOackaVHGqMt29CL1lHa4BDtj+9KDz6cf2R2yvt30a3Z/rd20X+b+77XuB\nw5LO7L11LrB/gCn1cyfwKkmre/8uzqXgAfxhN5AZ4nNNdx9ELhWcA7wT+IGkW3vvfbQ3IzaO3XuB\n63r/SR0C3j3gfOZk+2ZJ24E9dO/i7iWzxRuTGeIRUaTMEI+IIqU4RUSRUpwiokgpThFRpBSniChS\nilNEFCnFKSKKlOIUEUX6f5mIKgZZX6ziAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10bf64e10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn.metrics import confusion_matrix\n",
    "C = confusion_matrix(yts,yhat)\n",
    "\n",
    "# Normalize the confusion matrix\n",
    "Csum = np.sum(C,1)\n",
    "C = C / Csum[None,:]\n",
    "\n",
    "# Print the confusion matrix\n",
    "print(np.array_str(C, precision=3, suppress_small=True))\n",
    "plt.imshow(C, interpolation='none')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualizing the Weights\n",
    "\n",
    "To see the problem with the logistic classifier, it is useful to plot the weights for each digit.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAAE0CAYAAAAv9jd/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJztnXecVdX1xdeJvcbeYuwaCyoq9oId\nFbtRsSCoERsaLFHsXRE1sYsFGyqKgqLYERuIDVDBbgSjxobYwJLE3N8f4PmtvcI8mWFm3hvu+n4+\nfLIn+82b++4957zj2WvvnYqigDHGGGNMWflNtS/AGGOMMaaaeDNkjDHGmFLjzZAxxhhjSo03Q8YY\nY4wpNd4MGWOMMabUeDNkjDHGmFLjzZAxxhhjSk0pNkMppQVSSvemlCallD5IKe1b7WsyDSOl1DWl\n9HJK6aeU0s3Vvh7TcFJKs6WUek+Zk9+llEallLav9nWZhpFSui2l9ElK6duU0jsppT9V+5rM9JFS\nWjGl9GNK6bZqX0tTM3O1L6CZuArAvwAsCqA1gAdTSq8WRfF6dS/LNIB/AjgXQDsAc1T5Wsz0MTOA\nDwG0BfAPADsA6JdSWr0oinHVvDDTIC4AcHBRFD+llFYG8FRKaVRRFCOqfWGmwVwF4KVqX0RzMMOf\nDKWU5gKwB4DTiqKYWBTFUAD3A+hY3SszDaEoigFFUdwH4MtqX4uZPoqimFQUxZlFUYwriuK/RVEM\nAjAWwDrVvjZTf4qieL0oip9++XHKv+WreElmOkgpdQDwNYAnqn0tzcEMvxkCsBKAn4uieIf+v1cB\nrFal6zHGTIWU0qKYPF99YttCSSldnVL6HsBbAD4B8FCVL8k0gJTSvADOBnBcta+luSjDZmhuAN/I\n//cNgHmqcC3GmKmQUpoFwO0AbimK4q1qX49pGEVRHIHJa+umAAYA+Knyb5ga5RwAvYui+LDaF9Jc\nlGEzNBHAvPL/zQvguypcizFGSCn9BkAfTNb1da3y5ZjppCiKn6fIEZYEcHi1r8fUj5RSawBbA/hb\nta+lOSmDgPodADOnlFYsiuLdKf/fmvBRvDFVJ6WUAPTG5OSGHYqi+HeVL8k0HjPDmqGWyOYAlgHw\nj8nTE3MDmCmltGpRFGtX8bqalBn+ZKgoikmYfFx7dkpprpTSxgB2weT/EjUtjJTSzCml2QHMhMkT\ndPaUUhk29TMq1wBYBcBORVH8UO2LMQ0jpbRISqlDSmnulNJMKaV2APYBMKTa12bqzXWYvIltPeVf\nLwAPYnIG7wzLDL8ZmsIRmJyG/TmAvgAOd1p9i+VUAD8A6A5g/yn2qVW9ItMgUkpLAzgUkxfcT1NK\nE6f826/Kl2bqT4HJIbGPAHwF4GIA3YqiGFjVqzL1piiK74ui+PSXf5gsNfmxKIovqn1tTUkqiqLa\n12CMMcYYUzXKcjJkjDHGGDNVvBkyxhhjTKnxZsgYY4wxpcabIWOMMcaUGm+GjDHGGFNq6lWfZd55\n5y0WXnjhproW8yt88cUX+Pbbb1NjvNfcc89dLLjggo3xVqYBfPnll5g4cWKjPEsAmGuuuYoFFlig\nsd7O1JMJEyZg0qRJnpszAI05N+eZZ55ioYUWaoy3Mg1k3Lhx44ui+NWNS702QwsvvDB69uzZ8Ksy\n08UJJ5zQaO+14IIL4sQTT2y09zP148ILL2zU91tggQVw7LHHNup7mmnnr3/9a6O914ILLoju3bs3\n2vuZ+tGjR49Ge6+FFloIZ511VqO9n6k/nTp1+mBaXufKvcYYAMDMM8fl4Le//W22f/op9tucddZZ\n6/RNmjSp3n+7Ur2z2WabLfz8r3/9q87fm9I+wJhSwvPhhx9iQfc55pgj2w2dJ/z+Ovd+85uWrbpp\n2VdvjDHGGDOdeDNkjDHGmFLjzZAxxhhjSk2L1QxV0hj8+9//zvZ///vfbLPOAYgxTs3EmTBhQp1/\ny7qEpoOfFwCMGzcu26xpWXXVVcPrWrVqle2JEycG31dffZXtDz/8MPh+/PHHbKs2ZaaZZprGq265\n8L1ijRAQnwXPKQD4+uuvs61zZ/XVV8/2J598Enw851jT8I9//CO8jp97mzZtgm+ppZbKNuuHAOCz\nzz6b6vUDwM8//5xtncNleNZ1UUkHwj+rJqSudVDvO/+sv8M+1ax5na0/n3/+eba///774GPNkN5r\nnkecyajzgseDrgn8LGeZZZbg42ep38Nzzz036kLX8qbEJ0PGGGOMKTXeDBljjDGm1FQ9TMZHbWPH\njg0+PkbXo3j+WY/k+OfZZ58923rkxz4OpQDAO++8k+211lor+PjoTo/xll566Wzzsbz+DT0C5vtQ\n6YixpaDH7RyO0nvGR7R6z9jHKdwrrLBCeN2wYcOyzSFOZeTIkeFnDo3xswOA3/3ud9nWo90ZhTnn\nnDPblcarPs+PP/4429ttt13wcajqyy+/DL4NN9ww2xxC47AbAGy00UbZ1nvP1/z+++8HH38Gvg4A\nmG+++ab6HkD8fC1xvimVPsN//vOf8PMaa6yR7Ztvvjn42rdvn+3NNtss+F555ZVsr7POOtnWlO6+\nfftmW8sw7LHHHtl+6KGHgo/fp23btsH39ttvZ5vXFqDlp3j/GhzuWnnllYPvo48+yvZzzz0XfBx6\n1jHA5TD4/i255JLhdXPNNVe29XuKv181hMbzctFFFw2+SmU6+Pf0uTb2c56xR40xxhhjzK/gzZAx\nxhhjSo03Q8YYY4wpNVXXDLEmZO211w4+1ix89913wbfMMstkWzUoHJNcaaWVsj1+/PjwOk6zfuml\nl4KP03nnmWeeOq9f4+PffPNNtjVVm3/m1wEx1XFG0CxoHP+LL77Itmq3ll122WxvueWWwffII49k\ne4MNNsj2KqusEl7Hep8///nPwde6detsL7LIIsF30kknZfu0006r8z2//fbb4GOtio5bTS2uZTju\nrzoC1tGonojH/QcfxNY/PB/13nC6L+uJtt566/A61j6wfg+IWgjVTLCOQDVgr732WrZV78BjUH0t\nRYPCz4t1QECcj/p55p133mxvtdVWwcfaPNaJAcCYMWOyzbqP0aNHh9fxs1Vt5quvvpptLskATO6F\n+QtPPPFE8PF41PWfr2X++ecPvpbyLBn9PmAdj641PId1DLzxxhvZ1nR21n99+umn2eaxAcS1VBvQ\n/vOf/8y2an/4GekzYd+KK64YfKxRevPNN4OP9w6N8Z3Z8kaGMcYYY0wj4s2QMcYYY0pNs4fJNOWO\nw0N8NK4+DVM8++yz2dawy+abb55tDkcNGDAgvI7Ta3fdddfg4zCcVtHlv6fHjXzEqGE5DtlpKrFW\nBG2J8JGtHoVymExTKzlcMnjw4ODjY21O0dVqyfz+nBIPANtuu2229UidK7ZuvPHGwcfp3e+9917w\n8bGs/j0OG7WkY3kND/G90SNqvqeLLbZY8D3//PPZ1hIInN7L82/EiBHhdRyyvPPOO4Nvxx13zPbj\njz8efJwirCGZQw45JNuPPvpo8HFVa62GzSEZTcmvJhoa4NIBHH4CYohewye8tr7wwgvBx2Hrrl27\nBt8FF1yQbS5Zsdxyy4XX8Vq69957B98xxxyTbS2FwGEYDZ9zKI7nPhDnrcoY+L5waZVqU6kyt5aW\n4J9VPsLhKa3O/oc//CHbXLEfiN9bfC3rr79+eB2vg7vttlvwXXvttdneZJNNgu+pp56a6vsDce7p\n9yB/Vt0f8Hcof7dO7W9MCy1npTbGGGOMaQK8GTLGGGNMqfFmyBhjjDGlplmEKpzyqSnXHNPVOF+l\nLtOsPdIYOOs0OJ58wAEHhNf17t072/fee2/wcVonl5oHou6E00sBYPHFF8/2euutF3ysl3j33XeD\njzUwqjPRFOFaQVs0cLxa06053qtp8RzH5/gxEFPhK+l7ODWUXwdEjdLAgQODj/UNqiPhkgqaks+l\nFzTFlFNYK5VlqDVUh8XPjHV4QLynt99+e/Bxeq/qb1gDwNoO1QOceeaZ2dZWPKwhUi0Xa1C6desW\nfHwtOqf4mWlJDNYhqc6kmpowHVt///vfs82fBwA6duyYbdU9seZN11Ju/aCtT1hjxi1wVHvWp0+f\nbF966aXBx6UzjjrqqOC78MILs63zj3Vqqgt6+eWXs60lG4YMGZJt1frVEvxsdYzxdw7rgICob+VU\ndwB48skns73ffvsF3x//+Mdsczr99ddfH17Ha37Pnj2Db4kllsi2atZ4beHvViCWtWnXrl3w8R6A\ny18AUd/m1HpjjDHGmOnEmyFjjDHGlJpmz+fmipJADJvxsTkQU641vX3NNdfMtqbW8xExH2trlWmu\ntsrHwUBMNdT0eQ7RaYd0/nya2supo1rZk49C1VeraMiT74XeT65ey5W/gfi8uEs2EO8n/z1N2ebQ\nmB71n3rqqdnWLuocEunevXvw8XHuTTfdFHxcOkDLRXBKq4YMaimdF4hjW691tdVWy7aGL/m1WiaC\nj+A1ZZ6rSXM16gUXXDC8jkMYJ5xwQvCtuuqq2db15K677sq2ViIfNGhQtjXFu1evXtnW0g8cJtNK\n3M0dJuMwxe9///vg42epawiHgy+55JLg43nFXeSBGLri+w7E0ASHZDhcAsS1dJdddgm+6667Lts6\n/vg5aLjr1ltvzfbYsWODr0uXLtnmchxALD9Qy/D6ydILIEoQtNo3p7u//vrrwcdhYn3PffbZJ9tX\nX311tjVEzWFpLanBa8RDDz0UfLyWaBiaZQb33Xdf8HH4brPNNqvz77399tvBpxKOacEnQ8YYY4wp\nNd4MGWOMMabUeDNkjDHGmFLTLJohTrVVvQh3ndbu1Bz305LkrBniLtZAjEFyil/nzp3D6ziuqeW8\nWYfE2gb9WdPnOWVetUasi9B4P7cw0O7htRTn5mei6YzcJkFLwbO2Yvnllw8+7kbPqdH6c//+/bOt\n92/YsGHZZj0BEDUSmk7Lab8cNweAhx9+ONuqS+P4u94HjrNrfLzaaCydNTdaImD48OHZ5nRlIGo/\nNOWVyxyofoTvGz/Dm2++ObyO/95ee+0VfHy/NU2cdV7PPPNM8G244YbZ5nY+QEw71rR7vi86Dpob\nnldvvfVW8PGaqO1TuGwEp9kDUduo7RC4i7im6/M6xd3MVUPGa6K2KGJdl3amZ00I644A4NBDD822\n6lP473EqPRDbeujYaYz07MaCx6B+/9x///1TfR0QtUbaFmXnnXfOtpYj4HZUPK70u/Wggw7K9t13\n3x18PHZUs7bWWmtl+5RTTgk+Hju6XrAuSb8H+TpVA6VlBaYFnwwZY4wxptR4M2SMMcaYUtMsYTI9\n5mM4LVmPfbmisR7v85GwHq3xsS+n0OrxLYfetGoxH/1z53u9Fk2t5fQ/rY7L4SGt9svhPA2v1dLx\nLacWa7doTtHVDtTchblSV/Cjjz46+LjjNVcQ12q4Z599drZPPvnk4OOjck3h5rGp3dE5hKZHzvw+\nGibi6rG19OyA/51HXA5Bw5f8Wk4vB2IKu94bPiLXlFdOz37ggQeyrSnXHArQI31eM/h1QFwX9Kie\nn4V2++aO9nocz+NTy2U0N/xMNBTQoUOHOn/vhhtuyDaHvoC4FmlYhDva6/jg9+HyFSx9AGIIVMub\nDB48ONv8DADgwQcfzLaWpOA1ksPzQFxDeP0AYiXrWpqbOi95/dTQIocrdaxyRW+dN9yBnmUZQJyX\n/Fw1zM9VzrXkBX9/a0V5HmM6jli+wmsCECUO+r3B5U04BNhQfDJkjDHGmFLjzZAxxhhjSo03Q8YY\nY4wpNU2iGdL4J8d7OcUOiCmgGmseMGBAttddd93gY42NxlRZC8RpztpJnWPGqnFhtKw5x141HZlT\nb7V8OMc4Of0YiCnB2mKgmmjLjYkTJ2ZbtVQc063UZVrTfrlFBqfMAlEfwvF/1YPws9WuyNye47zz\nzgs+Ho/atZ7HnL4n64K0jQCnP9da+w0toc/jVcc5v3bgwIHBx92sOS0diGnx2lWbf2ZNlo6XVVZZ\nJds6b7fbbrts69zfdNNNs81dz4HYakJLafCz19YjrC/STvGait7U8N9T3eEtt9yS7W233Tb4WBNy\n2mmnBR/fC10HOWVZy2WwFo+fpWp4evTokW0ticE6p6FDhwYf6wJVE8IaT9UP8nPnkhtAbemEGF1D\nuMWIPhNuUaHd4Vnfpu/JGk/VUnE5El4HVOP1t7/9LdtanoLXUi5jAQBXXXVVtnV+XXTRRdlWbSJ/\nT2r6PK/rqnXTvzEt+GTIGGOMMaXGmyFjjDHGlJpmOePlsJWmNnOanaatclhJu6Bz2uDmm28efBy+\n4SN87W7M4RkN7a299trZ5o67iqbW85GmHivzMZ+GBThNVjtjVxMOiwH/e4xZF9rJfdSoUdnWZ8mp\novfcc0/wceo7p8xqRXIOv2iYkVPy99xzz+DjVFEdm/yzVh7mUgF69K6Vbafma+5u57+g956rvms6\n7AorrJDtfffdN/i46q9WB+Yw8h133BF8V155Zba5u71Wi+aUYR1LfM06j/iIX7usb7LJJtk+7LDD\ngo9DbzrmK1XJ11BVU8P3Qktw8LqhoRXtbl7Xe2p16o8//niqNhDDoXwfdB1fcskls62p2iNHjsy2\nlj7hEJCGYkeMGJHt4447Lvg4PVv/XrXm3dQYN25ctrU8B/s++eST4OP7yRWhgbgW6e9xiFfHMVff\n5/CXyhZ4XGkobNCgQdnWzgyHH354trUkBJf3aNOmTfDxOqCyDA776VxvCLUzMowxxhhjqoA3Q8YY\nY4wpNd4MGWOMMabUNIlmSDUUlcrgc+lvjpMCMW1ddTSVujezZoRTm1UTwroTTXXn91h00UWDr2/f\nvtnW9HnWxmj6N3dM17RHjo1qm4tqxrlVY8LaKtXmsLZK9QzceVzvC8eezz333ODbaKONss3l+VXb\nwKm2L774YvC1atUq26xTAYDLLrss26ov4M+uz4s/n37WSs+r2l3PWSMFxHII+vm5+7eWxGC9nc7N\nnj17Zrtt27bBx9ou1pxo923WRaivf//+2dbnwuNCu6xzl23VmTz11FPZ1hIfrFXgTu3VgD+vpvVz\nirlqsHgs33jjjcHHqeg6PjkVntcvIOq1eA6o/oV1ho899ljw7bTTTtnmsghAXHd1HHFpB27FA8R1\nXls1VRPVRPHc0zWRNV5aJoF1lFr6hPWKmhbfvn37bGurDv6u5Wd5/vnnh9ddcMEF2dbPw5pY/dus\nDVbtbvfu3bPN6woQ5/7jjz8efPxanZcNKXnhkyFjjDHGlBpvhowxxhhTapokTMZpx0A8wtIO9hxS\n0CrT3PFaQ0edOnXKNqfrArGq6ZZbbpntc845J7yOu5Rr+ITDLpzGCcTjXK2mzCl+etzOKbsakvjp\np5+yXUvpnxoK47Cj+rhKM4e0gJjmqVVo+UhYw3J8L/jIlJ/Pr8HdsDXVm5+DplTz0bV+Vj6GrU9V\n21/CxDpHmgs9TuY0V60yzSFgDTdwN3MNU3NqPYdHgdhVW4//GT4C1w7bnO577733Bh8/J/08Rxxx\nRLZ1/vE1azj9vffeq/P3mhueDxoiYQmChi45JHn77bcHH5cc0HWWw1oa/uJyC5yer1XX+ZlfffXV\nwcddBri6MQB06dIl21qmhOePljfgysQ8ToH/LaHSnOjaxuVbtEQEV8rXMPGYMWOyrWUnmN133z38\nzKUKWDoAAC+88EK2uTp169atw+v4vmsJES5xo2G4wYMHZ1vn/RNPPDG1ywcQQ3E6NrlUgI65hoSz\na+db1xhjjDGmCngzZIwxxphS0yxhMq5Qq9ksq6++era1wiQ3pNOKqvyemu3FDSb5CJ+zY4CYOaFH\nfnyEqWEArrDK2XBAzJbhCtpAPCbVMA8f82kosZbge6bhLlb7a9YIH9tr5kSvXr2yvf/++wffiSee\nmG0+JtVquMstt1y2tZkmZ0TokTOHu/TImY9oNWuJwxXVzhCrDzo3ee5ccsklwccVqZ9++ung48yi\niy++OPjuv//+bOvc4XvFITqtcM33/uGHHw4+rp6rc59DWnvvvXfwbbzxxtnW6rxcuZrDRkDMZOKQ\nTzXgcafry/bbb59trQretWvXbFdq7snSBCCGFrUZJl8Lr59ciRgANthgg2zvtttuwXf88cdnW9c9\nDqGpHIGfCVecVp/O6WqurRpi5Z/1ujgrWdceDodylXggzm+d68OHD8+2NjLl7yOWO+j8YvmDhuG4\nabJWs+eQuzZX5mbqOv547mvTZ85eboyuDT4ZMsYYY0yp8WbIGGOMMaXGmyFjjDHGlJom0QyphoL1\nPlrJllMmNc7N8UpNS+S0Uk6zB2Lc8U9/+lO2Nc2SY6qaBs+xStWxcOVqTffjeK6morLWQWPZHMPl\nLr61TCUtlcZ3uRyBVmzeaqutsq1p6qwj47gwV6YGYrqpaimGDh2aba76CsT4OMeugZiyyxVagZjC\nrVXVNQ24llA9FVd61o7RO+64Y7ZVM8TaP02XZm3eLrvsEnyccrvzzjtnW7tos65FtRasr2PdABDH\nlur5hgwZkm0dB6wl07/HOiRd2xpDq9BQtNwDj0MthXDddddlW8cnlzTR92T9pVYHvueee7LNY0d1\nJqwL1efMawhXAQeA++67L9uqM+R1Qv8eV9/mNHsgrvmVOiU0N1otnat2c8o6EEuA6Hcmr4P6nclz\nkUstADHVniuZa+Xv3r17Z1t1tvxdqGs8z2etDr3SSitlW58X7x10XeexqZrihuCTIWOMMcaUGm+G\njDHGGFNqmiRMpnAaroZPONyl6a581KbH0Zy6x01OgXiszc3jNH2Rq+Nq+uIBBxyQbW0gy2EyTV/k\nz6NVOD/44INsa8VMPiqsNnxcrBVb+fNWChFoNWc+Dtcjbw5RcHowENNB2aelCXr06JHt008/Pfj4\nyFnDHBwa4+rhQDyG1THGYTMeY0BsUqghl2pUF+fnpOUQeBxyRXYghje0VAJ/Zr3fnBavVZL59/jY\nW9Oju3Xrlu0zzjgj+DidWD/PpZdemm2tTMyVb7lxJBDHK6fZA3HcaXituSuJ8zqlaxavL/r5+Dnz\nuAbimjxs2LDg45C9/h6XPuFSFxyqAWKVcE3pZukCl2QAYuhIS2Jw1WJdC7g0gj4fDvXp/atmmEw/\nA4e/VIrBMhAN7XO5Cg5DA8BFF12UbQ15cqVpvhZupAsAnTt3zrbKJLi8ia6zHGZnaQIQyybwdzcQ\n5RAqsdHyKtOLT4aMMcYYU2q8GTLGGGNMqfFmyBhjjDGlplk0Qxz3U10Qp59rKmylFg6sJ9IYNcf/\nOWVe06NZS8K6JiBqGzSOzqXFtUQ9p6lqyjV30lUtjqb9V5Np7eys8XhOg9dnwvoQvkdATHnWWDbr\nPH766adsjxgxIryOtSPcegCIacY77LAD6kJ1XDzG9LOylkJ1MaxV0fGuHdGbAy7joOX7edyrBoB1\nBKoBY02AltDn1HrVYXCqLsf8tVQCd7nmjvJAvPd9+/at85pVf8Cfb5111gk+1o9oaxAed5oyzOUk\nmgNOZ2ZtGhDHL69RAHDXXXdlW9dB1pmoJopLKNx4443Bd+655071WlTHwmuwpsGzLkhTrrkMhM4/\nfkYPPfRQ8K233nrZZs0hEOeC6vmqiWqpeD17/vnng4+fn5ZJ4DIlrCEDgFVWWSXbXFIDiJqsV155\nJdujR48Or+Nr0e/rnj17ZvvAAw8Mvn322Sfb+py5pAJ/JwPA2LFjs63lAPh7aptttsH04pMhY4wx\nxpQab4aMMcYYU2qaJEymIR8+mtRjZQ5HaViJU+GffPLJ4OOO1P379w8+Dqnxcb6m9PHruEu9+hZZ\nZJHg4+N2TicE4tGuHvnx36+Url/L8NGkVtF+7bXXss2hBSAet2tXcD6m1xRu7nC84YYbZlsrKfMx\nrB4d873VMbDmmmtm+8MPPww+TiXm0A8Q06/5PYD4WXVMVyNMxiEgfS5cpVbHOYcX9Ria56Mex2+x\nxRbZPuuss4KPxw8f1es84jGhXdb5mo888sjg41C0hsn4OF4rV/Nn4DABEENjmsrf3NXieWxrmQb+\nvBqO4rGsac8cgtKQJFf3147z3NGeq0X369cvvI6r+WsYjtP6TzvttODj8LOGW7ljeqW09FtuuSX4\ntHwGoyHQ5kTDVlyJWcexfq8wXBVcOzNwWFPvNae3s8RBuyhwJ4FHH300+HgNVukAh9c0RZ7Xcq2A\nzqFfLeXS2N+ZPhkyxhhjTKnxZsgYY4wxpcabIWOMMcaUmibRDGkMl+PzmnbJ8UJNdeS0SE335r+h\nKX4cZ+Q2F6ptYDRtlFMNtew364k0pqopyAxrH1Q7UsudzhnWKXDbBSC2b9Duw6z/Ut3FIYcckm3V\nfzGs9+E0fiC2Z9H0cS5nz5oSIOppdGxyerd2dGetm94Hvk79rNWANSKq2WBtiaZjsx5G4/VcKkFb\nJbC+Y//99w++UaNGZZvvDbfRAOI9VI0Za/Y4ZRyIc3zHHXcMPk7xHjhwYPBx2rF2++bSGjrumhte\nJ7SjOJehUH0PtzTRNkGcpt6lS5fgY82Xplnz/eWO6Nr26Lrrrsu2plxzF3tuF6Hvo9osTrlWnQnP\nYx3T/NpqtMZheF3Scgc8zvR5sX5Ry2HwOL799tuDj+cpa/KAuH7y95SWoOCSJvrdx1phLTnDJUYq\ntUjRli+8DnBZCeB/S6FMLz4ZMsYYY0yp8WbIGGOMMaWmScJkmrrJqYGa0sqpvhqm4GNMTaPjNEhO\nzwTicRof/+lxI4foNAzC16W/x9VCNXWaj/n4b+t7asXRSuG1WoKP6TV9ksNT2mmZQw9aNZjTqvVZ\ncoVaTvnUVPfDDz8825rqzmOOU0MB4Nhjj822hkPbt2+fbT2+5VCQpsVyan21wyqKhjC4W7ym3fNR\n+gEHHBB8rVq1yrZ2G+dq0pzmDMRO2hzq1iP3SiUqtttuu2zr5+G5qsfoHJo+++yzg2+ZZZbJtoZS\n+ahe1yH9fM2JliPgkKHOTb4vPOaBOAY4nAbEMCqXtgCAoUOHZpvH+R577BFed84552RbS6vwuqfj\niO+7hpG4DIQ+Zw7LaShRS5pUEx6PWjqAw19c0R6I90LDXRzG0nWQ7xmn4AOxpEKHDh2yrd9THB5X\neQBXoteyKyxX0dAbh+71OXMF8wCQAAAgAElEQVRpFYXHdGOURfDJkDHGGGNKjTdDxhhjjCk13gwZ\nY4wxptQ0iWZIO8BzDF5TJDmdXuP/rPXgFGggxlg1tfLFF1/MNsccVffBcUaNf3KsXFs4cOny+eab\nL/g4dq7aKS4Fr/qMloimVLPOQzvTs55B4/j8nDllXX/m8graooE1StpahUvBc0drIOq6tJuytudg\nll9++Wxrqiint2rH6Wqjc4DHoeqb+HmqVoZLUWy//fbB9+CDD2ZbNWDPPvtstlkXsd9++4XXvfDC\nC9nW+Tdo0KBsq86rY8eO2dY1g1P5NcWbP59qBHltqCUNGKeXA1G/xBooAOjevXu2OQUfAG666aZs\naysSbnvEWjggrmf8nHv37h1ex6UzeO0EYskG1V8dfPDB2eb0fwDYa6+9ss1tVoD42fV58XeRzvfm\nplJqP38vqs6Jv0NVK8OlBLRExDPPPJNtnVO8hrH2Z8CAAeF1rC3S7zDWY+rc48+jZXJ47HCaPRC/\nX1VX29ilEXwyZIwxxphS482QMcYYY0pNk5wT6pE6p9lpyjWHKTTtmdND9fiMj1S5uzYQU5258iWn\nBQIx9VCPhznt9/333w++zz77LNta/ZpTCDV8xyE7PYpvKXD5g4033jj4OAVaQ6UjR47MNlerBWK6\nvoZL+Mj7ggsuyLambvJRuYYPOP1UQ1+cPq/jg8OteiTM1Va1OrOmydYSWoGa0211vPKR+yOPPBJ8\nnNatx9ccQtxggw2Cj+cVh+Euu+yy8DoOU2tok1OpO3fuHHw8zjQMwuuSznc+juc5DMRnrSnrzQ2H\nM7SsB98XLUfAKfMaWuS19Pzzzw8+/uw77LBD8J100knZ5irhWj6Fw8bagV2/DxgupdGtW7fge+ml\nl7Kt1e45rKnlTTR0VKtwWEnTzXl+PfbYY8HH8oTzzjsv+Pg5aJkElhZcddVV2dZQFK/rWp6CZROa\nPs+hbe22wGEyXTs53KZdG3TdnV58MmSMMcaYUuPNkDHGGGNKjTdDxhhjjCk1TaIZ0rgfx+o1hsvx\nUE3349/T+DinTGpMlXUhe+65Z7ZZ6wMAL7/8crY1ds0xW/08HKNWbQzH9LVUOsfxNfWwpcDaGU3r\nZC2CllBgHdfrr78efPz8NE2bY9unn356trXkP7do0E7pnM6t+pbWrVtnu1evXsHH40X1RKxLqmWN\nkMIaECBqr1SXxz7VzXEHeG1HwunLOgf4nrIGQPWCPLZUA7bFFltkm9PlgTiP9bqefvrpbKtmgjVK\n3NKg1uB1Q3VsrCXRDt9c7oC1WkBMd9dyGfz3brnlluDjMg3897TcCOvUdM3g58Up+ABwxRVXZPuD\nDz4IPta4aDsO/nwtdZ3lNUU/Hz+jSvq2rl27Bh9rbHTsnHzyydlmbZjq7ljHxW1WAKBPnz7Zvuii\ni4KP10/VDPGarGUEWMfY1M/SJ0PGGGOMKTXeDBljjDGm1DRLCU4+atOwEh/dDR8+PPj4GFbT6Lga\n6muvvRZ8u+++e7bfeeedbGuIjlNKOeUXiOm8GobjcJem6HKoQTv+6vFxS4RDIM8991zwcUiEu80D\n8Rhdj8r5+F1TgjmFe5tttsm2pojz+NBqtZwOran7/fr1y7aGwrhkg4YdWioa0uPjcq0Az+ETDSNz\niLldu3bBx3NnyJAhdV4Lz1ud32uttVa2NXzHc45LAwAxfKrlKzgEq6E3roZdy3CoQLuZs3RAQwoc\nktR1iKu8a1iTw6Nnnnlm8B166KHZ5nmr5UY4rV9DYfwZdB3nZ8IhQCCuQxp20dIuLR1de7iCs37H\n8DzlsiEAcO2112Zb17pdd90127zuaSkSLnWi4TtOkdf1n8PQGvJ89913s73CCisEHz/bpi5H45Mh\nY4wxxpQab4aMMcYYU2q8GTLGGGNMqWkWzRCX9Nby3hzb1m62HK9UrQPHDzUdj7VHlToY88+ahssl\n3TUdmfUN+nk4FqvX3NhddquNdpnmz6cl+RmOEQOxxYdqyriMOz8v1ZFwOq22Z2Gf/m0ef5pyrNqm\nGQEdkzzn9PPyOFeNAesDtL3DNddck+155503+FifwmUo9LpYV6CtV1gXoRqzW2+9Ndv6PHl8qkaJ\n30evuZZg/ZvqnFhHw+0qFO48DkQ9h651vC7qe1588cVTvZZhw4aF1z3wwAPZ5mcORH0Pvx8QU8O1\ntENLKmcxvWgLGC5bos+Lta/qYy0lr6tAnPtcgkJLmIwaNSrbo0ePDj7Wden787XsuOOOwcdzj793\ngeZ9zjPWt7MxxhhjTD3xZsgYY4wxpaZZwmTTih6JcfVJPdbmtHutQltXN9uZZpop/MzH7/q3+bWa\nvsuv5aNpIB5h6pEwhyRm9GPeSp+Pj+WBeH+1MipXHueuy+PHj6/zdRoe4b+nz4vRaqtlgJ9TQ8sH\naHqv/lwXmhrO8Nz56quvgo+fp4b2pjW0qWtBLYfGmEqduvkzaTV/5sknn6zTp+FnLivCoVEglqVg\n2QLPUyCWYdC1lMNf+v4ttcp7Y6PyCg4l6Xjn78K+ffsGH3//tGrVKvj4e4u/d0888cTwOh5jGkIb\nOXJktm+++ebg47Id+pwrrRcqnWlKfDJkjDHGmFLjzZAxxhhjSo03Q8YYY4wpNTUtkuA4cX30HNP6\n2kpdcKc1Rq3vwaXuNU2Q9SottZtyY6D3tlKH5jnmmGOq76G6I/25LsqoC2qJsE6CW0mYpkXvNWs2\nNJW/Lg2klq+oBD/nMuuC6kOl+8R6ydatW0/ze/K6yDqkgw46aJr/dps2bbK99tpr1/l7tfqcfTJk\njDHGmFLjzZAxxhhjSk2qT+paSukLAB/86gtNU7F0URR158zWAz/LqtNozxLw86wBPDdnHPwsZyym\n6XnWazNkjDHGGDOj4TCZMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PG\nGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482Q\nMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1Hgz\nZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXe\nDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qN\nN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS\n482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm\n1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOM\nKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYY\nY0qNN0PGGGOMKTXeDBljjDGm1HgzZIwxxphS482QMcYYY0pNKTZDKaWnUko/ppQmTvn3drWvyTSc\nlFKHlNKbKaVJKaW/p5Q2rfY1mfpD8/GXfz+nlK6o9nWZhpFSWial9FBK6auU0qcppStTSjNX+7pM\nw0gprZJSGpJS+ial9F5KabdqX1NTUorN0BS6FkUx95R/f6j2xZiGkVLaBsCFAA4EMA+AzQC8X9WL\nMg2C5uPcABYF8AOAu6t8WabhXA3gcwCLA2gNoC2AI6p6RaZBTNnEDgQwCMACALoAuC2ltFJVL6wJ\nKdNmyMwYnAXg7KIoni+K4r9FUXxcFMXH1b4oM938EZO/SJ+t9oWYBrMsgH5FUfxYFMWnAB4BsFqV\nr8k0jJUBLAHgb0VR/FwUxRAAwwB0rO5lNR1l2gxdkFIan1IallLavNoXY+pPSmkmAG0ALDzl2Paj\nKUfxc1T72sx00wnArUVRFNW+ENNgLgPQIaU0Z0rpdwC2x+QNkWl5pDr+v1bNfSHNRVk2QycCWA7A\n7wBcB+CBlNLy1b0k0wAWBTALJp8ibIrJR/FrATi1mhdlpo+U0lKYHFK5pdrXYqaLpzH5JOhbAB8B\neBnAfVW9ItNQ3sLkk9q/pJRmSSlti8lzdM7qXlbTUYrNUFEULxRF8V1RFD8VRXELJh/37VDt6zL1\n5ocp/3tFURSfFEUxHsBf4WfZ0jkAwNCiKMZW+0JMw0gp/QbAowAGAJgLwEIA5sdkfZ9pYRRF8W8A\nuwJoD+BTAMcB6IfJm9wZklJshqZCgakfA5oapiiKrzB5MjqUMmNxAHwq1NJZAMDvAVw55T86vwRw\nE/wfKi2WoiheK4qibVEUCxZF0Q6ToysvVvu6mooZfjOUUpovpdQupTR7SmnmlNJ+mJyB9Gi1r800\niJsAHJVSWiSlND+Abpic8WBaICmljTA5fO0sshbMlFPasQAOn7LOzofJOrBXq3tlpqGklNaY8r05\nZ0rpeEzOEry5ypfVZMzwmyFM1picC+ALAOMBHAVg16IoXGuoZXIOgJcAvAPgTQCjAJxX1Ssy00Mn\nAAOKoviu2hdippvdAWyHyWvtewD+A+CYql6RmR46AvgEk7VDWwHYpiiKn6p7SU1HcvKGMcYYY8pM\nGU6GjDHGGGPqxJshY4wxxpQab4aMMcYYU2q8GTLGGGNMqfFmyBhjjDGlZub6vHjuuecuFlhggaa6\nFvMrTJgwARMnTmyUYpHzzDNPsdBCCzXGW5kGMH78eHz33XeNVvhznnnmKRZeeOHGejtTT7744otG\ne55zzjlnMd988zXGW5kG8PXXX+P7779vlGc511xzFfPPP39jvJVpIB9//PH4oih+dXGs12ZogQUW\nwPHHH9/wqzLTxcUXX9xo77XQQgvhzDPPbLT3M/Wjse/9wgsvjHPOOadR39NMO6eddlqjvdd8882H\nLl26NNr7mfpx3XXXNdp7zT///Dj66KMb7f1M/TnxxBM/mJbX1WszZAzz73//O/w8fvz4bM8222x1\nvlZPpGaZZZZsT2vdq5Tif7jNNNNM2f7++++Dj69F3991towxxlgzZIwxxphS482QMcYYY0qNN0PG\nGGOMKTUzpGbov//97zS9jvUiP/0U+8/NPvvs0/R7P//8c52+mWeOt5d1LS2VBRdcMNuff/558PHn\nnXvuuYNvrbXWyrZqel588cVsL7PMMsG32mqrZfuLL76Y6t8CgDnnnDPbTz31VPCxlmmuueYKvhVX\nXDHbqkP64Ycfsv3ll18GH+ucFOuQTEtiWserrqv/+te/ss26PJ1H0/q3K/3ef/7zn/Azr7uzzjpr\n8P3mN/5v/PpS13NQXWgl+JmoZrQ+Y6JaeNQYY4wxptR4M2SMMcaYUtMkYbJKYSo9wqz0Wj6GXWGF\nFYLvn//8Z52/t9lmm2WbQzBvvPFGeN0666wz1b8FAM8++2y2O3ToEHwcytHw2tNPP53tb7/9Nvj4\nKFLDcHPMMcdUXze1n5sTPd7kz6ShKi4Up0XjODz17rvvBh9/vrFjxwYfh6ratGlT59/m57fGGmsE\nH9/bzz77LPj4iH3kyJHBx3+vUoHKr776qs6fZ5SQmYZEeRzoWF5iiSWyrZ+fj9157jc01OGQyLTB\nz0FDTr/73e+yrWFqDmm//fbbwffaa69lm5/l+++/H163wQYbZJvXXAD49NNPs73SSisF3+jRo7P9\nwQexVAyH4XWt4ZC50hLCNU1FpbWIZQATJkzINs9lRe8lr8mVvsP092rlmXglMcYYY0yp8WbIGGOM\nMaXGmyFjjDHGlJpmT61XjRCnS3/88cfBN++882Z7ySWXDL4333wz24svvnjwsfbjnXfeyfbf//73\n8Lqtttoq23/4wx+Cj+PVzz33XPBxTPz3v/998HH6vKYXcky8UiuLH3/8MfhqSRfB8V2NC7Meh+PO\nQNSYfPTRR8G39tprZ/v1118Pvpdffnmq76+6gBNOOCHbPXr0CL5DDjkk27179w4+1p5xmj0QNUNa\nFmHYsGHZVm0Y68hUC1NrcDpspdIPlVJsWdcFxM+v+pTvvvtuqraWLmCNmc4x1iboNbNPdWXLLrts\ntnWt0RIZLZ0xY8aEn5deeulsc7kKIGp6HnnkkeDjdfbBBx8MPp6brC1SPeedd96ZbdWe8d/+7W9/\nG3ynnHJKtg899NDgY/2S/r1JkyZlW7V+Sy21VLZrSZvZWPBnqPT5Jk6cGHyffPJJtr/++uts63cP\nazP1u5y1mbo+c4N3nZf89/g7H6hcwqSxqZ1vWWOMMcaYKuDNkDHGGGNKTZOEySqFdbQCMKflcqgI\nAOaff/5scygFAN56661s6zEppwPuvffe2R41alR4HR/zDR8+PPj4M3Tt2jX4BgwYkG0+4gPi0bGm\nkfLxrR7vr7zyytnmEAwAzDPPPKgWetTKz4RDnEA8QtXjTT5u5yNTIB7ha5rs5Zdfnm2uYs3vBwBH\nHHFEtvV+cchAU3Q33XTTbH/zzTfB98ADD2Sbnx0QQzxaXoHDOhomqjU4JKtjko+611tvveB74okn\nsq0VxTlspiE0fhYcQjzssMPC6zhsdffddwffo48+mm0NRa+yyirZXn311YOPP5+GcXkdqpSaXctw\nSQcN+XLleJUV8OfV8frYY49le+DAgcF3xRVXZPvJJ5/MdufOncPr+N4ut9xywXfXXXdle9111w2+\ne++9N9t33HFH8HE5gHbt2gXfLrvsku0uXboEH48XDusA/zsmqgl/N+kz4XVKQ048nzWMxd+hQ4cO\nDb6DDjoo2xyquvjii8PrjjzyyDp9q666arZ79eoVfPwZNt544+Dj7xENa/LY1HA5l1dojPR8nwwZ\nY4wxptR4M2SMMcaYUuPNkDHGGGNKTbO04+D4v6bQctyvY8eOwcetLRZeeOHg43RmjQufc8452eb4\np2qL1lxzzWxrqi1rHTS1nuPvnK6rf0/1DKwh6tOnT53vybocoLq6E9WD8LPVdGSOGaumZ8MNN8w2\np3ECMZ1eWzsstthi2eby/zxuAOC9997LdqtWrYLvlVdeyfa4ceOCj3UCRx11VPDdeuut2daWGxyj\nZr0XELUwtaYZ0rm5yCKLZFs1Njwn+B4CcSzr3GRNimrHWIPCmhNt9cC6P9aUAcDxxx+fbdYuAVEz\n1LZt2+BjveKgQYOCj0sH6DU3Z3pvfdBnyanpqiVhXYamyHOaNbcvAqI+65prrgm+7bbbLts8BzSt\nn+fVZZddFnxbbrlltnWO8fzbaaedgu/mm2/OtpbL2HbbbbOtpS123XXXbOvYqSVY+6OlVvh7RfWl\nrInUe8YlTbS9Fd9rLimieqxrr70227vttlvw8Rho3759ndfFOjEAOOuss7I9ZMiQ4ONSD6pFXn75\n5bPdGHPUJ0PGGGOMKTXeDBljjDGm1DRLBep//OMf2db0Za4C/dRTTwUfH5Wff/75wcehJA6fAMDW\nW2+dbU7P1K7I/fr1y7am+3GoRbusczo9p9ID8fhRKxNzuG333XcPvmeffTbbmkLI96+50c/O4QRN\nrefjW60kzUevgwcPDr4ddtgh29o1m/8ep8xrGQYOpWjKLF+Lhvb22WefbN9+++3Bx2EADS1wmrZW\nTdXnzmj4oimo1BVaQyschtRr4xCDpkRfeeWVU7WBGMIcMWJE8HGl93PPPTfbmtbM3cw1tM7zm8cO\nEKvMayjsgAMOyLZW1N5oo42yreUX9J4xtVQdnu+Zlgfgsf3qq68GH4dDdR3s2bNntjWEsckmm0z1\nPVUewGVLuLI/EMeKlm/gZ6mVsTm1vlOnTsHHJVQuvfTS4OPwrkojqhkO1THG41PLerCUQMsksPRj\n++23r9OnMgYua8DfP2+88UZ43dFHH51tLm8DAGPHjs32PffcE3x77rlntnWM8fcBl08B4lqqawTf\nMy0B05BU+9qZycYYY4wxVcCbIWOMMcaUGm+GjDHGGFNqmkUzxLFZ1VdwzJPTo4Goo2GdABBTpFVn\nwnBMmtO7gRj/1DRE1hqphoc/D2sNgKhr0c7AnIKsaaSsXXnmmWdQK3D6IhDj1foZWHOiLTE45s9l\nBIB4X0aOHBl83GGbU4eff/758DpO6+d0YCBqX7RNwUMPPZTtJZdcMvg4FVVj+txShNuQAFFPVKmr\nenPBGg59LtxKRFvS8NxkzRcQy1loyQpOnWWdDhB1PK1bt862pvqylmmNNdYIPtaxaXo0awW22Wab\n4LvuuuuyrRoo1r/pGOQWMZrmX03NkOoiuCP8XnvtFXz8nFW3s//++2f7jDPOCD5uc8Md5gHgvvvu\nyzbrh1jrA0QN1n777Rd83MaD5zAQ9SOsNQOiLkjbLLHeRjVlvG5oqxgtadLUsM5FS3CwVkbHaiVt\nJref+uMf/xh8XFZAS5Pw/L7qqquyrfrLG2+8Mduq/eQyCbxuA1H3d//99wffFltskW0dOzzGWUMM\nRI2cfkdbM2SMMcYYU0+8GTLGGGNMqWmWM3s+gubQAxArxmqqHqfda3d4DptpWImrJnNqIFfPBIAe\nPXpkWztj83G4hhb22GOPbGu4ho8w9eifP6se0XKohY+7gf8P0Wj6YHPAVXuBWLFY05M5HKUpknxM\nyuFPIFYa16N4Ph7nSqz8O0AM8eh18XEuH8kCsYSChnCvv/76bGuYlsNd2rWeU3T1uPaXa2vOZ8mh\nHO0KzWE8LQnAr9USEksttVSdPk7r1mNvPurm9F4N61x00UXZ1vnH82qzzTYLPi6doKEAfvZaBZ3D\nwVqlnD9PLaXSK5yCrWsIp1Xrcz7ppJOyreGG008/Pdsa+mBpAYfM+G8BwHnnnZftP//5z8HXtWvX\nbGvqPnc+P/HEE4OP10sthcBjWtdg/f5hfpEAVOMZ63rAa52G77kaPH8XAcANN9yQ7Ycffjj4uFo1\nz18gSk24u712n7/tttuyrSn/HH7Vqtk831R6wb+n5RU4fKjfRZWoVA6jLmp3ZhtjjDHGNAPeDBlj\njDGm1HgzZIwxxphS02iaIdZpaLzu4IMPzrZ2v+YWFZrezpoC1WVw6W/9vaFDh2b7L3/5S7a13QK3\nEdAS4axD0q7k/BmOO+644OPWGVrunTsDa5oga4hYWwT8ryajOVHNC2t/vvvuuzpfq+nQ/Pm0HAHH\n/C+//PLgO+igg7LNcW3tVH3hhRdm+/HHHw++zp07Z1vT+vlZ6hjjUgybbrpp8LGOhdODgfhZtWv2\nL79XLc0Q6+mAqK/QtHGO82uZCB7nHTt2DL6rr74626rD4nY8rEF55513wutYC7HssssGH2sftO0E\npxqrLqJSCxVea7jUg/5eLaE6Sm5RwWnHQNQCaQo567q0Xc2ZZ56ZbdbQAcCHH36YbZ4PXH4DiOuu\n6jZZ18J6RCC2WeLvECBqCXXerr/++tnWecvrF+ucgP/XpFRDM6Sp9bxOqSaKxy7re4CopbrpppuC\nr3v37tnmMhNAfLZcNkTvO5fNUP0lt2TZdtttg491w/p5eI3U/QHPYdUA8nPS++fUemOMMcaYeuLN\nkDHGGGNKTZOk1usxI6eba0djDilwyAyIx6TagZer0D722GPBx8dpnDKvx9/HHHNMtrXy6vHHH5/t\n0047Lfi4kqemCXLVWz2i5XBNpc70/fv3D75fUrUbcvTX2HB4R0M9/LN2u+fq0ZoiyT/vvPPOwccp\nmhwC0eqqXJaB07KBGFJ7//33g4/TgDXsx2G58ePHBx8f6WsYhe+Dvueiiy4KoHmP4nmM6ufgcaiV\nfPm5aCVaPvZ+4IEHgo+P8XlcA7H67CmnnJJt7S5+5JFHZlur7PI81urXfMSvFdK5wr2mK/PRfaXy\nA7WEpi/z2qOV1jn8oKUtOnTokO0LLrgg+LhitMoReC3nCsZaNoRLmGgphF/mAwCstNJKwcdhHk0h\n33zzzbPN3xN6nRoS5DnNKf/A/1cz19IcTQWH9zTUN2bMmGwfdthhwcdrh3aV5/CUft/169cv21x1\nHIjhUK5izfMJiGE5fc5t27bNtoaoOaR66qmnBh9/32lVaw4Fa0ic74NW9m/I+uqTIWOMMcaUGm+G\njDHGGFNqvBkyxhhjTKlpNM0QtyDQ1HpOe9MUdo55cmdzAFh33XWzrWnqHIPcc88967wubpehLSL+\n9Kc/ZVtj2awT6NSpU/CxBoPL1QPAhAkTsq1xU9YsvPTSS8HH9081C6p9qCYcT9c2CZy2rdfMGg3V\nn/Bn1/fkNOpDDz0021pqnsfYU089FXycuqnaM45RH3300cHHegYu1wDEGPXqq68efJy++8ILL6C5\n0bINnKar9551BapJ43R61cZxR3vVIXFqtWoEDznkkGyzdoC7ngNRB6gl+lmT9Pbbbwcfa0R4XAFx\nXul18efRNN1aRTUUPP+4hAEQdUKq9eByFnrPLr744myvs846wcdp+HxvtcM8lzD5RZfzCzxXBg0a\nVOfvabkR1v5o6ZEjjjhiqtcIxPID+pybutyFlrVgnZBqC7k9jbaOad++fbZVz8r3XnU03OqInysA\n3HvvvdnmNfGKK66o87p0LWF9j2rPeO3WEg2s1+PvTyBqDnX88XdKY2gwfTJkjDHGmFLjzZAxxhhj\nSk2zdK3n9FcNAfHxuKaKcvqkph5yGian6wLxuJ9DJJq6z++hR5icLqzHjVxxmI9rgXj83rdv3+Bb\nbrnlsq3d2eebb75sa1fpWoKPQjmdFoghoT59+gQfh6e0mjFXMtWwUl3Hvlp9l5+Jpmcutthi2dbw\nHYd0NBWVn6Uet2uYheHX6ns2Bzp+OE1XQwEcDtb5x6Ex9XFojKveArGatJaQ4K7hd955Z7b1PnHq\nL1e2BYAHH3ww21rplt9fPyuHM7mCPRDXF322fI/q0zm7qdGUfw5hc0VtIIaftTo8v1ZLjPA81vAX\nzyUec1pWg8tZqBSCnxemxXQAABPxSURBVAmHUoBYPkPXbu4QsNVWWwUfrw36eyxd+Pzzz4Pvl7Wg\nqUqYaAiSQzs6rlgWoj6u4q0heh7XI0aMCD7+ztluu+2Cj8vRcMjshBNOCK/jcCuHOIFYUb5Lly7B\nxxIRlh8AscODlkJ48cUXs60lDzgUrL/XEHwyZIwxxphS482QMcYYY0pNk4TJ9Nicj+mHDBkSfJzF\npSE0VplrJVtW32uGGh8dcgVqbdTKR/F87ArE0J4e+fFRpB5V8+fTyqFcSVeP6floXrNEagmuPKxh\niM8++yzb3BQQiMf0GvbgY/tx48YFH9+Lbt26ZVsr0vLxrWZTcciTQ69APGrVTBcODXF11aldJ8NN\nMjWjsDnQzBT+uVLVYm7YCcSsPA41AvEz6vzje6NZoAMGDMg2V5nWyuMcgtEj8Ntuu22q1wjE+ajh\nQp5/HJYGgH333TfbmkHD4Y1aCpPptfCaomsph261kfHWW2+dba3szuNes4c4+4szdbUJLofNuPEy\nEL8bBg4cGHy8BnPmLxDnqobkOeyiY3P48OHZ1vWrqSv8a5iMf9ZKzzx2da3hz9eqVavg46rdKv04\n8MADs60NczlbljPNNOOOs7i1IjnLVTR8zfeWv0OA2J1A4TGn+wqHyYwxxhhjGhFvhowxxhhTarwZ\nMsYYY0ypaRLNkHYf5vR27gYPAMcee2y2ubMyAEyaNCnbmtbJMXBNuevZs2e2W7dunW3VfXDclK8D\nAO66665sazzy1ltvzbaml3M3c9UzcAVNrdqrWo5ahXUKGpPmGPhGG20UfA899FC2tUI5/8z3D4ha\no0suuSTbWkGVY8uqx2LdgN53jrlrB3Qex3rNXFJBf4+rqmsaabVRDRPraDTVmO+96h1Yx6Dzj3VC\nXP0aAJZeeulss05OO6mzhueOO+4IPi7ToKnFfC2qTeMxqboW1iPo+OE1pJbQ+cdrCOstgfictawH\nr6W9evUKPk53V/3IhRdemG3WnOh8WHXVVbPNGk4gapu0MjZ/Hl2f+W+odoqr0x933HGoi6bWCNWH\n1157LfzMejrV2fI81ft5wAEHZPuee+4JPtZZaVV3Xne5QjSXgwFiBeqnn346+FhvpvOZyyuoVotL\nAOgzYT2brsFcDqAx8MmQMcYYY0qNN0PGGGOMKTVNEibjZnhAbLLHx7UAsOmmm2ab026BmP6qYSQ+\n2tVj7RVXXDHbfJTMjd0A4Kijjsp27969g49DWpyOCcTGsJxiDAB/+ctfss1Hg0AMA+rv8fFgLaXv\nKnxUqU0y+UhTj0nXXHPNbGvFcD4m1Wa6nGJ96aWXZlvv0ejRo7PNx/JADFvpNXOqrx7D8hGxNoPk\n8I8ecfNxf3On7/4aGibk+8jzBoip9ssuu2zwcahKwxs8z/SecpVaXgu0ESenza699trBx/Nd1wyu\naK8hSg6haco/f1a+RuB/K9DXCvoseQ5oyJPDgNp0mCUIXD1cfdrkmFOrb7jhhmxzhXAgViHXVHBe\nJ3jtBOIz4YbNQGx0zSnjQBw7GkrkEi1aHb2ac1M7LPD3m1YM5wr7GiriEhj6ncmv1RAaV5vn99DK\n8Nw895RTTgk+Xme1Efnuu++ebb3P/N2rEhu+D/p72tR7evHJkDHGGGNKjTdDxhhjjCk13gwZY4wx\nptQ0WjCcY3tcshuInaQ1Tjtq1Khsn3rqqcF30UUXZZtbPQDAbrvtlm1OJwRiR15OodXYK6cCcsov\nEOPhGl9lTYq28eAUb21vwDoT1VJwB+NaQjUvrOF56623go/jwmeffXbwcXo764eAWLpdtQ5cup31\nINr5nvUnmlLN768lGrhVh5aFZy0Mj2Egah/4PYA45mrtOWunaWbdddcNP7PGRjUNfK8mTJgQfNym\nQX08rzp06JBt1eWwDkvfg3UMXOIAiG0NdIwwqj/grtrqY+2iPs9qohpILnegz4s1UtragnWbW2yx\nRfCxhog7iAMxBZu1W9rqhEuT7LLLLsF3+eWXZ1t1m08++WS2Tz755ODj1Hpdu7mMh45pbZ9UK8w+\n++zhZ9ZdqfaNf9YyBtzSRNPi+bX8/anweqalOLjNkrZP4bHDulogfoe2b98++FjPpiU8+PtH2ww1\ntpavNr+BjTHGGGOaCW+GjDHGGFNqmiRnVEMKfJyllVHffPPNbHOFYSBW4eTjVCB2ktcUv4MPPjjb\nG264YbY1rMNH8VxpFQCee+65bGuKPB/Na8oxH8NyuuKMgnYOZvj+7rPPPnX6ttpqq+DjSuOakslH\n7hwq1XBau3btsq3Hylwpdf311w8+Dr1pmiqHezW8xuNYw7R8HF2phEI10HRUDpu98MILwcep2hpe\n4G73WnqCw1rnnXde8HHoiksqPPbYY+F1XCH30UcfDT4eg3POOWfwcUVjDVPzdXH6NRBDbzrGufJy\ntcOcHKbTyt+8Th199NHBxyUUNFTF4SmtTs0Vvrt16xZ8HMLgccVhNwC48cYbs83hSCCmxS+xxBLB\nx+E1lh8AwA477JBtHX/8fbP55psHn3ZMrxV07eGSAGPGjAk+Tj/XEHKlED2vPVxiAIhhfy6xweVS\ngBh21FDY888/n2297/w9r6V3ONyr6wzPN60az89Zq6M3ZJ31yZAxxhhjSo03Q8YYY4wpNd4MGWOM\nMabUNJpmiGN0GsfneJ7GIDkV9pBDDgk+1hdoKfhjjjkm29odnkv7czxZy/pzOmOfPn2Cj7UH2oGd\ntSuqC+LSARrP5ZRv1SFpCmitwnopTd/lMv+cgg/EdF5NheW2Alp6gd+H4+ia8s+aiK+//jr4WPfA\nJeOBmPLP2iUgxupVo8R/X0sF8HOvtkYIAFZYYYVsa+sQ/hyapsufi9O2gagB0/RbTrvXccB6LtaH\naSsG1jvo3OR5xM8PiGnV2kKFdRH6WVknoZohHnfV1gzxeFL95TnnnJNtLfnBmo3zzz8/+K644ops\na9o9r6Vdu3YNPi5nwTpR1eXwPdOWL6xRUq0pt0Tq27dv8PEYU80otwbRzur8bDWdvblhvZ62F+Ln\nvPjiiwdfv379ss1lZICoUdQSN7yederUKfi41RGvz3feeWd4HZe74U73v/a3eS1RvRLrNnmcAjHV\nXktJ8D379ttvMb34ZMgYY4wxpcabIWOMMcaUmkYLk/FRqB4lc5VprQLNqc69evUKvv79+2dbU28P\nO+ywbOux7913351tPhbldHkgHhvqUWSbNm2yrSESDuXMNttswXfVVVdlm9P/gZg6ymmIQEyT1fes\nJbiKtqZI8mfQas78/EaOHBl8nCqqIVau9s0VVbW6Kj9bDpUAMRTbtm3bOq9LKxZz+FVTX7mkgoYL\nOfRUC13rOcSl18NzVX1cgVvDlxye0nT9q6++OttadoCfE6fKDhkyJLyO1wWez0AMtWs6NldM1uvi\n0Ld29ObnwuMFiJVva6mDvVbl52fEJUWAGBJ65JFHgu/www/P9v777x98nD7dunXr4OPQJoc+tKTB\nZpttlm0u1wDEUKaWPmFZg4beODyk6+Vaa62V7TPOOCP4uKyHVsrmta054NIBGiLksJKWgeCxqxXl\nOQzJXSGAOFc0xMprOctCtDQI31sNqd5yyy1TvUYgfl9rKYS3334727p34Ovi+6Wv1TITWsl6WvDJ\nkDHGGGNKjTdDxhhjjCk13gwZY4wxptQ0WgCcU2E1Xsixe9WZfPLJJ3X+HmtuNEW3c+fO2b7rrruC\nj/UGnO692mqrhddxzFHjzpy+qLHR+++/P9vcqR2ImhSNcfLnU20Mp5frPaomqnHhtHX9DBxP1lRR\n1oSoj1shaNo96wg4Hq76AoZTcoHYHkD1PQxr24AYq9f08Uqp/EwtpNZzaQstd89p4zrHOCVf05D5\nfW677bbg23XXXbM9ePDg4ONnzZoN1W9wzF/TdLnVi2otOBVcdQSssdHPyu+j65BqiKoJjydtT8Aa\nr1lnnTX4WLPBbY6AqO979913g4/H/fHHHx98nL7P2iJtucFzTp/llVdemW1u2wHE1G1e7/W1rDkB\ngFdeeSXbO+20U/DxuFUdYHPD16L6Hh5z2q2d56Xqs1hLqXo61jk++eSTwcfj6sEHH8y26qjOPffc\nbHPLHCCOD9UG85jjOap/W3+P57CuEbzuNoaWzydDxhhjjCk13gwZY4wxptQ0SZ6oHjPzUZt2Jecj\nzi222CL4uHKwHtOPHTs223rcycdufDys1Uj322+/bOtxO/89TskFYiVMrabJR5N6FM+ViVdeeeXg\nq/aR7bTCx+98LA/E42mtzM2hUg2TXXLJJdnmbujA/6YP/8K9994bfuYwraYV89/WauV8vKoprJw+\nriEXfpZaRoDDr7UQJmM05ZSvlY/fgdjFXufmwIEDs61p9+3bt8+2PicePxxm0XDztddem22tTt2h\nQ4ds77XXXsF36KGHZvuyyy4LPi6foWEJvi5Nya921em60HA6hz50jvH6xqFRIIbbdLxec8012X7x\nxReDj+fqsGHDsr3zzjuH1w0aNKjO6+IKx/vuu2/wcVkUDevcfPPN2T722GODr2PHjtkeMWJE8PF3\nU7VLmHApCw1H8fry8ccfB98cc8xR53tySFKlBHwPtZQMj3++LzyfAKBHjx7Z5jUXiN+FOod4/GnJ\nmbquH4ifXbs2MI2xztbmLDfGGGOMaSa8GTLGGGNMqfFmyBhjjDGlpknacWhZf07P1FYMnC6tJcL5\n9zTOzZoFjWuy/ub666/PtsZJWc+gHea5FLx2a+aYKpeyB2I3bNUacFxY4621VOa/EhwLVh0X62q4\nZIL+rKm3nDLJuiMgaju4JL/qlXjMaYz9nnvuyTanl+o16xhjbZhqbWpNC1QJTq1nG4ifUTUhXCaC\nyyYA8Tlp6i+n32oLB9ZXcRkF1Uzw3GHNCRBTellXCEQNouoWOBVXSziwNm2ZZZZBS0BLd9xxxx3Z\nPvDAA4OPU9q1+zdrUHQOsPaPy38AURPJGsgHHnggvO6oo47KNus0gTjm9LuBW4PofPvrX/+abW1J\nwe+jpSSqrROqC11L+V7rmsVrn2oqWRP18MMPBx+Xi3n22WeDjzvH33DDDdnu06dPeF2rVq2yrWNl\njz32yLZ+n3bq1Cnb3GoEiN8NqvGaNGlSthvSYqM++GTIGGOMMaXGmyFjjDHGlJpmic1wyrKGqjiV\nTkMffHSoR2Q77rjjVN8DiCEoruTJ3XgB4PXXX8+2ptpy92vudKxUCgfp0SeH75r6yK+p4ONq/Qx6\nJM1wZVRNU69UiZWPSTl1k7tWAzF1U1O9+TifnysQj5krVUFWOCxX6yEzvd8MV/nV+cFH9fw6IIbQ\nHn/88eDjarN6bwYMGJBtrhSvYWN+Dw2fLLrootnm8gcAsNJKK2Vbw0HcgVtDe1rNuSWg696qq66a\n7Upp8Bp25HBimzZtgo+r7d99993Bx+FQLnXCNlA5fMffB5o+X6na/bhx47Kt3xtc3qRSxflaQj8D\nr4M6hziEzGUsAODNN9/MtoZ7OUSo9/O+++7LNoesNTTKZQt0LeXvPq0kzZXNtWr2uuuuO9W/DcTy\nCmwDjb/u+mTIGGOMMaXGmyFjjDHGlBpvhowxxhhTappEM6TaEY71aQyXtRfq4/TJbbfdNvhYj6Od\nlhnWFimclqtxdI6pqs6JO29PnDgx+L799ts6/15L1QnVhcZs+blreQUuHaBlBCqVZdDY9i9oaXZO\n5dQ4t7bLYLgMvaadM7WuC6oEawW03MPnn3+ebW2hwnH//v37Bx9rubQ9AOuLtJUNa/i4rMZVV10V\nXsed6RV+Ftr6geejasdYh6HjSvUILQFtbcHrGWtqgFheQlObV1tttWzrfXjuueeyrTo91n6wrks7\nxfO44jYaQNSTcno3EMuU6BjjeautjPi+tBTNkK5RvEZqmxy+Z6q/4We5+eabBx+319ExwLpC/m7d\ne++9w+u4rIV+71544YXZbteuXfDxe2q6PpcOUO0gj8emXoN9MmSMMcaYUuPNkDHGGGNKTbOk1vPR\nvFaFZTS9lavEaviEj8807MLp0qNHj852pZCW/m0+etVKrxxO0bAYX4seOfMRYK12wm4s9EhzWtPU\nGwqHZupznMrH7xrebcmhMabSWOOKtVxqAohhJQ1FcMXhxRZbLPg4TMEV4IEYsuN7r8fq/Dw5ZRwA\nhg4dmm1dT7jitVbZ5evUEBqng+tRfUuBQ0dcYgCIoRZdl4YPH55tDSvxPdP7wmEsLp2hFeY5nKLh\nLl4vtQMB/z1dP/h9NLyt1Y9bIvyZtPs8h5x0bo8ZMybb3bp1C75rr722zr/H4eazzz4721q64pFH\nHsm2hjy7d++ebZ1fPJ91vXj11VezrSUvmnMNnrG/kY0xxhhjfgVvhowxxhhTarwZMsYYY0ypqelW\n6ZVi9xwz5u72QIxDc5lzLXnO8VbtZsx6CS7jr2jMlv826zH0PWd0zVB9qKY2p5KGrQxU0mtVaq/C\nujzuiA7EruVchgKIKcRcvl9fx+/JqcRA1CPonOYUb51jrMOolI6tuoWWQqV5xM9SdTvaFX160THF\n16WlFir9Xl3vUZ/fmxHh7xX9juQxPnjw4ODjchLPP/988LHOiucNd7AHon5Wx9FNN92UbS37wHo2\nbbnEaf3ff/998PE8rbQeNQb+RjbGGGNMqfFmyBhjjDGlJtXniDGl9AWAD371haapWLooioV//WW/\njp9l1Wm0Zwn4edYAnpszDn6WMxbT9DzrtRkyxhhjjJnRcJjMGGOMMaXGmyFjjDHGlBpvhowxxhhT\narwZMsYYY0yp8WbIGGOMMaXGmyFjjDHGlBpvhowxxhhTarwZMsYYY0yp8WbIGGOMMaXm/wDvsASC\nzvYGoQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10bf675c0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "W = logreg.coef_\n",
    "nlabel = W.shape[0]\n",
    "plt.figure(figsize=(10, 6))\n",
    "for i in range(nlabel):\n",
    "    plt.subplot(2,5,i+1)\n",
    "    plt_digit(W[i,:])\n",
    "    plt.title('{0:d}'.format(i))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What you see is that each weight is a very blurry version of the digit.  The blurriness is due to the fact that weight must correlate with all shifts, rotations and other variations of the digits.  As a result, the weights begin to correlate with other incorrect digits leading to poor performance."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Run an SVM classifier\n",
    "\n",
    "We now try an SVM classifier.  The parameters are given by \n",
    "\n",
    "https://martin-thoma.com/svm-with-sklearn/\n",
    "\n",
    "This website has a nice summary of the main equations for SVM as well.  That site trained on 40000 samples and tested on 20000.  But, to make this run faster, we will train on 10000 and test on 10000.  If you increase to 40000 training samples, you can get past 99% accuracy.\n",
    "\n",
    "First, we import the SVM package and construct the SVC with the parameters."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## First let us try the linear kernel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn import svm\n",
    "\n",
    "# Create a classifier: a support vector classifier\n",
    "# svc = svm.SVC(probability=False,  kernel=\"rbf\", C=2.8, gamma=.0073,verbose=10)\n",
    "svc = svm.SVC(probability=False,  kernel=\"linear\", C=2.8, gamma=.0073,verbose=10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we get the training and test data.  The features are re-scaled to be between -1 and 1."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We fit the training data.  Again, this will take several minutes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Xs = X/255.0*2 - 1\n",
    "ntr = 10000\n",
    "nts = 10000\n",
    "Xtr = Xs[Iperm[:ntr],:]\n",
    "ytr = y[Iperm[:ntr]]\n",
    "Xts = Xs[Iperm[ntr:ntr+nts],:]\n",
    "yts = y[Iperm[ntr:ntr+nts]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[LibSVM]"
     ]
    },
    {
     "data": {
      "text/plain": [
       "SVC(C=2.8, cache_size=200, class_weight=None, coef0=0.0,\n",
       "  decision_function_shape='ovr', degree=3, gamma=0.0073, kernel='linear',\n",
       "  max_iter=-1, probability=False, random_state=None, shrinking=True,\n",
       "  tol=0.001, verbose=10)"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "svc.fit(Xtr,ytr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Save the results in case you want them without re-running the above the code."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pickle\n",
    "\n",
    "with open( \"mnist_svc.p\", \"wb\" ) as fp:\n",
    "    pickle.dump( [svc, Xs, y, Iperm, ntr, nts], fp)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can reload the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Can skip this step if you run the classifier directly\n",
    "import pickle\n",
    "with open( \"mnist_svc.p\", \"rb\" ) as fp:\n",
    "    svc, Xs, y, Iperm, ntr, nts = pickle.load(fp)\n",
    "Xtr = Xs[Iperm[:ntr],:]\n",
    "ytr = y[Iperm[:ntr]]\n",
    "Xts = Xs[Iperm[ntr:ntr+nts],:]\n",
    "yts = y[Iperm[ntr:ntr+nts]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Measure the accuracy on the test data.  The prediction can take several minutes too -- SVMs are *very* slow!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "yhat_ts = svc.predict(Xts)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since even the prediction (sometimes called inference) is slow with SVMs, we will save the results in `pickle` file. Instead of running the prediction again, you can recapture the data with the following comamnd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "with open(\"mnist_svc_test.p\", \"wb\") as fp:\n",
    "    pickle.dump([yts,yhat_ts,Xts], fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "with open(\"mnist_svc_test.p\", \"rb\") as fp:\n",
    "    yts,yhat1,Xts = pickle.load(fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuaracy = 0.905000\n"
     ]
    }
   ],
   "source": [
    "acc = np.mean(yhat_ts == yts)\n",
    "print('Accuaracy = {0:f}'.format(acc))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should get an accuracy of around 91%.  Again, had you trained on all 50,000 samples, it would have been much better -- close to 98.5%.  But, even this result is much better than logistic regression."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Looking at the Support Vectors\n",
    "\n",
    "Let's take a look at the support vectors.  We see there about over 2500 support vectors.  So, about quarter the training samples were used as SVs.  This is partly why the prediction was so slow."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(2685, 784)"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "S = svc.support_vectors_\n",
    "S.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us plot some support vectors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAACOCAYAAAAy0AzYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAADyBJREFUeJzt3XtQVVUbx/F1MEUJYTJNAdPAMi21\ni5cuU42pzWSlXS2ZIayxicLRNGtizCbJQiO7CJlaJiroNEr94XRxypy0Rs1iki5m2UVmorS0BAEp\nkP3+9+az1o5zNuxzOJz1/fz3W+21z5qX824fN49rBRzHUQAAALaK6+gFAAAAdCSKIQAAYDWKIQAA\nYDWKIQAAYDWKIQAAYDWKIQAAYDWKIQAAYDWKIQAAYDWKIQAAYLXTvFwcCATYrtoCjuMEwnVvvkPW\nOOI4Tp9w3ZzvkR14FsEHIT2LeDMEIByqOnoBAKBCfBZRDAEAAKtRDAEAAKtRDAEAAKtRDAEAAKtR\nDAEAAKtRDAEAAKtRDAEAAKtRDAEAAKtRDAEAAKtRDAEAAKt5OpsMgDRx4kRjrE8feQzOli1bRC4o\nKAh636uuukrkTz75JOicr7/+WuSVK1ca15w4cSLofQDANrwZAgAAVqMYAgAAVqMYAgAAVqNnCGiH\nK664whibP3++yIFAQGTHcYLet6WlReTU1FTjmh49eojcpUsXkR977DFjTmFhochLly4N+tkAol9R\nUZExlp2dLXJiYqLIcXH+vA9pbm4WOT09XeTq6mpfPieceDMEAACsRjEEAACsRjEEAACsRjEEAACs\nFgilmfP/FwcCoV+MTstxnEDwq9om1r5DI0eONMbefvvtVucsWbIk6H2PHj0q8po1a4xrsrKyRF60\naJHIaWlpQT8nKSnJGKurqws6LwQVjuOM8uNGbmLtewR3PIv+W0ZGhshz5swxrvn4449FPnjwoMhu\n/zAjmDFjxhhjeXl5IhcXF4v80EMPef4cH4X0LOLNEAAAsBrFEAAAsBrFEAAAsFpUbrqob1IXHx9v\nXNPU1CSyW+9Tt27dPH/2kCFDRM7NzRV5+PDhxhx9g6n8/HyRV6xYYczx0quF6FVRUWGMpaSkROSz\ny8rKRO7Zs6fIy5Yti8g60DGmTJlijM2dO1dkt/4O3V9//SVyeXm5yDk5OW1YHcLtp59+EnnmzJkR\n+dxevXoFvebnn3+OwEr8xZshAABgNYohAABgNYohAABgtajYZ0g/PE7vhZg8ebIxp7KyUuQ//vjD\nuGbChAk+rM47/felbnvRHDt2LFLL8Yy9PTqn3bt3i+zWL3L8+HGR3fqbGhoa/FgO+wy1YuDAgSK7\nHap7xhlniDxx4kSR9YN6lVKqa9euItfU1Iis91oqZT5/9R7NqVOnGnM2bdpkjIUDz6LoU1VVZYzp\nzxH9z7yvvvoqrGsKgn2GAAAAgqEYAgAAVqMYAgAAVqMYAgAAVouKTRdbWlpE1g+ga2xsNOacf/75\nreb/mneqI0eOGGPr1q0TecaMGSInJye3ek+llLr55ptFjuZmaXRes2bNEvniiy8OOmfp0qUi+9Qs\njVPExZl/x3z55ZdFvueee0Tu3r27MefkyZMi6xt8Hj582JhTWloq8gcffCCy3lCtlFJXX321yFu2\nbBHZbXPHSDVQI/L0xn394On+/fsbc/SDWDu4YbpNeDMEAACsRjEEAACsRjEEAACsFhU9Q3rfwogR\nIyLyuW6/px8/frzIs2fPFln/Pb5SSn3//fciV1dX+7A64F96f5BSShUWFoqsH0xcUlJizFm4cKG/\nC4PB7X/3u+++u9U527ZtM8bmzZsn8p49ezyvZdCgQSIXFxcb1+zYsUPkgwcPiux20DRiw5lnnmmM\nffbZZyIPGDBAZL0/SCmzJ64z4s0QAACwGsUQAACwGsUQAACwGsUQAACwWlQ0UHcUt9OYV69e3eoc\ntybGyy+/3Lc1AUopVV5eLvKkSZOMa/QTyu+9916RN2zYYMxxO7Uc7aOf0J2ZmWlco28sq2+6uH79\nemOO47R+qLrbBrCvvfaayJMnTxZZb7JXyvxHI4MHDxa5vr6+1XWg89Abpt3+PNM3XbzppptE1jfl\njBW8GQIAAFajGAIAAFajGAIAAFazqmfouuuuE9ltAzLd/v37RdY3ugO8SkxMNMZ27twp8pAhQ0Q+\n7TTz/6pr164VWe8Roj8oMhYsWCCy289q69atIpeVlXn+HP1A1WXLlhnXDBs2zPN9CwoKRKZHKHbo\nfbH6z7q2ttaYc9ttt4lcWVnp/8KiEG+GAACA1SiGAACA1SiGAACA1WK6ZyghIUFk/Xfsp59+ujFH\n7xEaO3asyL///rs/i4M19L1dKioqjGvOO++8Vu9RWlpqjN13330iux0ijPDbvXu3yDfeeKNxzfDh\nw0XevHmzyGlpacacAwcOiKzvNdWjRw9P61RKqePHjxtj69at83wftF2vXr2MMf3PmVBccMEFIrvt\nb6X3Hubn54v83HPPGXNOnDjheS2xgDdDAADAahRDAADAahRDAADAahRDAADAajHTQK03SytlNpie\ne+65Ih8+fNiYc+2114qcmpoqst6QFgq9wVIppf755x/P90Hn0LNnT5H37t0rcnp6ujFH/z5s3LhR\n5GnTpvm0Ovht8eLFIs+aNcu4pm/fviLrh1+6GTRokMiffvqpyG6N+DNnzhRZb97ftWuXMcetqRrh\n89JLLxljWVlZYfmsuro6kfXDf5OSkow5NFADAABYiGIIAABYjWIIAABYLaD/DrHViwOB0C+OsDFj\nxhhjbr06p2psbDTGvvvuO5H1HqH4+HjPa9Pv6fbZ+qZbNTU1nj/HL47jBMJ172j+DoVC36jzjjvu\nMK55/fXXRY6Lk3/n+PXXX405/fv392F1UaXCcZxR4bp5NH+P+vXrZ4z17t3b8330DV4bGhpE/vzz\nz405gwcPFlnv/9B7IJXq2GdNMLH4LLrsssuMsUsuuUTkJ554wrhGHxs9erTIKSkpxhz9GrcNH3Xz\n588X+cUXXxS5paUl6D2iTEjPIt4MAQAAq1EMAQAAq1EMAQAAq8XMPkNt0b17d2PsoosuEvmXX34R\n2W0Phm+//VZk/Xe3GRkZxhy9h6Br166tLxYdIjk5WeQPP/xQ5EsvvTToPWpra0WOwf4gnOLQoUMh\njXm1fPlykfX+IDf6vkPR3B9kiy+//NIYe/zxx0WeO3eucc0bb7wh8urVq9u9lkceecQY0w9vHTp0\nqMj333+/MacT9hEZeDMEAACsRjEEAACsRjEEAACsRjEEAACsFjMN1Pv27TPGbr311nbfd/v27SIf\nO3bM8z1WrVpljE2fPr3Na4I/9ANVV65caVyjb5DmdsiqTt9UcdKkSW1YHWw3Y8YMkXNyckQOBMz9\nCPUNXv1osoW/pk6daozph4jrzdLhUlRUZIzpzzi9Ydptg+MRI0b4u7AOwJshAABgNYohAABgNYoh\nAABgtZg5qDWa7Nq1S+RRo8wz4urq6kQ+55xzROag1vZJSEgQOSsry7gmLy9PZP1n4Ka+vl7khx9+\n2LimpKRE5Obm5qD3jUHWHtTaFvr3VSmlfvzxR5H79u0rclNTkzFH3zR2//79Pqyu48TCs0j3zTff\nGGM7duwQ+cEHH4zUcoJ6//33RZ4wYYJxjf7Mi7KeWA5qBQAACIZiCAAAWI1iCAAAWC1m9hmKlAED\nBhhjjz76qMh6j5DeH6SUUuPGjROZAxTb56yzzhJ5z549Irv93ELxww8/iJybmyvy1q1b23Rf4FTP\nP/+8Mab3COnee+89Y6yz9wjZIDU11Rg7cOBAB6wkNDfccIPI+mHVSimVnZ0tcnFxsch79+71f2E+\n480QAACwGsUQAACwGsUQAACwGsUQAACwGpsueqQ3hillHqioW7x4sTE2b94839bkt2jf6Gzs2LHG\n2IYNG0Tu16+fyI2Njcact956S+Tly5cb13zxxRciNzQ0hLpM27HpogeHDh0yxvR/FFBbWyvy6NGj\njTnR3IjbFtH+LGqLRYsWGWNDhw4V+ZZbbonUcnyh1xH6QbOZmZmRXI6OTRcBAACCoRgCAABWoxgC\nAABWY9PFIBYuXCjyAw88EHTOihUrRH7yySd9XZPt3A5d1XuEdH///bcxVllZKfKVV15pXOM25lV8\nfLzIOTk5xjX69+zVV19t9+ciet1+++0i9+7dO+icd999V+RY6w+yRWlpqTG2fft2kd0OgH7hhRfC\ntiYv3DYD1XuG3DYajna8GQIAAFajGAIAAFajGAIAAFajZ0jz1FNPiZyXlyey2+9C9f1p9B6hpqYm\nn1aHtkpOTjbGnn322Q5YibtXXnlFZL3Hae3atZFcDnwWFyf/3qk/I/T/rpS5p9XTTz/t/8IQcfv2\n7TPGioqKRM7PzzeuaWlpEXnTpk0iV1dX+7A608iRI0Xetm2bcU19fb3IhYWFYVlLOPFmCAAAWI1i\nCAAAWI1iCAAAWI1iCAAAWM3qBurrr7/eGNMPUO3SpYvIa9asCToH4TVnzhxjLDExUeQ777xTZL3B\nTynzMMEpU6YY1yQlJYm8efNmkTMyMow5w4YNE1lvfCwpKTHmHD16VOT169cb16Dzuuaaa0TWvyNu\n9O+AW+MtYoO+6WpNTY1xjX7gd0FBgchu34+PPvpI5J07d4qckJBgzMnOzhZ53LhxIrsdej1+/HiR\nO+OGoLwZAgAAVqMYAgAAVqMYAgAAVgvoB6y1enEgEPrFUWjBggUi6xsqKqVUt27dRM7MzBT5zTff\nNOY0Nze3f3FRxHGcQLju3dm/QwhZheM4o8J18872PaqqqhL57LPPDjpn2rRpIrsd8BnreBb9KyUl\nRWR9E84LL7zQmNOnTx+R09PTPX9uWVmZyLNnzzau+fPPPz3fN4JCehbxZggAAFiNYggAAFiNYggA\nAFiNYggAAFgtpjddTEtLEzk3N1dkt82j9M2vysvLRT558qRPqwNgC33z1lC4begJe/32228iT58+\nvYNWEpt4MwQAAKxGMQQAAKxGMQQAAKwW0z1D+u/p4+Jk7bdx40ZjzjPPPBPWNQGwzzvvvCPyXXfd\nJbJ+UK9SSq1atSqsawLwL94MAQAAq1EMAQAAq1EMAQAAq1l1UCtCw+GI8AEHtaLdeBbBBxzUCgAA\nEAzFEAAAsBrFEAAAsBrFEAAAsBrFEAAAsBrFEAAAsBrFEAAAsBrFEAAAsJrXg1qPKKWqwrEQRI2B\nYb4/3yE78D1Ce/Edgh9C+h552oEaAAAg1vBrMgAAYDWKIQAAYDWKIQAAYDWKIQAAYDWKIQAAYDWK\nIQAAYDWKIQAAYDWKIQAAYDWKIQAAYLX/AQI1BVkwGhD6AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1a24457860>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "nplt = 4\n",
    "nsv = S.shape[0]\n",
    "Iperms = np.random.permutation(nsv)\n",
    "plt.figure(figsize=(10, 4))\n",
    "for i in range(nplt):        \n",
    "    plt.subplot(1,nplt,i+1)        \n",
    "    ind = Iperms[i]\n",
    "    plt_digit(S[ind,:])   "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that the support vectors look like digits we want to recognize. They are like the 'match filters'."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Now let us see the performanc of SVM using the RBF kernel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# from sklearn import svm\n",
    "\n",
    "# Create a classifier: a support vector classifier\n",
    "svcrbf = svm.SVC(probability=False,  kernel=\"rbf\", C=2.8, gamma=.0073,verbose=10)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[LibSVM]"
     ]
    }
   ],
   "source": [
    "# Can skip this step if loading from previous result \n",
    "svcrbf.fit(Xtr,ytr)\n",
    "import pickle\n",
    "with open( \"mnist_svcrbf.p\", \"wb\" ) as fp:\n",
    "    pickle.dump( [svcrbf, Xs, y, Iperm, ntr, nts], fp)\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Can skip if you run the previous cell\n",
    "import pickle\n",
    "with open( \"mnist_svcrbf.p\", \"rb\" ) as fp:\n",
    "    svcrbf, Xs, y, Iperm, ntr, nts = pickle.load(fp)\n",
    "Xtr = Xs[Iperm[:ntr],:]\n",
    "ytr = y[Iperm[:ntr]]\n",
    "Xts = Xs[Iperm[ntr:ntr+nts],:]\n",
    "yts = y[Iperm[ntr:ntr+nts]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# skip if load from saved result\n",
    "yhat_ts = svcrbf.predict(Xts)\n",
    "\n",
    "with open(\"mnist_svcrbf_test.p\", \"wb\") as fp:\n",
    "    pickle.dump([yts,yhat_ts,Xts], fp)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# skip if you run the previouse cell\n",
    "with open(\"mnist_svcrbf_test.p\", \"rb\") as fp:\n",
    "    yts,yhat_ts,Xts = pickle.load(fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuaracy = 0.972100\n"
     ]
    }
   ],
   "source": [
    "acc = np.mean(yhat_ts == yts)\n",
    "print('Accuaracy = {0:f}'.format(acc))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that the RBF kernel provides more accurate results. \n",
    "Now let us plot some errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAACbCAYAAABhy5/wAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAEGtJREFUeJzt3W+QFVV6x/Hfg0Akm6T8g39SS9yV\nZGUiVnghaHALFDK+EozrWGXxYoUKW1agNpBYqZRLFW9iZI2JJi8oIJDUlkQLy8oQWFxccZwyy7Ib\nU1EjIWFAQCVrAN1VEiBBFE5edE/qnu6e233v3Hu77z3fT1VXzek53X3u9HN7nun7zGlzzgkAACBU\nE8oeAAAAQJlIhgAAQNBIhgAAQNBIhgAAQNBIhgAAQNBIhgAAQNBIhgAAQNBKTYbM7D0z6y9zDPE4\nbjKznWb2kZl9bGYvm9mMNhxn3K/XzJaamTOzb7RqXN2uKnEkSWa22MwOmNlZM/uRmd3chmM09XrN\nbF48rtrFmdlAq8fYbaoSQ2Y21cz2mdnPzOy0mf3YzL7ahuM0/XrjmDlXE0N/3erxdasKxVFH3uvj\njKPNZnbIzC6Z2bJWjqsZlb4zZGYTO3SoKyR9V9IMSddJ+idJOzt07MLM7EpJ35L0b2WPpZt0Ko7M\n7CuSnpP0u4piapek73Ywjutyzu11zv3C6CJpkaSzkr5f8tAqr4Pn8Kyk35F0jaQrJf2ppF1ViaEa\ns2piiT/MCurUeeyS9/rbklZKerPsgUiSnHOlLJL+VtIlSf+r6CT9kaQvS3KSlks6LukHku6S9JPE\ntu9J6o+/niDpUUlHJf1M0guSrhrn2K6Kx3H1GN//OUl/Ho/xlKRNkqbE35sq6UVJpyV9LGlvPMbU\n621iXJsUBc9rkr5R1rmr0lKlOJL0TUnfq2lPiMf1W1WKo5rjf0fSd8o+h2UvVYqhxL4nSFocj+Pa\nqsRQPJ5fK/u8VW2pahzF+6z7Xi/zWiTph5KWlX7+Sg6e/w+AuD0aOFslfUHSlAKB8/uS/lHStPiE\n/pWkbTV9T9dZHh1jXPdJOlFn3H+p6E7SVZJ+UdEdgG/H3/t2HEiT4mWeJMt6vY2MT9Jtkv45DsLX\nRDJUuTiS9HuSdtdsc5mk85JWVyWOavr/vKQzku4q+/xVYalKDNX03S/pQjyGLXXGXca1yEn6T0kn\nJW2X9OWyz19VlqrFUdw/971eRhzV9CcZqhM402vW5QXOQdX85S3plyV9Jmlik2OaJukDSUvG+L5J\nOifpV2vWzZX0bvz1Hyv6iC31l1NW4BQc02WKEqG5cfs1kQxVLo4k9cWxcZekyZLWKvrL6VtViKPE\nPr4u6d3Ri1roS1ViKLHvyyUtkbR0jO+XEkOS5sfxfYWk9ZIONPsae22paBzVfa9X4FpUiWSoqjVD\n/9FA3y9J+vu42PC0okC6qKj2pyFmdo2kPZI2OOe2jdHtGkWZ9hs1x/x+vF6S/kzSEUl7zOyYmT3a\n6DgyrJS03zn34xbsKyQdjSPn3IikpYp+QZxQdHv53yX9JKN7GXFUa6mkrS6+GmFMpVyLJMk5dz6+\nDj1qZrMyupQSQ865HzjnLjjnTktaLelGSb/ein33sNLiSPnv9bKvRdVQchb9rrKz6Ik16+ZI+rim\nfZmiLHY0iz4k6at1jnG2zrKmpt+Vkt6S9ETOmCdI+h9JXyzw+mZK+lBxlp98vUXHJ2mHpE8U3ZY+\nqej2+X9JWl92Nl2FpUpxlNjmCkW3p/uqEEc1fX9F0ueq+Usw9KWqMRRvd0TS16oUQ4mfwVlJv1H2\nOazCUrU4KvJeLzuOVJE7Q2X/h8IpSdNz+hyWdLmZ3aPors0aRZ+jjtok6XEzW+qcez++u3OHc26n\nJLmokr4uM/slSS9L2uecq5v1OucumdkWSX9hZt90zn1oZl+UdItz7mUzWyRpRFHx238ryugvjvV6\ni4xP0jJFt8xHbZf0d5L+psC2IahEHEmSmd0q6V8Uffa+XtIuF90x8pQUR6O+LulHzrmjDWzT6yoR\nQ2b2m5ImKvqP1sskrVJ0R+D1ZN8yYsjMZiqqG/lXRfUvf6KorOBg3raBqEQc1ch9r5d1LTKzyYoS\nMZM0ycwul3TBOXepgdfXOmVmYpJ+W1H1+mlJf6iMLDrut0zRxw4fxv3ek195/4iibPqMohO2rsFx\nLI2Pe05+FnvDGP0vl7RO0jFFwXFQ0qr4e38Qj++coo9H1o71esfxc3tN1AxVLo7i/fww3v5jRYWP\nX6jTt5Q4UnRhW172eavSUpUYknSnon85Ho2hf5A0vyoxJGlh/PrOxT+DHZK+Uvb5q8pSlTiqOU6h\n93oZ1yJFv8dcYrmrrHM3WhEOAAAQpKoWUAMAAHQEyRAAAAgayRAAAAgayRAAAAgayRAAAAhaQ/MM\nmRn/ehYA55y1a9/EUDB+6py7Jr9bc4ijMHAtQgsUuhZxZwhAO7xf9gAAQAWvRSRDAAAgaCRDAAAg\naCRDAAAgaCRDAAAgaCRDAAAgaCRDAAAgaCRDAAAgaCRDAAAgaCRDAAAgaCRDAAAgaCRDAAAgaCRD\nAAAgaCRDAAAgaCRDAAAgaBPLHkAnbd682Wvfe++9qT7XXnut1zYzr33gwIHUNvfdd5/XPnr0aLND\nBAAAHcadIQAAEDSSIQAAEDSSIQAAEDSSIQAAELSeLqBevXq11162bJnXvnjxYmqbN99802vPmjXL\na8+cOTO1zSuvvFJ3mzNnzuSOFUBv2LBhQ8PbrFixog0jkQ4dOuS1+/r62nIcFHf//fen1g0ODrbl\nWMnzPzw87LWHhoZS22zfvr0tY6k67gwBAICgkQwBAICgkQwBAICgmXOueGez4p0r4Nlnn/Xac+bM\n8dp33313apvjx4977dmzZ3vt119/PbVNcmLGxx9/3GuvXbs2f7AV4pyz/F7N6bYYQtPecM7Nzu/W\nnCrF0cjIiNeeMWNGSSPJl6whkapdRxTKtaiZOrMiFi5c6LWLxObAwIDX7oEaokLXIu4MAQCAoJEM\nAQCAoJEMAQCAoPV0zVA7vPDCC6l1DzzwgNdO1hDcfPPNbR1Tq4XyOX0zmvlsP/m5fZbk/B8rV65s\n+DgV05M1Q1nnv5k5gjZu3Oi1s2IkGRNFNDOWZM1jlXAtaq+s3//J2AzlWsSdIQAAEDSSIQAAEDSS\nIQAAEDSSIQAAELSeflBrO2zZsiW1LllAfcMNN3RqOGihVhXHNiM5GVryAYo9MPFZ0DpVpNzMJHsI\nR9ZDYhHhzhAAAAgayRAAAAgayRAAAAgaNUMIQtZn5YODg+Peb3KCslZZt26d1+7v78/dhjqj9kv+\njKV0XVnyQZetkozhZIxI+TVC7RobukOR60iouDMEAACCRjIEAACCRjIEAACCRjIEAACCRgE1gpBV\nbJp06NCh1Lq+vr5xH7sVha9FJs9LFvJW+Wnk3SqrKL3IzzlrQs9anZrck6J6IBt3hgAAQNBIhgAA\nQNBIhgAAQNC6tmbo1ltv9drPP/987jarVq3y2i+99FJLx4TuNjw8nFqXV+uRpayHZWbVPKHznHNl\nD2FMWfG8cuXKEkYCVAt3hgAAQNBIhgAAQNBIhgAAQNC6omZo0qRJqXXJep+rr7461eezzz7z2jt2\n7PDaa9asSW3z1FNPee0pU6Z47U2bNtUfrKjdqKKsc518UGu75npJxkNWbVKybmNkZMRrF6k7ynqN\naK2sB/4WmcOqHYpcZ5JxUyTGqSHqXZ2az6obcWcIAAAEjWQIAAAEjWQIAAAEjWQIAAAErSsKqJcs\nWZJaN3XqVK998ODBVJ/Fixd77QMHDnjtxx57LLXNO++847WPHz/utW+88cb6g5W0b9++3D7orCIP\n2Mwqju3v7/faQ0NDdb8vtaYAtZmJGnkIZ/tlne9OTaq5ceNGr50VZ3kPBc4aa7KoNhnjxBVCwJ0h\nAAAQNJIhAAAQNJIhAAAQtK6oGXrxxRdT686ePeu1n3vuuVSfY8eOee3k5+m7d+9ObbNt2zav/fnn\nn+eO7+233/ba69evz90G1ZNVG5FXL9GqeoqseqU8AwMDLTk2ypc1gWJyEs0isZbXJznRaJZknRE1\nQwgBd4YAAEDQSIYAAEDQSIYAAEDQzDlXvLNZ8c4tdNNNN6XWvfXWW177k08+SfWZM2eO107OK5P1\n+fntt99edyyXLl1KrVuwYIHX3rt3b919VJ1zzvJ7NaesGKq6vAezJueYkSr/QM03nHOz27XzKsVR\nkYfqJmuCmqkHaoVGrvejktfNTuJa1FpFzn+R+ay6TKFrEXeGAABA0EiGAABA0EiGAABA0EiGAABA\n0Lpi0sXDhw+n1p04ccJrT58+PdXngw8+aPlYPvroo9S6/fv3t/w46F1ZEyx26mGfaL2+vr6yh1BY\n1uSOxB7AnSEAABA4kiEAABA0kiEAABC0rqgZynLnnXd67awJFG+77ba6+zhy5Ehq3cWLF7128vP0\n6667LrXN8uXLvfbTTz9d97gIW39/f8PbDA0NtWEkCM3w8HBqXV7N0IYNG1LremAiPsDDnSEAABA0\nkiEAABA0kiEAABC0rq0ZSs4hdMcdd6T6TJkype4+Pv3009zjvPrqq1573rx5qT5PPvmk1z5//rzX\nfuaZZ1LbnDt3LvfY6E0rVqzI7ZOcD6ZTD/LE+GXNI8X5A6qNO0MAACBoJEMAACBoJEMAACBoJEMA\nACBo5pwr3tmseOceMXnyZK+9Z8+eVJ/58+fX3UfWRGfNTLzXKc45a9e+Q4yhkZERr13kwZhmbTsF\nnfKGc252u3ZepThKnt8sVX6YayO/A0Z1Kj65Fo1PcsLMIv+8sXHjRq/dAxNsFroWcWcIAAAEjWQI\nAAAEjWQIAAAErWsnXeyUCxcueO177rkn1Wf37t1eOzkxY9ZEjQ8//LDX3rx5c7NDRMUVqRFKTrKI\n7lHk/CZrN8qqw8h66CoA7gwBAIDAkQwBAICgkQwBAICgUTPUoKwHrC5YsMBrJx/uOnfu3NQ2ybkc\nHnnkEa+9du3a1DY7d+702sl6JlRD1oM686xZs6YNI0EnDAwMeO3BwcFUn+T8LgsXLvTaWXORtaKu\nqJl5ZpKS1yqgF3FnCAAABI1kCAAABI1kCAAABI1kCAAABI0HtXbAQw89lFr3xBNPeO3rr78+dz+H\nDx/22qdOnUr1efDBB732yZMniwzRw8MRxyfvwaxZEyxW+UGeTQrmQa1JzTz4tIgiE3MWmQCy0eOU\nGZtci8aHB7VK4kGtAAAA+UiGAABA0EiGAABA0Jh0sQO2bt2aWrdr1y6vfcstt3jtVatW5e530aJF\nqXXTpk3z2s3UDGF88uo2sibYQ+8wS5e55NWRFdGKeqAsVaoRAsrCnSEAABA0kiEAABA0kiEAABA0\nkiEAABA0Jl1EChOdjU/eeypr8rzkU+u3b9/e0jGVINhJF5uRnBxPSj/ZvpkC6uQEekNDQ6k+VY41\nrkWtVeT3fdY/AHQ5Jl0EAADIQzIEAACCRjIEAACCRs0QUvicfnyaeVBnqJ/TNyuEOALXIrQENUMA\nAAB5SIYAAEDQSIYAAEDQeFAr0GLJeYSS88MMDAx0cjgAgBzcGQIAAEEjGQIAAEEjGQIAAEEjGQIA\nAEGjgBposb6+vrKHAABoAHeGAABA0EiGAABA0EiGAABA0BqtGfqppPfbMRBUxpfavH9iKAzEEcaL\nGEIrFIqjhp5aDwAA0Gv4mAwAAASNZAgAAASNZAgAAASNZAgAAASNZAgAAASNZAgAAASNZAgAAASN\nZAgAAASNZAgAAATt/wA15A7UTq5LDAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1a2c0fb470>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "Ierr = np.where((yhat_ts != yts))[0]\n",
    "nplt = 4\n",
    "plt.figure(figsize=(10, 4))\n",
    "for i in range(nplt):        \n",
    "    plt.subplot(1,nplt,i+1)        \n",
    "    ind = Ierr[i]    \n",
    "    plt_digit(Xts[ind,:])        \n",
    "    title = 'true={0:d} est={1:d}'.format(yts[ind].astype(int), yhat1[ind].astype(int))\n",
    "    plt.title(title)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that again a human would not have made these errors, but the digits in error are much less clear."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Looking at the Support Vectors\n",
    "\n",
    "Let's take a look at the support vectors.  We see there about 5000 support vectors.  So, about half the training samples were used as SVs, more than that for the linear kernel.  This is partly why the prediction was so slow."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4815, 784)"
      ]
     },
     "execution_count": 57,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "S = svcrbf.support_vectors_\n",
    "S.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "We can plot some of the vectors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAACOCAYAAAAy0AzYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAADztJREFUeJzt3X9MlWUfx/Hr4K9s/qAfxEjCympW\nVlJYypxWayunacmWrVX0g0ypWRG1TMVcfxhZW2UpLhdGhqazLZuFWpvTxYyoxXCTYek0C3EaoWEQ\n6v389zzP93vdcp8D5xzO4Xq//vvcu677XM883Xyfmy/XFfI8zwAAALgqpbcXAAAA0JsohgAAgNMo\nhgAAgNMohgAAgNMohgAAgNMohgAAgNMohgAAgNMohgAAgNMohgAAgNP6RzI4FAqxXbUDPM8Lxere\nfIeccczzvLRY3ZzvkRt4FiEKwnoW8WYIQCwc7O0FAIAJ81lEMQQAAJwW0a/JEJ7W1laRBw4caI0Z\nPHhwvJYDAAC6wJshAADgNIohAADgNIohAADgNIohAADgNBqoI5Sbm2tdq6ysFLm5uVnkgoKCmK4J\nAAB0H2+GAACA0yiGAACA0yiGAACA0+gZClBeXi7yfffdZ40ZPny4yMXFxSLv3Lkz+gsDgDjo39/+\nMVFWViby448/bo1JSZH/X7u0tFTkV155JQqrA6KDN0MAAMBpFEMAAMBpFEMAAMBp9AwpK1asEDk/\nP1/k06dPW3OeeeYZkfXv04H/p3vMjDGmtrZW5MzMTJGLioqsOStXrozuwuCkAQMGiKx7eUpKSqw5\nnueJXFNTY40ZN26cyDt27OjmChFPGRkZIt9zzz0iz50715qTk5MT8efU1dWJPG3aNGvM77//HvF9\nu4s3QwAAwGkUQwAAwGkUQwAAwGkUQwAAwGlON1CvW7fOujZlyhSROzs7RV6yZIk1h0ZWRMKvUXDU\nqFFdzlm+fLl1ra2tTeSKioqeLQx93pAhQ6xre/bsETkrK0vklpYWa44eo/+IxBhjtm3bJnJVVVXY\n60R86D8YMsaYvLw8kdPS0mLy2TfeeKPImzdvtsboZ2VTU1NM1mIMb4YAAIDjKIYAAIDTKIYAAIDT\nnOoZys3NFVn3BxljzLBhw0R+4YUXRH733XejvzD0afqgy+zs7IjvoQ+9NMbeDI2eIWhTp04V2W9D\n2BEjRoisD5aeMWOGNUf3q7355pvdXSJiaOTIkSJXV1eLrDdYDIdfD9n+/ftF1v0/qamp1hz9s1X3\nEBljzNChQ0WmZwgAACBGKIYAAIDTKIYAAIDT+kzPkN/vJL/66iuR9WFyP/30kzVH92GcOHEiCquD\nSy644AKR6+vrRU5PT4/K50ycODEq90Hfoff7KS0tFfm8886z5lRWVopcUFAgcnt7e5RWh1jSPbHG\n2HvghdMj1NHRIfLWrVtFnjNnjjXnyJEjXd6zX79+1jXdvzZ58uTAtcQSb4YAAIDTKIYAAIDTKIYA\nAIDTKIYAAIDTkraBWm/GdNddd1ljxo8fL/K+fftEnj59ujXnr7/+isLq0FdddNFFIn/55ZfWmFtu\nuUXkAQMGBN53165dIusm7DFjxlhzYrkBGRLfqlWrrGuPPfaYyLpxdf78+dacZcuWRXVdiA3d/K6b\n4+fOnWvN0Ru+agcOHLCuPfvssyJ//fXX4S7xnM6cOWNd0z+PdY433gwBAACnUQwBAACnUQwBAACn\nJW3P0Pbt20W+9dZbA+csX75c5KNHj0Z1Tej7JkyYILLuSwtHc3OzdW3mzJkiL168WGS/nqGff/45\n4s9G8tAbyRYVFYms+4OMMebs2bMiL1q0SGT6g5KD33Nl/fr1ImdlZQXeR29aqHsTZ82aZc3xO4jV\nBbwZAgAATqMYAgAATqMYAgAATqMYAgAATkvaBmq9EVQ4DdRr166N1XLgCH2idzg2bdoksj5Z3Bh7\nY8a8vLzA+x46dCjitSB5rFu3TuS7775b5NbWVmvOHXfcITJN9snhmmuuEfnbb7+1xgwePLjLe7S3\nt1vXiouLRV6xYkU3VhcfoVDIuqb/N586dSpmn8+bIQAA4DSKIQAA4DSKIQAA4LSk7Rny+/2itm3b\nNpHb2tpitRw4wvO8wDH//vuvyHpjM314pjHGfP/99yJnZGSIfPLkSWvOli1bAteC5LBw4ULr2p13\n3imy7gl54IEHrDn0CCWnNWvWiBzUH2SMMYcPHxZ59uzZ1piqqqoerStc+kDY22+/XeTc3FxrzujR\no0W++eabrTF6Y9Hdu3d3b4Fh4M0QAABwGsUQAABwGsUQAABwWtL2DC1YsCBwTF1dncidnZ2xWo5w\n7bXXilxYWBjxPfz2RNJ9JYg/fRCmn+PHj4usf39eUlJizRk7dqzIuj9E7zFjDP0hyUz3COkDVY0x\n5syZMyLrvgv+/fuOlJTI30v88ccfIuvDy7tr4MCBImdmZors15uk90VLT08X+fzzz7fmbNiwQWS/\n/wZi2SOk8WYIAAA4jWIIAAA4jWIIAAA4jWIIAAA4LWkbqONl5syZIuvmMj/vvfeeyBdffHHEn+vX\ncEYDde+rqKgQ+d5777XG6A0Tv/nmm8D76s0c9YGb/Nsnt+zsbJFffPFFkfVBvcYYU1lZKbJumL7h\nhhusObp5XzfZIjE1NDSIHM7B43rM3r17rTFLly4Vuby8XGS/Q6Nff/11kVNTUwPX0tHRIXJNTY3I\nb7/9tjVn8+bNgfeNJ94MAQAAp1EMAQAAp1EMAQAAp4XCOXjyv4NDofAHx5jeQNHv8Mtly5aJrDd1\nGjVqlDXno48+EjknJyfwc2JBb7hmjH9fQSx4nhd8Cm43JdJ3KBry8/Ota/r38vpQYb9/2+uuu07k\nxsbGKKyuV/3oeV5O8LDuSeTvkd5wzhhjvvjiC5F1v8eePXusOWVlZSLrg1tnzJhhzdGHBNfX13d5\nT2PsQ0ITiavPonnz5lnXXn31VZEvueSSHn+O34HnQTXBrl27rGu6B662trZnC4uusJ5FvBkCAABO\noxgCAABOoxgCAABO69M9Q3pfjhMnTog8adKk6C8sSugZSg5+/RaPPvpol3P8/purqqoSWR98qA9u\nTQLO9AzpZ8+qVausMU888USX9zh16pR1zW+vsf/X1NRkXdP7oA0fPlxkvx6Rp59+WmTdN9mbeBb9\njz6s+cMPPxRZH6gaLTt37hRZ74FmTHCfUS+jZwgAACAIxRAAAHAaxRAAAHAaxRAAAHBaUhzUOmbM\nGOuaXyOgNnbs2B5/9oEDB0Q+ffq0NUZvmKY3Pps1a1bg5+iD7ubPnx/uEhFH1dXVIodzoOKxY8dE\n9jvsd8qUKSKvXr1a5IcffjjcJSLOdANyULO0n0GDBlnX9AG/H3/8scgbN2605ug/sigpKRH5pZde\nsub4PdPQu4YMGWJdmzNnjsjhNEzrxuZ9+/YF3kM37k+cOFHk/fv3W3Ouv/56kf3+ICDR8WYIAAA4\njWIIAAA4jWIIAAA4LSl6hpYsWWJdS0mJTR2ne5GOHj0qsu7tMcaY8ePHi5yRkRH4Ofp3ubpH6J13\n3gm8B2JPH4Z52223iezXu9bS0iKyPlDR74Dguro6kSdMmBDROhEffhufLl68OHDe33//LbLeMO+t\nt96y5vhtqhhEP3sefPBBkVtbW605v/76a8Sfg+i68sorRfbrB8vOzhZZb8y7Y8cOa47uX/vtt99E\n9ut5XLhwocjTpk0TeeTIkdYc3fO4adMma0yi480QAABwGsUQAABwGsUQAABwWlL0DMWT7uXRPSLd\n4XfoaiIfjugqv/2g1q5dK7LuETp+/Lg1Z/LkyV1+zlVXXWVd03tTHTx4sMt7oHeUlZVZ19LS0kT+\n5ZdfrDHTp08XuaGhocdreeSRR6xrun/psssuE9nvEOHvvvuux2vBuenDcgsKCqwx+t/Nb58hvXdP\nRUWFyIWFhRGvrbGx0br2+eefi6z7gc6ePWvNCadPNtHxZggAADiNYggAADiNYggAADiNYggAADgt\nKRqoi4qKrGs33XSTyP369Qu8z9ChQ0W+8MILe7awc1i0aJHIW7dutcbU1tbG5LMRPt1c6rfRpf5e\n6c3z9AGFxtjNj2+88YbI8+bNs+bozfw++eQTnxUj3p566imR8/PzrTH6oNPi4mJrTHcapvVzTx/U\nefnll1tz9KaKetPFZNwML9npRnf9PDDGfs4cOnTIGvPyyy+LvGHDhh6v7bnnnrOu6Z9feoNjv407\n9SHCyYg3QwAAwGkUQwAAwGkUQwAAwGlJ0TPktwGd38Z1QUaPHi2yPoTTjz4kVm/KaIwxr732msiH\nDx8Wmf6gxFRTUyNyenp64Jwnn3xS5JUrV1pj9MGGuh9I95gYY8zs2bNFLi8vD1wLYk/3afgdEK0P\nyNS9aMYYs337dpH1JnUjRoyw5gwbNkxk/b0pLS215rz//vsiNzc3W2MQW3qj3qVLl4rs19+qN+p8\n6KGHrDH654i+T2ZmpjXn/vvvF/n5558XOSsry5qj6e+d/nlnjDEnT54MvE+i480QAABwGsUQAABw\nGsUQAABwWsivB+acg0Oh8AcjaXmeFwoe1T2J9B36559/RB40aFDE9/A7tLCjo0Nk3bcxdepUa87e\nvXsj/uwE96PneTmxunm8vkcbN24UOS8vL3CO3zNVH/Cr/fnnn9a19evXi6z3L2pvbw9cS7JLxmdR\nSUmJyH49Ntq4ceNEnjRpkjXm0ksvFVl/F/32nQrid4i4PjRa98198MEHEX9OLwvrWcSbIQAA4DSK\nIQAA4DSKIQAA4DSKIQAA4LSk2HQR6Cm/DRXDOdxX6+zsFFk3YRtjzBVXXCFyS0tLxJ+DxLBmzRqR\nc3LsPszU1FSR/ZrhGxsbRV6wYIHITU1N1hy/5nwkvnA2MtR++OGHHn+u3/elvr5eZP3HHJ999pk1\nx9UNX3kzBAAAnEYxBAAAnEYxBAAAnEbPEJzQ1tZmXTty5IjIfgcd6s3wCgsLRd69e7c1hx6hvmPL\nli1dZkD79NNPRb766qtFzs3NteaE079YXV0t8urVq0VuaGiw5vg9n+CPN0MAAMBpFEMAAMBpFEMA\nAMBpFEMAAMBpnFoPSzKeFI2E0ydOrUfv4lmEKODUegAAgCAUQwAAwGkUQwAAwGkUQwAAwGkUQwAA\nwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkU\nQwAAwGkUQwAAwGn9Ixx/zBhzMBYLQcIYGeP78x1yA98j9BTfIURDWN+jkOd5sV4IAABAwuLXZAAA\nwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGkUQwAAwGn/AeI7C6Fa\n/MWZAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1a2c0e08d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "nplt = 4\n",
    "nsv = S.shape[0]\n",
    "Iperms = np.random.permutation(nsv)\n",
    "plt.figure(figsize=(10, 4))\n",
    "for i in range(nplt):        \n",
    "    plt.subplot(1,nplt,i+1)        \n",
    "    ind = Iperms[i]\n",
    "    plt_digit(S[ind,:])        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "In this demo, we specified the parameters for the SVC. In the lab, you will be asked to find the optimal parameters through cross validation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}