Matrix Calculus.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Derivatives with respect to matrices\n",
    "\n",
    "What is the gradient of a scalar function with respect to a matrix argument?\n",
    "\n",
    "Remember that the gradient of a scalar valued function with respect to a vector argument is a vector of the same size, similarly the derivative of a scalar with respect to a matrix is a matrix of the same size.\n",
    "\n",
    "The trace (sum of diagonal elements) appears often in scalar valued functions of matrices, so we start with its definition and an example:\n",
    "$\\newcommand{\\trace}{\\mathop{\\text{Tr}}}$\n",
    "\\begin{eqnarray}\n",
    "\\trace X = \\sum_i X(i,i)\n",
    "\\end{eqnarray}\n",
    "\n",
    "\\begin{eqnarray}\n",
    "f(X) & = & \\trace \\left( \\begin{array}{cc} X_{1,1} & X_{1,2} \\\\ X_{2,1} & X_{2,2} \\end{array} \\right) = X_{1,1} + X_{2,2} \\\\\n",
    "\\frac{df}{dX} & = & \\left( \\begin{array}{cc} \\partial{f}/\\partial X_{1,1} & \\partial{f}/\\partial X_{1,2} \\\\ \\partial{f}/\\partial X_{2,1} & \\partial{f}/\\partial X_{2,2} \\end{array} \\right)  =  \\left( \\begin{array}{cc} 1 & 0 \\\\ 0 & 1 \\end{array} \\right)\n",
    "\\end{eqnarray}\n",
    "The trace is like the inner product of the matrix $X$ (viewed as a vector) and the identity matrix (viewed as a vector). In the following, we will introduce a systematic way for computing these defivatives where\n",
    "we will make use of the Kronecker delta symbol that is defined as\n",
    "$$\n",
    "\\delta(i, j) = \\delta(j, i) = \\left\\{ \\begin{array}{cc} 1 & i = j \\\\0 & i\\neq j \\end{array} \\right.\n",
    "$$\n",
    "\n",
    "Before delving into the business of computing the derivatives, let us warm up a little. An important property of the trace is that the matrices can be rotated in a trace\n",
    "\\begin{eqnarray}\n",
    "\\trace A X = \\trace X A\n",
    "\\end{eqnarray}\n",
    "\n",
    "Using the Kronecker delta, the trace operator can be written as\n",
    "\\begin{eqnarray}\n",
    "\\trace X = \\sum_i \\sum_j X(i,j) \\delta(i, j)\n",
    "\\end{eqnarray}\n",
    "\n",
    "To prove this result, we will write the product explicitly using the index notation\n",
    "\\begin{eqnarray}\n",
    " (A X)(i,j) &=& \\sum_k A(i, k) X(k, j)\n",
    "\\end{eqnarray}\n",
    "\\begin{eqnarray}\n",
    "\\trace A X & = & \\sum_i \\sum_j \\left(\\sum_k A(i, k) X(k, j) \\right) \\delta(i, j) \\\\\n",
    "& = & \\sum_k \\sum_i  A(i, k) X(k, i)   \\\\\n",
    "& = & \\sum_k \\sum_i   X(k, i)  A(i, k) \\\\\n",
    "& = & \\sum_k \\sum_r \\sum_i   X(k, i)  A(i, r) \\delta(r,k) \\\\\n",
    "& = & \\sum_k \\sum_r (X A)(k, r) \\delta(r,k) = \\trace X A\n",
    "\\end{eqnarray}\n",
    "\n",
    "Let's now focus on the derivative of the trace of a product. We define\n",
    "\\begin{eqnarray}\n",
    "f(X) = \\trace A X & = & \\sum_i \\sum_j \\left(\\sum_k A(i, k) X(k, j) \\right) \\delta(i, j) \n",
    "\\end{eqnarray}\n",
    "where $A$ is $I \\times K$ and $X$ is $K \\times I$. Note that the result matrix in the trace must be a square matrix as otherwise the result is rectangular matrix and the trace does not make any sense.\n",
    "\n",
    "The gradient with respect to $X$ will be a object of the same size as $X$, with each entry defined by\n",
    "$$\n",
    "\\frac{d f}{d X} = [\\frac{\\partial f(X)}{\\partial X(u, r)}]\n",
    "$$\n",
    "Formally, by the chain rule we have \n",
    "\\begin{eqnarray}\n",
    "\\frac{\\partial f(X)}{\\partial X(u, r)} & = & \\frac{\\partial f(X)}{\\partial X(k, j)}\\frac{\\partial X(k, j)}{\\partial X(u, r)} = \\frac{\\partial f(X)}{\\partial X(k, j)} \\delta(u, k) \\delta(r, j)\n",
    "\\end{eqnarray}\n",
    "\n",
    "We can now sum over the indices to remove the Kronecker delta symbols\n",
    "\\begin{eqnarray}\n",
    "\\frac{\\partial f(X)}{\\partial X(u, r)} & = & \\sum_i \\sum_j \\sum_k A(i, k) \\delta(u, k) \\delta(r, j)  \\delta(i, j) \\\\\n",
    " & = & \\sum_i \\sum_j  A(i, u)  \\delta(r, j)  \\delta(i, j) \\\\\n",
    " & = & \\sum_j  A(j, u)  \\delta(r, j) =  A(r, u) = (A^\\top)(u, r) \\\\\n",
    "\\end{eqnarray}\n",
    "\n",
    "We have just shown that \n",
    "\\begin{eqnarray}\n",
    "\\frac{d f}{d X} = A^\\top\n",
    "\\end{eqnarray}\n",
    "\n",
    "\n",
    "## Example: Derivative of $\\trace A^\\top X A$ \n",
    "We will use the results above \n",
    "\\begin{eqnarray}\n",
    "f(X) & = &  \\trace A^\\top X A = \\trace A A^\\top X  \\\\\n",
    "\\frac{d f}{d X} &= & (A A^\\top)^\\top = A A^\\top \n",
    "\\end{eqnarray}\n",
    "\n",
    "Now, this quantity isnew.\n",
    "\\begin{eqnarray}\n",
    "f(X) & = &  \\trace A^\\top X A = \\trace A A^\\top X  \\\\\n",
    "\\frac{d f}{d A} &= & ?\n",
    "\\end{eqnarray}\n",
    "\n",
    "\\begin{eqnarray}\n",
    "f(A) = \\trace A^\\top X A & = & \\sum_i \\sum_j \\left(\\sum_k \\sum_l (A^\\top)(i, k) X(k, l) A(l, j) \\right) \\delta(i, j) \\\\\n",
    "\\frac{d f}{d A} &= &[\\frac{\\partial f(A)}{\\partial A(u, r)}] \\\\\n",
    "\\frac{\\partial f}{\\partial A(u,r)} & = & \\frac{\\partial}{\\partial A(u,r)} \\sum_i \\sum_j \\left(\\sum_k \\sum_l (A^\\top)(i, k) X(k, l) A(l, j) \\right) \\delta(i, j) \\\\\n",
    " & = & \\sum_i \\sum_j \\left(\\sum_k \\sum_l \\delta(k, u) \\delta(i, r) X(k, l) A(l, j) \\right) \\delta(i, j) \\\\\n",
    " & & + \\sum_i \\sum_j \\left(\\sum_k \\sum_l A(k, i) X(k, l) \\delta(l, u) \\delta(j, r) \\right) \\delta(i, j) \\\\\n",
    " & = & \\sum_i \\sum_j \\sum_k \\sum_l \\delta(k, u) \\delta(i, r) X(k, l) A(l, j)  \\delta(i, j) \\\\\n",
    " & & + \\sum_i \\sum_j \\sum_k \\sum_l A(k, i) X(k, l) \\delta(l, u) \\delta(j, r)  \\delta(i, j) \\\\\n",
    " & = & \\sum_i \\sum_j \\sum_l  \\delta(i, r) X(u, l) A(l, j)  \\delta(i, j) \\\\\n",
    " & & + \\sum_i \\sum_j \\sum_k  A(k, i) X(k, u) \\delta(j, r)  \\delta(i, j) \\\\\n",
    " & = & \\sum_j \\sum_l  X(u, l) A(l, j)  \\delta(r, j) + \\sum_j \\sum_k  A(k, j) X(k, u) \\delta(j, r)  \\\\\n",
    " & = & \\sum_l  X(u, l) A(l, r)  + \\sum_k  A(k, r) X(k, u)  \\\\\n",
    " & = & \\sum_l  X(u, l) A(l, r)  + \\sum_k  (X^\\top)(u,k) A(k, r)   \\\\\n",
    " & = & (X A)(u, r)   +  (X^\\top A)(u,r)    \\\\\n",
    " & = & (X+X^\\top) A\n",
    "\\end{eqnarray}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example: Derivative of $\\trace A X A^\\top$ \n",
    "This looks quite similar\n",
    "\n",
    "\\begin{eqnarray}\n",
    "f(A) & = & \\trace A X A^\\top  =  \\sum_i \\sum_j \\left(\\sum_k \\sum_l A (i, k) X(k, l) (A^\\top)(l, j) \\right) \\delta(i, j)  \\\\\n",
    "& = & \\sum_i \\sum_j \\sum_k \\sum_l A (i, k) X(k, l) A(j, l) \\delta(i, j)  \\\\\n",
    "\\frac{\\partial f}{\\partial A(u,r)} & = & \\sum_i \\sum_j \\sum_k \\sum_l \\delta(u,i) \\delta(r,k) X(k, l) A(j, l) \\delta(i, j) + \\sum_i \\sum_j \\sum_k \\sum_l A (i, k) X(k, l) \\delta(u,j) \\delta(r,l) \\delta(i, j) \\\\\n",
    "& = &  \\sum_l  X(r, l) A(u, l)  +  \\sum_k  A (u, k) X(k, r)   \\\\\n",
    "& = & \\sum_l  A(u, l) (X^\\top)(l, r)  +  \\sum_k  A (u, k) X(k, r)   \\\\\n",
    "& = & (A X^\\top)(u,r) + (A X)(u, r) \n",
    "\\end{eqnarray}\n",
    "\n",
    "So the derivative is\n",
    "$$\n",
    "\\frac{d f(A)}{d A} = A (X^\\top + X)\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Matrix Factorization by Alternating Least Squares (ALS)\n",
    "\n",
    "We define the Frobenious norm\n",
    "$$\n",
    "\\|E\\|_F = \\sqrt{\\trace E^\\top E}\n",
    "$$\n",
    "\n",
    "By using the results that we have just derived we obtain we can derive the derivative of the \n",
    "Frobenious norm of the error matrix $E$ as\n",
    "\n",
    "\\begin{eqnarray}\n",
    "E^\\top E &= & (X - M C)^\\top   (X - M C) \\\\\n",
    "& = & X^\\top X + C^\\top M^\\top M C - X^\\top M C - C^\\top M^\\top X \\\\\n",
    "\\|E\\|_F^2 & = & \\trace \\left( X^\\top X + C^\\top M^\\top M C - 2 X^\\top M C \\right)\n",
    "\\end{eqnarray}\n",
    "\n",
    "By the linearity of the trace\n",
    "\\begin{eqnarray}\n",
    "\\|E\\|_F^2 & = & \\trace \\left( X^\\top X + C^\\top M^\\top M C - 2 X^\\top M C \\right) \\\\\n",
    "\\frac{d\\|E\\|_F^2}{d C} & = & 2 M^\\top M C - 2 M^\\top X \\\\\n",
    "C & = &  (M^\\top M)^{-1} M^\\top X\n",
    "\\end{eqnarray} \n",
    "\n",
    "\\begin{eqnarray}\n",
    "\\|E\\|_F^2 & = & \\trace X^\\top X + \\trace C^\\top M^\\top M C - 2 \\trace X^\\top M C \\\\\n",
    " & = & \\trace X^\\top X + \\trace  M C C^\\top M^\\top  - 2 \\trace C X^\\top M  \\\\ \n",
    "\\frac{d\\|E\\|_F^2}{d M} & = & 2 M C C^\\top - 2 X C^\\top \\\\\n",
    "M & = & X C^\\top (C C^\\top)^{-1}\n",
    "\\end{eqnarray} \n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$B = (A^\\top A)^{-1} A^\\top X$\n",
    "\n",
    "$A = X B^\\top (B B^\\top)^{-1}$"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [default]",
   "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.6.5"
  },
  "toc": {
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": "block",
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}