Chapter 18 - The Singular Value Decomposition.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy as sp\n",
    "import scipy.linalg\n",
    "import sympy as sy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.set_printoptions(precision=3)\n",
    "np.set_printoptions(suppress=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# <font face=\"gotham\" color=\"purple\"> The Singular Values </font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Technically only square matrices have eigenvalues and vectors, however we want to extend a similar concept for any $m \\times n$ matrices. \n",
    "\n",
    "If $A$ is an $m \\times n$ matrix, then $A^TA$ and $AA^T$ are both symmetric and orthogonally diagonalized."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 4, -3],\n",
       "       [-3,  2],\n",
       "       [ 2, -1],\n",
       "       [ 1,  0]])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = np.array([[4, -3], [-3, 2], [2, -1], [1, 0]]); A"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute eigenvalues and eigenvectors of $A^TA$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 30, -20],\n",
       "       [-20,  14]])"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ATA = A.T@A; ATA"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [],
   "source": [
    "D1, P1 = np.linalg.eig(ATA)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check if $P$ is an orthonormal matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1., 0.],\n",
       "       [0., 1.]])"
      ]
     },
     "execution_count": 57,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P1@P1.T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([43.541,  0.459])"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "D1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The square roots of eigenvalues of $A^TA$ are called <font face=\"gotham\" color=\"red\"> singular values</font>  of $A$, denoted by $\\sigma_1, ..., \\sigma_n$ in decreasing order. \n",
    "\n",
    "We can also show that singular values of $A$ are the lengths of vectors $A\\mathbf{v}_1,..., A\\mathbf{v}_n$, where $\\mathbf{v}_i$ is the eigenvalue of $A^TA$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The length of $A\\mathbf{v}_i$ is $\\|A\\mathbf{v}_i\\|$\n",
    "\n",
    "$$\n",
    "\\|A\\mathbf{v}_i\\| = \\sqrt{(A\\mathbf{v}_i)^TA\\mathbf{v}_i} = \\sqrt{\\mathbf{v}_i^TA^T A\\mathbf{v}_i}=\\sqrt{\\mathbf{v}_i^T(\\lambda_i\\mathbf{v}_i)} = \\sqrt{\\lambda_i}=\\sigma_1\n",
    "$$\n",
    "\n",
    "where $\\sqrt{\\mathbf{v}_i^T\\mathbf{v}_i} = 1$ and $\\lambda_i$'s are eigenvalues of $A^TA$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    " # <font face=\"gotham\" color=\"purple\"> Singular Value Decomposition</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$A$ is a $m\\times n$ matrix. However $AA^T$ and $A^TA$ are symmetric matrices,then both are orthogonally diagonalizable.\n",
    "\n",
    "$$\n",
    "AA^T = U\\Sigma\\Sigma^T U^T=(U\\Sigma V^T)(V\\Sigma U^T)\\\\\n",
    "A^TA = V\\Sigma^T \\Sigma V^T = (V\\Sigma^T U^T)(U \\Sigma V^T)\n",
    "$$\n",
    "\n",
    "where $\\Sigma\\Sigma^T$ is a diagonal matrix with all eigenvalues of $AA^T$ and $\\Sigma^T \\Sigma$ is a diagonal matrix with all eigenvalues of $VV^T$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because both $AA^T$ and $A^TA$ are symmetric, then $UU^T= U^TU=I_{m\\times m}$ and $VV^T= V^TV=I_{n\\times n}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have implicitly shown the singular value decompositions above, one of the most important concept in linear algebra.\n",
    "\n",
    "<font face=\"gotham\" color=\"red\">\n",
    "$$\n",
    "\\Large\n",
    "SVD:\\quad\n",
    "A_{m\\times n} = U_{m\\times m}\\Sigma_{m \\times n} V^T_{n \\times n}\n",
    "$$\n",
    "</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So next question: what is $\\Sigma$? \n",
    "\n",
    "It is an $m\\times n$ main diagonal matrix, with all singular values on the main diagonal. Rewrite"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "A^TA = V\\Sigma^T \\Sigma V^T = V\\Sigma^2 V^T\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Post-multiply both sides by $V$\n",
    "\n",
    "$$\n",
    "A^TAV = V\\Sigma^2\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is the matrix version of $A\\mathbf{v}_i = \\lambda_i \\mathbf{v}_i$, but here the matrix of interest is $A^TA$ rather than $A$. Similarly it can be written with singular values\n",
    "\n",
    "$$\n",
    "A^TA\\mathbf{v}_i = \\sigma_i^2\\mathbf{v}_i\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because $U$ and $V$ are not unique, we tend to standardize the solution by arranging $\\sigma_1 \\geq \\sigma_2 \\geq \\sigma_3\\geq ... \\geq\\sigma_r$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Why we only arrange $r$ singular values? Because it is the rank of $A$, so is the rank of $A^TA$. Explicitly $\\Sigma$ looks like\n",
    "\n",
    "\n",
    "$$\\Sigma =\\left[\\begin{array}{cccccc}\n",
    "\\sqrt{\\lambda_{1}} & & & & &\\\\\n",
    "& \\sqrt{\\lambda_{2}} & & & &\\\\\n",
    "& & \\ddots & & &\\\\\n",
    "& & & \\sqrt{\\lambda}_{\\mathrm{r}} & &\\\\\n",
    "& & & & \\ddots &\\\\\n",
    "& & & & & 0 \n",
    "\\end{array}\\right]\n",
    "=\\left[\\begin{array}{cccccc}\n",
    "\\sigma_1 & & & & &\\\\\n",
    "& \\sigma_2 & & & &\\\\\n",
    "& & \\ddots & & &\\\\\n",
    "& & & \\sigma_r & &\\\\\n",
    "& & & & \\ddots &\\\\\n",
    "& & & & & 0 \n",
    "\\end{array}\\right]\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can do the same for $AA^T$ and get\n",
    "\n",
    "$$\n",
    "AA^TU = U \\Sigma^2\n",
    "$$\n",
    "\n",
    "or \n",
    "\n",
    "$$\n",
    "AA^T\\mathbf{u}_i = \\sigma_i^2\\mathbf{u}_i\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have shown why $A_{m\\times n} = U_{m\\times m}\\Sigma_{m \\times n} V^T_{n \\times n}$ holds."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To perfomr a SVD on $A$, we just need two equations and this is also a mannual procedure to decompose any matrix.\n",
    "\n",
    "$$\n",
    "A^TA = V\\Sigma^T \\Sigma V^T\\\\\n",
    "AV = U\\Sigma\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## <font face=\"gotham\" color=\"purple\"> Reformulate SVD</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rewrite $SVD$\n",
    "\n",
    "$$\n",
    "AV = U\\Sigma\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "vector version is \n",
    "\n",
    "$$\n",
    "A\\mathbf{v}_i = \\sigma_i \\mathbf{u}_i\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There two implications from the equation above: $(a)$ $A$ can be decomposed into\n",
    "\n",
    "$$\n",
    "A = \\sum_{i}^r\\sigma_i\\mathbf{u}_i \\mathbf{v}_i^T\n",
    "$$\n",
    "\n",
    "$(b)$ We can compute $\\mathbf{u}_i$ by using \n",
    "\n",
    "$$\\mathbf{u}_i = \\frac{A\\mathbf{v}_i}{\\sigma_i}$$"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}