II_02_Projection_onto_subspaces.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "+ This notebook is part of lecture 15 *Projections onto subspaces* in the OCW MIT course 18.06 by Prof Gilbert Strang [1]\n",
    "+ Created by me, Dr Juan H Klopper\n",
    "    + Head of Acute Care Surgery\n",
    "    + Groote Schuur Hospital\n",
    "    + University Cape Town\n",
    "    + <a href=\"mailto:juan.klopper@uct.ac.za\">Email me with your thoughts, comments, suggestions and corrections</a> \n",
    "<a rel=\"license\" href=\"http://creativecommons.org/licenses/by-nc/4.0/\"><img alt=\"Creative Commons Licence\" style=\"border-width:0\" src=\"https://i.creativecommons.org/l/by-nc/4.0/88x31.png\" /></a><br /><span xmlns:dct=\"http://purl.org/dc/terms/\" href=\"http://purl.org/dc/dcmitype/InteractiveResource\" property=\"dct:title\" rel=\"dct:type\">Linear Algebra OCW MIT18.06</span> <span xmlns:cc=\"http://creativecommons.org/ns#\" property=\"cc:attributionName\">IPython notebook [2] study notes by Dr Juan H Klopper</span> is licensed under a <a rel=\"license\" href=\"http://creativecommons.org/licenses/by-nc/4.0/\">Creative Commons Attribution-NonCommercial 4.0 International License</a>.\n",
    "\n",
    "+ [1] <a href=\"http://ocw.mit.edu/courses/mathematics/18-06sc-linear-algebra-fall-2011/index.htm\">OCW MIT 18.06</a>\n",
    "+ [2] Fernando Pérez, Brian E. Granger, IPython: A System for Interactive Scientific Computing, Computing in Science and Engineering, vol. 9, no. 3, pp. 21-29, May/June 2007, doi:10.1109/MCSE.2007.53. URL: http://ipython.org"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<link href='http://fonts.googleapis.com/css?family=Alegreya+Sans:100,300,400,500,700,800,900,100italic,300italic,400italic,500italic,700italic,800italic,900italic' rel='stylesheet' type='text/css'>\n",
       "<link href='http://fonts.googleapis.com/css?family=Arvo:400,700,400italic' rel='stylesheet' type='text/css'>\n",
       "<link href='http://fonts.googleapis.com/css?family=PT+Mono' rel='stylesheet' type='text/css'>\n",
       "<link href='http://fonts.googleapis.com/css?family=Shadows+Into+Light' rel='stylesheet' type='text/css'>\n",
       "<link href='http://fonts.googleapis.com/css?family=Philosopher:400,700,400italic,700italic' rel='stylesheet' type='text/css'>\n",
       "\n",
       "<style>\n",
       "\n",
       "@font-face {\n",
       "    font-family: \"Computer Modern\";\n",
       "    src: url('http://mirrors.ctan.org/fonts/cm-unicode/fonts/otf/cmunss.otf');\n",
       "}\n",
       "\n",
       "#notebook_panel { /* main background */\n",
       "    background: #ddd;\n",
       "    color: #000000;\n",
       "}\n",
       "\n",
       "\n",
       "\n",
       "/* Formatting for header cells */\n",
       ".text_cell_render h1 {\n",
       "    font-family: 'Philosopher', sans-serif;\n",
       "    font-weight: 400;\n",
       "    font-size: 2.2em;\n",
       "    line-height: 100%;\n",
       "    color: rgb(0, 80, 120);\n",
       "    margin-bottom: 0.1em;\n",
       "    margin-top: 0.1em;\n",
       "    display: block;\n",
       "}\t\n",
       ".text_cell_render h2 {\n",
       "    font-family: 'Philosopher', serif;\n",
       "    font-weight: 400;\n",
       "    font-size: 1.9em;\n",
       "    line-height: 100%;\n",
       "    color: rgb(200,100,0);\n",
       "    margin-bottom: 0.1em;\n",
       "    margin-top: 0.1em;\n",
       "    display: block;\n",
       "}\t\n",
       "\n",
       ".text_cell_render h3 {\n",
       "    font-family: 'Philosopher', serif;\n",
       "    margin-top:12px;\n",
       "    margin-bottom: 3px;\n",
       "    font-style: italic;\n",
       "    color: rgb(94,127,192);\n",
       "}\n",
       "\n",
       ".text_cell_render h4 {\n",
       "    font-family: 'Philosopher', serif;\n",
       "}\n",
       "\n",
       ".text_cell_render h5 {\n",
       "    font-family: 'Alegreya Sans', sans-serif;\n",
       "    font-weight: 300;\n",
       "    font-size: 16pt;\n",
       "    color: grey;\n",
       "    font-style: italic;\n",
       "    margin-bottom: .1em;\n",
       "    margin-top: 0.1em;\n",
       "    display: block;\n",
       "}\n",
       "\n",
       ".text_cell_render h6 {\n",
       "    font-family: 'PT Mono', sans-serif;\n",
       "    font-weight: 300;\n",
       "    font-size: 10pt;\n",
       "    color: grey;\n",
       "    margin-bottom: 1px;\n",
       "    margin-top: 1px;\n",
       "}\n",
       "\n",
       ".CodeMirror{\n",
       "        font-family: \"PT Mono\";\n",
       "        font-size: 100%;\n",
       "}\n",
       "\n",
       "</style>\n",
       "\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from IPython.core.display import HTML, Image\n",
    "css_file = 'style.css'\n",
    "HTML(open(css_file, 'r').read())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true,
    "format": "row"
   },
   "outputs": [],
   "source": [
    "from sympy import init_printing, Matrix, symbols\n",
    "from IPython.display import Image\n",
    "from warnings import filterwarnings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true,
    "format": "row"
   },
   "outputs": [],
   "source": [
    "init_printing(use_latex = 'mathjax')\n",
    "filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Projections onto subspaces"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Geometry in the plane"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Projection of a vector onto another (in the plane)\n",
    "* Consider the orthogonal projection of **b** onto **a**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "format": "row"
   },
   "outputs": [
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAhIAAAEOCAYAAADcwf7CAAAACXBIWXMAAAsTAAALEwEAmpwYAAAA\n",
       "GXRFWHRTb2Z0d2FyZQBBZG9iZSBJbWFnZVJlYWR5ccllPAAAIuxJREFUeNrs3QtwXNV9x/GzkrUr\n",
       "ydiS7Elj0gle2yQYEmoxNqVpAK0HEaWmjEWCgQkdECS1p2FaljYpTprihaZgT9sgJ2Smdplabppp\n",
       "wUwtDwkEbMfr4hBeDnJIAwHbrMzLgLElP2Q9bG3Pf++51vX1PuV93F19PzOHfWi1e3V38f3tOf97\n",
       "ji8ejytkx+fzNcql3md97A0AAJSqYhfkJKxbO7sBAADzJZseiRx2ls8Xk0u9z4LsDQAA6JHIJUR0\n",
       "6IuZ0vT1EHsEAAB6JHIJEjETJMR2vd8IEwCACY8eiexCRIcjRIgWfR+1EgAAjpH0SGQMEXKmRky3\n",
       "BtePenVr5gwOAMBERo9EZhFHiBgyTUgPRZjdgwkQpkPUBQFI+W8EPRJp/wGV4YuNum3SrUu3Zt1q\n",
       "dXtNWaeBLtZtod6HUfYWKuyzHzSf8bAJzbP05zzGngHgNoldkFGTPXyh/3GVIDGob0uo6DLDHnxT\n",
       "Q6UEh2bzeZY2z/Hj7YQIAASJcdD/eHZn+LkEjG72FMo4QEgY7tTt1jQP62JPAUiFGglgYoflPt06\n",
       "9NVZuq1P8TDCMgCCBIC0gSJmAsUduo04frSJM5MAECQAZGTmS7lfnb6eTJQ9A4AgASCbECG1EqF4\n",
       "PP6EGhvm6GHvACBIAMg2RNjBISL/4dRmAAQJALmGCGVO97yXPQSAIAEgpxDhCBMR9hIAggSAnEME\n",
       "ABAkABAiABAkABAiABAkABAiABAkABAiAIAgARAiAIAgAYAQAYAgAYAQAYAgAYAQAYAgAYAQAQAE\n",
       "CYAQAQAECQCECAAECQCECAAECQCECAAgSACECAAgSACECEIEAIIEAEIEAIIEAEIEAIIEAEIEABAk\n",
       "AEIEABAkAEIEIQIAQQIAIQIAQQIAIQIACBIAIQIACBIAIQIACBIACBEACBIACBEACBIACBEAQJAA\n",
       "CBEAQJAACBEAQJAAQIgAQJAAQIgAAIIEQIgAAIIEQIgAAIIEQIgAAIIEAEIEAIIEAEIEABAkAEIE\n",
       "ABAkAEIEABAkAEIEABAkAJRDiNDbEtKtmXcFAEECIETksh1B3br11W26dfHOACBIAISIbLcjoi/k\n",
       "9Rebu6K8OwAKYRK7AKicECHDGMrqfZhp7tquW5gaDQCFQo8EUAEhwjWMISGiX7fb9LZQ6AmAIAEQ\n",
       "ItK+flidPoyxWreg3pYu3h0AhcbQBlDGIcKI6dagJtgwhglQfebvTyYk/9H7I8InFSBIAISIFPRr\n",
       "d+ttWagvo2m2db6+aDU3t+jH7qyAt6DbtHlpHnMbn1SAIAEQIjKHiWiK7WzSF486QoR9vwSJZeUc\n",
       "KPS2x0xxqYSJliQP2cTwDlB41EgAZR4i0mynhIjNJkRIYFiu2wbzY+mh2Gx6KsqW3v99UlCqr+5y\n",
       "/ahXtw4+rQBBAiBEjC9AiKUmMEjPwwLdVul2g769wASLRG+F4/Fl+X7oFlNWrYRTh4QMPrEAQQKo\n",
       "6BChn7Mxz9so4WGPvrzb9ESs1du61vUtXkKEBIpDus02gaPc3guZ9juqr0akuXolNqWrFwFAkAAq\n",
       "JUTIWQcxc5kve01Pw0oTEpLWQOi/QR633Ny8u4zeB5kvo0tZdRFR/Xc4T3ONmcswn1iAIAFUbIgw\n",
       "36bleR5U1mmb7fnaVr19W6QXwtycneHhdr1Ek96eVo+/B41m2u83zV3BJKd1yj5dL0WYfGoBggRQ\n",
       "cSHCHAzl27PMPimnLMrsk3eZbvl8kp6GQ+b6/DShQx6zpUzeAwkHsp8u0dudqv6hx7xXAAgSQMWF\n",
       "iLA5GN5q7lpvvlXn/cBnAoI9bLG0XM/McPTcSM9De6bpvmU+DaYDBwgSQEWFiCTDGFIQuDDNt+p8\n",
       "hYm1jt6GlSm2TYY+7CGNnR7a90FTSCl1EF2mDiLKpxIgSAATridCWXMZOIcxmgt1UJRgoNtmc+aG\n",
       "cvRKtOr71rgea09UJVaZXoxS73cZ+pH9LnUQso8K0mMDIL+Y2RIoXIgQEXMZLlQPhAkFSx09D1tM\n",
       "r8RO/bNVyjorY6kpqJQCS3n8EnO5QT9uuQf2u+ynsAkQsyiYBOiRAAgR1sE85hzGyPfkT/r5JBC8\n",
       "ZEKEBIjExFOOh8h1Z29Dqwkd9hTZN5R4n7ebCaVCyqqDaCdEAPRIAISIM59fCh7XmIP61XnsibCH\n",
       "LJa5J54yQeaQftxy87jZJmjs9MD+luAgvRBBZU0o1cWnEKBHAiBEpJdYgdNRw3BG2JAZKU1bkunJ\n",
       "TF2D9CjMSRYiHI/LWHhZxH192oRSujUTIgCCBECIyHzQl14Ae8hhpXOIQ65LkaQaG6KQJmtg7Ml0\n",
       "6qZMQJVloaSz8HJpCfazPaGUvY8TE0qxHgZAkAAIEdmz6xUkRDinpU61QqcMReRlhc50QaZI+zmm\n",
       "rDqIUKFPfQVAkAAqMUS4J4q62x7KUMVboTNVkCnUPs5pQikABAmAEJE5TLjrFYq2QqcjyBxSp5/J\n",
       "ke/9y4RSAEECIEQU0Kl6BdO2pDjw532FThNY3KeI5mvfMqEUQJAACBFF6JVw1itkkvcVOk1Ayfe+\n",
       "jSirDiKorAmlKKQECBIAIaKAnBNFtaY56Ht6hU4mlAJAkAAhogTKfYVOU0gZNfs1Ygopo3zaAIIE\n",
       "QIgoXpgouxU6mVAKAEEChAhvyXaFzuWlXKEzyYRSzdRBABCstYGMBxA5aGT5DTtKiDhtG+RMi7Xp\n",
       "AkCaFTrtngh7hc5VJd6XEiJiHg9lAAgS8CAJEduyPeYQIhKvLwHAXiRLQsTaDL8iIWGpCQ17TYCQ\n",
       "momdKskcE0X8O0JmP0qYlGXQu/nfAQBBArmSruvt5bKxpQwRpp5Bah2WmEBwtayFken3XCt0yu0F\n",
       "Jd6HQXsfKquQkrkgABAkMD7mYBwiRGR8bRmasCeNWp7rUIT0OpjFtLaUcP9Jz0NEtzt1W62sCaWo\n",
       "gQBAkMBZHVykm32zOdj5PLydJQkRrmGMDSZE7HU9RoYsZpuJqNKFiQUl3H9hEyKiyppQKsanH0A2\n",
       "OGsDlRB2ih4iZBhDt0ftkKWsYYwbkoQI6WXYo8bOvvDavrMnlJJ9yIRSAAgSIETk+PtyWmPYdOtn\n",
       "+zsyhPGSsooipQdijrsWwqzsudn0VtgLcHlpvzW7JpRqZkIpAAQJECJy/DaurHkRHlRWt342vyPD\n",
       "FHYtwxmLX8nPdVtpgoa9PPjVmYY1irjP7AmloooJpQDkATUSmHAhwpyVIAfPFnNXrzmoZmTOsFiQ\n",
       "bG4IM4whIULChgSMVaWcRMq1bYlTOE3rNgEixicJAEEChIjcD6YrHHffK8+Vy9kJ7nBg1smQAGEv\n",
       "Cb7cKz0Qjv0VUUwoBaAAGNpArgelJVJkqFvctJfMN3Gvhwh7GMMOETI3xlktd10GwxiysFaPCRFh\n",
       "s7AWIQIAQQIlOzBJkaGcfTDbfPPeaw6ga9zrRBTm9S/4epOa9f1xhIiovtio20xlDWNcJwdV+fxL\n",
       "EJCiSNPuznGTJEDJ78gwxpxSzUCZ5O+VOggZvpDWpbcryKyUAAgS8AI5aMo37gXmm/cc+RZuH1QL\n",
       "1TPh86lmn2/XDqVefOgD9UZVXMVzfYqouZRhDKkN6Haclil/U6tpK00PS7ZLektwkH2x3Au1EOYM\n",
       "FOmteVNZwxhBZqUEQJCAl6xyn+ZovoWvdQSNfIUHfVBUHbpJz8PLSs37vFLvrp2kqncngoEvMe9B\n",
       "tuRgemoYQ4Zn1NhpmXPMRFtzzG0JEY+aszPSkvDgoWGMsB0ezN8aZlZKAAQJeM3aDPfPNutNnG14\n",
       "kG54+Ya/Trd5Y4+YK/UIcnBsSPzMp7L6ti0HVNcZCvI89loYe81jErdNmLDXzPA8Uwchf5sEKyaU\n",
       "AlB0nLWBrKXqvjdLYds3Z5uDdLbhQVYXDclBUI2djpnMpnhcf+P2Kee37Dv1bet34yqW3eslhi1k\n",
       "G08bjnDMD2EPa+z18nshE0opq6dFLsPMBQGAIIGKZ0JD0Bz8QuayIctftw+UMtSx2HH/vMR9PhXR\n",
       "YSKbHoqmJAflJaYHwi4iXeae6tpDAaLRBIh2+5IhDAAECZR5QBgrTky2bLYOEPaCUA3jfIn+eFyl\n",
       "O+tAnvdBHSbk4NqRZe/EfDMMI7USraYHQtbK2ODhAMGEUgA8hxoJ5HIwS3VWRqu5THoQjls9BdIT\n",
       "IWdN9I/jpbM9dbHlVO9EanZxpPRC7DHbvlxZZ19scP29TR7Z7x3K6okJmR6IDkIEAIIEytFKs2y2\n",
       "8yAnZ2rYhYkp51HQYaJPt4hu8s36NmXN5zCeIJHpACq9Eyt0mIiZ+gnXdiTqIuzAINflrI1UU1mv\n",
       "LMb8GGkCRMjMgSHByJ5QKsrHEABBAuVoi/k2v9nMtSCXexwhYlmyYY0UoaJLt6AJFLsyPNw9rJHt\n",
       "N3GZfGqbDhPRJIHCXmirSY31ppzWE2GWCM96XozDN38qpFvwbHaweY72XYuD9zkmlIoyoRQAL6NG\n",
       "AtnaKRMvmR4I59kN8u1+bbYhwh0orHkiTjynP4qBLHojxqPFBAqZErtLJUJM4iyT5SYErTG9LBtM\n",
       "D0Wr+fuazN+1LIfX6tFBQHpEel2BR4ohZWgiaJpT0ISeU77zqwNysVpZE0pRSAnA03z6Hyr2QrY7\n",
       "y+eLWAfAeIS9kY/9KWdtDD+rlL9OqaF3lQp8IsnDbpPAMfZLid6FbTm+VL85kFvNPJ+p+Ug2dCGF\n",
       "l8tzLbzUIaLZBJ+Z490na17v7/3blz4MUQMBgCBBkEDWIeLY60pNvswc6N0H4Vlx51kY6YPEQd1e\n",
       "UdaU2LFEi6dfHtycuSGFl7NNj8Tes1kzQ4eJRvP683L93Y+GTu6a9djeZj4dAMoJQxvwRIiQYkyf\n",
       "dbbFOsdDe+NnnsoZTPPU05ScHhlXWS/oZeaLWJWvv23qj9+QoYhmHSi69OWtOfxq//RAdYhPB4By\n",
       "Q7ElPBEirIN6YsjBeTZHst6EZEFCijbXm+ueWKRKB4oOZZ3umq12E0IAgCAB5BoiHJxBIJuehdWm\n",
       "5iFibrckO+2zRGFCtumObP4G/dgonxAABAng7EKEkFDQnyZINDqur9chIpy4Zg2B3Ot4jpI7fPOn\n",
       "JEjcn+FhvY4QBAAECeAsQoQy93emCRJ2MeJ2HR46XD/rNCFkZi7LjPt8vmCeA0S7bhJsVujW8MKB\n",
       "wcE0D+9gSAMAQQLIQ4hwB4IUj5MeCZnEqv3MFHJaCInoMNGYIUA06talr74ps0jmIUDIhFJRfXWj\n",
       "hJl3Bk6Mfv/Vvk2tT71Vp5LP5LmJIQ0A5Y6zNuC1EKHMGRypehSCiZbqeeKJANGhrNNI7cXCkoWI\n",
       "iBpbAGu1GluSezwBImheJ3GWxuGRUbVp39Fd7w+eXHxfzwE7QLhPbe03r39asMlhG3qYrAoAQQKE\n",
       "iNRh4swZLX2JEBFSmZ9HDurrEgdqnw4IjsebngfphYgpawGsqLk/Jotj6dtdOQQI54qciZVNf/r2\n",
       "sfffODx8/T0vH9jhPvCr05c/75z64zdirsdIiMh2sq2FSil6MwAQJECIyFo8y3U25CwOq1eixYSK\n",
       "sKmD6DIH63CSwCCP75Y1LrL5pq9DRId57kQvg9RBvHxw6FvffPGDVKef9rh6I5I9Tl53e5Z7g94I\n",
       "AAQJECIKKGK+3d95qe/Sqcqaa+Je0wvRd2YPSDyqQ0SPSjMcYgJEyPxcQoqSOoifvzew7o7n3v9a\n",
       "a/rtibl6I5Jtg71UOACUz7/7TJGdy0GSKbLLJESYDUp0/bdsUVv6r1ZXN2dav8L0WrypEtNyn/7Y\n",
       "ZHUQOkBsv3XH/vZsaxX0c8RNb0SQMzUAVArO2kBlhghLh/ynVbU2xFU8mOnBJjxIr0WX4+DfaOaD\n",
       "eNMOETpA7LnrhQ/n3fLMe6EcCx53peqNSBJqlspS5rod1C1umtyez6cNAEEChIhiiKvYcXX8v+Tq\n",
       "frX/oSx/K3H2hhRlmjoICRcrEing4FDffbs++vP2n79z/obY4V+PY4tiKovJskxYWKOsJc136iZL\n",
       "tMuaILK42Etm2XMAIEiAEFFodaru2yNqZGCGmvGZXl/vzVn0SvTddVHTvz5/zXlPKevMj4Z3B06c\n",
       "+MGrfT+44sl9Tf/8m4MPn8XmhJOcqZFsGyQ8LNOX03S72rQ5+j57VdKVfPIAeOY4QI1ELgdNaiTK\n",
       "KUSMbWiitmHFIXWovynelHKSKqmDGBqN/zBQ5Vskt4+MjMa3vjfwvzKE4ZHPnyx1vsfclJBxiE8h\n",
       "gFLjrA1UdoiwdJ5UJ/+mSTU17PHtuW9OfM49rgBhzwexQoeIxH3/03t0996jI62OCaVKGR6kuWsj\n",
       "5PYWPokACBIgRBSa3rZqX/Vf6mvrzlXnLlc+9T17kiqpgzgRjz80yeebLLd3fjR49Jn3j/9Jkgml\n",
       "ih0gpA7ibmXVSQCAZ1EjgcoOEWNhomtUjb5Vr+prdqqdnTIfxMGvnP+2hAsJETIfxJrX++9a+LO3\n",
       "pnggRCzVF5tNr8Ny3RbE43EZhvTxCQRAjwQIESVLzVV/pS82zpt0wS2Dx37/1kmT30nMB7H53YHH\n",
       "btvx3pILvbOpdjHlDTo8MHwBgCABQsRZvnpQWVNbx5SK94zj9xtH3/VfPvyfU78aeK9d+V6Yq7/d\n",
       "H1E/efvYK6/1D19b6joIV2+EDGU0yXV3iDD1EgBAkAAhIjcyUZRPlg3fKKdgGNmsSdFyqjfiE8Oq\n",
       "9uYD6uS+H6lfvfixD5984/iXSj2EkcJeR3BYosPEBnNdwsWjfBoBECRAiBhfmOjUWyOvKxNGNThD\n",
       "QrZODlTFN/xj3XdveuStey716P7WwWGvDg0yX4TUSchMlnavxHwTMnaqM8/gAICSodgSZRAiTh1m\n",
       "u5S1qFX/eH67un70SzetOXSP1/e7DhPLlFVkKcFBhjpkSGOVvn+BvrTnjmCYA4A3jhVMSJXLgXXi\n",
       "TEjl7cJKX1D/p1u3eTn80mr9zoX5FANAftEjgTILEYkoF1NWz8T2LH9hFyECAAgSIEQ4w4Tepngo\n",
       "3jfplSwe3ME7CwAECRAiThl+ounGk7HaQ77GExdneOi94ztlFABAkEDFhYihx6ddfuKVybv9i/r+\n",
       "uzo41Dg6UKU+6PF/kOLhvTpERHh3AaBwOP0TZREihrqnzaw6d2S9v3XgSl/daGIyifefqTv47q9r\n",
       "br7kjsOD+ua2JL/WwbsLAAQJTPAQMbK94WF/y8CtvqYTic/rkd8GRt55wb9ybsfhez5+xam/xP1r\n",
       "65WKR3mHAYAggQkaIoafagxXXzC4oqZlSJb5VicOTIrvfrr2yblfOXLN3IvOeLhMcz3TcTvCOwwA\n",
       "BAlMwBAhdRDVwaEuf9vAHLktdRC9W+t+94mT/ra5XzmYal2MHkeQuNecIgoAIEhgooSIVHUQzzxZ\n",
       "Hbr+/iOZTvOUILFYWbNedvIuAwBBAhMoREgdRM2VA7dWTTN1EK8Ghg+8GPiLWbf0//v1V2T1FPYp\n",
       "nt3WHBMAAIIEKj5ESB3EpD84fn9Ny3Cd3JY6iNefqOu66JbDt0+5MKenipnLCO80ABAkUOEhQuog\n",
       "Jn16cIO/7fgMuZ2YD2JH3a4ZXzjafNEt43lGmXTKt4naCAAgSKCCQ4RdBxG49uipZcAPPFt3YMoH\n",
       "tQtmtKcspMwW62kAAEEClRoiEvNBtB293S6kPPpqYOjDFwNflzqI/LwCvREAUGxMkU2IKDipgxh9\n",
       "LzBQ03L4qxIipA7ita4p/3DOhYO1+QsRAIBSoEeCEFEwiTqIiwce97dZE0rFB6rUvqfrd8xsP3LF\n",
       "3A7eIwAgSIAQkSxAdE+bWX3ecHfg2mPN9n3vba3fN+1I4MqZZ18HAQAgSKBSQ8TIjqkb/W3HFtt1\n",
       "EMdeCwy+vjVw3SV39P+MdwcAKg81EoSIvEjUQRysGam5/Ei7sw5i8tzBulQhwufzhXSLZ9lCvMMA\n",
       "QI8EKixEJKuD6H2ivjt4/ZHrqIMAAIIECBHJA0SSOoi3npj8f7837L8meH12dRDxeGKZbx/vIAAQ\n",
       "JDCBQsTI81Oi/raxhbX6X649XPt23TWfvPbgDt4RAJhYqJEgRGRteGvjA4k6iMuOtkiIGNrnP7H/\n",
       "sanfbrjkeEOAEAEABAkQIpIGiCeabjwZqz3kv6p/uazOKXUQ+x+f8kjgvKGaGdf3PzC+7fe1miLK\n",
       "zeb2Et0edRRXvqTbUt5pAPA2hjYIESkl6iDmDG31LxqYY9+3/+nJu5oG/Itn5HE+CB0Y1ugLCQ07\n",
       "ddui22zd5uu2Rv9sfjweX8a7DgAECZRRiHDXQRz8Zd0Hkw/UfnlG/ocw5pt2tQ4MWxzhQoJFImDo\n",
       "6zv1z9by7gOA9zC0QYg4PUBsb3h49GDNCbsOQiaUkjqIaZ8b+HiB6iCadFvuDBHCBAc7PNzNuw8A\n",
       "BAl4OETIhFIn3qztk4W1qqadqJYJpaQOQiaUGm8dRLbS9DbY98/2+Xyz+RQAgPcwtDHBQ4RMKFUd\n",
       "HOryt1l1EFJIefCXdTvPOeL/8owSr4uhA8ZOHSDsmxIk9vJpAACCBDwQIqSQsurckfX+1rE6iL7n\n",
       "6/bXfVC7ZDqncgIACBKEiFSkDqLmyoEOGcKQ21IHUdNb+63Gtr5Ob/39vvn2dXcNBQDAG6iRmEAh\n",
       "IkkdxOhHj099UOog/CUMETKHRIof2fNIbOATAQDeRI/EBAgRiYW1Pj24wd92fEbi2/1AlTryYv3z\n",
       "U1uO/NH0az3xJ8t8EYdcp3/e7QgSnPoJAAQJFDtE2HUQgWuPttj3HX6x7t3AO7V/PLXEhZQuEiA2\n",
       "y3wR+vKQsgor7bM0ljGsAQAECRQ5REgdhL/t6O12IeXx1wMD1bvrbp+66NAj6lJv/d06KNzg6IGw\n",
       "6yJkOGMtIQIACBIoYoiQOojqiwceqGkZqZXbJz+aNNr3bP3q6df2/7X6tHf/fh0YVumLVXwSAIAg\n",
       "gRKECKmDqDp/8DF/2/GPJw7MMh/Ejvqnp3/hSJtH6iAAAAQJeC1EJKuDOPhM/e7JHwVap3urDgIA\n",
       "QJCAl0LEyI6pG/1txxbbdRCDbwSO+V6r/+I0JpQCABQJ80iUYYiQOojEwlqXH2mXECF1EMNPNdxV\n",
       "+6nBcwKECABAEdEjUUYhIrEuxmcHfuJvG2qQ21IHcez5+p+es/DIn1a3ldc+MGdj+Pg0AABBAgUO\n",
       "EVIHUX3ecHfg2mPN9n1Hn5v825r9/kXnUAcBACBIECJScddBDLxSe7A6Vrf4HIYwAAAeQI2ER0PE\n",
       "8NbGB5x1EMNv+UekDqL+4uPTqYMAAHgFPRIeCxHDTzTdWHXh8TX+q8bqIA4/M/lHDW2Hb1GfZF8C\n",
       "AAgShIgkISJRBzFnaKt/0cAc+77D2895PnCo5sYG6iAAAAQJQkSqEDHy/JSov23gSncdxFSGMAAA\n",
       "HkeNRAlDhNRBnJQ6iMuOtkiIkIW1qIMAAJQTeiRKECKkDsI39/ha/1VDU+W2TCg1+pv6dXUt/V/z\n",
       "8sJaAAAQJEoYImRCKd95Q//hXzQwS25LIeXxl+t/Uf1hzc0B6iAAAAQJQkSyEGEvrOVvHauDOP5y\n",
       "XW/V27V/Vs8QBgCAIIFUIWJke8PDk64c6KiedqJaHje0O3DUt6f27+va+jrVJew3AABBghCRJEQM\n",
       "/ayxY2Tv4L01LafXQQSkDuJ89hkAgCCBJCHiue8vuvOSK1/6lX8edRAAAIIEsgwRH5u+b+/ujQs+\n",
       "nDL/oyd99aOJnyfqIHprW+oJEACAClUR80j4fL6QbuFShYi7vv4v/fv3zZ499YoPPy8hQuoghp9o\n",
       "vKnukoEgvRAAgDI9tnbo1l7xQUL/kRF9sU23aLFDRCj0bN1rv5kb/94Pv9FQVX8yUQcxvLVhZeD8\n",
       "wSn+RYce4WMIAChjclzdqI+znbo1VlyQkJSkW0xfXaHbrng83lOcV75uxpzZv3uuu/uGum3bFqoL\n",
       "PvM7n9RBjDx/zvbq6SPV/qv6vsVnDwBQ7vRxVY6xm3S7U7ce6aGoiCBhAkQiJek209zdVYzX/kzw\n",
       "i+d/974Llu7eMzewePGmxH3Du+rfHH66IVhz2ZEQHzsAQIXpNpdyvF0nX+DdgcKnE0c5hAfpUpEN\n",
       "DzvCg9Nq3foKvR3/9Hc3fOcb3300MR9E/2/rh++888SP128Z3sfnDABQoWrNsTfgul/q/zrli3xZ\n",
       "nLWhw06fDhMSFGIpgkShBXW74Qf/9lh1KPTZ0V9s+mhb+KH3mJESADARBJLc12d/gS+LHgknOUND\n",
       "X0R0a3HcfZf+OzoL8FoSILp0+0O5qdt63fbr14rwuQIAVDozjLHOcdd2OQbr42DUvqPsaiRk43WT\n",
       "MLFQt13m7nCed1yjVKnqq2+axDWi2+ckRPCxAgBMIB3mUoYyFsrx1xkiyjJIuAJFs756l24z9YG/\n",
       "OU8hIqKsIZSgbt/UTUJLS/HOCgEAoPRMr7z0/t+rj4FBd4Cwlf3MljKkYc7ikCDRcxY7TCbd6DQh\n",
       "ot0ECbkdIkQAACYgOa4uTBUgKiZImDDRM94Q4ai5kOAg4z5dZkyIEAEAmLD08a87m8dVTdQdJF02\n",
       "unUp6xxZSVvNhAgAAHIz4RbtMnNShE2TECHjPn3mZ4QIAAAIEilDhB0UetxhgRABAABBIlWACJmQ\n",
       "IL0R7e7CEUIEAADjU9E1EqYOQkKDDGF0JTt9hRABAABBwh0gnBNKSXAIJpv5khABAMDZqbihDTOh\n",
       "VNgEiFlmGVRFiAAAgCCRLkCcNqFUugk0CBEAABAk7FAQUq4JpTI8nhABAECelG2NRKoJpQgRAAAU\n",
       "T9n1SKSbUIoQAQAAQSJTGIgoqw4i60BAiAAAYAIHCdeEUuFsFxIhRAAAUFierpEwdRASGpwTShUt\n",
       "RJgAk+kxzXyMAAAECW8FCOeEUjGVYkKpQoYIo1k/T0y2xRkY5LoUeuomtRkECQDAhOXTB1mvhQgp\n",
       "oowo60yMcKoJpYoQIuznkh6QxeZmv7lsMJeb9PO38zECABAkSh8g7Aml+kyAiI7zefIWIszzBZW1\n",
       "WmiD60cSKrI6YwQAgEpV8qENM0wQNQd/mVCq2SshQpgekWTDKmFCBABgoitZj4T5ph/Rze6J6Dyb\n",
       "A3Ohz86QWgl9MdPc7JXCTz4+AICJrug9EqaQUgKEfbCXHoiIl0OEEUlxHQCACauoPRKuCaXCeaph\n",
       "KEaIsGfUlO3uozcCAABLUSakOpsJpbwQIoT0mJi1PWJ8bAAAMMfiQvZImDqIxIFeWYWUnXl87qKF\n",
       "CNff00eRJQAAloLUSORjQimvhQghZ3AQIgAAKGCQMBNKJcKDbrP0gTevp0mydgYAAN6RtxoJ14RS\n",
       "7eOdC4IQAQDABAoSZg0KObjLpfQ+dBViQwkRAAB4z7iHNkwdhISGqGlBQgQAABM8SJilu5szBIiI\n",
       "GjsN8qwnlCJEAABQnpINbXSYy54UB3U7RBSkDoIQAQBAmQYJM3ujnHXR7bo/ZAJEUOVxQilCBAAA\n",
       "5c09tCEhosEEBnuYo9sEi6hMDU2IAAAAZwQJM2vjCnOzxdRBOCeUihRjgwgRAACUD+fQRpfrZzcp\n",
       "a0KpWLE2hhABAEB5qTIH8Ii+aHH97AI1VnhJiAAAAGcev5W1oNa2NI9Zrw/sBQ0UhAgAAMo3SEjx\n",
       "pJytEU3xmKCyVrwMEyIAAMBpx/FCLiNOiAAAoLJVleqFCREAABAkCBEAABAkCBEAAMDjQYIQAQAA\n",
       "QYIQAQAAihckCBEAABAkCBEAAKB4QYIQAQAAQYIQAQAAihckCBEAABAkCBEAAKB4QYIQAQAAQYIQ\n",
       "AQAAihckCBEAABAkCBEAAKB4QYIQAQAAQYIQAQAAihckCBEAACBtkNBhodkEBkIEAADILUhojbqF\n",
       "CREAAGA8QSKk2zwdHkKECAAAkGuQaDaXEUIEAABIxqeDQfIf+Hx9+qLB3BzQ7fOECAAA4FSVIkS0\n",
       "O0KEqNetM1nxJQAAIEi4tSe5r0W3dTpMxEzQAAAAE1zSoQ3XsIbo1S1qWkz/TpRdBwAAJiUJER0m\n",
       "REh46NKtm9oIAACQVZDQgrotpNcBAABk8v8CDADQOAM3nI09hgAAAABJRU5ErkJggg==\n"
      ],
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Image(filename = 'Orthogonal projection in the plane.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Note that **p** falls on a line, which is a subspace of the plane &#8477;<sup>2</sup>\n",
    "* Remember from the previous lecture that orthogonal subspaces have A**x** = **0**\n",
    "* Note that **p** is some scalar multiple of **a**\n",
    "* With **a** perpendicular to **e** and **e** = **b** - x**a**\n",
    "* Thus we have **a**<sup>T</sup>(**b** - x**a**) = 0 and x**a**<sup>T</sup>**a** = **a**<sup>T</sup>**b**\n",
    "* Since **a**<sup>T</sup>**a** is a number we can simplify\n",
    "$$ x=\\frac { { \\mathbf { a }  }^{ T }\\mathbf { b }  }{ { \\mathbf { a }  }^{ T }\\mathbf { a }  }  $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* We also have **p** = **a**x\n",
    "$$ \\mathbf { p } =\\mathbf { a } x=\\mathbf { a } \\frac { { \\mathbf { a }  }^{ T }\\mathbf { b }  }{ { \\mathbf { a }  }^{ T }\\mathbf { a }  }  $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* This equation is helpful\n",
    "    * Doubling (or any other scalar multiple of) **b** doubles (or scalar multiplies) **p**\n",
    "    * Doubling (or scalar multiple of) **a** has no effect"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Eventually we are looking for proj<sub>**p**</sub> = P**b**, where P is the projection matrix\n",
    "$$ \\mathbf { p } =P\\mathbf { b } \\\\ P=\\frac { 1 }{ { \\mathbf { a }  }^{ T }\\mathbf { a }  } \\mathbf { a } { \\mathbf { a }  }^{ T } $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Properties of the projection matrix P\n",
    "    * The columnspace of P (C(P)) is the line which contains **a**\n",
    "    * The rank is 1, rank(P) = 1\n",
    "    * P is symmetrix, i.e. P<sup>T</sup> = P\n",
    "    * Applying the projection matrix a second time (i.e. P<sup>2</sup>) nothing changes, thus P<sup>2</sup> = P"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Why project?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(projecting onto more than a one-dimensional line)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Because A**x** = **b** may not have a solution\n",
    "    * **b** may not be in the columnspace\n",
    "    * May have more equations than unknowns\n",
    "* Solve for the closest vector in the columnspace\n",
    "    * This is done by solving for **p** instead, where **p** is the projection of **b** onto the columnsapce of A\n",
    "$$ A\\hat { x } =\\mathbf { p }  $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Now we have to get **b** orthogonally project (as **p**) onto the column(sub)space\n",
    "* This is done by calculating two bases vectors for the plane that contains **p**, i.e. **a**<sub>1</sub> and **a**<sub>2</sub>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Going way back to the graph up top we note that **e** is perpendicular to the plane\n",
    "* So, we have:\n",
    "$$ A\\hat { x } =\\mathbf { p } $$\n",
    "* We know that both **a**<sub>1</sub> and **a**<sub>2</sub> is perpendicular to **e**, so:\n",
    "$$ { a }_{ 1 }^{ T }\\mathbf { e } =0;\\quad { a }_{ 2 }^{ T }\\mathbf { e } =0\\\\ \\because \\quad \\mathbf { e } =\\mathbf { b } -\\mathbf { p } \\\\ \\because \\quad \\mathbf { p } =A\\hat { x } \\\\ { a }_{ 1 }^{ T }\\left( \\mathbf { b } -A\\hat { x }  \\right) =0;\\quad { a }_{ 2 }^{ T }\\left( \\mathbf { b } -A\\hat { x }  \\right) =0 $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* We know that from ...\n",
    "$$ \\begin{bmatrix} { a }_{ 1 }^{ T } \\\\ { a }_{ 2 }^{ T } \\end{bmatrix}\\left( \\mathbf { b } -A\\hat { x }  \\right) =\\begin{bmatrix} 0 \\\\ 0 \\end{bmatrix}\\\\ { A }^{ T }\\left( \\mathbf { b } -A\\hat { x }  \\right) =0 $$\n",
    "* ... **e** must be in the nullspace of A<sup>T</sup>\n",
    "* Which is right because from the previous lecture the nullspace of A<sup>T</sup> is orthogonal to the columnspace of A"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Simplifying the last equations we have\n",
    "$$ {A}^{T}{A} \\hat{x} = {A}^{T}{b} $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Just look back at the plane example in &#8477;<sup>2</sup> example we started with\n",
    "* Simplifying things back to a column vector **a** instead of a matrix subspace A in this last equation does give us what we had in &#8477;<sup>2</sup>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Solving this we have\n",
    "$$ \\hat { x } ={ \\left( { A }^{ T }A \\right)  }^{ -1 }{ A }^{ T }\\mathbf { b }  $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Which leaves us with\n",
    "$$ \\mathbf { p } =A\\hat { x } \\\\ \\mathbf { p } =A{ \\left( { A }^{ T }A \\right)  }^{ -1 }{ A }^{ T }\\mathbf { b }  $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Making the projection matrix P\n",
    "$$ P=A{ \\left( { A }^{ T }A \\right)  }^{ -1 }{ A }^{ T } $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Just note that for a square invertible matrix A, P is the identity matrix\n",
    "* Most of the time A is not square (and thus invertible) so we have to leave the equation as it is\n",
    "* Also, note that P<sup>T</sup> = P and P<sup>2</sup> = P"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## Applications"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Least squares"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Given a set of data points in two dimensions, i.e. with variables (*t*,*b*)\n",
    "* We need to fit them onto the best line\n",
    "* So, as an example consider the points (1,1), (2,2), (3,2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* A best line in this instance means a straight line in the form\n",
    "$$ {b}={C}+{D}{t} $$\n",
    "* Using the three points above we get three equations\n",
    "$$ {C}+{D}=1 \\\\ {C}+{2D} = 2 \\\\ {C}+{3D}=2 $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* If the line goes through all points, we would give a solution\n",
    "* Instead we have the following\n",
    "$$ \\begin{bmatrix} 1 & 1 \\\\ 1 & 2 \\\\ 1 & 3 \\end{bmatrix}\\begin{bmatrix} C \\\\ D \\end{bmatrix}=\\begin{bmatrix} 1 \\\\ 2 \\\\ 2 \\end{bmatrix} $$\n",
    "* Three equation, two unknowns, no solution, **so** solve ...\n",
    "$$ { A }^{ T }A\\hat { x } ={ A }^{ T }b $$\n",
    "* ... which for the solution is\n",
    "$$ \\hat { x } ={ \\left( { A }^{ T }A \\right)  }^{ -1 }{ A }^{ T }b $$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1 & 1\\\\1 & 2\\\\1 & 3\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡1  1⎤\n",
       "⎢    ⎥\n",
       "⎢1  2⎥\n",
       "⎢    ⎥\n",
       "⎣1  3⎦"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = Matrix([[1, 1], [1, 2], [1, 3]])\n",
    "A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1\\\\2\\\\2\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡1⎤\n",
       "⎢ ⎥\n",
       "⎢2⎥\n",
       "⎢ ⎥\n",
       "⎣2⎦"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b = Matrix([1, 2, 2])\n",
    "b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\frac{2}{3}\\\\\\frac{1}{2}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡2/3⎤\n",
       "⎢   ⎥\n",
       "⎣1/2⎦"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(A.transpose() * A).inv() * A.transpose() * b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Thus, the solution is:\n",
    "$$ b=\\frac { 2 }{ 3 } +\\frac { 1 }{ 2 } t $$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.4.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}