Module1-PreliminaryPython/M1D_NumPy_and_Linear_Algebra.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "4nmji1trStg4"
   },
   "source": [
    "<html>\n",
    "    <summary></summary>\n",
    "    <p float=\"left\">\n",
    "         <div> <p></p> </div>\n",
    "         <div style=\"font-size: 20px; width: 800px;\">\n",
    "              <h1>\n",
    "               <left>Intro to Python: NumPy and Linear Algebra Review </left>\n",
    "              </h1>\n",
    "              <p><left>============================================================================</left> </p>              \n",
    "             <pre>Course: BIOM 421, Spring 2024\n",
    "Instructor: Dr. Brian Munsky\n",
    "Contact Info: munsky@colostate.edu\n",
    "Authors: Will Raymond, Dr. Luis Aguilera, Dr. Brian Munsky\n",
    "</pre>\n",
    "         </div>\n",
    "    </p>\n",
    "\n",
    "</html>\n",
    "\n",
    "\n",
    "\n",
    "<details>\n",
    "  <summary>Copyright info</summary>\n",
    "\n",
    "```\n",
    "Copyright 2023 Brian Munsky\n",
    "\n",
    "Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:\n",
    "\n",
    "1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.\n",
    "\n",
    "2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.\n",
    "\n",
    "3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.\n",
    "\n",
    "THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS \"AS IS\" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\n",
    "```\n",
    "<details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "hKUM83zzMfrg"
   },
   "source": [
    "# 1. Basic NumPy\n",
    "\n",
    "<html>\n",
    "<summary></summary>\n",
    "<p float=\"left\">\n",
    "<img align=\"left\" src=\"https://numpy.org/images/logo.svg\" width=\"150\" height=\"150\">\n",
    "</p>\n",
    "</html>\n",
    "\n",
    "In this notebook we will take a look at one of the most widely used data manipulation packages within Python, Numpy!\n",
    "\n",
    "numpy and scipy are ubiquitous in the world of scientific computing in Python. **NumPy**, NumericPython, is a library provides access to a lot of generalized array functions for linear algebra or manipulation that is not available with the standard list typing in Python. In addition to providing its powerful ```ndarray``` type, ```numpy``` also includes a myriad of packages for random simulation, basic statistics, sorting, and various matrix transforms.\n",
    "\n",
    "The core of numpy is the array, ```ndarray```, so let's take a look at that first.\n",
    "\n",
    "[Numpy Documentation](https://numpy.org/doc/stable/user/whatisnumpy.html)\n",
    "\n",
    "\n",
    "<img align=\"right\" src = \"https://media.springernature.com/full/springer-static/image/art%3A10.1038%2Fs41586-020-2649-2/MediaObjects/41586_2020_2649_Fig1_HTML.png?as=webp\" >\n",
    "\n",
    "\n",
    "------\n",
    "```Reading: Kinder, Nelson Section 2```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "faGG0LXFw8bQ"
   },
   "outputs": [],
   "source": [
    "  # Plotting style - here we set some options for pltting later in this notebook.\n",
    "  import matplotlib.pyplot as plt\n",
    "  from cycler import cycler\n",
    "  import matplotlib.cm as cm\n",
    "  from matplotlib.lines import Line2D\n",
    "  colors = ['#ef476f', '#06d6a0','#7400b8','#073b4c', '#118ab2',]\n",
    "  #colors = ['#fa8174', '#b3de69', '#bc82bd','#ccebc4','#ffed6f','#81b1d2']\n",
    "  font = {'family' : 'monospace',\n",
    "          'weight' : 'bold',\n",
    "          'size'   : 12}\n",
    "  plt.rcParams.update({'font.size': 12, 'font.weight':'bold', 'font.family':'monospace'  }   )\n",
    "  plt.rcParams.update({'axes.prop_cycle':cycler(color=colors)})\n",
    "  plt.rcParams.update({'axes.prop_cycle':cycler(color=colors)})\n",
    "  plt.rcParams.update({'axes.prop_cycle':cycler(color=colors)})\n",
    "  plt.rcParams.update({'xtick.major.width'   : 2.8 })\n",
    "  plt.rcParams.update({'xtick.labelsize'   : 12 })\n",
    "  plt.rcParams.update({'ytick.major.width'   : 2.8 })\n",
    "  plt.rcParams.update({'ytick.labelsize'   : 12})\n",
    "  plt.rcParams.update({'axes.titleweight'   : 'bold'})\n",
    "  plt.rcParams.update({'axes.titlesize'   : 10})\n",
    "  plt.rcParams.update({'axes.labelweight'   : 'bold'})\n",
    "  plt.rcParams.update({'axes.labelsize'   : 12})\n",
    "  plt.rcParams.update({'axes.linewidth':2.8})\n",
    "  plt.rcParams.update({'axes.labelpad':8})\n",
    "  plt.rcParams.update({'axes.titlepad':10})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "F62JOW7VKQpc"
   },
   "source": [
    "## 1.A.  Making your first array\n",
    "\n",
    "There are multiple ways of making an array, but one of the easiest is passing ```np.array``` a list or list of lists, or even a list of lists of lists!\n",
    "\n",
    "We will look at some other ways in a little bit.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "7enN7MQdc_6a",
    "outputId": "1da912d6-4dd9-4a9a-b6dd-4934f5e27bee"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "my_list = [0,1,2,3,4]\n",
    "my_np_array = np.array(my_list)\n",
    "\n",
    "print(my_np_array)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "AayRme4tU96Q"
   },
   "source": [
    "Now that we have this object, let's compare its methods to the methods of a regular list:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "QvBE6-bdUmv7",
    "outputId": "134dc79f-4b34-43b0-e1c6-bdec8fe0d145"
   },
   "outputs": [],
   "source": [
    "print('The functions of np_arrays are:')\n",
    "print(dir(my_np_array))\n",
    "print('\\nThe functions of lists are:')\n",
    "print(dir(my_list))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "DQ9yAPwTVIT2"
   },
   "source": [
    "Numpy arrays have way more functions associated with them than a regular list! Let's take a look at some of the more commonly used ones:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "9DKQJ6PQPC1z"
   },
   "source": [
    "## 1.B.  Array information methods\n",
    "\n",
    "Now that we have an array in memory, how do we look at some of its features? How big is it? What is its shape? Thankfully, ```numpy``` has some useful methods to get information such as shape, dimensions, and information type.\n",
    "\n",
    "Here's some common ones:\n",
    "\n",
    "given ```arr = np.array([[0,1,2,3,4]])```\n",
    "\n",
    "| Method    | Description  | Result\n",
    "| ----------- | ----------- | ------------ |\n",
    "| arr.shape   |  returns a tuple of the size of each dimensions | ```(1,5) ```   |\n",
    "| arr.ndim | number of dimensions |```2 ```   |\n",
    "| arr.dtype   |  the datatype of the info inside the array | ```dtype('int64') ```    \n",
    "| arr.size | integer total entries of this array | ```5```  |  \n",
    "| arr.nbytes | the number of bytes in RAM consumed by storing only the information inside the array | ``` 40 ```|\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "vPnY4LuYO8bM",
    "outputId": "7ead1820-c5fc-4e49-f895-dcad3eed6e61"
   },
   "outputs": [],
   "source": [
    "arr = np.array([[0,1,2,3,4]])\n",
    "\n",
    "desc = {'shape': '',\n",
    "    'ndim': 'dimensions',\n",
    "    'dtype': '',\n",
    "    'size': 'elements',\n",
    "    'nbytes':'bytes'}\n",
    "\n",
    "arr_methods = {'shape': arr.shape,\n",
    "    'ndim':arr.ndim,\n",
    "    'dtype':arr.dtype,\n",
    "    'size':arr.size,\n",
    "    'nbytes':arr.nbytes}\n",
    "\n",
    "for x in [\"shape\", \"ndim\", \"dtype\",\"size\",\"nbytes\"]:\n",
    "    dropdown = x\n",
    "\n",
    "    print('The %s of the array is %s %s' % (dropdown, str(arr_methods[dropdown]), desc[dropdown] ))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "6PKZyosKS8nD"
   },
   "source": [
    "## 1.C.  Array allocation methods\n",
    "\n",
    "```numpy``` also has some nice methods to preallocate arrays of any size or shape you wish, that you can then fill later with some values.\n",
    "\n",
    "| Method    | Description  | Examples  | Result\n",
    "| ----------- | ----------- | ------------ | ------------ |\n",
    "| np.zeros([n1,n2...]) |  Returns an ndim array of zeros  | ```np.zeros([10,10,10]) ```   | a (10 x 10 x 10) 3d array of zeros\n",
    "| np.ones([n1,n2...]) | Returns an ndim array of ones  |```np.ones([3,5]) ```   | a (3 x 5) 2d array of ones |\n",
    "| np.array([n1,n2...]) | converts an iterable to an array | ``np.array([1,2,3])``  |  a 1d array of [1,2,3] |\n",
    "| np.empty([n1,n2...]) | an empty ndim array | ``` np.empty([3,3]) ```| a (3 x 3) empty array - note values will be from whatever memory it allocated.\n",
    "\n",
    "Specific cases for easily allocating ranges/identity matrix\n",
    "\n",
    "| Method    | Description  | Examples  \n",
    "| ----------- | ----------- | ------------ |\n",
    "| np.linspace(start, stop, n_steps) | returns a linearly spaced array from start to stop with n_steps | ```np.linspace(0,10,11) ```  |  \n",
    "| np.logspace(start, stop, n_steps) | returns a linearly spaced array from 10^start to 10^stop with n_step | ``` np.logspace(0,10,11)```|\n",
    "| np.arange(stop,) | returns a linearly spaced array from a given start or 0 to stop with integer steps (or optional step sizes) | ``` np.logspace(0,10,11)```|\n",
    "| np.eye(N)   | gives a 2d identify matrix of NxN|  ```np.eye(3) ```   | a (3 x 3) identity matrix (diagonals are 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "fs7v6BGRbMYg",
    "outputId": "fe1d70ce-5385-47ad-f15c-724fc7838d7a"
   },
   "outputs": [],
   "source": [
    "print('np.zeros:')\n",
    "print(np.zeros([3,3]))\n",
    "\n",
    "print('')\n",
    "print('np.ones:')\n",
    "print(np.ones([3,3]))\n",
    "\n",
    "print('')\n",
    "list_of_lists = [ [y for y in range(3*x,3*x+3)] for x in range(0,3)]\n",
    "print('np.array() of a list of lists')\n",
    "print(np.array(list_of_lists ))\n",
    "\n",
    "print('')\n",
    "print('np.empty:')\n",
    "print(np.empty([3,3]))\n",
    "print('')\n",
    "print('np.eye(3):')\n",
    "print(np.eye(3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "_7DYx1D5cxLl",
    "outputId": "aee5e426-84aa-4960-a485-fb77f1c0a65c"
   },
   "outputs": [],
   "source": [
    "# Using linspace \n",
    "start = 1 \n",
    "stop = 14\n",
    "n_steps = 49\n",
    "\n",
    "print('np.linspace(start, stop, n_steps):')\n",
    "print(np.linspace(start, stop, n_steps))\n",
    "\n",
    "print('\\nnp.logspace(log10(start), log10(stop), n_steps -- in logspace):')\n",
    "print(np.logspace(np.log10(start), np.log10(stop), n_steps))\n",
    "\n",
    "print('\\nnp.arange(stop):')\n",
    "print(np.arange(stop))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "dlES1vmKbOeS"
   },
   "source": [
    "## 1.D.  Array manipulation methods\n",
    "\n",
    "### Slicing/ Indexing rules\n",
    "\n",
    "**array[ row, column, ...., n ]**\n",
    "\n",
    "Each dimension in a ```numpy``` array can be sliced or indexed just like its a python list:\n",
    "\n",
    "For example, consider this array:\n",
    "\n",
    "```\n",
    "arr = np.array([[1,2,3], \\\n",
    "                [4,5,6], \\\n",
    "                [7,8,9]])\n",
    "```\n",
    "\n",
    "| Example    | Description  | Result  \n",
    "| ----------- | ----------- | ------------ |\n",
    "| arr[0] | row 0 | ``` [1,2,3] ```  |  \n",
    "| arr[0,0]| element in row 0, column 0 | ``` [1]```|\n",
    "| arr[:,0]| column 0 | ``` [1,4,7]```|\n",
    "| arr[0:2,0:2]| top left block of size 2 x 2 | ``` array([[1,2],[4,5]])```|\n",
    "| arr[::-1]| reverses the rows | ```[[7,8,9],[4,5,6],[1,2,3]] ```|\n",
    "| arr[::-1, ::-1]| reverses the rows AND columns | ```[[9,8,7],[6,5,4],[3,2,1]] ```|\n",
    "\n",
    "\n",
    "### Array logicals\n",
    "\n",
    "```numpy``` arrays will accept logical indexing or any indexing as long as the shape and maximum indices are correct.\n",
    "\n",
    "| Example    | Description  | Result  \n",
    "| ----------- | ----------- | ------------ |\n",
    "| arr[arr > 5 ]| returns all elements greater than 5| ```array([6, 7, 8, 9]) ```|  \n",
    "| arr[1][arr[1] > 5 ]| returns all elements in row 1 greater than 5| ```array([6]) ```|  \n",
    "\n",
    "### Other common manipulation methods:\n",
    "\n",
    "| Method    | Description  | Examples  \n",
    "| ----------- | ----------- | ------------ |\n",
    "| arr.reshape() | will resize an array, the new size must have the same amount of elements!| ```arr.resize(1,9)```  |  \n",
    "| arr.T | transposes the array | ``` arr.T```|\n",
    "| np.expand_dims(array, axis) | adds an axis to a given array| ``` np.expand_dims(arr,axis=0)```|\n",
    "| arr.flatten()   | flattens an array to 1 dimension |  ```arr.flatten()```   |\n",
    "| np.diag(array,offset)   | gets the diagonal of an array based on offset from middle diag |  ```np.diag(arr,0)```   |\n",
    "| np.where(statement)   | returns indexes matching the statement |  ```np.where(arr < 5)```   |\n",
    "| np.at_least2d(arr)   | adds an axis to a 1d array if needed, leaves any ndim >= 2 alone|  ```np.atleast_2d(arr)```   |\n",
    "| np.tile(iterable, N) | repeat an iterable object N times along an axis | ``` np.tile(arr, N)```|\n",
    "| np.hstack(iterable of iterables) |Horizontally stack the arrays provided in an iterable (list of arrays) | ```np.hstack([arr1,arr2])``` |\n",
    "| np.vstack(iterable of iterables) | Vertically stack the arrays provided in an iterable (list of arrays) | ```np.vstack([arr1,arr2])```|\n",
    "| ___.astype(dtype) | convert a given array to a particular data type | ```np.array([1.2,1.4,1.5]).astype(int)```|\n",
    "| ___.tolist() | convert a given array to a python list or lists of lists| ```np.array([1.2,1.4,1.5]).tolist()```|"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "8NYrjNUngD-1",
    "outputId": "5fda1706-12f8-44dc-c621-4ffc62253f46"
   },
   "outputs": [],
   "source": [
    "# Let's try some of those out.\n",
    "arr = np.array([[1,2,3], \\\n",
    "                [4,5,6], \\\n",
    "                [7,8,9]])\n",
    "\n",
    "print('\\nOriginal:')\n",
    "print(arr)\n",
    "\n",
    "print('\\nReshaped:')\n",
    "print(arr.reshape(1,9))\n",
    "\n",
    "print('\\nFlattened:')\n",
    "print(arr.flatten())\n",
    "\n",
    "print('\\nTranspose:')\n",
    "print(arr.T)\n",
    "\n",
    "print('\\nadd an extra dimension to the front:')\n",
    "print(np.expand_dims(arr, axis=0))\n",
    "\n",
    "print('\\nmiddle diagonal:')\n",
    "print(np.diag(arr, 0))\n",
    "\n",
    "print('\\nindexes of array where the element is less than 5:')\n",
    "row,col = np.where(arr < 5)\n",
    "print(f'rows:    {row}')\n",
    "print(f'columns: {col}')\n",
    "\n",
    "print('\\nconvert the array to the floating point datatype')\n",
    "print(arr.astype(np.float64))\n",
    "\n",
    "print('\\nconvert the array to a python list')\n",
    "print(arr.tolist())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.E.  Using Logical Indexing to Manipulate Arrays\n",
    "When procesing large amounts of data, we often want to perform the same task again and a gain to similar entries in a maxtrix or vector. Logical indexing is a fast and convenient way to accomplish such tasks."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Let's start with a large matrix:\n",
    "A = np.array([ [y for y in range(3*x,3*x+15)] for x in range(0,12)])\n",
    "print(f'Our Original matrix is:\\n {A}')\n",
    "\n",
    "# We will construct a new matrix B that is saturated at a value 15.  \n",
    "# That is, all numbers greater than 15 are set at 15.\n",
    "\n",
    "# Because we do not want to change A itself, we will start by making a copy.\n",
    "B = np.copy(A)\n",
    "\n",
    "# Now we can apply our logic:\n",
    "B[B>15]=15\n",
    "print(f'\\nOur saturated matrix is:\\n {B}')\n",
    "\n",
    "# Let's try another.  This time, we want a matrix 'C' where all values are between 0 and 15, but\n",
    "# when the numbers exceed those limits, we reflect their values (e.g., 17 which is 2 greater than 15\n",
    "# becomes 13 which is two less than 15, or -6 which is 6 less than zero becomes 6).  For very large \n",
    "# or very small numbers, we may need to run multiple rounds of reflection, so we will use a while loop.\n",
    "C = np.copy(A)\n",
    "while (np.max(np.abs(C)))>15:\n",
    "    C[C>15] = 30 - C[C>15]\n",
    "    C[C<0] = -C[C<0]\n",
    "print(f'\\nOur reflected matrix is:\\n {C}')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "R6LBk0wFfD4-"
   },
   "source": [
    "## Questions and practice\n",
    "Answer the following questions to solidify your understanding of this material. \n",
    "You should now also be able to answer questions Q6 - Q9 in Homework 2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "4kkYmBPkp_hc"
   },
   "outputs": [],
   "source": [
    "## Generate an array of all zeros of shape 3,3,3,3. Double check that the shape is that using\n",
    "## the internal shape function. Convert this array to an integer array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "GWhQdm2kq0he"
   },
   "outputs": [],
   "source": [
    "## Generate a 4x4 Identity matrix and then remove the final row."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "_P0fdehbrEeP"
   },
   "outputs": [],
   "source": [
    "## Repeat the following array vertically 5 times to end up with a 5x3 matrix\n",
    "arr1 = np.array([[1,2,3]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "w0y4HAMgeula"
   },
   "outputs": [],
   "source": [
    "## transpose the following numpy array\n",
    "arr2 = np.array([[1,2,3],\n",
    "                 [4,5,6]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "hxBIIS_urcM3"
   },
   "outputs": [],
   "source": [
    "## Use array logicals to return a flattened array of values greater than 0 in the following array\n",
    "## Take the sum of all the negative values obtained from the logical operation.\n",
    "arr3 = np.array([[-1, 2, 3, -31, -0.3 ],\n",
    "                 [4, 0, 5, 39.3, -37,]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "J5N9AFcir26e"
   },
   "outputs": [],
   "source": [
    "## use np.where to find any indexes in the following array that is equal to zero\n",
    "arr4 = np.array([[-1, 2, 3, -31, 0 ],\n",
    "                 [4, 0, 5, 39.3,-37,],\n",
    "                 [5, 4, 0, 0, 1] ])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "fnOEzvm6sIUf"
   },
   "outputs": [],
   "source": [
    "## use np.sum with the correct axis=? argument to take the sum of the following array\n",
    "## along each column\n",
    "\n",
    "arr5 = np.array([[-1, 2, 3, -31, 13 ],\n",
    "                 [14, 0, 5, 39.3,-37,],\n",
    "                 [5, -3, 0, 0, 1] ])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "BMUi9_rxeu_F"
   },
   "source": [
    "# 2.  Linear Algebra in Numpy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Hs_77nq8eNKk"
   },
   "source": [
    "## 2.A.  Matrix math with Numpy\n",
    "\n",
    "NumPy provides a fast implementation of many linear algebra functions stored in ```np.linalg```.\n",
    "\n",
    "```numpy``` arrays also can be treated as matrices for arithmetic operations!\n",
    "Heres some examples of some basic arithmetic operations with two 2x2 matrices:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "fL8ejwYwzGMp",
    "outputId": "5d24a550-005d-43ec-ff2a-ce2e6642fde9"
   },
   "outputs": [],
   "source": [
    "A = np.array([[1, 2],\n",
    "              [3, 4]]) #matrix 1\n",
    "\n",
    "B = np.array([[5, 6],\n",
    "              [7, 8]]) #matrix 2\n",
    "n = 5\n",
    "\n",
    "examples = {'matrix multiplication': A @ B,\n",
    "            'element wise multiplication': A * B,\n",
    "            'matrix division': A @ np.linalg.pinv(B),\n",
    "            'element wise division': A/B,\n",
    "            'scalar multiplication': n*A,\n",
    "            'matrix subtraction': A - B,\n",
    "            'matrix addition': A + B,\n",
    "            'scalar addition': n + A,\n",
    "            }\n",
    "examples_str = {'matrix multiplication': 'A @ B:',\n",
    "                'element wise multiplication': 'A * B:',\n",
    "                'element wise division': 'A / B',\n",
    "            'matrix division': 'A @ np.linalg.pinv(B):',\n",
    "            'scalar multiplication': 'n*A',\n",
    "            'matrix subtraction': 'A - B:',\n",
    "            'matrix addition': 'A + B:',\n",
    "            'scalar addition': 'n + A',\n",
    "            }\n",
    "\n",
    "print(f'A:\\n{A}\\n')\n",
    "print(f'B:\\n{B}\\n')\n",
    "\n",
    "for dropdown in [\"matrix multiplication\",\"element wise multiplication\", \"matrix division\",\n",
    "                 \"element wise division\",\"scalar multiplication\", \"matrix subtraction\",\n",
    "                 \"matrix addition\",\"scalar addition\"]:\n",
    "\n",
    "    print(examples_str[dropdown])\n",
    "    print(f'{examples[dropdown]}\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "1H5V6SCT2Eee"
   },
   "source": [
    "Short list of linear algebra functionality in numpy\n",
    "\n",
    "| Method    | Description  |\n",
    "| ----------- | ----------- |\n",
    "| np.cross | cross product of vectors|\n",
    "| np.dot | dot product of arrays |\n",
    "| np.linalg.det | determinant of a matrix |\n",
    "| np.linalg.inner | inner product |\n",
    "| np.linalg.outer | outer product|\n",
    "| np.matmul  | matrix product |\n",
    "| np.trace | trace of a matrix  |\n",
    "| np.linalg.norm | norm of a vector or matrix, takes arguments for nuclear/frobenius |\n",
    "| np.linalg.matrix_rank | rank of a matrix |\n",
    "| np.linalg.cond | condition of a matrix |\n",
    "| np.linalg.inv  | matrix inverse (multiplicative) |  \n",
    "| np.linalg.pinv  | matrix psuedo-inverse (Moore-Penrose) |\n",
    "| np.linalg.svd | Single Value Decomposition |\n",
    "| np.linalg.null* | Null space |\n",
    "\n",
    "For more specific commands check out the linalg module in numpy:\n",
    "https://numpy.org/doc/stable/reference/routines.linalg.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "PI-mbvbl3U8-"
   },
   "source": [
    "### Eigen vectors and values\n",
    "\n",
    "The eigenvector is a characteristic vector of a data set such that any linear transformation applied to the data set only scales the eigenvector linearly. The eigenvectors scale with an associated eigenvalue. Another way of wording that is if I have some set of data, any linear transformation will also linearly transform its eigenvectors, dictated by its eigenvalues.\n",
    "\n",
    "| Method    | Description  | Examples  \n",
    "| ----------- | ----------- | ------------ |\n",
    "| np.linalg.eig | returns the eigen vectors and values of a matrix| |  \n",
    "| np.linalg.eigvals | returns the eigenvalues of a matrix |\n",
    "\n",
    "For a very in depth explanation  check out the following youtube video: [\"Eigenvectors and eigenvalues | Chapter 14, Essence of linear algebra\"](https://www.youtube.com/watch?v=PFDu9oVAE-g)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "cp5oSUx8BPmI",
    "outputId": "6ce91724-630c-49c1-eb89-e38c49cdfb21"
   },
   "outputs": [],
   "source": [
    "# Create a random symmetrix matrix\n",
    "random_matrix = np.random.randn(5,5)\n",
    "random_matrix = random_matrix+random_matrix.T\n",
    "\n",
    "# compute eigenvalue/eigenvector pairs\n",
    "eigenvals, eigenvecs = np.linalg.eig(random_matrix)\n",
    "\n",
    "print('Random symmetric matrix:')\n",
    "print(random_matrix)\n",
    "\n",
    "print('\\nEigen vectors:')\n",
    "print(eigenvecs)\n",
    "print('\\nEigen values:')\n",
    "print(eigenvals)\n",
    "\n",
    "# Note that for symmetric matrices, the eigenvalues and eigenvectors are always real.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "y3v-7TCNkujR"
   },
   "source": [
    "## 2.B.  Matrix Transformations\n",
    "\n",
    "A linear transformation allows for one to transform a vector in one othonormal basis vector space (i.e., where all basis vectors are orthogonal to one another) to another, potentially more convenient cordinate definition.  \n",
    "\n",
    "To achieve such a transformation, we can multiple a vector in the original coordinate system by a unitary matrix (a matrix whose collumns are orthogonal to one another).\n",
    "\n",
    "For example, suppose that we examine a point in 2D denoted by [x,y] in the regular cartesian coordinate system. We can rotate this coordinate systme about the origin [0,0] by an angle $\\theta$ (defined as positive in the counter clockwise direction), using the linear transformation:\n",
    "\n",
    "$\\begin{bmatrix}\n",
    "x'\\\\\n",
    "y'\n",
    "\\end{bmatrix} =\n",
    "\\begin{bmatrix}\n",
    "\\cos(\\theta) & -\\sin(\\theta)\\\\\n",
    "\\sin(\\theta) & \\cos(\\theta) \\\\\n",
    "\\end{bmatrix} * \\begin{bmatrix}\n",
    "x\\\\\n",
    "y\n",
    "\\end{bmatrix}\n",
    "$\n",
    "\n",
    "You chould confirm on your own that this transformation matrix is indeed a unitary matrix (i.e., that all columns are orthonormal) by showing that if $\\mathbf{v}_i$ and $\\mathbf{v}_j$ are both columns of the transformation matrix, then the dot product $\\mathbf{v}_i\\cdot \\mathbf{v}_j$ is one if $i=j$ and zero if $i\\ne j$. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 476
    },
    "id": "TdSpbrTItep6",
    "outputId": "d6a1b217-7d28-4a10-a6ab-01ac9e3142da"
   },
   "outputs": [],
   "source": [
    "xy = np.random.randint(0,100,size=(2,100)) # a random constellation of points at integer values\n",
    "center = np.array([[30],[15]]) # axis of rotation for our image.\n",
    "\n",
    "print(xy - center)\n",
    "# perform a rotation:\n",
    "thetas = [.02, .1, .2]\n",
    "transformation_mat = lambda theta: np.array([[np.cos(theta), -np.sin(theta)],\n",
    "                                            [np.sin(theta),  np.cos(theta)]])\n",
    "\n",
    "# we are going to translate the points to the new center point, then apply the transformation matrix\n",
    "# application of the rotation:  A^-1 (2x2) @ (xy (100x2)- rotation point (2x1)).T (2 x 100)\n",
    "new_xy1 = transformation_mat(thetas[0]) @ (xy - center) + center\n",
    "new_xy2 = transformation_mat(thetas[1]) @ (xy - center) + center\n",
    "new_xy3 = transformation_mat(thetas[2]) @ (xy - center) + center\n",
    "\n",
    "plt.plot(center[0], center[1], 'x')\n",
    "plt.scatter(*xy)\n",
    "plt.scatter(*new_xy1)\n",
    "plt.scatter(*new_xy2)\n",
    "plt.scatter(*new_xy3)\n",
    "plt.legend(['center', 'original','.02 rad', '.1 rad', '.2 rad'], bbox_to_anchor=(1.05, 1.05))\n",
    "plt.title('rotating a 100 2D points w/ rotation matrix')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "mGtqR_dSxefY"
   },
   "source": [
    "Notice how we had to translate our constellation to a new center then rotate, then translate it back?\n",
    "\n",
    "With clever application of linear algebra, we can do all of these in one step using an Affine Transformation Matrix which can be constructed as the product of multiple transformation matrixes. Note how above we (1) translated the matrix to its new center, (2) then rotated, and then (3) we translated back to the original origin. We have to keep that order here too when we apply these transformation matrices, so the order is $MT = M3\\cdot M2\\cdot M1$, where\n",
    "\n",
    "$\n",
    "M1 = \\begin{bmatrix}\n",
    "1 & 0 & -\\Delta x\\\\\n",
    "0 & 1 & -\\Delta y\\\\\n",
    "0 & 0 & 1\n",
    "\\end{bmatrix}$\n",
    "\n",
    "$M2 = \\begin{bmatrix}\n",
    "\\cos(\\theta) & -\\sin(\\theta) & 0\\\\\n",
    "\\sin(\\theta) & \\cos(\\theta) & 0\\\\\n",
    "0 & 0 & 1\n",
    "\\end{bmatrix}$\n",
    "\n",
    "$\n",
    "M3 = \\begin{bmatrix}\n",
    "1 & 0 & \\Delta x\\\\\n",
    "0 & 1 & \\Delta y\\\\\n",
    "0 & 0 & 1\n",
    "\\end{bmatrix}$\n",
    "\n",
    "---\n",
    "Applying the transformation to a set of points:\n",
    "\n",
    "$\\begin{bmatrix}\n",
    "x'\\\\\n",
    "y'\\\\\n",
    "1\n",
    "\\end{bmatrix} =\n",
    "M3\\cdot M2 \\cdot M1 \\cdot \\begin{bmatrix}\n",
    "x\\\\\n",
    "y\\\\\n",
    "1\n",
    "\\end{bmatrix}\n",
    "$\n",
    "\n",
    "You can try other affine transformation matrices as well -- but for now we are just demonstrating the common linear algebra transformation in ```numpy```."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 476
    },
    "id": "lNz-gIN9wyLa",
    "outputId": "4616ed26-7911-49be-ed81-d8f80c7bded8"
   },
   "outputs": [],
   "source": [
    "# Let's try it:  First, we add in a column of 1s to our x,y constellation so now its 3x100:\n",
    "print(len(xy[0,:]))\n",
    "xy_appended = np.vstack([xy, np.ones([1,len(xy[0,:])])])\n",
    "\n",
    "M1 = lambda rx, ry, theta: np.array([[1, 0, -rx],\n",
    "                                    [0,  1, -ry],\n",
    "                                    [0,  0,  1]])\n",
    "M2 = lambda rx, ry, theta: np.array([[np.cos(theta), -np.sin(theta), 0],\n",
    "                                     [np.sin(theta), np.cos(theta), 0],\n",
    "                                     [0,              0,             1]])\n",
    "M3 = lambda rx, ry, theta: np.array([[1, 0, rx],\n",
    "                                    [0,  1, ry],\n",
    "                                    [0,  0,  1]])\n",
    "\n",
    "\n",
    "MT = lambda rx, ry, theta: M3(rx,ry,theta)@M2(rx,ry,theta)@M1(rx,ry,theta)\n",
    "#Set the center point about which to rotate: \n",
    "rx = center[0,0]\n",
    "ry = center[1,0]\n",
    "\n",
    "#Set the angle to rotate counter-clockwise (radians)\n",
    "theta = 0.2\n",
    "\n",
    "xy_translated = MT(rx,ry,theta) @ xy_appended\n",
    "\n",
    "plt.plot(rx, ry, 'x')\n",
    "plt.scatter(*xy_appended[:2]) # only plot along the first two dimensions\n",
    "plt.scatter(*xy_translated[:2])\n",
    "plt.scatter(new_xy3[0,:],new_xy3[1,:],marker='.')\n",
    "\n",
    "plt.legend(['center', 'original',f'{theta} rad (new approach)',f'{theta} rad (old approach)'], bbox_to_anchor=(1.05, 1.05))\n",
    "plt.title('rotating a 100 2D points w/ affine rotation matrix')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "tocGXL9wkype"
   },
   "source": [
    "# 3. SVD and PCA\n",
    "\n",
    "Principal component analysis (PCA) is a specific linear transformation that chooses a set of ordered basis vectors (i.e., a new coordinate system) where each orthoganal axis is a \"Principal Component,\" usually arranged in an order such that each new vector captures the next greatest percentage of variation in the given dataset.  In other words, the first principal component is the direction that has the most variation, the second has the next most, etc.\n",
    "\n",
    "Singular value decomposition is an efficient linear algebra method that generalizes the eigenvalue-eigenvector decomposition to allow a similar decompositon to be performed on rectangular matrices. SVD decomposes a $m\\times n$ matrix $\\mathbf{X}$ as \n",
    "\n",
    "$\\mathbf{X} = \\mathbf{U}\\mathbf{\\Sigma}\\mathbf{V}^*,$\n",
    "\n",
    "where $U$ is an orthonormal matrix with m rows and $p=\\min(m,n)$ columns; $\\Sigma$ is a $p\\times p$ diagonal matrix with non-negative real values called the singular values and are usually ordered from largest to smallest; and V is a second orthonormal matrix with $n$ rows and $p=\\min(m,n)$ columns ($\\mathbf{V}^*$ is the conjugate transpose of $\\mathbf{V}$).\n",
    "\n",
    "We can use SVD to obtain the largest orthonormal eigenbasis for any m x n matrix, which is precisely the operation needed to perform PCA! We will show two implementations below, but first we need some data!\n",
    "\n",
    "A 2D multivariate gaussian is one of the best ways to visualize this particular transformation (without dimensionality reduction)! So let's generate one:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 458
    },
    "id": "tqXfuo5p63A-",
    "outputId": "d998c5fc-d5a4-4e8e-cf68-20f09ad7f454"
   },
   "outputs": [],
   "source": [
    "# Specify some parameters to define our distribution:\n",
    "mu = np.array([3.6,2.7])\n",
    "sig = np.array([[2,-1],[-1,2]])\n",
    "\n",
    "# Generate random data from multivariate gaussian:\n",
    "random_gaussian = np.random.multivariate_normal(mu,sig, size=(1000))\n",
    "\n",
    "# Plot the data:\n",
    "plt.scatter(random_gaussian[:,0], random_gaussian[:,1])\n",
    "plt.xlim([-10,10])\n",
    "plt.ylim([-10,10])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "olxKfu1Ji0PO"
   },
   "outputs": [],
   "source": [
    "# Let's get some help on how to use SVD.\n",
    "np.linalg.svd?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "_5xOoEGiO5_g",
    "outputId": "aad93346-c7a7-43b3-d4db-9b3367b45c96"
   },
   "outputs": [],
   "source": [
    "# PCA using  SVD\n",
    "\n",
    "# PCA requires CENTERRED data (zero mean), so let's subtract the mean:\n",
    "centered_gaussian = random_gaussian - random_gaussian.mean(axis=0)\n",
    "\n",
    "# Now to compute the SVD transformation of the zero-mean data:\n",
    "U, S, V = np.linalg.svd(centered_gaussian,full_matrices=False) #single valuar decomp with numpy\n",
    "\n",
    "# This results in:\n",
    "# U = unitary matrix (possibly complex) for the range (columns space) of the matrix. \n",
    "# S = singular values (non-negative real values)\n",
    "# V = unitary matrix (possibly complex) for the \n",
    "\n",
    "print(U.shape)\n",
    "print(S.shape)\n",
    "print(V.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 458
    },
    "id": "AjWhg_BcjEJy",
    "outputId": "a08caf8f-2c6e-4fd6-9815-3ada92ebe589"
   },
   "outputs": [],
   "source": [
    "X_new = centered_gaussian @ V # transform our original data after being centered\n",
    "\n",
    "plt.scatter(random_gaussian[:,0], random_gaussian[:,1])\n",
    "plt.scatter(X_new[:,0], X_new[:,1])\n",
    "plt.xlim([-10,10])\n",
    "plt.ylim([-10,10])\n",
    "\n",
    "#plot the eigen vectors by their eigen values!\n",
    "plt.quiver(*mu, *V[:,0], color=['k'], scale=np.sqrt(S[1]))\n",
    "plt.quiver(*mu, *V[:,1], color=['k'], scale=np.sqrt(S[0]))\n",
    "plt.quiver(*np.mean(X_new,axis=0), *(np.eye(2))[:,0], color=['k'], scale=np.sqrt(S[1]))\n",
    "plt.quiver(*np.mean(X_new,axis=0), *(np.eye(2))[:,1], color=['k'], scale=np.sqrt(S[0]))\n",
    "\n",
    "plt.legend(['Original','Transformed'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "f3ZUpBRlnPZ5"
   },
   "source": [
    "Note how the original gaussian has been placed on its 2 new axes perfectly by rotating the original data. Now what happens if we apply PCA and start removing axis / dimensions? Well for one, our reconstruction will be imperfect. You can calculate the exact amount of variance captured by including each axis, which we will do next."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Let's start by creating a data set that contains a mixture of two different types of \n",
    "# multivariate random variables:\n",
    "\n",
    "# We start by creating a random mean and covariance matrix for the two types of random variables.\n",
    "nVariables = 20\n",
    "\n",
    "# Let's generate a pair of random means for two types of data\n",
    "mn1 = np.random.rand(nVariables)\n",
    "mn2 = np.random.rand(nVariables)\n",
    "\n",
    "# Let's also generate a pair of random covariance matrices.\n",
    "sig1 = np.cov(np.random.rand(nVariables,6))\n",
    "sig2 = np.cov(np.random.rand(nVariables,6))\n",
    "\n",
    "# Next, we generate many samples from each distribution:\n",
    "nSamples = 1000\n",
    "samples1 = np.random.multivariate_normal(mn1,sig1,size=nSamples)\n",
    "samples2 = np.random.multivariate_normal(mn2,sig2,size=nSamples)\n",
    "\n",
    "# Now to combine the two data sets into a single normalized data set:\n",
    "X = np.vstack([samples1,samples2])\n",
    "std = np.sqrt(np.var(X,axis=0))\n",
    "X = (X-np.mean(X,axis=0))/std"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now, let's plot these two sample sets on a couple of their axes:\n",
    "ax0 = 0\n",
    "ax1 = 3\n",
    "plt.scatter(X[:1000,ax0], X[:1000,ax1])\n",
    "plt.scatter(X[1000:,ax0], X[1000:,ax1])\n",
    "plt.legend(['Set 1','Set 2'])\n",
    "\n",
    "# Note that the two data types are similar but not exactly the same.\n",
    "# Try different combinations to see if you can separate out the two colors."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "AOBWfEeoug2p"
   },
   "source": [
    "Now we can \"Select\" components by discarding axis along our matrices. Using only the first two components (first 2 indexes) of each matrix will let us construct some $\\hat{X}$ that is a lossy representation of our original data, $X$. Play with the number of components included, and see how the error of $\\hat{X}$ vs $X$ changes in relation to number of components used!\n",
    "\n",
    "We can also calculate the explained variance by squaring the sum of the singular values kept divided by the total singular values squared!\n",
    "\n",
    "$ Explained\\: Variance = S[:N_{components}]^2 / S^2$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "e1L4ZemAt9w3",
    "outputId": "cba68f3c-ed91-4204-8e75-aac0366b235b"
   },
   "outputs": [],
   "source": [
    "# Compute the SVD transformation of the zero-mean data:\n",
    "U, S, V = np.linalg.svd(X) #single value decomp with numpy\n",
    "\n",
    "components = [0,1,2]\n",
    "X_hat = U[:,components] @ (np.diag(S[components]) @ V[components,:]) # apply the transformation to the data\n",
    "explained_variance = np.sum(S[components]**2/np.sum(S**2))\n",
    "X_reduced = (X @ V[components,:].T)\n",
    "\n",
    "print('MSE transformed data to original data:')\n",
    "print(np.sum(np.subtract(X_hat, X)**2))\n",
    "print('')\n",
    "print(f'Total explained variance of components: {components}')\n",
    "print('{} %'.format(explained_variance*100))\n",
    "print('')\n",
    "print('reduced dimension:')\n",
    "print(X_reduced.shape)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "gwTdHE102ek4"
   },
   "source": [
    "So by keeping just a few principal components, we are able to reconstruct most of the variations in our data! Another way of saying this is that we can effectively reduce the dimension of our data set and keep track of how much information we lose.  \n",
    "\n",
    "But how does that affect the differences between the two different populations?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 463
    },
    "id": "L688HZdqq1ns",
    "outputId": "10b33e7f-55a1-49ef-c9f1-2d2b0de96932"
   },
   "outputs": [],
   "source": [
    "# Plotting our reduced representation based on the first 2 components\n",
    "components = [0,1]  # Pick which components to include.\n",
    "X_reduced = X @ V[components,:].T # Project data onto the\n",
    "\n",
    "plt.scatter( *X_reduced[:1000,:].T, )\n",
    "plt.scatter( *X_reduced[1000:,:].T, )\n",
    "plt.xlabel('PC 1')\n",
    "plt.ylabel('PC 2')\n",
    "\n",
    "plt.legend(['Set 1','Set 2'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "zXtV_vM6yZE1"
   },
   "source": [
    "Notice any trends? Using the first few PCA vectors, it is much easier to distinguish between the two data types. Later in the course we will explore this for classification of different types of data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "LlwS6h_e3FJv"
   },
   "source": [
    "## Questions and practice\n",
    "You should now be able to complete problems 12 and 13 in Homework 2, and also try the following on your own."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "yNHquo6D3Ebz"
   },
   "outputs": [],
   "source": [
    "## Apply a scaling affine transformation matrix to a constellation of points.  This transformation should\n",
    "## scale (i.e. strectch or compress) the points in the x and y directions relative to some spatially fixed \n",
    "## origin point. \n",
    "## Hint: [ Scaling x   0       -rx ]\n",
    "##       [ 0       Scaling y   -ry ]\n",
    "##       [ 0           0         1 ]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "fyVm_kXO27C0"
   },
   "outputs": [],
   "source": [
    "## Apply SVD to an image, then reconstruct the image with varying amounts of singular values\n",
    "## What do you notice about using more singular values for reconstruction vs less?\n",
    "## How does the \"quality\" of the reconstructed image relate to the explained variance?\n",
    "## Hint: You may have to apply SVD to each channel if you are using a H x W x Color image."
   ]
  }
 ],
 "metadata": {
  "colab": {
   "provenance": [],
   "toc_visible": true
  },
  "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.11.5"
  },
  "widgets": {
   "application/vnd.jupyter.widget-state+json": {
    "054083162da645f9b247be27173efe1e": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "FloatProgressModel",
     "state": {
      "_dom_classes": [],
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "FloatProgressModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/controls",
      "_view_module_version": "1.5.0",
      "_view_name": "ProgressView",
      "bar_style": "success",
      "description": "",
      "description_tooltip": null,
      "layout": "IPY_MODEL_3658851275694ee1a15f9d5a305c1e22",
      "max": 35624,
      "min": 0,
      "orientation": "horizontal",
      "style": "IPY_MODEL_3e7a113c1fa34f8b9810b782447b3e8c",
      "value": 35624
     }
    },
    "0fb04728349c4dde94aaf9897492743f": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "HTMLModel",
     "state": {
      "_dom_classes": [],
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "HTMLModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/controls",
      "_view_module_version": "1.5.0",
      "_view_name": "HTMLView",
      "description": "",
      "description_tooltip": null,
      "layout": "IPY_MODEL_80c5f442014946df8957cc48dcca541e",
      "placeholder": "​",
      "style": "IPY_MODEL_c6017ca79bae42bfa03eb72781ab031f",
      "value": "100%"
     }
    },
    "1b8264c7ea5a421584e4871c2d8062ea": {
     "model_module": "@jupyter-widgets/base",
     "model_module_version": "1.2.0",
     "model_name": "LayoutModel",
     "state": {
      "_model_module": "@jupyter-widgets/base",
      "_model_module_version": "1.2.0",
      "_model_name": "LayoutModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "LayoutView",
      "align_content": null,
      "align_items": null,
      "align_self": null,
      "border": null,
      "bottom": null,
      "display": null,
      "flex": null,
      "flex_flow": null,
      "grid_area": null,
      "grid_auto_columns": null,
      "grid_auto_flow": null,
      "grid_auto_rows": null,
      "grid_column": null,
      "grid_gap": null,
      "grid_row": null,
      "grid_template_areas": null,
      "grid_template_columns": null,
      "grid_template_rows": null,
      "height": null,
      "justify_content": null,
      "justify_items": null,
      "left": null,
      "margin": null,
      "max_height": null,
      "max_width": null,
      "min_height": null,
      "min_width": null,
      "object_fit": null,
      "object_position": null,
      "order": null,
      "overflow": null,
      "overflow_x": null,
      "overflow_y": null,
      "padding": null,
      "right": null,
      "top": null,
      "visibility": null,
      "width": null
     }
    },
    "1c63e0a213d2422fad1013a4c7dc40c6": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "HTMLModel",
     "state": {
      "_dom_classes": [],
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "HTMLModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/controls",
      "_view_module_version": "1.5.0",
      "_view_name": "HTMLView",
      "description": "",
      "description_tooltip": null,
      "layout": "IPY_MODEL_cee19a8eb4454b238b640ea1f7669a2f",
      "placeholder": "​",
      "style": "IPY_MODEL_7abc07cf6e164036be85a891356d0a79",
      "value": " 35624/35624 [00:40&lt;00:00, 870.16it/s]"
     }
    },
    "3658851275694ee1a15f9d5a305c1e22": {
     "model_module": "@jupyter-widgets/base",
     "model_module_version": "1.2.0",
     "model_name": "LayoutModel",
     "state": {
      "_model_module": "@jupyter-widgets/base",
      "_model_module_version": "1.2.0",
      "_model_name": "LayoutModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "LayoutView",
      "align_content": null,
      "align_items": null,
      "align_self": null,
      "border": null,
      "bottom": null,
      "display": null,
      "flex": null,
      "flex_flow": null,
      "grid_area": null,
      "grid_auto_columns": null,
      "grid_auto_flow": null,
      "grid_auto_rows": null,
      "grid_column": null,
      "grid_gap": null,
      "grid_row": null,
      "grid_template_areas": null,
      "grid_template_columns": null,
      "grid_template_rows": null,
      "height": null,
      "justify_content": null,
      "justify_items": null,
      "left": null,
      "margin": null,
      "max_height": null,
      "max_width": null,
      "min_height": null,
      "min_width": null,
      "object_fit": null,
      "object_position": null,
      "order": null,
      "overflow": null,
      "overflow_x": null,
      "overflow_y": null,
      "padding": null,
      "right": null,
      "top": null,
      "visibility": null,
      "width": null
     }
    },
    "3e7a113c1fa34f8b9810b782447b3e8c": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "ProgressStyleModel",
     "state": {
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "ProgressStyleModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "StyleView",
      "bar_color": null,
      "description_width": ""
     }
    },
    "7abc07cf6e164036be85a891356d0a79": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "DescriptionStyleModel",
     "state": {
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "DescriptionStyleModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "StyleView",
      "description_width": ""
     }
    },
    "80c5f442014946df8957cc48dcca541e": {
     "model_module": "@jupyter-widgets/base",
     "model_module_version": "1.2.0",
     "model_name": "LayoutModel",
     "state": {
      "_model_module": "@jupyter-widgets/base",
      "_model_module_version": "1.2.0",
      "_model_name": "LayoutModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "LayoutView",
      "align_content": null,
      "align_items": null,
      "align_self": null,
      "border": null,
      "bottom": null,
      "display": null,
      "flex": null,
      "flex_flow": null,
      "grid_area": null,
      "grid_auto_columns": null,
      "grid_auto_flow": null,
      "grid_auto_rows": null,
      "grid_column": null,
      "grid_gap": null,
      "grid_row": null,
      "grid_template_areas": null,
      "grid_template_columns": null,
      "grid_template_rows": null,
      "height": null,
      "justify_content": null,
      "justify_items": null,
      "left": null,
      "margin": null,
      "max_height": null,
      "max_width": null,
      "min_height": null,
      "min_width": null,
      "object_fit": null,
      "object_position": null,
      "order": null,
      "overflow": null,
      "overflow_x": null,
      "overflow_y": null,
      "padding": null,
      "right": null,
      "top": null,
      "visibility": null,
      "width": null
     }
    },
    "8fa5a4cd197c4a19b703b9e3178d546e": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "HBoxModel",
     "state": {
      "_dom_classes": [],
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "HBoxModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/controls",
      "_view_module_version": "1.5.0",
      "_view_name": "HBoxView",
      "box_style": "",
      "children": [
       "IPY_MODEL_0fb04728349c4dde94aaf9897492743f",
       "IPY_MODEL_054083162da645f9b247be27173efe1e",
       "IPY_MODEL_1c63e0a213d2422fad1013a4c7dc40c6"
      ],
      "layout": "IPY_MODEL_1b8264c7ea5a421584e4871c2d8062ea"
     }
    },
    "c6017ca79bae42bfa03eb72781ab031f": {
     "model_module": "@jupyter-widgets/controls",
     "model_module_version": "1.5.0",
     "model_name": "DescriptionStyleModel",
     "state": {
      "_model_module": "@jupyter-widgets/controls",
      "_model_module_version": "1.5.0",
      "_model_name": "DescriptionStyleModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "StyleView",
      "description_width": ""
     }
    },
    "cee19a8eb4454b238b640ea1f7669a2f": {
     "model_module": "@jupyter-widgets/base",
     "model_module_version": "1.2.0",
     "model_name": "LayoutModel",
     "state": {
      "_model_module": "@jupyter-widgets/base",
      "_model_module_version": "1.2.0",
      "_model_name": "LayoutModel",
      "_view_count": null,
      "_view_module": "@jupyter-widgets/base",
      "_view_module_version": "1.2.0",
      "_view_name": "LayoutView",
      "align_content": null,
      "align_items": null,
      "align_self": null,
      "border": null,
      "bottom": null,
      "display": null,
      "flex": null,
      "flex_flow": null,
      "grid_area": null,
      "grid_auto_columns": null,
      "grid_auto_flow": null,
      "grid_auto_rows": null,
      "grid_column": null,
      "grid_gap": null,
      "grid_row": null,
      "grid_template_areas": null,
      "grid_template_columns": null,
      "grid_template_rows": null,
      "height": null,
      "justify_content": null,
      "justify_items": null,
      "left": null,
      "margin": null,
      "max_height": null,
      "max_width": null,
      "min_height": null,
      "min_width": null,
      "object_fit": null,
      "object_position": null,
      "order": null,
      "overflow": null,
      "overflow_x": null,
      "overflow_y": null,
      "padding": null,
      "right": null,
      "top": null,
      "visibility": null,
      "width": null
     }
    }
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}