pca_reprojection_visualization.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# numerical\n",
    "import numpy as np\n",
    "from sklearn.datasets import load_digits \n",
    "from sklearn.decomposition import PCA\n",
    "\n",
    "# plotting\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython import display\n",
    "import matplotlib.gridspec as gridspec\n",
    "\n",
    "# To make matplotlib figures appear inline in the notebook\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Digits data\n",
    "digitsX, digitsY = load_digits(return_X_y=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAM0CAYAAACbDZ8eAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X90nmd93/HPN4QfKQHLJcDGj1pJWPnVHctLuhbOmOUt\nGYXBIq91mpay2KXDLaedTbc1KYPFPpQdswFW1nZturNGhrasCT2Rx88Sd7FZaDtITuTT9QzaLpaB\nhvAjkURCKRvk2h/3rYOiPPb1VXxLzy193q9zdOJIX9339VzSR9etR9/nuqOUIgBezhv2AACsPYIP\nGCL4gCGCDxgi+IAhgg8YWvfBj4hfj4i3dV0LLBcRr4uIj6+X455VKaW3b5JmJX1D0kOS5iX9oaSf\nlnReB8cel/SFYT/G5ON/WNKXJE1JunDJx18p6RPt/HxF0glJ/2TA4yySrkuc7+mSJiV9rj3n/2n/\n/6LHMfZrJP1vSV9vj/OKATX/th3bFck5uH/5HKy3N0mj7WM+f5jjWA8r/mtLKU+TtEXSIUnXSfov\nwx3SmnptKeVCSX9H0uWS3ipJEfEjkm6V9F5Jz5P0bDVBeu2yz79W0oOS/tnZThIRT5L0B5JeKumH\n1PwQeJmkByT93ZUMOCKulPROSXskPU3S35d077KaSyXtkvTFxCEX52BM0jZJv7iS8WRFxBNW47i9\nNOyfgJWfjrNathqo+SZ8RNL3tf8/JemXlnz8F9R8M90n6afU/HR9wdJaSU9Vs4o8omYleVjSc9pj\n3yXpa2pW2Pf06fFL+g+SPiQp1KzK/7ry+U9VczVwjaT/K+nys9T+VPuYz3k1VXNl9oZKzcckvXrQ\n17gyB/9e0oeX/P+TJb2rnY8vSfp1SRe0HxuX9AVJb5H01fZYr1vyuVOSfk3SR9RcmVwhaZOaH6Zf\nkXRazQ/a89r63ZLuXPL5L5J0u5ofrJ+VdPWSj10g6d3tMRYk3dm+73Pt9+Ti993LBhz35ZI+3X7e\npyW9fMnHjkt6u6RPtl/bj+txXJGthxX/UUopn1LzxXzF8o9FxA9J+nk1X8AXqPnCDzrG1yW9StJ9\npZQL27f7JN0o6cZSytMlXSrpllV5EI9DRDxfTVDukfRCSc+X9IHKp/1TNd9ct0r6fTWr/5lcIelj\npZSHzzKGD0XE/BnePtTWPEHNlckzI+IvIuILEfErEXHBkuPskvTNUspHqg/80ed/npqv218sefch\nSd+r5mrgBZKeq+bKZ9HfkHRR+/5rJf1GRLxwycd/XNI71FyZ3Cnpl9WE/xJJ29VcKe0ZMJanqgn9\n70h6lpofrv8pIl7SlrxL0mVqQvzdahakR9Rc/UjSSPt990fLjvvdkj4s6T9Keoak90j6cEQ8Y9mY\n97TnfZKkfzVovs5m3QW/dZ+ayVzuakk3l1L+tJTyV5IOrPC4/0/SCyLiolLKw6WUPz7HcXZhOiLm\n1XxTnpD079R8Q0j1y+RrJf1uKeXbar5Br4mIJ56h9hm145VSXlNKGTnD22vasmdLeqKkH1Hzw3nx\n8nzxV5SntY9hX2XsS01HxEOSPi/py5JuaI8Vkt4o6c2llAdLKQ+1x75m2ee/rZTyzVLKCTWhunrJ\nx46WUj5ZSnlEzdf/Gkm/WEp5qJQyq2bVfv2AMb1G0mwp5eZSyrdKKfdI+j1JuyLiPEk/KWlfKeUv\nSynfLqX8YSnlm4nH+o8l/Xkp5X3tcd8v6TN69K9wN5dS/qyU8g01i9NY4riPsl6D/1w1l1fLPUfN\nN8eizw+oOZs3qFk9PhMRn46I19Q+YQ1MtMHaUkp5U/vFfqD92N880ye1Vwg7JP12+66jkp6i5htr\nkAfOdrwV+Eb7318upXyxlPJVNavWq9v3H5D0vjZUWROleZ5nXM3l9UXt+58p6bsk3b145aHmV4hn\nLvncufYKb9FpNd8ni5Z+j1yk5ofW6WX1zx0wpi2SfmDpVY+k1+k7VxhPUfOk5ko9Z9n5B43h/iX/\n/itJF670JOsu+BHx/Wom4c4BH/6imie6Fj3/LId6zMsSSyl/Xkr5MTWXUO+U9IH2kq5vPqvmG/aH\nz1LzejVf3w9GxP1qnlx7is58uX9M0ivP9ngj4qMR8fAZ3j4qSaWUOTW/ii2d36X//oeS/kVE3N+O\n6/mSbomI6876iJtjn1Dze/m72nd9Vc0PmpcuufLYVJonAhdtXvaYvkfNFeOgsX1Vzaq/ZVn9Xw4Y\nzuclnVh21XNhKeVn2uP8tZpfFx/zMCoP875l5z/bGB63dRP8iHh6uwL/V0m/VUr5kwFlt0jaExEv\njojvknS2v9l/SdIzImLTknP8REQ8s73sm2/f/UhHD6EzpXmW5+clvS0i9rRzc15E/L2I+I227FpJ\nB9VcBi6+/bCkVy/7fXHR+9R8M/9eRLyoPd4zIuItEfHq9ryvWvKcyPK3Vy051s2Sfi4inhURmyW9\nWc2TklIT/O9bMqb7JO2V9KvJhz8p6cqI2Np+nf6zpMMR8SxJiojnRsQrl33OwYh4UkS8Qs0l+q2D\nDtz+SnSLpHdExNMiYouaef6tAeUfkvS9EfH6iHhi+/b9EfHidly/Kek9EfGciHhCRLwsIp6s5knD\nR9Q8hzDIR9rj/nhEnB8RPyrpJfrO/HViPQT/g0t+v/s3ai4bH/NkiySVUj6q5kmRO9Q8AbT4O/pj\nfrcqpXxG0vsl3dteqj1HzZ+x/jQiHlbzRN817aV175RSPiDpR9X8Lnmfmh9kvyTpaET8oJpV41dL\nKfcveftvaublxwYc75tqnuD7jJonrb4m6VNqLlv/5wqH93Y1z0b/mZq/5d+j5gk0lVIeWDomSd9W\nczl+xicVl43zK2qedV98Au+69jH9cUR8Tc2Vy9In7+6XNKdmjn5b0k+3X/sz+Tk1z/Dfq+aq8nfU\nhHj5OB6S9I/UPCdwX3ued6r5K4PUPOH2J+08PNh+7Lz2uad3SPpk+333g8uO+4CaH07/Us2vX78g\n6TXtr0ydifZPBBtSRLxY0v+S9ORSyreGPR6srYgYV3N1+LxabfJ4PynpJ0op/6CL4w3TeljxVyQi\ndkbEk9tLzHdK+iChR0deKunUsAfRhQ0XfDW/L35ZzTOq35b0M8MdDjaCiJhW86vgu4c9li5s6Et9\nAINtxBUfQAXBx8Z6uSlSuNQfooiYVdPi+m01PfUfk/Sz2T9t9UlEjKp54uuJ6+3J1GVfh0VTpZSf\nHc6IVh8r/vCt+ktOrV5u+vi9dllD0mNCHxHnZ953NiutXy0EvyfaZpbfV/uCi/ZPku+KiM9FxJei\n2T3ogvZj4+2r3t4SEV+NiNmIeN3isSJiKiJ+LSI+EhFfl7QjIjZFxHsj4isRcToi3tq+mEQRsTsi\n7lzy+S+KiNsj4sGI+GxEXL3kYxdExLvbYyxExJ3tuD7Rlsy3LbwvG3Dcl0fzGoiF9r8vX/Kx4xHx\n9oj4ZEQ8FBEfj4jFnvyhaMf/yYg4HBEPSDpwhved187n6Yj4cjvPm9pjjEZEiYg3RMTnJP33YT6m\nRQS/J+KxLznl5aaP4+Wmq+AH1HTxPVtt9+GA9+1u33aomd8LJf3KsuNsl/RiNbsmDd+ZXqjP2+q/\nqdkY4mE1GyoUNTvgjKjZaOPrki5dUvsySafKdzaY+Jakpy75+C1qXn4qNS9kee+Sjz1BzUYcL1ny\nvr2Sjrf/3q12Iwg1bcD/Y9k4b1LzUtjz1LwoZuuAxzKqZVtKLTvu6yV9atnn/JGk3e2/j0t665KP\nvUnN/gBr+XWYX/L2z9vxf25Z7aD3/YGkNy35/xeqebHP+Uvm5ZJhf78tfevF7xvmJkopxyJiu5pV\n9iI1q93iy00X60JNgBet+stNl7zvfDUv4un9y03PwUQp5djSd0TEbg1+affy9y1/bKfVzNmzz/I5\nQ8Wlfk+UR7/klJeb9segx7X8fcsf2/eouSL7UuU4Q0Pw+2VS0pWS/rZ4uel68n5Jb46IiyPiQjW7\nAP1u6fGfNQl+j5RHv+SUl5uurQ/GozcWuW0Fn/uban4V+oSaXoa/VjPfvUUDzzoUvNwU54gVH9IG\nerkpcnhW31z7ctO/pebmFjDBpT5giEt9wNCqXOpHRCeXEbt25a4+Dx06VK05duxYteb666+v1szN\nzaXGlFFKCeYqp8u52uhKKVGrYcUHDBF8wBDBBwwRfMAQwQcMEXzAEMEHDBF8wFCve/UzzSaSdMkl\nZ3oJ+Hds3ry5WvPggw9Wa66++upqjSTdeuvAl8avGuYKK8GKDxgi+IAhgg8YIviAIYIPGCL4gCGC\nDxgi+IChoTXwXHbZZdWaTLOJJF166aAbuzzavffeW625/fbbqzWZcUvdNqUwV+gaKz5giOADhgg+\nYIjgA4YIPmCI4AOGCD5giOADhgg+YGhonXuZ7Z3uvvvu1LEynWYZ2fOtNeYKXWPFBwwRfMAQwQcM\nEXzAEMEHDBF8wBDBBwwRfMBQrxt4jh07tgYj+Y7MmObm5tZgJI/GXKFrrPiAIYIPGCL4gCGCDxgi\n+IAhgg8YIviAIYIPGBpaA0+muSN777WMTMNJ5nzDuM8bc4WuseIDhgg+YIjgA4YIPmCI4AOGCD5g\niOADhgg+YChKKd0fNKJ60EsuuaR6nLvuuit1vr1791Zrdu3aVa3JjOnyyy9PjSmjlBLMVU52rtDM\nVa2GFR8wRPABQwQfMETwAUMEHzBE8AFDBB8wRPABQwQfMDS0zr2MN77xjam66667rlpz9913V2uu\nvvrq1Pm60mU3GnOFRXTuARiI4AOGCD5giOADhgg+YIjgA4YIPmCI4AOGVqWBB0C/seIDhlblbrlr\n3Vp5/Pjxas3IyEi15sCBA9Wa6enpxIhyhtGGmpmrTZs2VWsOHjxYrVnvc7Ve0bILYCCCDxgi+IAh\ngg8YIviAIYIPGCL4gCGCDxhalQaetTY/P1+t2b59e7VmfHy8WtNlU8owzM3NVWsyc5WpWe9ztZGx\n4gOGCD5giOADhgg+YIjgA4YIPmCI4AOGCD5gqNcNPGNjY6m6TONNxszMTCfHGYatW7em6pgrSKz4\ngCWCDxgi+IAhgg8YIviAIYIPGCL4gCGCDxgi+IChoXXu7d+/v1qTuZedlLvXW0bmvnLDsG/fvmpN\ndq4y9xDM6OtcIYcVHzBE8AFDBB8wRPABQwQfMETwAUMEHzBE8AFDUUrp/qARnRw022ySuR9cxrZt\n26o1XW45VUqJruYq28SUuc9gRmZbtJMnT3ZyLqnbudroSilRq2HFBwwRfMAQwQcMEXzAEMEHDBF8\nwBDBBwwRfMBQr++dt9YyTSl9vWdcRLVno1Nr3cCDbrHiA4YIPmCI4AOGCD5giOADhgg+YIjgA4YI\nPmCIBh48LmvdMIRuseIDhgg+YIjgA4YIPmCI4AOGCD5giOADhgg+YIjgA4Z63bmXvc/b0aNHqzVX\nXXVVtWZ8fLxaMzU1lRjR2lvrudq+fXu1pq9zBVZ8wBLBBwwRfMAQwQcMEXzAEMEHDBF8wBDBBwxF\nKWXYYwCwxljxAUMEHzC0Kr36EbGmvz9keuynp6erNZlbYGfOlVVKiT7O1W233VatyczVjh07MkNK\nGcZcrVellOoWyKz4gCGCDxgi+IAhgg8YIviAIYIPGCL4gKFe77k3NjaWqrvjjjuqNQsLC9Wa0dHR\n1Pn6qMu5yuzft57nCqz4gCWCDxgi+IAhgg8YIviAIYIPGCL4gCGCDxjqdQPPxMREqu7kyZPVmsxG\nHDfccEPqfH2UnavMBhqZG2uu57kCKz5gieADhgg+YIjgA4YIPmCI4AOGCD5giOADhgg+YKjXnXuT\nk5OputnZ2U6OlelY66su5+rw4cPVmkwnJPqLFR8wRPABQwQfMETwAUMEHzBE8AFDBB8wRPABQ1FK\n6f6gEdWDjoyMVI+zf//+1Pky205l7vWWqcncVy6rlBJdzdW+fftS51zLucrcrzArO1do5qpWw4oP\nGCL4gCGCDxgi+IAhgg8YIviAIYIPGCL4gKGh7cBz4MCBak22KSVj586d1Zoum3O61Me56rI5B2uP\nFR8wRPABQwQfMETwAUMEHzBE8AFDBB8wRPABQ0PbgWdsbKx6nKmpqdT5tm7dmqqrydxCKzumzC2m\nsrvK9HGuMo/vyJEjnR2LHXjy2IEHwEAEHzBE8AFDBB8wRPABQwQfMETwAUMEHzBE8AFDQ9t6a2Zm\nplqT6VjL1mW2r7rqqquqNbOzs4kR5brRstZ6rm644YZqTeYefKdPn06Nqcu5Qg4rPmCI4AOGCD5g\niOADhgg+YIjgA4YIPmCI4AOGVmXrLQD9xooPGCL4gKFV6dXvajfU6667LlV36NChas29995brbn8\n8surNXNzc6kxZXS5cyxzhUXssgtgIIIPGCL4gCGCDxgi+IAhgg8YIviAoaHtuZf5e/KuXbtSx9q7\nd2+15qabbqrWXHbZZdWaY8eOpcbUJeYKXWPFBwwRfMAQwQcMEXzAEMEHDBF8wBDBBwwRfMDQquy5\nl9kw4ZJLLqkeJ7uRw1133ZWqq7n00ks7OU5WdnMJ5oqNOFaCjTgADETwAUMEHzBE8AFDBB8wRPAB\nQwQfMETwAUMEHzA0tK23MrdpynSsZesy20Bt3ry5WtPlbaGymCt0jRUfMETwAUMEHzBE8AFDBB8w\nRPABQwQfMETwAUND23qrS5lmkttvv72Tc1155ZWpukzzyjC2k2KuNj623gIwEMEHDBF8wBDBBwwR\nfMAQwQcMEXzAEMEHDG2IBp6MTOPKTTfdVK3J7IYjSddff321pq9NKczV+kYDD4CBCD5giOADhgg+\nYIjgA4YIPmCI4AOGCD5gaGi30Mo4dOhQqq6rWz5dccUV1Zpbb701Naa1xlxhJVjxAUMEHzBE8AFD\nBB8wRPABQwQfMETwAUMEHzBE8AFDve7cy9xTTcptA5WR6TTbu3dvJ+fqGnOFlWDFBwwRfMAQwQcM\nEXzAEMEHDBF8wBDBBwwRfMDQqtw7D0C/seIDhgg+YGhVevXX+nbGIyMj1ZqpqalqzcTERAejyRvG\nrZ8zc3XzzTdXa3bu3NnFcNK4TXYet8kGMBDBBwwRfMAQwQcMEXzAEMEHDBF8wFCv99zL2r17d7Vm\nZmZm9QeyDmTm6uTJk6s/EAwVKz5giOADhgg+YIjgA4YIPmCI4AOGCD5giOADhnrdwJPZNELKNaVM\nTk5Wa0ZHR1Pny5idne3sWBnMFVaCFR8wRPABQwQfMETwAUMEHzBE8AFDBB8wRPABQwQfMNTrzr1M\nl5mU6yLL3EIr07E2Pz+fGJF04MCBVF1XsnO1ZcuWak1mrg4fPlytWVhYyAxpzecKrPiAJYIPGCL4\ngCGCDxgi+IAhgg8YIviAIYIPGBpaA8/ExES1JtMkIklHjhw51+FIkvbt21et2bNnTyfnWok+ztX+\n/furNcOYK+Sw4gOGCD5giOADhgg+YIjgA4YIPmCI4AOGCD5gaGgNPJmdbLI7uFx77bXVmrGxsdSx\naqanpzs5zkqs9Vxt3bo1dayaYcwVcljxAUMEHzBE8AFDBB8wRPABQwQfMETwAUMEHzAUpZTuDxrR\n/UHPItOcc/z48WpNpuEke6uqjFJK9HGu7rjjjmrN0aNHqzXrfa7Wq1JK1GpY8QFDBB8wRPABQwQf\nMETwAUMEHzBE8AFDBB8wRPABQ0PbeqtLma2pNm3aVK2ZmprqYDT9lpmrkZGRao3DXG1krPiAIYIP\nGCL4gCGCDxgi+IAhgg8YIviAIYIPGFqVrbcA9BsrPmCI4AOGVqVXv6vdUDM740rS7OxstabLHV+7\n0uXOsdm5OnXqVLVmz5495zia7rHLbh677AIYiOADhgg+YIjgA4YIPmCI4AOGCD5gqNd3y838fV6S\ntmzZ0sXpdPr06WrN6OhoJ+eSuv3bdObv81J34898bS6++OJOziXxd/yV4O/4AAYi+IAhgg8YIviA\nIYIPGCL4gCGCDxgi+IChXt80M3ODRynXwLOwsFCtyWxmkbmhpJQfe1e6PF/mWOt5rsCKD1gi+IAh\ngg8YIviAIYIPGCL4gCGCDxgi+IAhgg8Y6nXnXnbrra1bt1ZrNm3aVK2ZmZmp1vS1yyw7V2NjY9Wa\nTMfdyZMnqzV9nSuw4gOWCD5giOADhgg+YIjgA4YIPmCI4AOGCD5gqNcNPBMTE6m68fHxak2mceXw\n4cOp82VMTk52dqyMnTt3puqYK0is+IAlgg8YIviAIYIPGCL4gCGCDxgi+IAhgg8Y6nUDT1bmPm5d\nGR0dXbNzrQbmChIrPmCJ4AOGCD5giOADhgg+YIjgA4YIPmCI4AOGet3Ak92BJ3OrpgMHDpzjaBrT\n09OdHKdrXc7VDTfccK7DkdTfuQIrPmCJ4AOGCD5giOADhgg+YIjgA4YIPmCI4AOGCD5gqNede5n7\nvEnSvn37OjnfkSNHqjVruXXVSqz1XE1NTVVr+jpXYMUHLBF8wBDBBwwRfMAQwQcMEXzAEMEHDBF8\nwFCUUoY9BgBrjBUfMETwAUOr0qsfEWv6+0OmJ3xkZKRaMzY21sFo8kop0ce52rRpU7Vm27ZtHYwm\nbxhztV6VUqJWw4oPGCL4gCGCDxgi+IAhgg8YIviAIYIPGOr1nnvZO8Bu3769WnPw4MFzHU6vdTlX\nXd1ZGP3Fig8YIviAIYIPGCL4gCGCDxgi+IAhgg8YIviAoVXZc6+rDRNmZmZSdVu3bq3WZDaOyJ6v\nK11uLnHPPfek6jKbjWz0udro2IgDwEAEHzBE8AFDBB8wRPABQwQfMETwAUMEHzBE8AFDve7cm52d\nTdXNz89Xa9b69lgZXXajnTp1KlWXmau1vj1WBp17eXTuARiI4AOGCD5giOADhgg+YIjgA4YIPmCI\n4AOGet3Ak2k2kXLbQE1PT3dSk20qyuiyKWVubi5V19VcHT16tFrT17na6GjgATAQwQcMEXzAEMEH\nDBF8wBDBBwwRfMAQwQcM9bqBp8t75508ebKT42R3p8mMva/3zsuMvat78GXPRwNPHg08AAYi+IAh\ngg8YIviAIYIPGCL4gCGCDxgi+ICh84c9gLOZmppK1R0+fLhak9kNZnR0tFozMTGRGFG++agrR44c\nSdVlGm+6mqurrroqMaK1nyuw4gOWCD5giOADhgg+YIjgA4YIPmCI4AOGCD5giOADhjZE516mi2z3\n7t3VmuPHj1drMveVG4Yu5+raa6+t1mTmKnN/PQwHKz5giOADhgg+YIjgA4YIPmCI4AOGCD5giOAD\nhlbl3nkA+o0VHzBE8AFDq9Kr39XtjPfv35+qGxkZqdZkdsfN3CZ7YWEhNaZMT/zc3Fxnt35mrrCI\n22QDGIjgA4YIPmCI4AOGCD5giOADhgg+YIjgA4Z6vdlm1vz8fLUm0+CSqck0wGTHNAzMFSRWfMAS\nwQcMEXzAEMEHDBF8wBDBBwwRfMAQwQcMrcpmm33cKeXAgQPVmszOM+Pj46nzZZpSSim93FWGuVrf\n2IEHwEAEHzBE8AFDBB8wRPABQwQfMETwAUMEHzBE8AFDQ9t6K9PVle38ysjeW64m07EmSVNTU52c\nT5K2b99ercnOVUS1qWtdzxVyWPEBQwQfMETwAUMEHzBE8AFDBB8wRPABQwQfMDS0Bp7Z2dlqzdjY\nWOpYXTX6ZBpOjh8/3sm5VuL06dPVmm3btqWOtdHnCjms+IAhgg8YIviAIYIPGCL4gCGCDxgi+IAh\ngg8Y2hD3zss8hp07d1ZrpqenuxhO2jDuB8dcbXzcOw/AQAQfMETwAUMEHzBE8AFDBB8wRPABQwQf\nMDS0HXgyJicnU3ULCwvVmo2+GwxzhZVgxQcMEXzAEMEHDBF8wBDBBwwRfMAQwQcMEXzAEMEHDPW6\ncy97n7fdu3dXa+bn589tMD3HXGElWPEBQwQfMETwAUMEHzBE8AFDBB8wRPABQwQfMLQq984D0G+s\n+IChVWnZ7equptkNJCcmJqo1U1NTnZyvy3bWLu8Ay1xhEXfLBTAQwQcMEXzAEMEHDBF8wBDBBwwR\nfMAQwQcMrUrLbleNFtm7to6OjnZxOs3OzlZrsnvbZXTZlJKdqy1btlRrIqr9H+t6rjY6GngADETw\nAUMEHzBE8AFDBB8wRPABQwQfMETwAUO9vmnmzMxMqi7TTNLVzSKzTSnZhpquZOfq1KlT1Zo9e/ZU\nazJztX379tSYTpw4kapDd1jxAUMEHzBE8AFDBB8wRPABQwQfMETwAUMEHzBE8AFDve7cy9zDTZLu\nueeeak1me65MN1qmS3AY1nqu5ubmqjWnT5/ODAlDwIoPGCL4gCGCDxgi+IAhgg8YIviAIYIPGCL4\ngKFeN/CMjIx0dqzMNlAXX3xxtaavDTxdzlVme7FMkw8NPP3Fig8YIviAIYIPGCL4gCGCDxgi+IAh\ngg8YIviAoSildH/QiOpBx8bGqsfJ7BYjSQcPHqzWZBpOMmOamJjIDCnV6FNKicxcbd26tXqs7L3z\nDhw4UK3paq527tyZGFG3c4Vmrmo1rPiAIYIPGCL4gCGCDxgi+IAhgg8YIviAIYIPGBpaA09mx5js\nbjeZhpNMTaZhKNMsJOUaZbJNKZm5OnXqVGpcmV2G1vNcgQYeAGdA8AFDBB8wRPABQwQfMETwAUME\nHzBE8AFDBB8wNLR7583Pz1drjh8/njrW3NxctWZhYaFac/To0WrN5ORkakxdyszViRMnUsfKzFXm\nfH2dK+Sw4gOGCD5giOADhgg+YIjgA4YIPmCI4AOGCD5gaFW23gLQb6z4gCGCDxhalV79rnZDnZ6e\nTtVldqEdHx8/x9F0r8udY2+77bZUXWauduzYca7D6Ry77Oaxyy6AgQg+YIjgA4YIPmCI4AOGCD5g\niOADhoZ2t9zMHVmzd4DtysmTJ6s1Y2NjnZ0v+7fpPs7VzMxMtWbbtm2dnY+/4+fxd3wAAxF8wBDB\nBwwRfMAQwQcMEXzAEMEHDBF8wNDQbpqZ2RAiK3PDyNnZ2WpNHzfrkLqdq8yNSNfzXCGHFR8wRPAB\nQwQfMETwAUMEHzBE8AFDBB8wRPABQwQfMDS0rbcy3Whzc3Op823evLlak7kdV2ZbrS676LLbSa31\nXGVux5WBsVFKAAABlklEQVSZq8y5sth6K4+ttwAMRPABQwQfMETwAUMEHzBE8AFDBB8wRPABQ0Pb\nemt+fr5ak7mXnZRrXrnxxhurNZmmlMx97KTc9lVZmbnK3MtO2vhzhRxWfMAQwQcMEXzAEMEHDBF8\nwBDBBwwRfMAQwQcMDW0Hni5lmkkyDS6Tk5PVmmxTysTERLVmGLvKdDVXhw8frtZk52rnzp3VGnbg\nyWMHHgADEXzAEMEHDBF8wBDBBwwRfMAQwQcMEXzA0NB24OlSV805u3fvrtZkGnP6rKvmnMxcZRpz\nMBys+IAhgg8YIviAIYIPGCL4gCGCDxgi+IAhgg8YIviAoV537mW67aTcdlIjIyPVmvHx8WpN9h51\nay3TbSd1N1c7duyo1vR1rsCKD1gi+IAhgg8YIviAIYIPGCL4gCGCDxgi+IChVbl3HoB+Y8UHDBF8\nwBDBBwwRfMAQwQcMEXzAEMEHDBF8wBDBBwwRfMAQwQcMEXzAEMEHDBF8wBDBBwwRfMAQwQcMEXzA\nEMEHDBF8wBDBBwwRfMAQwQcM/X+Fuss0/H3QjwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10910ce50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# select interesting digits\n",
    "dig_idxs = [np.where(digitsY==dig)[0][0] for dig in range(0,10)] \n",
    "\n",
    "# iterate over no. of pca components\n",
    "for comp in [1,2,3,4,5,6,7,8,9,10,20,30,40,50,60,64]:\n",
    "\n",
    "    # pca init\n",
    "    pca = PCA(n_components=comp)\n",
    "    # projection\n",
    "    digitsX_r = pca.fit_transform(digitsX)\n",
    "    # reprojection\n",
    "    digitsX_reproj = pca.inverse_transform(digitsX_r)\n",
    "    n_rows, n_cols = 10, 3\n",
    "    \n",
    "    # visualization\n",
    "    fig = plt.figure(figsize=(int(n_cols*1.4), int(1.4*n_rows)))\n",
    "    \n",
    "    plt.gca().cla()\n",
    "    \n",
    "    gs1 = gridspec.GridSpec(n_rows, n_cols)\n",
    "    gs1.update(wspace=0.025, hspace=0.005)\n",
    "    for i in range(n_rows):\n",
    "\n",
    "        plt.subplot(gs1[3*i])\n",
    "        if i == 0: plt.title(\"Digits \\n\")\n",
    "        img = digitsX[dig_idxs[i],:].reshape(8,8)\n",
    "        plt.imshow((img - img.min())/(img.max()-img.min()), cmap=\"gray\",vmin=0,vmax=1)\n",
    "        plt.axis('off')\n",
    "                   \n",
    "        plt.subplot(gs1[3*i+1])\n",
    "        if i == 0: plt.title(\"PCA C={}\\nReprojection\".format(comp))\n",
    "        img_r = digitsX_reproj[dig_idxs[i],:].reshape(8,8)\n",
    "        plt.imshow((img_r - img_r.min())/(img_r.max()-img_r.min()), cmap=\"gray\",vmin=0,vmax=1)\n",
    "        plt.axis('off')\n",
    "                   \n",
    "        plt.subplot(gs1[3*i+2])\n",
    "        if i == 0: plt.title(\"Reprojection \\nError\")\n",
    "        img = np.abs(img - img_r)\n",
    "        plt.imshow(img, cmap=\"gray\",vmin=0,vmax=1)\n",
    "        plt.axis('off')\n",
    "        \n",
    "    display.clear_output(wait=True)\n",
    "    display.display(plt.gcf())\n",
    "    plt.savefig(\"PCA_Reproj_Comp{}\".format(comp))\n",
    "    plt.clf()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"https://lh3.googleusercontent.com/-hE1hYMSjgJQ/WNwAkVbIoHI/AAAAAAAA0qI/RzqhXRPoz6QVhbem0z5yKDW8GpCHB2QiwCJoC/w576-h1024-rw/Digits_PCA_reprojection.gif\"/>"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}