Exams/Midterm-Revised-Solutions.ipynb

{
 "metadata": {
  "name": "",
  "signature": "sha256:d764361672de78f0d813fed685b18cac3bbe7553d270ea72bfebebf81be2dd6d"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Instructions\n",
      "----\n",
      "\n",
      "This is a \"closed book\" examination - in particular, you are not to use any resources outside of this notebook (except possibly pen and paper). You may consult help from within the notebook using ? but not any online references. Violation of any of these exam rules is a violation of the University's honor code and will result in penalty under that code.\n",
      "\n",
      "\n",
      "You have 1 hour and 45 minutes to complete the exam.\n",
      "\n",
      "- <font color=red>Use a new cell for different parts of a question</font>\n",
      "- <font color=red>The maximum score for the mid-terms is 60</font>"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import os\n",
      "import sys\n",
      "import glob\n",
      "import matplotlib.pyplot as plt\n",
      "import numpy as np\n",
      "import pandas as pd\n",
      "import scipy.linalg as la\n",
      "%matplotlib inline\n",
      "plt.style.use('ggplot')"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.set_printoptions(formatter={'float': '{: 0.3f}'.format})"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%load_ext rpy2.ipython"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "**Question 1 (10 points)**. \n",
      "\n",
      "Given the 2 matrices\n",
      "```\n",
      "A = np.array([[1,2,3],[4,5,6]])\n",
      "B = np.array([[1,2,3,4],[5,6,7,8],[9,10,11,12]])\n",
      "```\n",
      "\n",
      "Perform matrix multiplication of `a` and `b` using the following methods:\n",
      "\n",
      "1. Using nested `for` loops without the `dot` function (4 points)\n",
      "2. Using numpy (2 points)\n",
      "3. Using R (start the first line of a new cell with %%R) and passing in the A and B arguments defined in Python (4 points)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[1,2,3],[4,5,6]])\n",
      "B = np.array([[1,2,3,4],[5,6,7,8],[9,10,11,12]])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 37
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "C = np.zeros((A.shape[0], B.shape[1]))\n",
      "for i in range(len(A)):\n",
      "    for j in range(len(B)):\n",
      "        for k in range(len(A[i])):\n",
      "            C[i,j] += A[i,k] * B[k,j]\n",
      "C"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 38,
       "text": [
        "array([[ 38.000,  44.000,  50.000,  0.000],\n",
        "       [ 83.000,  98.000,  113.000,  0.000]])"
       ]
      }
     ],
     "prompt_number": 38
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A.dot(B)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 39,
       "text": [
        "array([[ 38,  44,  50,  56],\n",
        "       [ 83,  98, 113, 128]])"
       ]
      }
     ],
     "prompt_number": 39
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%R -i A,B\n",
      "print(A%*%B)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "text": [
        "     [,1] [,2] [,3] [,4]\n",
        "[1,]   38   44   50   56\n",
        "[2,]   83   98  113  128\n"
       ]
      }
     ],
     "prompt_number": 40
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%R -i A,B -o C\n",
      "C <- A%*%B"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 41
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "C"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 42,
       "text": [
        "array([[ 38.000,  44.000,  50.000,  56.000],\n",
        "       [ 83.000,  98.000,  113.000,  128.000]])"
       ]
      }
     ],
     "prompt_number": 42
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "**Question 2 (10 points)**. \n",
      "\n",
      "1. Read the `iris.csv` data set into a Pandas DataFrame. Dispaly the first 4 lines of the DataFrame. (2 points)\n",
      "2. Create a new DataFrame showing the mean `SepalLength`, `SepalWidth`, `PetalLength` and `PetalWidth` for the 3 different types of irises. (4 points)\n",
      "3. Make a scatter plot of `SepalLength` against `PetalLength` where each species is assigned a different color. (4 points)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "! head -n6 iris.csv"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "SepalLength,SepalWidth,PetalLength,PetalWidth,Name\r\n",
        "5.1,3.5,1.4,0.2,Iris-setosa\r\n",
        "4.9,3.0,1.4,0.2,Iris-setosa\r\n",
        "4.7,3.2,1.3,0.2,Iris-setosa\r\n",
        "4.6,3.1,1.5,0.2,Iris-setosa\r\n",
        "5.0,3.6,1.4,0.2,Iris-setosa\r\n"
       ]
      }
     ],
     "prompt_number": 8
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "df = pd.read_csv('iris.csv')\n",
      "df.head(4)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "html": [
        "<div style=\"max-height:1000px;max-width:1500px;overflow:auto;\">\n",
        "<table border=\"1\" class=\"dataframe\">\n",
        "  <thead>\n",
        "    <tr style=\"text-align: right;\">\n",
        "      <th></th>\n",
        "      <th>SepalLength</th>\n",
        "      <th>SepalWidth</th>\n",
        "      <th>PetalLength</th>\n",
        "      <th>PetalWidth</th>\n",
        "      <th>Name</th>\n",
        "    </tr>\n",
        "  </thead>\n",
        "  <tbody>\n",
        "    <tr>\n",
        "      <th>0</th>\n",
        "      <td> 5.1</td>\n",
        "      <td> 3.5</td>\n",
        "      <td> 1.4</td>\n",
        "      <td> 0.2</td>\n",
        "      <td> Iris-setosa</td>\n",
        "    </tr>\n",
        "    <tr>\n",
        "      <th>1</th>\n",
        "      <td> 4.9</td>\n",
        "      <td> 3.0</td>\n",
        "      <td> 1.4</td>\n",
        "      <td> 0.2</td>\n",
        "      <td> Iris-setosa</td>\n",
        "    </tr>\n",
        "    <tr>\n",
        "      <th>2</th>\n",
        "      <td> 4.7</td>\n",
        "      <td> 3.2</td>\n",
        "      <td> 1.3</td>\n",
        "      <td> 0.2</td>\n",
        "      <td> Iris-setosa</td>\n",
        "    </tr>\n",
        "    <tr>\n",
        "      <th>3</th>\n",
        "      <td> 4.6</td>\n",
        "      <td> 3.1</td>\n",
        "      <td> 1.5</td>\n",
        "      <td> 0.2</td>\n",
        "      <td> Iris-setosa</td>\n",
        "    </tr>\n",
        "  </tbody>\n",
        "</table>\n",
        "</div>"
       ],
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 9,
       "text": [
        "   SepalLength  SepalWidth  PetalLength  PetalWidth         Name\n",
        "0          5.1         3.5          1.4         0.2  Iris-setosa\n",
        "1          4.9         3.0          1.4         0.2  Iris-setosa\n",
        "2          4.7         3.2          1.3         0.2  Iris-setosa\n",
        "3          4.6         3.1          1.5         0.2  Iris-setosa"
       ]
      }
     ],
     "prompt_number": 9
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "df1 = df.groupby('Name').mean()\n",
      "df1"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "html": [
        "<div style=\"max-height:1000px;max-width:1500px;overflow:auto;\">\n",
        "<table border=\"1\" class=\"dataframe\">\n",
        "  <thead>\n",
        "    <tr style=\"text-align: right;\">\n",
        "      <th></th>\n",
        "      <th>SepalLength</th>\n",
        "      <th>SepalWidth</th>\n",
        "      <th>PetalLength</th>\n",
        "      <th>PetalWidth</th>\n",
        "    </tr>\n",
        "    <tr>\n",
        "      <th>Name</th>\n",
        "      <th></th>\n",
        "      <th></th>\n",
        "      <th></th>\n",
        "      <th></th>\n",
        "    </tr>\n",
        "  </thead>\n",
        "  <tbody>\n",
        "    <tr>\n",
        "      <th>Iris-setosa</th>\n",
        "      <td> 5.006</td>\n",
        "      <td> 3.418</td>\n",
        "      <td> 1.464</td>\n",
        "      <td> 0.244</td>\n",
        "    </tr>\n",
        "    <tr>\n",
        "      <th>Iris-versicolor</th>\n",
        "      <td> 5.936</td>\n",
        "      <td> 2.770</td>\n",
        "      <td> 4.260</td>\n",
        "      <td> 1.326</td>\n",
        "    </tr>\n",
        "    <tr>\n",
        "      <th>Iris-virginica</th>\n",
        "      <td> 6.588</td>\n",
        "      <td> 2.974</td>\n",
        "      <td> 5.552</td>\n",
        "      <td> 2.026</td>\n",
        "    </tr>\n",
        "  </tbody>\n",
        "</table>\n",
        "</div>"
       ],
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 10,
       "text": [
        "                 SepalLength  SepalWidth  PetalLength  PetalWidth\n",
        "Name                                                             \n",
        "Iris-setosa            5.006       3.418        1.464       0.244\n",
        "Iris-versicolor        5.936       2.770        4.260       1.326\n",
        "Iris-virginica         6.588       2.974        5.552       2.026"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "iris_map = {k: i for (i, k) in enumerate(df.Name.unique())}\n",
      "colors = [iris_map[name] for name in df.Name]\n",
      "plt.scatter(df.SepalLength, df.PetalLength, c=colors);"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXEAAAEECAYAAADeaATWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VOW9+PHPObOvWSGQhCVsAmFVZJMqICpqBREaWm6X\nq62tCm1ta2+xvbfa66+LeLV6W6393dtW+2ttRRQF9w0rYt0QBMO+gwFC9klmn3N+fwxMMiRkJmEy\nmYHv+/Xy9crJWeabh/E7Z57zPN9H0XVdRwghRFZSezsAIYQQ3SdJXAghspgkcSGEyGKSxIUQIotJ\nEhdCiCwmSVwIIbKYMdEBq1evZv369SiKwsCBA7ntttswmUzpiE0IIUQCnd6JV1dX88Ybb3Dvvfdy\n//33o2kaGzZs6PSClZWVKQ0wVTIxLokpORJT8jIxLokpOd2NqdMkbrfbMRgMBAIBIpEIgUCA/Pz8\nHgmkp2ViXBJTciSm5GViXBJTcrobU6fdKU6nk+uuu47bbrsNs9nM+PHjGTduXLdeSAghROp1eid+\n7NgxXnjhBR5++GF+//vf4/f7Wb9+fbpiE0IIkYDSWe2Ud999ly1btnDLLbcA8Pbbb7Nr1y6+8Y1v\nxI6prKyM+xpQUVHRg+EKIcS5a+XKlbGfy8vLKS8vT3hOp90pxcXFPP300wSDQUwmE1u2bGHYsGFx\nx3T0QlVVVV2JOy1cLhcej6e3w4gjMSVHYkpeJsYlMSWnuLi4WzfBnSbxwYMHc+mll7J8+XIURaGs\nrIw5c+Z0O0ghhBCplXCc+Pz585k/f346YhFCCNFFMmNTCCGymCRxIYTIYpLEhRAii0kSF0KILCZJ\nXAghspgkcSGEyGKSxIUQIotJEhdCiCwmSVwIIbKYJHEhhMhiksSFECKLSRIXQogsJklcCCGymCRx\nIYTIYpLEhRAii0kSF0KILCZJXAghspgkcSGEyGKSxIUQIoslXGOzqqqKBx98MLZ9/PhxFi9ezDXX\nXNOjgQkhhEgsYRIvLi5mxYoVAGiaxi233MLkyZN7PDAhhBCJdak7ZevWrRQVFVFYWNhT8QghhOiC\nLiXxDRs2MGPGjJ6KRQghRBclncTD4TAbN25k2rRpPRmPEEIAoChKb4eQFRL2iZ+yadMmhgwZgtvt\njvt9ZWUllZWVse2KigpcLlfqIkwRs9mccXFJTMmRmJKXiXF1NSZd16n68EO2PfEEuUOHMuL668kZ\nMKBXY0qXlStXxn4uLy+nvLw84TlJJ/ENGzZwySWXtPt9Ry/k8XiSvWzauFyujItLYkqOxJS8TIyr\nqzG1HDrE6gULCJ48J9DUxITvfhdd13stpnRwuVxUVFR0+bykulP8fj9bt25lypQpXX4BIYToiqDH\nE0vgACc++US6VjqRVBK3Wq384Q9/wGaz9XQ8QojznKO4mCHXXguAajIx9utfR9O0Xo4qcyXdnSKE\nEOlgzstjxq9+xfhbbsHkdOIeNqy3Q8poMu1eCJFxFEXBaDZjMBpRVElTnZHWEUJklGBtLZ/89rc8\nddVVrJ43j5qNG3s7pIwmSVwIkVEa9+5l86OPAuCvr2fzI4+gyt34GUmfuBAi/XSdpt278dXV4R40\nCFv//rFdRqsV1WhEC4cBsEmZj05JEhdCpF3tJ5/w3A03EAkEyL/gAq7+y1+wFxcDkFNezpWPPsrG\nhx7CNWiQjE5JQJK4ECKtFEXh4GuvEQkEAKjbuRPPoUOxJK4aDAy4+mqKP/c5DFYrGCVNdUZaRwiR\nVrquU9hmlrfBYsGWn9/uOIPTmc6wspYkcSFE2vWfMYOr/ud/qNm2jcFXXIF7xIjeDilrSRIXQqSd\nye1m4DXXMOjaa1NaE+UURVF65LqZSMbtCCF6TaoTrR6JcHT9et745jfZ9sc/EqitTen1M5HciQsh\nzhmNO3fywpIl6JrG3uefx5qXx5AFC3o7rB4ld+JCiE6FGhoINTT0dhhJCTQ1obcZjth08GAvRpMe\nksSFEGdU98knrL72Wp6eOzcrpr/nDB1K8fTpAFhyciibO7eXI+p50p0ihOhQsL6eV7/1LTyHDwPw\nyje+waLXX8dSUNDLkZ2ZtU8f5vzud7RUVWHJzcUxcGBvh9TjJIkLITqm6+iRSGxTC4chC2ZOWgoL\nsZxHU/WlO0UI0SFzfj5zHnkEa14elpwcrvy//xdLnz69HZY4jdyJCyHOqM/FF1Oxbh2AJPAMJUlc\niPNc8/790dol/fuT08HMyc6Sd9Pu3bRUVeEaOBBnWVnKYvIdO0bDrl1YcnLIHTUK1WxO2bXPNQmT\neEtLC48++ihHjhwB4NZbb2WETJEV4pzQvH8/z15/Pb6aGgxWKwvWrCGvTV2TzjRs387qefMIe71Y\n8/K4/rnncA0detYxBWpqeP2WWzj24YcAXP3YY5ReccVZX/dclTCJ/+lPf2LixIn84Ac/IBKJEDhZ\neUwIkf0a9+3DV1OD2eUi6PFQs2VL0km8bscOwl4vEF28oXHv3pQk8ZZjx2IJHGD7E08w8KqrpBzt\nGXSaxL1eLzt27GDZsmUAGAwG7HZ7WgIToifVf/opxz78kLzhwym88EKM5+n72llaysz77qPx4EFc\nJSXkjxyZ9LmuAQNiPyuqiqPNwg5nw5qXh61PH3wnTgBQeuml500dlO7oNIlXV1fjdrt55JFHOHjw\nIGVlZdx4441YLJZ0xSdEyjXu3MnqefNi9ayve/JJ+s2Y0ctR9Q4tHObtO++MraIz/+mnkz43f+xY\n5j/zDDVbt1J00UXkjh6dkpjsJSXMX7WKI+vX4+zfn6IpUySJd6LTJB6JRNi/fz833XQTw4YN47HH\nHuPZZ59l8eLFsWMqKyuprKyMbVdUVOByuXou4m4ym80ZF5fElJxUx1RdUxNL4AD1O3YwbO5cFEXp\ntZggWgxKi0RQDQaALsXTWVynEuCpa7e97tETJ2IJHKDpwAGGzJmT3Gu7XORcfjlDLr+803i701au\niRPpP2FCp9c9G5n4PgdYuXJl7Ofy8nLKk+ja6jSJFxQUkJ+fz7BhwwCYOnUqzz77bNwxHb2Qx+NJ\nOuh0cblcGReXxJScVMfkGDAAe1ER3uPHUU0mii6+mObm5l6NCaDm44/56L/+C/egQUxYtgx7SUmX\nr9FRXC2HD7PpoYdoPnaMyf/2b+SPGxfb5xw0CFtBAb7aWgwWCwVjxnS5LboTU2/L1JgqKiq6fF6n\nSTw3N5fCwkKqqqooLi5my5YtlJaWdjtIITKBY+BA5q9eTdP+/TiKisjpQj9wT2k5dIi1ixfHHhTq\nus6Me+89+24ETeP9X/yCvWvWAHD8o4+oWLcutjCxq6yM+c89h+fgQRz9+5NzwQVn93oi7RKOTrnx\nxhv5zW9+QzgcpqioiNtuuy0dcQnRo5yDBuEcNKhb54ZbWmioq0O3Wrv8QDTc3Ey4pQVTTk50/chT\nv/f7Ywkcot0auqbBWXYlaOEwnkOHYttBj4eI3x93jL1fP0x2O0ZZDi0rJZx2P3jwYH75y19y3333\ncccdd8joFHFe81dX887y5Tw+aRIbli/Hf3IERTK8R47w6te/zhPTpvHx/fcTamyM7XOUlDDx5Cgw\no83GpDvuOOsEDqCazUz9yU8wnJwsM+XOO7H16xfbH6yr44Nf/IInpk3jzaVL8R09etavKdJLaqcI\n0QU1mzez+5lniAQC7Hr6aWo2b0763MNvvcVn77xDJBBg8yOPUL99e2yf0eFgwne+Q8Wbb1Lx5pv0\nmTQpZTEXTZtGxVtvUbFuHeU33YTBZovtq62s5NM//pFIIMDB117j6Hvvpex1RXrItHshukKNv+9R\n1OTvg8xOJ9c89hiRYJCw39/uWkaHo9M+6eYDBwg1N+MoLcWcm5t8zIpyxq6j00d+dOXvEZlB/sWE\n6ILCCRMo/+pXseTmUv61r1E4fnzS5+YMHsyb3/8+r3zzm2z905+w5ecnfW7t5s08dcUVrLrqKt79\n938nWF/fnfDbyR8zhonLlmHJzWX4ggX0nzIlJdcV6SN34kJ0gbWwkCl33820O+9Es1gwdGHi26G3\n3sJfVwdA9aZN1O/ejevk8N3OKIrCp489FnvwuXv1asbefDMFeXnd+yPaMOfmcuEddzD25psxOhxx\nXS0iO0gSF1nLe+QIdTt2YM3PJ6+8vEsJtf7TT2nYswdn//4UTJwYVyXPf+wYJz75BC0UonDcuHar\nwxgsFlyFhV0eZ+xuM00dRcHahSSc36bonMFiwZzCkSSqyXReLaJwrpEkLrKS7+hRXvjyl2nYvRuA\na/78Z0pOzhxMpGHbNl5YsgRfbS2KwcC8v/+dvifXZdT8frb+z/+w+dFHARh67bXM+NWvMHeh6+NM\ncocN48Jvf5va7dsZMHNm0rVGdF1n2MKFaKEQNZ9+yribb05JoSlxbpA+cZGVvMePxxI4wJ61a5Oe\nnt20fz++2loA9EiEqg8+iO0L1NWxZ80aLDk5WPPz2fviiwQ7WOm9O5Nwqv75T3Y98wyBxkY+WLGC\npjbjtxOxFRUx4fbbmfO//0ufyZO7/Nri3CV34iIr2QoL4yrdDbjssqQTq7OkBKPdHutj7ttmGrol\nL48pd95J/a5daKEQ+SNHYmpTY0MLBjnx4Yds/uAD+k+ZQp9Jk5JesKBw7FiaP/uM5s8+w2i3t7sT\nb9q9mwMvvYSjf39KZs7EetpiDFIESnREkrjISvbSUuY//TTVmzZhLyqi8GSxpGTkT5jAdX//OzWV\nleQMHEifiy5q3amqHHjlFfY+/zwAfcaOZcCVV8Z2123dypo29S1uWLuWggsvTOp1+1x0EfOfeQbP\n4cMUjB6Nu81DTd/Ro6ytqMBbXQ3ApDvu4MLvf18St0hIkrjIWq6hQzvvG9Y09EgExWRqt6vwooso\nbJu8Twp7vRz/+OPYdk1lJZGWFsjJAaILFrTVfPQoBUnGa7Ba6TtlCn07GMYXbGqKJXCAqg0bmPjd\n78aPJdd1tGAQVUpBizakT1yck1oOH+bNpUtZs2AB1V2YhWhyuxl/662x7XE334y5zSiS/AsuwHJy\noo01L4+8FBWMshUVMfiqq6IbisK4b3wjLoF7P/uM9T/8IWsWLOCzN95Aj0RS8roi+yl6D3xfq6qq\nSvUlz1qmlp6UmBLrakyKovDWd77DrlWrgOgd8OJ169oNFTyTsNdLw44daKEQeSNHYjp5F35K8/79\n+I4fx9avH87Bg5OOK5FATQ0Nu3djdrnIGTEi1teuKAob77uPjb/+dXTbYKDi9ddxd7DW7bnw75cO\nmRhTcXFxt86T7hRxztEjkbhuj4jfH7cIRCJGux1bfj66pmFwONrtd5aV0X/cuJQnAUthIUVnGK/d\n9u/RI5HotH0hkO4UcS5SVabceSfmk6NKpv7kJzi6UAf/6D/+waq5c1k5Zw6HXnwxbuWb3qDrOuO+\n/nWsJ8eqj//Wt3CXlfVqTCJzyJ24OCcVTJhAxZtvEgkEsPXvH1e7uzP+6mrW/+QnBE/eZb/5ve9R\nMWYMziFDejLchHJGjeILr71GyOvFVlSEsYNvCOL8JElcnLNs3ehjVA0GjG0SvtFqRTm55uUp3iNH\nqK+txVJQgL2LK101HziAv64OZ0kJ1qKiLp1r7deP5D6KxPlEkrgQbZgLCph53338Y/lyIoEAn/vF\nL3C0KePavH8/axYtouXYMRzFxcxbuRJnkl0bjTt28Oz11xP0eCgcM4arHnsMe5JT74U4E0niQpwm\nf+JEPv/kk+jhMObTHjTWVlbGHjK2VFVRt3170kn88D/+Eeumqfn0U5r27pUkLs6aJHHRo8ItLdRu\n3ozn8GHyhg+n4LQJNr6qKk5s2YLJbqdw/Pi44XwRNUKtrYbDyiFybXm4fK7TL99tTbt3U79rF86S\nEvJGj243dd50hkUXbKdNhbcWJDvVJ1pP/BRFVbF0oYphxOejdssWfLW1FJaXx307EOe3pJL40qVL\nsdlsqKqKwWDgl7/8ZU/HJc4RtR9/zNolS9A1DZPTyXVPPBFL5IG6Ol679VaOf/QRANPvuosx3/pW\nbKr5UXsVz7vWgAIOzcEN2hdwBM7+gZ5n3z5Wz5tHsKkJFIXrn3026eXQzDk5TL/rLqo3b6bvhAlY\nu7DCTtHUqcx+6CGOffghw+bNI3fkyKTPPfLmm7z6zW8C4BowgHnPPIO9m+OKxbkl6Tvxu+++G6es\nhi3O4FQFwdPnjtVs2xZdtR0INTfTdOhQLIkH6+tjCRyiix2M+trXUC0WFEXhoPkAVt2KTbPRqDbi\nNbTg4OyTePORI0QCAfKGD6fl2DGOb9yYdBJv3LuXjx54AFdpKYfefBNnSQnuJGdtmnNyGLpoEcMr\nKtBOtkkyFEVh/0svxbY9hw/jq66WJC6ALowTl0I84kxaDh5k44oVbLr/fryHD8ftKxg9OrZqu9Fu\nx9Vm1qQlL4++bQpXDZs/PzYUUNd1hh0qZsJvod8tnzD7tRE4wqkZVucaOJAZ99xDn3HjuPDb36a4\nC0uSuQcNIuzzUbt9O2G/H3cXZmyGjCH2ufbyTu7bHHccR1eS+39K1/XWKflEqzDa+/ZN+nXFuS2p\naffLli3Dbrejqipz5sxhzpw5nR4v0+6Tcy7EFGpq4sUvfYnqk6u+l156KVf84Q8Y7XYAIs3NnNi8\nGc+hQ+QNH07hxRfHne/77DOqN2/G5HBQOGFC3ALAu5/8O299/wdAdOr8F159NSWLITRUVrLq6qtj\n9UeufvxxShO8p0/RNY2Gyko8hw7hGjSI3NGjk15c+LDrEC+41gKg6AoVDV8iz5dcv3i4pSWuT/xM\nD1PPhfdUOmRiTD067f6ee+4hLy+PpqYm7rnnHkpKShg1ahQAlZWVVFZWxo6tqKjA5UrdA6hUMZvN\nGRdXNsUUCYcJNDZidrkwtnkI2FBXR+2OHbHt2spKDJFI6zVcLnI7GYHhGjmSvh30Deu6Tm3lttbX\n9/sJezxxsWmahr+hAZPVitFmS3pRiOqGhrgCUp6DB3E6nUmf754+HfPMmQSDwaSOh+jf02hobN1W\ndEKmIC5jcv/+utOJY85MwkoYK1ZUpeMPjmx6T/WmTIwJYOXKlbGfy8vLKS8vT3hOUkk87+RTdLfb\nzeTJk9mzZ08siXf0Qpn2CQeZ+cmbLTGFGhvZ9thj7HzqKcquvppxt9yC5eSoDMXp5OLvf5/3fvEL\nACZ9//tgs6Xk7xp+/fVs+8tfiAQC9J8yBVtxcey6WijE4Vde4YMVK8gfOZJpP/1p0hNv7ANKyB0x\nnIZd0WJTRZMn0dzc3KXYuvNvV2IvxWQzEVJC5IcLcAZceALJXaPF2sIGx3rqjHVM9U5jUPNgFL39\nh062vKd6W6bGVNGmVn2yEibxQCCApmnYbDb8fj9btmxh0aJF3QpSZKeaLVv4YMUKADY/8gj9p0yJ\ndT+oZjOj/vVfKZ4xA0VRyBk+HMWYmpGrBRdeyKJXXyXc3Iy9uBhrm37gpr17efWWW0DXadi7l9yh\nQ5m0fHlSz248tccYctVcTAsXoSgK9Yf3kz92fEpi7kyeN48K7Yv4VT/OsAtbMLmV5RVFodL2Kfss\newF42fkiXwz/C7m+5EfGiHNXwv/bGhsbue+++4Do19cZM2YwfnzPv+FF5ji9AuDpFfSMDgcF3XxP\nRNQIPrMPg65iC9jb7bfk5GC22drVCtGCQWiTsAMdrIPpO3YMLRjEVlQUt5BC2O/j49/8JrY97Vf/\n2a3Yu8Pld+PC3eG+aFt4MeiGdm0RUNq0uQIaUk9cRCVM4n379o0lcXF+Khw3jkFXXsnBV1+l/5Qp\nFCW5HFkimqqxy7mTfzjWYcLE/KYbKGxpnSHZtGcPz3/pS7RUVTHmppuYdMcdsclA7qFDmfSDH/DR\n/ffjLClhzE03xd2F123ZwtrFiwk2NTH9rrsY+dWvxka+FIwuZ9gXF7Hn76vILx/FwNmzU/L3nI2w\nGma7q5IN9new6BbmNy0g3xvtstJ1nfG+CRw0HaBZbWaybyruUE6CK4rzhSwK0YuyKaZQUxPBxkbM\nbne7RRK6q9nq4S95f4aTXbv9wv2ZX7sAJaKgKArv/+d/8smjj8aOX/D88xROnBjbjvj9BGpqMFit\nWNpOj9d1XvrSlziyfn3sV4v/8Y+4NS0D3kYCjfWYHE5s7o5reHcm1f92TbZGnsj9S6wtBgYHcU3d\n56HNcPKA2U9YjWANWjFohg6vk03vqd6UiTHJohCiR5ncbkzujrsBuktFxYiRMNF63XbNHntYp+s6\ntraJWVEwnra2pC8vTH0fDZuuY/JqqFp0xIaiqnHT4w0WC+pp62xG8kx4C63YdRM2X0r/rE412Rrx\nKj7cETf2Nl0mqm7AgIHIyW4Sux5tC53WeyxL0MqZVtf0mX00GRtpphmn4gSZ1nHekCQueo3d72Be\n0/Wsd7yNQ3MwvWVG3J3nsPnzady7l+pPPuGi734X9/DhsX0ttmZed77KMdNRVF3lWuU6SjzR0Sm6\nrnPxHXcQ9nrxHDnC9Lvuwtmm1kiLpYW1Oc/SYGzAqJtYqHyBPG/ydUy6q95ex9M5TxFWwuSGc7mu\n8fpYGQFnwMnnPfN41/4O7kgOF7dMTnqCnc/s46Xc56k2VqPqKgu0hfRpkclA5wtJ4qJX9W0pYoF/\nIYquoGjxQ+ZsxcXMWLECExA6bUJNk6GJY6ajoIOmaOwwb2eAOjA2nd0xaBCX//736OEw6mkLQjQa\nG2gwRh+EhpUQn5kOk0fXknh3eiGPmI4QVqLfOhqMDTQaG1prwejQv7mYBb6FqLoh7sMsEY+xiWpj\nNRBtiz3m3fT1Fsks6/OEJHHR69RIJzMeVRWry0XotP5Le8TOZd6ZBJQgZt2MVbO0q0eiGI0dDne0\naXZUXUVTosfnRfKTjlVXdOpsdew37MNtzyHflx83Xttr8VJtPo6qq/QNFmENtn6A5EVaPyhUXcWm\ntR+No0Y67uvujFW3YdSNsQ+IPpG+ksDPI5LERVYyaiY+sn5Ii9oCwBXNVyU4o1WuP5cFTQs5YD5A\n33BfivzJr7BTb6tjVc6T6IqOoit8gcWxUSQhY4h1rtc5bI7WjxnnH8/Uhumxvvoifz+uVq6l2ljN\n4OBgcv2pGeft9rm5QVnEPvNeCvU+FHulMNb5RJJ4D1IUJePuiM4mps7OPVMVw7OlKEqHFf98Bl8s\ngQMcMh1kmDI87vXPFJOiK/Rp6dtpl8OZzm0yNMUKV+mKTpOhiXxOJnFDiMOm1gJg+8z7mGSYjFmL\nlikwhU0M8gxmsFKW8nbK9xZQ4CvE4XDQHO7a7FOR3SSJ94BwczOHXn2Vg6+/ztDrrqNk1qykF+rt\nKZqicdxxjG2WSvqHihnqG4YldKaxDvEiaoQjjsPsMe9mcLCMQb5BGMOtoz08Ng+f2DahojLOOwGn\nPzUli5utHvbY9nDCUM1odQwlLSWxvmJ7xE5OOJdGYwPoMCwwLC4xei1edti20WhoYoJvAnne9l0m\nZ0qkLdZmtti2EFZCjPdNxO1rHZWTE8nFoBuIKBEMupGccOvdtDliYlhwOHssuwEY5R+FKWxqd/2e\n+mDXdT3p+i/i3CFJvAec2LSJN779bQD2PPcc1z/3XNL1qntKo7WRNa5n0RWd3ZZdWHUrQ0LJVQSs\ntdXykvMFUGC3ZRfz9QX0b45+ZQ+ZQrzseoFaYy0AJwwnuDZ8Hcbw2b+1Dlj3857tXQD2m/axUPsC\nBScnA9kDduY1zafOWIdVt5Lva11hR1EUNts/Zovtk+h1TPuoiHwpqQUlIoYIb7nWcdh0CIAjxsPc\nEPkClmD0Ay/Pm8cX9MU0m5pxhlxxU9+NYRMzPJcyMjAKVTdQGCzssL6JEKkkSbwH+Gtr47YD9fW9\nFEmroBKIq1/dYGg/TT1kDNGst/8q7ld8sUkoAD6ldWB1RAnT0KY6X4OhgYgaxniWby1VVakztLab\npmgElPjp/w6/Ewft7/p1RafGWBPbDqgBwmqo3TFBUxBjxIihzcPEiBqhoc3rNhmaiCjxU9xzfXkM\nMA7E42s/WcQatFKkF4GudHgXHjFECBvCmENmSfAiJZJeFEIkr3DcOBz9+gGQO3w4eV1Yhqun5IRz\nKA5F754tmoUhwSFx+xvtDTyb/zR/cT7OUWf8jNuCcAG5J7sNHBEnfcKtY5DNIQuXeGdEJ5foMN17\nCeYku2k6o2kaIwMjMenRRFgaGkBOJMmZohpc7JuMqkff3mP947CHWu/CQ8YQW9ybeTL/Cd7KfQOv\nxRvbZwqbmNZySWyyzHTvJUl3OwEcdxznqfwneaZgFfX2+A9vr8XLutw3eDL/Cba4PyFsDJ3hKkIk\nT6bd9xDvZ5/hq6nBXlSE7WRC7+2YAmY/zYYWLLolrt9aM2g8l7+a46boKu4G3ci/1H0lbkahz+zF\na/Bi0+xxv4fo3WuTpQkFcAdyYqMxUqHBXk/AEMAVcmH3J7+yj67oeKweQoRwhVyYw6010Kvtx3km\nd1Vs+/LmKxjeNCK2rSkaTdYmdDTcwZy4O/VTOvq385t9PFHwF4JKtM54QaSQBbULY11Lu9w7edP5\neuz4GxoW0deb/MiYZGTidHKJKTky7T7D2EtKsJeU9HYYcaLTtts/YNXR0dp0GegdzDSxBe3YaD+u\nGcCgGZJeoaarcr150f/h/F37H07RFQy6AV3RMejxHyqnxofHtk+rCKhw8lzUdhOQOqMrOlqbttOI\nxE2bj5z2OskuzyZEZ6Q7RWCIGJjVfDl2zY5RNzG3+Zqka11nqhOOE/wt76/8Le8vfOLaTNgQju3L\nDxYw3jcBVVcpDpVQGhwYd+5njiP8Le8vPJH3/9jj2p10srUF7VztuRaTbsKm2bi8+Yq4fvGBwUEU\nh4pRdZXxvonkBZOfZCTEmciduACi44wXh7+E0WxEbTFk90M3Fd63v0tYifY5f2B/n6GB4eT4on3q\n5pCZyU1TmeC9EGPEGJdog6Ygbzveit2tv+V4k9LAgHZdSB3SoaSllH8JfgUFBUsw/luPw+/g6vC1\nhA0RzGFzh900QnSVJHERYwlacVlcePT2XRceaxMegwen5sTti3/A6LP4qDVHR4MUBvrETTVPJGD2\nU2eqw6RUByqGAAAgAElEQVSbyPPnn7HEalcouoJLax3bbdRNGPX46xoiBmyR9t82VF3FoTlpMjQB\nYNNtGPQuxKSDtZNvMaawGVP4jLuF6DJJ4iKhJlsTz+Q8hV/1Y9JNLKQiNj46ZAzxsWMjW63RMdkX\nGSdxYXhSUsk4ZAzyD/db7DPvBR3mGK5kePOIs54Mo+s6k1omA9CoNjLNOx1HkhOQjGEjs5sv55+O\ndwkR4hLvjNgYcSEykSRxkVC9oQ6/Gl0eLKSEqDPWkks0iQdMAbZZPj35IFCn0lLJaMsYHL7EI0l8\nRn80gQMosMX2CUO8QzsviJUkh9/BzMDsaH92FyoCArh8bq4MzI1udPFcIdJNkrhIyKk5owsUKDro\n0e1TTJqJz3kvw6f6UFGwara44XyaGuGE7QR1hjqKwkXRmZUnb7TNmpmccA6NxuhkodLQgJQOT9R1\nvfuLI0jyFlkiqSSuaRrLly8nPz+f5cuX93RMIsOouspM72waDPUna4e0vm0MEZV95r0cMh0EYHhg\nOMO01sUbTthOsNr9NCjR61ToXyT35HBEa9DKdZ75HDQfxKbZKA6UyIo0QnRRUkn8xRdfpLS0FJ8v\njetYibRTVZWw1v6pW5OhiXX2N7DqNvyKjyub51JwqnKfGuYz45HYsYdNRwgbwrGRF/WGetxaDjma\nm1pDLS1qC7ltFmBw+lyU+8acMaYzVTGM7VcV0HuuqJQQmS7hd9fa2lo2bdrE7Nmz5X+Uc1ijvYHX\n815ltW0VNc4TcfvyInlYdAt+1YcRY9wiCuaImTGBsbHtMf6xmEKtQ/aKwv0YFhyODkz0X0hOJPka\n2j6zj49yPmCl7e9UOT+LG6+tKzrHnEd5Pn8NH+S8h88sNxji/JTwTvzxxx/ny1/+styFn8M0o8YG\n+zscMke7RI4bj7MoUoHDF+37dvtyWEgFHrWp3RBDQ8TARZ5JDA6UoaBQECiM1SwBaDZ4+Nj2ERCt\nCNgv2B9nB0WrOrLftpcP7R9EzzUd5ouRJeScHBXjsTaxxvUsmqJx2HwIs25hQmii3GiI806nSXzj\nxo243W7KysqorKzs8JjKysq4fRUVFbhcrtRGmQJmsznj4upOTGE9TAvNGDFhx56S+tENWgPNhtax\n4T7Fh2bQ4mIza2YcigOrasXqOn0cuIsCTq5M32aItK7r+NU2H/4KhI3hpP5mXddpahOTruhEjJHY\nuQ3Ux02fbzQ24HA40lZPOxPfT5CZcUlMyVu5cmXs5/LycsrLyxOe02kBrCeeeIL169ejqiqhUAif\nz8eUKVNYtmxZpxeVAljJ6WpMmqqxx7mL9Y63cUfczPVcg6vNggVn46D7AK84XkJDY6p/GmM8Y2ML\nPwRMAT5wvcdOy07KgmVM98xIelp+o72Bta7naDY0UxQq4ormuTh9yd2JN9jrWe1+moAaYHhgBDOa\nLo1VFAyagrzreocd1u2YdTPXNy4kv4OFH3pKJr6fIDPjkpiS090CWElXMdy2bRtr1qxJanSKJPHk\ndDWmRlsDf8v9a6y298jAKGbVXZ6SLgRN1Wi0NaApGm5/TtxU9M+cR1jrfi62PddzDYM9ZUldd59z\nL3XmWky6mRa1mTLfEPp3YQ3IFksLmlnD4rdgDpnj9oWMIVpMzZh0M44uVDhMhUx8P0FmxiUxJSct\nVQxl6afMopDCfw89Wo3QZDShJqib0pXX1RWdj6wfxraH+JNbTegUR8CBy+zCE2r/P5wpbCI33DPV\nE4XIFkkn8dGjRzN69OiejEUk4PK7mdNyJevtb5MTyeFC76SUPcj7zPkZLzrXois609UZlHvGxKbO\nFwb6MNY/nl3mHZQFyygKdFwfvSPFgWKGmodxxHSY0YEx5AcLEp8khEiazNjMIqquMsQzlAH+Aaia\nocPlv7ojbAzzT/s7sSF879rfYWhgWKyLwhKyMLVhGhcZJ2GKmLpUfc8WtDOr4fLokmRhc0pnZAoh\nJIlnHVVX25U4TcU1h/lHYFbMRJQIIULtClgZNAOGYPcqDBojRowReasJ0RPk/yyBoikYFQPv2N8G\nYEhgaEpKwgohep58txWEDWG2WbbFtvdZ9hI0Bjo5QwiRKSSJC4wRI4ODrUMGi8JFmCKtw/l0Rafa\ncZxNuR9z1FmFpkY6uowQohdId4pA0RUmtEykf7g/ETVCUaBf3JjsOlstz7hXRcen63CDvoi+Lald\npV0I0T2SxHuJruodVgzsLdaglUHBwTidTpr9zXH7WtSW2AQjFGhWm+nLaUn81Hc6qcMtRFpJEu8F\ntY4a3rGvx6Zbma7NwOnLnBoOHU3oyovkYdNs+FQfZt1Mfjh+rHejvYH1jrdRdZUZLZ9rtwanEKLn\nSBJPM5/ZxxrXswTU6IPDEGGuCX4eJZK5s2FdPjcL9Qo8Bg8OzYG7Tb2WoCnIC661sYWFX1VbmB+6\nIWVj2IUQnZMknmaaGiGgtI788Bq8aIqGgcwe0uf0OzssIaspEXxtKhV6VS+aKn0qQqSLjE5JM1vA\nzsyWWdFaJbqBz7VciiGc2Qm8M5aQldnNc1B0BVVXmdVyeazSoBCi58mdeJqpusrw5gsoCZZiMVow\nt2R3wlN0hcHNZfxL6CsAOAJOWSdTiDSSJN5DWiwt+Axe7BEH9oA9bp+magSVIIqiYMKc2mqEvUDR\nFZz+zHk4K8T5RJJ4D/DYmlid8zRe1UtOOJd5TfNx+KP9yRFDmM3Oj9lo/wh0mEvytbmFEOJ00ife\nA04Yq/GqXiC6bFidsS62z2fys/HkmpMo8IH9fcLGzBkvLoTILpLEe4BNb7PKjA42rXUpM5NmxKG1\n7i8IF0ixKSFEt0l3Sg8o9BVylXI1h8wHGRIYSr6/de1HS9DK/KYFbLVtwabbucA7EkXL7j5xIUTv\nkSTeA0wRE2WeIQxRhna48o7bl8MM/6U4HA6aA80dXEEIIZIj3Sk9qLOl03RdlzVLhRBnLeGdeDAY\n5O677yYUChEOh7n44otZsmRJOmITQgiRQMIkbjabueuuu7BYLEQiEX7605+yY8cORo4cmY74zku6\nqhNRI9ElzWTijBCiE0n1iVss0VmF4XAYTdNwOtvX0BCp4TN7+cD1PseMx7jQdxFDmmWpNCHEmSWV\nxDVN40c/+hHHjx/nyiuvpLS0tKfjOm/tt+5n+8ml0t5wvEZ+OJ8Cb2EvRyWEyFRJJXFVVbnvvvvw\ner38/Oc/p7KykvLycgAqKyuprKyMHVtRUYHLlXlTsM1mc8bFdXpMuq4TMoRaD1AAA2mNOxvaKRNk\nYkyQmXFJTMlbuXJl7Ofy8vJYnu2Monc2hKIDq1atwmw2M2/evDMeU1VV1ZVLpoXL5cLj8fR2GHE6\niqnZ5mGtaw2NhgbGBMYy2TM1bqm03oipt0lMycvEuCSm5BQXF3frvIR34k1NTRgMBhwOB8FgkK1b\nt7Jo0aJuvZhIzOlzcUNkEWE1hCVsxRiWofxCiDNLmCEaGhp4+OGH0TQNXde59NJLGTt2bDpiO29Z\nghYsZHeJWiFEeiRM4gMHDuTee+9NRyxCCCG6SGZsCiFEFpMkLoQQWUySuBBCZDFJ4kIIkcUkiQsh\nRBaTJC6EEFlMkrgQQmQxSeJCCJHFJIkLIUQWkyQuhBBZTJK4EEJkMUniQgiRxSSJCyFEFpMkLoQQ\nWUySuBBCZDFJ4kIIkcUkiQshRBaTJC6EEFlMkrgQQmSxhGts1tTU8PDDD9PY2IiiKFx++eVcc801\n6YhNCCFEAgmTuNFo5Gtf+xqDBw/G7/fzox/9iHHjxlFaWpqO+IQQQnQiYXdKbm4ugwcPBsBqtVJS\nUkJ9fX1PxyWEECIJXeoTr66u5sCBAwwfPryn4hFCCNEFCbtTTvH7/TzwwAP867/+K1arNfb7yspK\nKisrY9sVFRW4XK7URpkCZrM5pXH5/WF27apF12HEiHxsNlNKYtJ1nb1762lo8FNa6qKoyImiKKkK\nu1sx9TaJKXmZGJfElLyVK1fGfi4vL6e8vDzhOYqu63qig8LhMPfeey8TJkzg2muvTXjRqqqqhMek\nm8vlwuPxpORakYjO2rUHWLbsdXQdHnhgFgsXDsFo7Npgn45i2rmzkYUL11Bf7+fyywdy//0z6dPH\nkpK4uxtTb5OYkpeJcUlMySkuLu7WeQmzjq7rPProo5SUlCSVwM8HDQ1B/uM/3uHUx99//Mc71NUF\nU3Lt55/fR329H4A33jjEvn2NKbmuEOLclLA7ZefOnaxfv56BAwfyb//2bwAsWbKECRMm9Hhwmcps\nNjBggIu6umiyLSlxYbEYUnLtAQNav+KpqoLTaU7JdYUQ56aESXzkyJE8+eST6Ygla7hcRh5+eA6/\n/vVGIhGN731vEjk5Xe8T78js2QO4884pvP/+UW68cSwXXJCTkusKIc5NSfWJd9W53id+yqkHjt1t\nwjPFdLbXPRuZ2FcoMSUvE+OSmJLT3T7xpEeniPZ6Isl6PCHef/84u3bVc+mlpYwZk5eya3/wwQne\nfvsII0fmM3VqPwoLU/PA9MgRL2++eQhVVZk9u5TiYntKriuESEySeIZ5552jfOMbrwDwwANGXn55\nIcOGuc/6ups31/GlL63F748A8PDDc7j++rKzvq7HE+J731vHu+9Gv31dfXUZ//3fs7DbU/OMQAjR\nOSmA1QlFUQiH9TOO045Eov91RyAAPl+g3ett3XqC0aMLuOKKwbhcZmpqfN17gdNUV3tjCRxg69aa\nlFzX643w0UfHY9vvv38UrzeckmsLIRKTJH4GXm+Ev/1tNwsWrOXBBzdTVxefcLdureerX32ZJUte\nYvPmui5d+5//rGbJkhf46ldf5NNPG2K/13WduXOHMHFiEXV1PpYuvTButMrZKCtzM3hw9CGpyaQy\ne/bAlFw3J8fEzTePi21/85vjUvaQVwiRmHSnnMGnn9bxgx+sA+Djj48zcmQ+V101AIATJwLcfvub\n7NgRTd5Ll3pYteo6+ve3Jbzu7t0ebrrpJZqaouPKfb4Q//u/V+JwRLsftm49wV//ui32ulOn9qOk\n5Oz7mIcPd/OnP81l794GioocTJhQcNbXBLBaDSxdOp7LLx+IqiqMGpWLyST3BkKky3mfxGtqAoTD\nOn36WDAYWrtNWlpCccc1NPhjPwcCEU6caO3mqK31EQhoccdXVfloaQnjdpsoKmotUxAMRmIJHODE\nCS9ebySWxKurvbF9ut4+Dq83zLFjfqxWA8XF7T806uuD+P0RCgst7ZJpUZENo1HF4TCipjDP5uSY\nGDzYjaKA0yl34UKk03l9y/Tpp/VcccUqpk17guefP0A43DraZPTofGbMKAFg6NAcpk5tHf5TXGzn\nZz+7BINBQVUVfvazSygpaU2oO3c28b3vrWPWrL/z7/++gb17m2P7Sksd/OQnUwGwWAz86EdT6NOn\ndULP7NmDKCqK3nnPnVsW+xnA4wnz2GM7mDXrSW644Tm2bo2vJrlrVxPXXPMMU6c+wV//uiuuD7y6\n2s+KFR9x2WV/5ytfeZGdO5vOqu3aevvto1x66d+57LInef/96pRdVwiR2HmbxMNhnX//93eorvYS\nDEZYtuwNPvusJba/qMjKo4/O4a23FvP00/MZNMgR26eqcN11g3j55UW8/PIirr9+SNxd76ZN1bzz\nzmfoOrz44j62bj0R2+fxhIhENJ5+ej5PPjmPd989gs/XmmzfeOMgixeP4s47p2CzGdmzp3Xa/a5d\nDfz85/8kHNY4fNjDww9vQj15S60oCg8+uJFDhzyEwxo/+cl6Dhxo/fCorKzjscc+RdN0KitreeGF\nvSlpxxMnAtxyy2s0N4doagry7W+/kbISBEKIxM6J7pR9+zx4PEEGDHCRn5/8NPXp04tZunQiwWCE\ngwebYgnxlLw8M3l5HV/P79dobg6h69HuFYul9dzTC2G1TfCqGh0ZoqoqqqpQVxdAVVu7cQwGlfvv\n/zC2vWjRBW3OVVEUYjVbrNb4fz673djmWCWue8hojB9hY7HEnxsK6ezd24SmNTFokB2HI7m3hqKA\n2dz695nNhrjXFUL0LMPdd999d6ovms6ZUB9/XMN1163m8ce3ceBAEzNmlGKztR+jbLFYCAZb7xCj\nCTTEN7/5KmvW7GXw4Bwuu6wkLhGeSTCo8de/7uSWW17jySd3kJNjZezYgliyVlUDPl+Yujo/ixeP\n5PLLB5GfH51Yo2k6O3c28Otff8T69Ue47baJjBmTH7u2y2Xh2LEWfL4w3/nORUye3I/c3OgHidls\nYODAHHbsqGX06AK+852L6Nevtb99xIgCduyoQ9N0/uu/ZjJpUt/YB4TbbSEnx8revQ1Mn17CjTeO\nicUE8PLLh6ioWMOf/1yJ0Whg4sS+ST2gtNuNTJlSwnvvVdGnj51HHpnDwIGOhOd1xen/dpkgE2OC\nzIxLYkpOd0vjZvW0e0VR+O533+Kpp3bFfvfiizcwfnz7kRenT7NtaYmwaNFatmxp7ep45ZVFSc2Q\nPHzYy9VXPx2rNmi3G3nttQoGD44mr2ef3U91tZfRowv44IOjjB3bhyuuiC5nt317I3PmtNYMvuCC\nfFatuo78fDOKovCzn73PyJH59O/v5KWX9jNr1oDYuZWV9Xz9668wa9ZA6ur8RCIaf/zjlWha60NV\nrzdCIBAhN9fM6cPbg0Gdmho/brcJp7P1w6q5OcznP/8su3e39rG/994SBgxIPhk3NUUfwLrdqX+w\nmYlTpDMxJsjMuCSm5Jy30+6vvXYoV189hNpaP/37O3C747s/DhxoZseOevr1czJqlDtWbdBqNTB0\naG4siRcU2HA44hPQrl1N7NvXSEmJk9Gjc2PdBHa7kbIydyyJDx6cE3f3X1hoY+nS12Pba9YsiP3s\ncJgoKrJz/Hh0FMqwYbmxrgtd1xkyJJfvfS86tNFgUFi4sHUVJafTRF2dnz//OboIxx13XNxu6r/d\nbjjjbMkjR5rZtauBfv3slJfnxe60LRYDixePxGBQCIU0wmENm61rb42eSN5CiMSyOonruk5TU5Dv\nfOcNAAYNcvP//l9rzfMjR1r4whfWUFUVfWD5xBPXctll0U87gwFuv/1CiorsnDjh48Ybx1BW5oyd\nu3t3E/PmrcbjCWIwKKxZsyA2trqgwMyKFTN5/PFKIhGNm24aGzeMcMKEQn7728tZt+4Q1147lPLy\n1rv7nBwTP/vZDF5+eT8ul5k5cwZhNrfeMl92WQk//el0Nm+uZtGiEYwb19rVMmiQk6eemsfjj1cy\ndGguixYNS7p+y8GDzSxY8Bw1NT4UBZ56ah7TphXF2iIQCHPffdG++FmzBnDjjaOTuq4QondldRJX\nVZVXXz0Q2z54sImjR5sZOjSajI8e9cYSOMCbbx5i5sySWOIbNszNXXdFh/u17ZIAOHQo+rAUoiv5\nbNtWGzdBZtSoHFasmNHhuU6nkQULhnDDDUPbJdmGhiBLl77GiBF5+P0R3nuviqlTr491bwwY4OBb\n3ypHVce2uy7A+PH5/PrXlwJdK8BVVdUSm8Kv6/Duu1VMn94PXddpaYnwzDO7Y8euW3eYxsZgSu+u\nFUXplaqMQpzrsnqIoaZpzJw5ILbdv78jblx1UZGdgoLW8duXXFLaLpFomtZhsiwtdWK1RrslFCXa\nd93R63d07ikdJa2cHBMTJvRl+/Y69u9vZO7cwR0+TE103a4mxH797HFdTRdf3C92DbvdyDXXtBbD\nmjy5X7tuqe5qaQnzyiuH+c//fJ/33qsmHD7z3yWE6LqsfrAJUFsb4OOPT1BT42PcuD6Ul+fG7d+z\np4lt22rp18/FmDG5SY0+gWii3L69kd276xk40M2YMfmYTGc/dK6mJsDq1XsJBMKYTAZsNiNLloxo\nNwSwJ+za1cT27bX07+9k3Lj82IcUQF1dkI8/riYY1JkwobDD2aDdsX79Ub74xeeB6IigF19cyNix\nXSuvm4kPoTIxJsjMuCSm5JzzDzZ9vggmk9ou2RUUWGKjNzoybJibYcPcXf5HUxSF0aNzGT06N/HB\nXdDSEubuuzfEtkeMyOOGG4bGjRbpKSNGuBkxouOytvn5ZubMKcXpdNLc3NzhMd1x6FBrm2uaTk2N\nF0hdjXQhzncZ352i69G7ueuue46bb36NgwdbEp+UwQoLLdx445jY9rJlE9OSwJN1prK73XXRRUW4\nXNGumYEDXQwdmtoPRSHOdwmzxyOPPMKmTZtwu93cf//96YgpzsGDzXzlKy8SCmls316L223hv/97\nZtY+JHM4jPzwh5NYsGAYdrsl9hD2XDVyZA4vvngDJ074KC11pqQioxCiVcIkPmvWLK6++mp++9vf\npiOedkIhjVCo9WFYba0PTdPbTWTprmPHfIRCOv36WdNWQjUnx8RFF/XJyH65njBkiIshQ1JTF10I\nES9h1ho1ahQOR2qnUXfFgAEO7r57OgC5uRaWL5+SsgS+eXMtM2euZNq0v7Jq1V6CQRk5IYTILpnT\nGXsGVquBr351FFddNRiLxRA3qeZshEI6P/3phthY8DvueIvp0/szaNC53b0hhDi3nHUSr6yspLKy\nMrZdUVHR7UIuZ+JyQWFhzlldw2w2x8UVDIbJzW39QLBaDVitlpTH3pWYMoHElJxMjAkyMy6JKXkr\nV7bWVSovL6e8vDzhOWedxDt6oUzs5+2o//lnP5tOKBThxAkv/+f/fI6iIlNaY8/EPnGJKTmZGBNk\nZlwSU3JcLhcVFRVdPi/ju1N6UlmZk8cfv4pwWD9j0SghhMhkCZP4gw8+yPbt2/F4PNx6661UVFQw\na9asdMSWFmazijk1M8yFECLtEibx22+/PR1xCCGE6IaMn7EphBDizCSJCyFEFpMkLoQQWUySuBBC\nZDFJ4kIIkcUkiQshRBaTJC6EEFlMkrgQQmQxSeJCCJHFJIkLIUQWkyQuhBBZTJK4EEJkMUniQgiR\nxSSJCyFEFpMkLoQQWUySuBBCZDFJ4kIIkcUkiQshRBZLuDzb5s2beeyxx9A0jdmzZ3P99denIy4h\nhBBJ6PROXNM0/vCHP/DjH/+YBx54gA0bNnDkyJF0xSaEECKBTpP4nj176NevH3379sVoNHLJJZfw\n0UcfpSs2IYQQCXSaxOvq6igoKIht5+fnU1dX1+NBCSGESI482BRCiCzW6YPN/Px8amtrY9u1tbXk\n5+fHHVNZWUllZWVsu6KiguLi4hSHmRoul6u3Q2hHYkqOxJS8TIxLYkrOypUrYz+Xl5dTXl6e+CS9\nE+FwWF+2bJl+/PhxPRQK6XfccYd++PDhzk7Rn3zyyU7395ZMjEtiSo7ElLxMjEtiSk53Y+r0Ttxg\nMHDTTTfx85//PDbEsLS09Kw/bYQQQqRGwnHiEydOZOLEiemIRQghRBcZ7r777rtTfdG+ffum+pIp\nkYlxSUzJkZiSl4lxSUzJ6U5Miq7reg/EIoQQIg1kiKEQQmQxSeJCCJHFEj7Y7IymaSxfvpz8/HyW\nL1/ebv8f//hHNm/ejMVi4bbbbqOsrOxsXu6sY6qsrGTFihUUFRUBMGXKFBYuXNjjMS1duhSbzYaq\nqhgMBn75y1+2OybdbZUopt5oq5aWFh599NFYfZ5bb72VESNGxB2T7nZKFFO626mqqooHH3wwtn38\n+HEWL17MNddcE3dcutspmbh64z21evVq1q9fj6IoDBw4kNtuuw2TyRR3TLrbKlFMXW6nsxnXuHbt\nWv2hhx7Sf/WrX7Xbt3HjRv0Xv/iFruu6vmvXLv3HP/7x2bxUSmL69NNPO/x9T7vtttt0j8dzxv29\n0VaJYuqNtvrNb36jv/HGG7quR+cotLS0xO3vjXZKFFNvvad0XdcjkYh+88036ydOnIj7fW/9v5co\nrnS31fHjx/WlS5fqwWBQ13Vdf+CBB/R169bFHZPutkompq62U7e7U2pra9m0aROzZ89G7+DZ6Ecf\nfcRll10GwPDhw2lpaaGhoaG7L5eSmIAz/r6ndfa6vdFWiWJKZn8qeb1eduzYwezZs4HoHAW73R53\nTLrbKZmYoPfeU1u3bqWoqIjCwsK43/fW+ylRXJDetrLb7RgMBgKBAJFIhEAg0G7GebrbKpmYoGvt\n1O3ulMcff5wvf/nL+Hy+DvefXjyroKCAuro6cnNzu/uSZx2Toijs2rWLH/7wh+Tn5/OVr3wlLZOX\nFEXhnnvuQVVV5syZw5w5c+L290ZbJYop3W1VXV2N2+3mkUce4eDBg5SVlXHjjTdisVhix6S7nZKJ\nqbfeUwAbNmxgxowZ7X7fG++nZOJKd1s5nU6uu+46brvtNsxmM+PHj2fcuHFxx6S7rZKJqavt1K07\n8Y0bN+J2uykrK+v0EyOdn7rJxFRWVsbvfvc77rvvPubOnct9992XltjuueceVqxYwY9//GNeeeUV\ntm/f3u6YdN/NJYop3W0ViUTYv38/V155Jffeey9Wq5Vnn3223XHpbKdkYuqt91Q4HGbjxo1Mmzat\nw/299e2gs7jS3VbHjh3jhRde4OGHH+b3v/89fr+f9evXtzsunW2VTExdbaduJfGdO3eyceNGli5d\nykMPPURlZSW//e1v445JpnhWKiUTk81mi91FTZw4kXA4THNzc4/FdEpeXh4AbrebyZMns2fPnrj9\n6W6rZGJKd1sVFBSQn5/PsGHDAJg6dSr79++POybd7ZRMTL31ntq0aRNDhgzB7Xa329cb76dk4kp3\nW+3bt48LLrgAl8uFwWBgypQp7Ny5M+6YdLdVMjF1tZ26lcSXLFnC7373Ox5++GFuv/12ysvLWbZs\nWdwxkyZN4u233wZg165dOByOHv06l0xMDQ0NsU/dU0nL6XT2WEwAgUAg1r3j9/vZsmULAwcOjDsm\n3W2VTEzpbqvc3FwKCwupqqoCYMuWLe2+Qqa7nZKJqTfeUxDtsrjkkks63Jfudko2rnS3VXFxMbt3\n7yYYDKLreka8p5KJqavtdFZDDE9RFAWA1157DYArrriCCy+8kE2bNvHtb38bq9XKrbfemoqXOquY\n3jDqJkMAAAC4SURBVHvvPV577TVUVcVisfDd7363x+NobGyMfR3SNI0ZM2Ywfvz4Xm2rZGLqjba6\n8cYb+c1vfkM4HKaoqIhbb721199TiWLqjXby+/1s3bqVb33rW7Hf9XY7JRNXuttq8ODBXHrppSxf\nvhxFUSgrK2POnDm92lbJxNTVdpJp90IIkcVkxqYQQmQxSeJCCJHFJIkLIUQWkyQuhBBZTJK4EEJk\nMUniQgiRxSSJCyFEFpMkLoQQWez/A1MeRO7ejAEkAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x116c52590>"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.scatter(df.SepalLength, df.SepalWidth)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 34,
       "text": [
        "<matplotlib.collections.PathCollection at 0x117e262d0>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAEECAYAAADJSpQfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X9sE+f9B/C3bWYnJA6pUy18SQS4WtYWq+2K0mYTGdEo\nqRiokQosUzOqtqraKAko8AdN1j/aSIhRiAo1bUpTMST+mbSs2iamVULRNAmYxB8JRmVuIUWjTAWR\njriBxCSExP7+kcbE9Y87++7xPc7zfkmVktyd753Hx6fnu0+es0Wj0SiIiGhBs1sdgIiIxGOxJyJS\nAIs9EZECWOyJiBTAYk9EpAAWeyIiBSzSs1JbWxsKCwtht9vhcDiwb9++hHWOHTuG8+fPw+VyobW1\nFV6v1/SwRESUHd1n9l1dXThw4EDSQn/u3DkMDw/j8OHDeP3113H06FHN1wsGg5klzQFm0k/GXMyk\nDzPpJ2OubDPpLvbp/vZqYGAAdXV1AICqqiqEw2GMjo6mfb2FNIgiyZgJkDMXM+nDTPrJmCvbTLou\n49hsNuzZswd2ux3r16/H+vXr45aHQiGUlZXFvi8rK0MoFEJpaWlWoYiIyFy6iv2ePXvwwAMP4Pbt\n29izZw8qKirw6KOPxq3DWReIiORly3RunD/96U8oKCjAc889F/vZxx9/DJ/PhzVr1gAAdu7cia6u\nrrgz+2AwGPfxo7Gx0Wh2IiIl9fX1xb72+Xzw+Xya22ie2d+9exeRSASFhYWYnJzEZ599hq1bt8at\nU11djZMnT2LNmjUYGhpCUVFRwiWcZIGuX7+uGTCX3G43xsbGrI4RR8ZMgJy5mEkfZtJPxlzLli3L\n6mRZs9jfunUL3d3dAIBIJILa2lo88cQT6O/vBwDU19dj9erVCAQC2LFjBwoKCtDS0pJxECIiEifj\nyzhm4pm9NhkzAXLmYiZ9mEk/GXMtW7Ysq+34F7RERApgsSciUgCLPRGRAljsiYgUwGJPRKQAFnsi\nIgWw2BMRKYDFnohIASz2REQKYLEnIlIAiz0RkQJY7ImIFMBiT0SkABZ7IiIFsNgTESmAxZ6ISAEs\n9kRECmCxJyJSAIs9EZECWOyJiBTAYk9EpAAWeyIiBbDYExEpYJGelSKRCDo7O+HxeNDZ2Rm3LBgM\n4sCBAygvLwcA1NTUYMuWLeYnJSKirOkq9p9++ikqKysxMTGRdPmqVavQ0dFhajAiIjKP5mWckZER\nBAIBrFu3DtFoNOk6qX5ORERy0DyzP378OLZt25byrN5ms2FoaAi7d++Gx+PBiy++iMrKStODEhFR\n9tIW+8HBQZSUlMDr9SIYDCZdx+v14siRI3C5XAgEAuju7obf709YLxgMxr1GY2Mj3G63wfjmcjqd\nzKSTjLmYSR9m0k/WXH19fbGvfT4ffD6f5ja2aJprMH/4wx9w+vRp2O123Lt3DxMTE6ipqcH27dtT\nvmBbWxv279+P4uJizZ1fv35dc51ccrvdGBsbszpGHBkzAXLmYiZ9mEk/GXMtW7Ysq+3Sntk3NTWh\nqakJAPD555/jxIkTCYV+dHQUS5Ysgc1mw+XLlwFAV6EnIqLc0dWNM8dmswEA+vv7AQD19fU4e/Ys\n+vv7Ybfb4XK50N7ebn5KIiIyJO1lHNF4GUebjJmA7HKFQk4AgMczJSKSlGPFTPrImAmQM1e2l3H4\nF7SUEwMDRaiv96C+3oOBgSKr4xAph8WehAuFnGhuduPGDTtu3LCjudkdO8snotxgsSciUgCLPQnn\n8Uyht3cMS5dGsHRpBL29Y8Ku2xNRchl14xBlq7o6jP7+ewDE3aAlotRY7ClnWOSJrMPLOERECmCx\nJyJSAIs9EZECWOyJiBTAYk9EpAAWeyIiBbD1kjIWCjkxPh4FZ7Imyh88s6eMzE1oVldXzAnNiPII\niz3pxgnNiPIXiz0RkQJY7Ek3TmhGlL94g5YyMjehmdPpRHFx2Oo4RKQTz+wXqFDIKex6usczhf/7\nP5uQ1yYiMVjsFyA+ApCIvo/FfoFhxwwRJcNiT0SkABb7BYYdM0SUjK5unEgkgs7OTng8HnR2diYs\nP3bsGM6fPw+Xy4XW1lZ4vV7Tg5J+fAQgEX2frjP7Tz/9FJWVlbDZEjswzp07h+HhYRw+fBivv/46\njh49anpIypzHM5V3hV5kBxGR6jSL/cjICAKBANatW4doNJqwfGBgAHV1dQCAqqoqhMNhjI6Omp+U\nFjR2EBGJpVnsjx8/jm3btsFuT75qKBRCWVlZ7PuysjKEQiHzEtKCxw4iIvHSXrMfHBxESUkJvF4v\ngsFgyvWSnfF/XzAYjHuNxsZGuN3uDKKK53Q6mUknM3ONjyceP7Ov77Isk1mYSR8ZMwHy5urr64t9\n7fP54PP5NLdJW+wvXbqEwcFBBAIB3Lt3DxMTE/jggw+wffv22DoejwcjIyOx70dGRuDxeBJeK1mg\nsbExzYC55Ha7mUknM3MVFwO9vVE0N8/+o+rtHUNxcRiZvryMY8VM+siYCZAzl9vtRmNjY8bbpS32\nTU1NaGpqAgB8/vnnOHHiRFyhB4Dq6mqcPHkSa9aswdDQEIqKilBaWppxEFIbO4iIxMpoIrS5bpz+\n/n4AQH19PVavXo1AIIAdO3agoKAALS0t5qckJbDIE4mju9ivWrUKq1atAjBb5Od79dVXzU1Flrp2\nrRAAUFExYXESIjIL/4KW4pw540ZDwxI0NCzBmTPy3Zgiouyw2FPMtWuFaG8virVAtrcXxc7yiSi/\nsdgTESmAxZ5iKiom4PeHY5Oo+f1hXrcnWiD4WEKKU1s7hhMnpgHwBi3RQsJin8euXl0MAFix4o6p\nryuqyM9NgcAWSzKLkWNKteORl3Hy1KlTbmzeXILNm0tw6pT8XTOc6IzMZuSYUvF4ZLHPQ1evLsau\nXfe7ZnbtKoqd5cuIE52R2YwcU6oejyz2REQKYLHPQytW3MGhQ/e7Zg4dCpt+3d5MfFQimc3IMaXq\n8WiL6pmfWJDr169bteukZJ3hLlUmUTdo9chmrETfEMu3988qCymT6Bu0Mo7VsmXLstqO3Th5TOaz\n+WRUOHui3DJyTKl2PPIyTh4z8szWdNvyWbBECw+LfZ4S1XamYksakQpY7POQqLYzVVvSiFTAYk9E\npAAW+zwkqu1M1ZY0IhWwGydPGXlma7pt+SxYooWJxV4nKyZNCoWcGB+Porg4+XJRbWcs8pSMahOH\nLTS8jKODFR0qc/usqytmVwxZjl1a+Y/FXoMVHSrsiiGZ8HhcGFjsiYgUwGKvwYoOFXbFkEx4PC4M\nmjdop6am0NXVhXv37mF6ehpPPfUUmpqa4tYJBoM4cOAAysvLAQA1NTXYsmWLmMQWsKJDZW6fTqcT\nxcXhnOyTKBV2aeU/zWLvdDrx9ttvw+VyYWZmBm+99RYuXryIRx55JG69VatWoaOjQ1hQq1lxgE9M\nODA1lbob59q1QgDJHyMosnNCq0uIFiYW+fym6zKOy+UCAExPTyMSiaA4yb9yC2dKXpDOnHGjoWEJ\nNm1y48yZxMcOzi1vaFiSsFxk5wS7hIjyk65iH4lEsHv3brz22mvw+XyorKyMW26z2TA0NITdu3dj\n3759+Prrr4WEVcW1a4Vob7//2MH29qLYWbzWcpGdE+zKIMpfuv6oym63o7u7G3fu3MHevXsRDAbh\n8/liy71eL44cOQKXy4VAIIDu7m74/f641wgGgwgGg7HvGxsb4XbL9aBsp9MpRSaHI5LkZ/ZYtnTL\nx8cTP2HN/l4uw7lEvrYZZHn/5mMmfWTMBMibq6+vL/a1z+eLq8epZPykqk8++QROpxMNDQ0p12lr\na8P+/fuTXu6Zj0+qSu3MGTfa22cvk/j9YdTWjulePjBQhObm2QO0t3cM1dXm3eAV+dpGyfT+zWEm\nfWTMBMiZS9iTqm7fvg2Hw4GioiJMTU3hwoUL2Lp1a9w6o6OjWLJkCWw2Gy5fvgwAmoWe0qutHcOJ\nE9NwOOxYujSxoM4tBxJv0IrsnGCXEFF+0iz2o6Oj6OnpQSQSQTQaxdq1a/HYY4+hv78fAFBfX4+z\nZ8+iv78fdrsdLpcL7e3twoOroKJi4rszi9TLUxHZOeHxTMHtdqXMRUTy4QPH50n3kc1IK2O6Fkmt\n7VKd2WtlEp1Xto+3oZDzu08b41ltC5g/joB84wQwUyZkzJXtZRz+Ba0ORloZ07VI6tkuVeulqEcL\nZpvXSkbaQfmIRlIFz+znSfZ/8VDIifp6D27cmP3/4tKlEfT3h3Sd6V27VoiGhiVx2544cUvzDF9r\nu3SZcpVXljMeI7+vqHGcT5Zxmo+Z9JMxF8/siYgoJRZ7DUYmgaqomIDfH45t6/eHdV2319pO1KMF\ns81rJT6ikUgfXsaZhzdo7+8X4A1arWV6yDZOADNlQsZcwvrsaZaRs7psz44LC2fgdDqy2taKvFYy\n0g6abqwmJrIbf6OyPUEgSoWXcSSl1WHCThHxrOpMyseOKJIfi72EtCYc44Rk4mlNRrfQ9ksLH4s9\nEZECWOwlpNUJwk4R8azqTMrHjijKD+zGmUe2O+9aHSYin0alRbaxAsRkMnqjNNtMIm/QqvLemUHG\nXOzG0cFIcUy3rVWdEyqdzVv1PzY9f+2caj0jj29Mt1+tsdBqJc02k5UnF2ScMpdxjHSvpNtWVOcE\nH/93n6ydR3oeDWn2+6c1Fnrm+jF7DiHKD0pcxtE7z0mmc+NkO/eNWXmtlKuPt5mMRS4/cqd770W9\nf1qvK2qun1wcjzJeLgHkzMW5cYiIKCUlir2o+VNEdU6w2+Y+Wcci3XsvKrORLi1R/wYofyhxGWeO\n1g2mbOfGEXWD1sh8L6Ll+uOtnpuDVnzk1rpBK+L9M3qDVsQcQkbJeLkEkDMXu3F0MHKQpttWVBcO\nH/93n6xnklqPhhTx/mmNRbrlouYQIvkpcRnHDKGQM+spCdJte/XqYly9uthINLKQkeMiHSPHBY8p\nSobFXgdRbZunTrmxeXMJNm8uwalTnPAq34hqRzRyXPCYolRY7DUYmXQs3bZXry7Grl33J7zatauI\nZ2N5RNRkdEaOCx5TlA6LPRGRAljsNYhqWVux4g4OHbrfunfoUBgrVtwR+auQiUS1Ixo5LnhMUTpp\nWy+npqbQ1dWFe/fuYXp6Gk899RSampoS1jt27BjOnz8Pl8uF1tZWeL1eXTvPp4nQRM2rM/cxO9U/\nShlbvwA5c1mRyUg7bzpax4WRbfne6SdjLiGtl06nE2+//TZcLhdmZmbw1ltv4eLFi3jkkUdi65w7\ndw7Dw8M4fPgwvvzySxw9ehR79+7NKoxRRvqPtRg5a0v3aDu3ezrlMiOTVhkh6vm1Ip+Lm26sRE5U\nl+53uXp1MWy2CJYvT1ymNRaLFqX+8xetbY2czYt8Hu9CypSPNC/juFwuAMD09DQikQiKv/evaWBg\nAHV1dQCAqqoqhMNhjI6OCoianpEJokTSM1mW2ZNWGWFknLL9XY287vzXTjZWVj3ib64r5vnn3Qld\nMVpjIWocteg5HnP970fGTPlK8y9oI5EIOjo6MDw8jGeffRbbtm2LW/7OO+/g+eefx8MPPwwA2LNn\nD37zm9/goYce0ty5WZdxjEwQNZ/ZH9mynSzLqonQjEw6JmpiMK3J5qyYqE7L1auLsXlzSdx+//zn\n21ix4o7mWORigrVMJ/yzaiI0qzOlymU1YX9Ba7fb0d3djTt37mDv3r0IBoPw+Xxx6+iZcSEYDCIY\nDMa+b2xshNttzpnW+Hji/p1OJ9xul67l8T8z7+zP4Ygk+Zkdbrc7bSa9ec2WyX6/P1bZ/q5a0r2u\nVmatbUWx2RL3a7PZdI2FqHFM3CZ+DKw+HmXMlCqXDPr6+mJf+3y+hJqcTEZz43zyySdwOp1oaGiI\n/ezjjz+Gz+fDmjVrAAA7d+5EV1cXSktLNV/PzBu0AwNFaG6efVN6e8dQXR3OaDkg5v/iZ8640d4+\n+xHT7w+jtvb+66fLpCevCHr3m2yssv1dtaR7Xa3X1tpWlFOn3Ni1a3a/hw6FsXat/rEQNY5zUh3n\nVh6PMmZKl8tK2Z7Zpy32t2/fhsPhQFFREaamprB3715s3boVjz32WGydc+fO4eTJk/jtb3+LoaEh\nHD9+XPcNWrO7cYzeoBX1xmZ709KqidCMTDpm5Q3aVGNl1ZPEZm/Q2rB8eWIR0hoLUeMIZN91ZtVE\naFbeoFWm2P/3v/9FT08PIpEIotEo1q5di4aGBvT39wMA6uvrAQC///3vcf78eRQUFKClpUXX9Xog\nv1ovrSJjJkDOXMykDzPpJ2MuIcVetHya4tgqMmYCUuey6iws3Zm90SmBs90WyG6cRBN1TBn5nfLt\nOLcSn1SlgW1a4lnVJpeu9dKMZ7aa3c67EI/Fhfg7LTRKnNlb1XppBhkzAYm5rGqTM7JfUdtmO065\nYvYxZcbvlC/HuQx4Zk9ERCkpUez5DE3xRD3/VOR+RW2bbd58tRB/p4VIics4c3iD1jy8QatvOW/Q\n6pNvx7mV+AxaHXi2Ya10E8IBxieqS/VsVa39GpHutdNNzmbkWJR18q90eURllnUsZKTEZRzKjXQd\nGXonMzO7m8PIfs2YsMzsiezysetFVOZ8HAsrKXUZR4uMH9lkzARk1mViZDIzI5mM7FeGCcu+T+RE\naKKI6oTjRGiZ45k9EZECWOzJFOk6MioqJuD3339cnt8fjpvzRVQ3h5H9amVK99qifp987HrhWMiD\nl3HmkfEjm4yZgOy6TPRMZpZqWyOZjOzX6IRlIiayEzkRmiiiOuE4EZp+SnXjqETGLgWtWSdv3Zo9\nHLN9zm+qzhet/RoZo8LCmbSvm6pDyAiZ3lO9RGXOx7GwCi/jLEAyPkJOy6lTbmzdWoKtW0sSHuOn\nd7+57nxhNwjlE17GmUfGj2z50qVgZL9XrizG1q3xj/H75JPb8Hq1H5xtVecL51syl4yZADlzsRuH\niIhSYrFfYKzqUjCyX6/3Dg4dut/ZcuhQWNdZvdH9GnlddoNQvuFlnHlk/MiWb10KRvZ75cpiANBd\n6L+/Xys6XzjfkjlkzATImYvdOBTHqrNMI/vNpsjP368VnS88m6d8wWJPGUvX5qi1HZD9WXK2ZGxD\nJXPwvdWP1+wpI9m2OVrVxsj2yIWL721mWOxJt1DIieZmN27csOPGDTuam92xMysj22X7uqLykvz4\n3maOxZ6ISAEs9qSbqEfxcbIsyhTf28xptl7evHkTPT09uHXrFmw2G5555hls3Lgxbp1gMIgDBw6g\nvLwcAFBTU4MtW7Zo7pytl9pkzJRtm6PoG7QqPQLQiIWUiROh6afZjbNo0SK89NJLWLlyJSYnJ9HR\n0YHHH38clZWVceutWrUKHR0dWYVQmZGDVcZOBKOPFhTByOtqzZgpiozvrYw4PvppXsYpLS3FypUr\nAQAFBQWoqKjAt99+m7CehX+blbeMdBNYPdlZsm4cI4/4k5HWIw1FycexIvlldM3+m2++wVdffYWq\nqqq4n9tsNgwNDWH37t3Yt28fvv76a1NDLkRGugms6kRIt99sl8nq2rVCtLcXxTK3txfFzvJFysex\novyg+4+qJicncfDgQbz88ssoKCiIW+b1enHkyBG4XC4EAgF0d3fD7/fHrRMMBhEMBmPfNzY2wu3O\n3dmSHk6nM2eZxscTPwnN7t+lmUnvtmZLt99sl5nJzPfP4Ygk+Zk949fPNFMuxiqXx7leMmYC5M3V\n19cX+9rn88Hn82luo2tunOnpaezfvx8/+clPsGnTJs0XbWtrw/79+1Gs8SeWqt+gHRgoQnPz7IHU\n2zuG6uqw7kx6thUh3X6zXWYWs9+/M2fcaG+fvYzi94dRW5v5a2eTSfRYyXjTUcZMgJy5sr1Bq1ns\no9Eoenp6UFxcjJdffjnpOqOjo1iyZAlsNhsuX76MQ4cOoaenR3Pnqhd7wNhEWlbdxEvXjWPkEX9G\niXj/jN6glbHLRMYCJmMmQM5cwrpxLl26hNOnT2P58uV44403AAAvvPACbt68CQCor6/H2bNn0d/f\nD7vdDpfLhfb29qzCUGas6kSYmHBgair53DjpMuVj50Suu3Dm5ONYkdw4xfE8+XQZxypmXNoQRbax\nAphJLxkzAXLm4pOq8kw+dl1Y1aFCRMax2BMRKYDF3iL5OLdHRcUE/P77jw/0+8OWXdMmoszw4SUW\nqq4Oo7//HoD8uSFXWzuGEyem4XDYsXRpbto9icg4FnuL5UuRn6+iYuK7G1dWJyEivVjsv5Pto/ZE\nkjETIG8uETghGS0UvGaP7B+1J5KMmQB5c4nACcloIVG+2MvYAiljJplziaDS70pqUL7YExGpQPli\nL2MLpIyZZM4lgkq/K6mB0yV8J9tH7YkkYyZA3lwi/rRd1KMSrcRM+smYS9hEaKrweKbgdrukaieU\nMRMgZ65r1wrhcESwdGnm26Yr6Dybp4VC+cs4lP/mHh+4aZM748cHsuOGVMFiT3nNyORs7LghlbDY\nExEpgMWe8pqRydnYcUMq4Q1ayntGJmfLx8noiLLBYk8LgpHJ2VjkSQW8jENEpAAWeyIiBbDYExEp\ngMWeiEgBLPZERArQ7Ma5efMmenp6cOvWLdhsNjzzzDPYuHFjwnrHjh3D+fPn4XK50NraCq/XKyQw\nERFlTrPYL1q0CC+99BJWrlyJyclJdHR04PHHH0dlZWVsnXPnzmF4eBiHDx/Gl19+iaNHj2Lv3r1C\ng6tApcf/EZFYmpdxSktLsXLlSgBAQUEBKioq8O2338atMzAwgLq6OgBAVVUVwuEwRkdHzU+rEJUe\n/0dE4mV0zf6bb77BV199haqqqrifh0IhlJWVxb4vKytDKBQyJ6GCOEEXEZlN91/QTk5O4uDBg3j5\n5ZdRUFCQsFzrGSjBYBDBYDD2fWNjI9zuzKajFc3pdEqRaXw8cSxns7ksSJOcLGM1HzPpw0z6yZqr\nr68v9rXP54PP59PcRlexn56exrvvvouf//znePrppxOWezwejIyMxL4fGRmBx+OJWydZINmeACPL\nU2mKi4He3iiam2cPst7eMRQXh6V6WIgsYzUfM+nDTPrJmMvtdqOxsTHj7TSLfTQaxUcffYSKigps\n2rQp6TrV1dU4efIk1qxZg6GhIRQVFaG0tDTjMHTf3ARds4//y2xyLyKi79Ms9pcuXcLp06exfPly\nvPHGGwCAF154ATdv3gQA1NfXY/Xq1QgEAtixYwcKCgrQ0tIiNrUiZHz8HxHlJ81i/8gjj+CPf/yj\n5gu9+uqrpgQiIiLz8S9oiYgUwGJPRKQAFnsiIgWw2BMRKYDFnohIASz2REQKYLEnIlIAiz0RkQJY\n7ImIFMBiT0SkABZ7IiIFsNgTESmAxZ6ISAEs9kRECmCxJyJSAIs9EZECWOyJiBTAYk9EpAAWeyIi\nBbDYExEpgMWeiEgBLPZERApYpLXChx9+iEAggJKSErz77rsJy4PBIA4cOIDy8nIAQE1NDbZs2WJ+\nUiIiyppmsf/FL36BX/7yl/jggw9SrrNq1Sp0dHSYGoyIiMyjeRnn0UcfRVFRUdp1otGoaYGIiMh8\nmmf2Wmw2G4aGhrB79254PB68+OKLqKysNCMbERGZxHCx93q9OHLkCFwuFwKBALq7u+H3+83IRkRE\nJjFc7AsLC2NfP/nkkzh69CjGx8dRXFwct14wGEQwGIx939jYiGXLlhndvencbrfVERLImAmQMxcz\n6cNM+smYq6+vL/a1z+eDz+fT3MZw6+Xo6Gjsmv3ly5cBIKHQzwVqbGyM/Tc/rCyYST8ZczGTPsyk\nn4y5+vr64mqpnkIP6Dizf++99/DFF1/g9u3baGlpwa9+9SvMzMwAAOrr63H27Fn09/fDbrfD5XKh\nvb3d2G9CRESm0yz2O3fuTLt8w4YN2LBhg2mBiIjIfI6urq4uq3b+wx/+0Kpdp8RM+smYi5n0YSb9\nZMyVTSZblE3yREQLHufGISJSAIs9EZECDPfZ6xGJRNDZ2QmPx4POzs6E5ceOHcP58+fhcrnQ2toK\nr9draSYrJndra2tDYWEh7HY7HA4H9u3bl7BOrsdJK5NVk+CFw2F89NFH+PrrrwEALS0t+PGPfxy3\nTq7HSitTrsfq+vXreO+992LfDw8P49e//jU2btwYt14ux0lPJiuOqb/85S84ffo0bDYbli9fjtbW\nVvzgBz+IW8eKGqWVK+OxiubA3/72t6jf74++8847CcsGBwejv/vd76LRaDQ6NDQUffPNN3MRKW2m\nf//730l/LlJra2t0bGws5XIrxkkrkxXjFI1Go++//370H//4RzQajUanp6ej4XA4brkVY6WVyaqx\nikaj0ZmZmehrr70W/d///hf3c6v+7aXLlOtxGh4ejra1tUWnpqai0Wg0evDgweg///nPuHWsGCc9\nuTIdK+GXcUZGRhAIBLBu3bqkE6YNDAygrq4OAFBVVYVwOIzR0VFLMwHWTO6Wbp9WjJNWJj3LzXbn\nzh1cvHgR69atAwA4HA4sXrw4bp1cj5WeTIB1EwZeuHAB5eXlePDBB+N+btUxlS4TkNtxWrx4MRwO\nB+7evYuZmRncvXsXHo8nbh0rxklPLiCzsRJ+Gef48ePYtm0bJiYmki4PhUIoKyuLfV9WVoZQKITS\n0lLLMlkxuZvNZsOePXtgt9uxfv16rF+/Pm65FeOklcmKcfrmm29QUlKCDz/8EFevXoXX68Urr7wC\nl8sVWyfXY6Unk5UTBv7rX/9CbW1tws+tOKa0MuV6nIqLi/Hcc8+htbUVTqcTTzzxBB5//PG4dawY\nJz25Mh0roWf2g4ODKCkpgdfrTft/oFz+n1xPprnJ3bq7u7FhwwZ0d3cLz7Vnzx4cOHAAb775Jk6e\nPIkvvvgiYZ1cnxlqZbJinGZmZnDlyhU8++yz2L9/PwoKCvDXv/41Yb1cjpWeTFaMFQBMT09jcHAQ\nP/vZz5Iut+LTRrpMuR6nGzdu4O9//zt6enrQ29uLyclJnD59OmG9XI+TnlyZjpXQYn/p0iUMDg6i\nra0Nfr8fwWAw4SEoHo8HIyMjse9HRkaSflzJZabCwsLYWdmTTz6J6elpjI+PC8sEAA888AAAoKSk\nBE8//XTMRQF/AAABt0lEQVRsnqE5uR4nPZmsGKeysjJ4PB786Ec/AgD89Kc/xZUrV+LWyfVY6clk\nxVgBQCAQwEMPPYSSkpKEZVYcU1qZcj1O//nPf/Dwww/D7XbD4XCgpqYGly5dilvHinHSkyvTsRJa\n7JuamnDkyBH09PRg586d8Pl82L59e9w61dXVOHXqFABgaGgIRUVFQj8e6cmkd3I3s9y9ezd2SWly\nchKfffYZli9fHrdOrsdJT6ZcjxMAlJaW4sEHH8T169cBAJ999lnCR9dcj5WeTFaMFTB7uWTNmjVJ\nl+V6nPRkyvU4LVu2DF9++SWmpqYQjUalOJ705sp0rHLSejnHZrMBAPr7+wHMTqS2evVqBAIB7Nix\nAwUFBWhpacllpKSZcj25261bt2IfwSKRCGpra/HEE09YOk56Mlk1Cd4rr7yC999/H9PT0ygvL0dL\nS4vlx5RWJivGanJyEhcuXEBzc3PsZ1aPk1amXI/TypUrsXbtWnR2dsJms8Hr9WL9+vWWj5OeXJmO\nFadLICJSAP+ClohIASz2REQKYLEnIlIAiz0RkQJY7ImIFMBiT0SkABZ7IiIFsNgTESng/wEeNOVh\n+xe3GAAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x118a1a090>"
       ]
      }
     ],
     "prompt_number": 34
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "**Question 3 (10 points)**.\n",
      "\n",
      "Given the following covariance matrix\n",
      "```python\n",
      "A = np.array([[2,1],[1,4]])\n",
      "```\n",
      "\n",
      "1. Show that the eigenvectors of $A$ are orthogonal. (2 points)\n",
      "2. What is the vector representing the first principal component direction? (2 points)\n",
      "3. Find $A^{-1}$ without performing a matrix inversion. (2 points)\n",
      "4. What are the coordinates of the data points (0, 1) and (1, 1) in the standard basis expressed as coordinates of the principal components? (2 points)\n",
      "5. What is the proportion of variance explained if we keep only the projection onto the first principal component? (2 points)."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import scipy.linalg as la"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 12
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[2,1],[1,4]])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 13
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "e, v = la.eig(A)\n",
      "v[0].dot(v[1])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 14,
       "text": [
        "0.0"
       ]
      }
     ],
     "prompt_number": 14
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "v[:, e.argmax()]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 15,
       "text": [
        "array([-0.383, -0.924])"
       ]
      }
     ],
     "prompt_number": 15
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "e = np.real_if_close(e)\n",
      "Ainv = v.dot(np.diag(1/e).dot(v.T))\n",
      "print Ainv\n",
      "print Ainv.dot(A)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 0.571 -0.143]\n",
        " [-0.143  0.286]]\n",
        "[[ 1.000  0.000]\n",
        " [ 0.000  1.000]]\n"
       ]
      }
     ],
     "prompt_number": 16
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "v.T.dot([0,1]), v.T.dot([1,1])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 17,
       "text": [
        "(array([ 0.383, -0.924]), array([-0.541, -1.307]))"
       ]
      }
     ],
     "prompt_number": 17
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "e.max()/e.sum()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 18,
       "text": [
        "0.73570226039551578"
       ]
      }
     ],
     "prompt_number": 18
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "**Question 4 (10 points)**.\n",
      "\n",
      "Find the minimum of the following quadratic function on $\\mathbb{R}^2$ \n",
      "\n",
      "$$f(x) = x^TAx +b^Tx +c$$\n",
      "where\n",
      "$$A = \\left(\\begin{matrix}13&5\\\\5&7\\end{matrix}\\right), b = \\left(\\begin{matrix}1\\\\1\\end{matrix}\\right) \\textrm {and } c = 2$$\n",
      "\n",
      "Under the constraints:\n",
      "$$g(x) = 2x_1-5x_2=2 \\;\\;\\;\\;\\;\\; \\textrm{ and } \\;\\;\\;\\;\\;\\; h(x) = x_1+x_2=1$$\n",
      "\n",
      "1. Use a matrix decomposition method to find the minimum of the *unconstrained* problem without using `scipy.optimize` (Use library functions - no need to code your own). Note: for full credit you should exploit matrix structure. (3 points)\n",
      "2. Find the solution using constrained optimization with the `scipy.optimize` package. (3 points)\n",
      "2. Use Lagrange multipliers and solving the resulting set of equations directly without using `scipy.optimize`. (4 points)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import scipy.optimize as opt"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 19
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[13,5],[5,7]])\n",
      "b = np.array([1.0,1.0]).T\n",
      "c = 2\n",
      "\n",
      "la.cho_solve(la.cho_factor(A), -b/2)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 20,
       "text": [
        "array([-0.015, -0.061])"
       ]
      }
     ],
     "prompt_number": 20
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def f(x, A, b, c):\n",
      "    return x.T.dot(A).dot(x) + b.T.dot(x) + c\n",
      "\n",
      "# check unconstrained solution\n",
      "usol = opt.minimize(f, [0,0], args=(A, b, c))\n",
      "usol.x\n",
      "\n",
      "cons = ({'type': 'eq', 'fun': lambda x: 2*x[0] - 5*x[1] - 2},\n",
      "        {'type': 'eq', 'fun': lambda x: x[0] + x[1] - 1})\n",
      "\n",
      "opt.minimize(f, [0,0], constraints=cons, args=(A, b, c))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 21,
       "text": [
        "  status: 0\n",
        " success: True\n",
        "    njev: 3\n",
        "    nfev: 14\n",
        "     fun: 15.999999999999996\n",
        "       x: array([ 1.000,  0.000])\n",
        " message: 'Optimization terminated successfully.'\n",
        "     jac: array([ 27.000,  11.000,  0.000])\n",
        "     nit: 3"
       ]
      }
     ],
     "prompt_number": 21
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "M = np.array([\n",
      "    [26, 10, 2, 1],\n",
      "    [10, 14, -5, 1],\n",
      "    [2, -5, 0, 0],\n",
      "    [1, 1, 0, 0]\n",
      "])\n",
      "\n",
      "y = np.array([-1, -1, 2, 1]).T\n",
      "\n",
      "la.solve(M, y)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 22,
       "text": [
        "array([ 1.000, -0.000, -2.286, -22.429])"
       ]
      }
     ],
     "prompt_number": 22
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "**Question 5 (10 points)**. \n",
      "\n",
      "Consider the linear transformation $f(x)$ on $\\mathbb{R}^3$ that takes the standard basis $\\left\\{e_1,e_2,e_3\\right\\}$ to $\\left\\{v_1,v_2,v_3\\right\\}$ where\n",
      "\n",
      "$$v_1=\\left(\\begin{matrix}10\\\\-10\\\\16\\end{matrix}\\right), v_2=\\left(\\begin{matrix}2\\\\-5\\\\20\\end{matrix}\\right) \\textrm {and } v_3=\\left(\\begin{matrix}1\\\\-4\\\\13\\end{matrix}\\right)$$\n",
      "\n",
      "1. Write a matrix $A$ that represents the same linear transformaton. (2 points)\n",
      "\n",
      "2. Compute the rank of $A$ (use any method you like). (2 points)\n",
      "\n",
      "3. Find the eigenvalues and eigenvectors of $A$. (2 points)\n",
      "\n",
      "4. What is the matrix representation of $f$ with respect to the eigenbasis? (4 points)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([\n",
      "    [10, 2, 1],\n",
      "    [-10, -5, -4],\n",
      "    [16, 20, 13],\n",
      "])\n",
      "\n",
      "A.dot(np.eye(3))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 23,
       "text": [
        "array([[ 10.000,  2.000,  1.000],\n",
        "       [-10.000, -5.000, -4.000],\n",
        "       [ 16.000,  20.000,  13.000]])"
       ]
      }
     ],
     "prompt_number": 23
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.linalg.matrix_rank(A)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 24,
       "text": [
        "3"
       ]
      }
     ],
     "prompt_number": 24
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "e, v = la.eig(A)\n",
      "e = np.real_if_close(e)\n",
      "print e\n",
      "print v"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 9.000  3.000  6.000]\n",
        "[[-0.577 -0.000 -0.120]\n",
        " [ 0.577  0.447 -0.241]\n",
        " [-0.577 -0.894  0.963]]\n"
       ]
      }
     ],
     "prompt_number": 25
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.diag(e)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 26,
       "text": [
        "array([[ 9.000,  0.000,  0.000],\n",
        "       [ 0.000,  3.000,  0.000],\n",
        "       [ 0.000,  0.000,  6.000]])"
       ]
      }
     ],
     "prompt_number": 26
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "**Question 6 (10 points)**\n",
      "\n",
      "Given the the function $f(x) = x^TAx + b^Tx$ and the Hessian matrix:\n",
      "\n",
      "$$A = \\left(\\begin{matrix}13&5&-3\\\\5&11&7\\\\-3&7&20\\end{matrix}\\right), b = \\left(\\begin{matrix}1\\\\1\\\\1\\end{matrix}\\right),$$\n",
      "\n",
      "complete the following code to find the first 3 vector directions for a conjugate gradient algorithm.\n",
      "\n",
      "Hint: Recall from linear algebra that the projection of $v$ on $u$ is \n",
      "\n",
      "$$\n",
      "\\frac{v \\cdot u}{u \\cdot u}u\n",
      "$$\n",
      "and that $u$ is conjugate to $v$ with respect to $A$ if \n",
      "$$Av\\cdot u = 0$$\n",
      "\n",
      "Also, recall that the gradient at $x_k$ of the quadratic function $f$ is given by $Ax_k + b$.\n",
      "\n",
      "\n",
      "<p>\n",
      "<font color=red>There are two lines of code to complete: the correction term in the function next_p, and the initial value of p = p0 in the main body of code.</font></p>"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def next_x (A,x,b,p):\n",
      "    return x + ((np.dot(p,-b))/(np.dot(p,A.dot(p))))*p # returns x_{k+1} given x_k, p_k, A and b"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 27
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def next_p(A,x,b,ps):\n",
      "    correction = np.zeros_like(x)\n",
      "    for i in range(len(ps)):\n",
      "        correction += ((np.dot(ps[i].dot(A),(-b-A.dot(x))))/(np.dot(ps[i].dot(A),ps[i])))*ps[i]\n",
      "    \n",
      "    return -b - A.dot(x) - correction \n",
      "                                       \n",
      "# Complete the indicated line above with the correction term so that this function returns \n",
      "# the gradient at x minus the correction term that makes the new p conjugate to all the others"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 28
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[13.,5.,-3.],[5.,11.,7.],[-3.,7.,20.]])   #  A\n",
      "b = np.array([1.,1.,1.])                           # Value of b\n",
      "x0 = np.array([0.,0.,0.])                          # Initial guess for xmin\n",
      "\n",
      "p0 = -b - A.dot(x0)                  # fill in p0 (5 points)\n",
      "\n",
      "x = x0 # initializes x for the for loop \n",
      "p = p0 # initializes p for the for loop\n",
      "ps = []  # start list to store the conjugate vectors\n",
      "xs = []\n",
      "ps.append(p0) # store first vector\n",
      "xs.append(x0) # store first vector\n",
      "\n",
      "for i in range(2):\n",
      "    x = next_x(A,x,b,p)   # gets x_{k+1}\n",
      "    p = next_p(A,x,b,ps)  # computes p_{k+1}\n",
      "    ps.append(p)          # stores p_{k+1} in list\n",
      "    xs.append(x)   \n",
      "\n",
      "ps  # prints the list of three vectors"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 29,
       "text": [
        "[array([-1.000, -1.000, -1.000]),\n",
        " array([-0.312,  0.075,  0.123]),\n",
        " array([-0.066,  0.118, -0.072])]"
       ]
      }
     ],
     "prompt_number": 29
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "xs"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 30,
       "text": [
        "[array([ 0.000,  0.000,  0.000]),\n",
        " array([-0.048, -0.048, -0.048]),\n",
        " array([-0.069, -0.044, -0.040])]"
       ]
      }
     ],
     "prompt_number": 30
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def f(x):\n",
      "    return 0.5*np.dot(np.dot(x,A),x) + np.dot(b,x)\n",
      "def fprime(x):\n",
      "    return np.dot(A,x) + b"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 31
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "opt.minimize(f,x0,method=\"newton-CG\",jac=fprime)\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 32,
       "text": [
        "  status: 0\n",
        " success: True\n",
        "    njev: 26\n",
        "    nfev: 6\n",
        "     fun: -0.07850070721357852\n",
        "       x: array([-0.083, -0.017, -0.057])\n",
        " message: 'Optimization terminated successfully.'\n",
        "    nhev: 0\n",
        "     jac: array([ 0.000, -0.000, -0.000])"
       ]
      }
     ],
     "prompt_number": 32
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 32
    }
   ],
   "metadata": {}
  }
 ]
}