chapter16_tensor_methods/tensor_basics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Tensors are higher order extensions of matrices that can encode multi-dimensional data\n",
    "\n",
    "![tensor_illustration](../img/tensor_cartoon.jpg)\n",
    "\n",
    "In this tutorial we will show how to manipulate tensors as NDArrays, and write from scratch functions to manipulate these as defined in [TensorLy](http://tensorly.github.io)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import mxnet.ndarray as nd"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. Creating a Tensor\n",
    "\n",
    "A tensor can be represented in multiple ways. The simplest is the slice representation through multiple matrices.\n",
    "\n",
    "Let's take for this example the tensor $\\tilde X$ defined by its frontal slices:\n",
    "\n",
    "$$\n",
    "   X_1 = \n",
    "   \\left[\n",
    "   \\begin{matrix}\n",
    "   0  & 2  & 4  & 6\\\\\n",
    "   8  & 10 & 12 & 14\\\\\n",
    "   16 & 18 & 20 & 22\n",
    "   \\end{matrix}\n",
    "   \\right]\n",
    "$$\n",
    "\n",
    "and \n",
    "\n",
    "$$\n",
    "   X_2 =\n",
    "   \\left[\n",
    "   \\begin{matrix}\n",
    "   1  & 3  & 5  & 7\\\\\n",
    "   9  & 11 & 13 & 15\\\\\n",
    "   17 & 19 & 21 & 23\n",
    "   \\end{matrix}\n",
    "   \\right]\n",
    "$$\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "In Python, this array can be expressed as a numpy array::"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X = nd.arange(24).reshape((3, 4, 2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[[  0.   1.]\n",
       "  [  2.   3.]\n",
       "  [  4.   5.]\n",
       "  [  6.   7.]]\n",
       "\n",
       " [[  8.   9.]\n",
       "  [ 10.  11.]\n",
       "  [ 12.  13.]\n",
       "  [ 14.  15.]]\n",
       "\n",
       " [[ 16.  17.]\n",
       "  [ 18.  19.]\n",
       "  [ 20.  21.]\n",
       "  [ 22.  23.]]]\n",
       "<NDArray 3x4x2 @cpu(0)>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can view the frontal slices by fixing the last axis:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[  0.   2.   4.   6.]\n",
       " [  8.  10.  12.  14.]\n",
       " [ 16.  18.  20.  22.]]\n",
       "<NDArray 3x4 @cpu(0)>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X[:, :, 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[  1.   3.   5.   7.]\n",
       " [  9.  11.  13.  15.]\n",
       " [ 17.  19.  21.  23.]]\n",
       "<NDArray 3x4 @cpu(0)>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X[:, :, 1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3. Basic Tensor Operations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3.1 Unfolding\n",
    "\n",
    "Also called **matrization**, **unfolding** a tensor is done by reading the element in a given way as to obtain a matrix instead of a tensor.\n",
    "\n",
    "It is done by stacking the **fibers** of the tensor into a matrix.\n",
    "\n",
    "![tensor_illustration](../img/tensor_fibers.png)\n",
    "Illustration: *Nonnegative Matrix and Tensor Factorizations*, Andrzej Cichocki, Rafal Zdunek, Anh Huy Phan, and Shun-ichi Amari, John Wiley & Sons, 2009.\n",
    "\n",
    "\n",
    "\n",
    "### Definition\n",
    "For a tensor of size $(I_1, I_2, \\cdots, I_N)$, the n-mode unfolding of this tensor will be of size $(I_n, I_1 \\times \\cdots \\times I_{n-1} \\times I_{n+1} \\cdots \\times I_N)$ and is obtained by reading the tensor as a matrix with the $n$-th dimension first. \n",
    "\n",
    "\n",
    "Specifically, given a tensor $\\tilde X \\in \\mathbb{R}^{I_1, I_2, \\cdots, I_N}$, the\n",
    "mode-n unfolding of $\\tilde X$ is a matrix $\\mathbf{X}_{[n]} \\in \\mathbb{R}^{I_n, I_M}$,\n",
    "with $M = \\prod\\limits_{\\substack{k=1,\\\\k \\neq n}}^N I_k$ and is defined by\n",
    "the mapping from element $(i_1, i_2, \\cdots, i_N)$ to $(i_n, j)$, with\n",
    "\n",
    "\n",
    "$$\n",
    "    j = \\sum\\limits_{\\substack{k=1,\\\\k \\neq n}}^N i_k \\times \\prod\\limits_{\\substack{m=k+1,\\\\m \\neq n}}^N I_m.\n",
    "$$\n",
    "\n",
    "### Convention\n",
    "\n",
    "   Traditionally, mode-1 unfolding denotes the unfolding along the first dimension.\n",
    "   However, to be consistent with the Python indexing that always starts at zero,\n",
    "   as done in tensorly, we will start indexing modes at zero!\n",
    "\n",
    "   Therefore ``unfold(tensor, 0)`` will unfold said tensor along its first dimension!\n",
    "   \n",
    "\n",
    "### Example\n",
    "\n",
    "For instance, using the $\\tilde X$ previously defined:\n",
    "$$\n",
    "   X_1 = \n",
    "   \\left[\n",
    "   \\begin{matrix}\n",
    "   0  & 2  & 4  & 6\\\\\n",
    "   8  & 10 & 12 & 14\\\\\n",
    "   16 & 18 & 20 & 22\n",
    "   \\end{matrix}\n",
    "   \\right]\n",
    "$$\n",
    "\n",
    "and \n",
    "\n",
    "$$\n",
    "   X_2 =\n",
    "   \\left[\n",
    "   \\begin{matrix}\n",
    "   1  & 3  & 5  & 7\\\\\n",
    "   9  & 11 & 13 & 15\\\\\n",
    "   17 & 19 & 21 & 23\n",
    "   \\end{matrix}\n",
    "   \\right]\n",
    "$$\n",
    "\n",
    "The 0-mode unfolding of $\\tilde X$:\n",
    "\n",
    "$$\n",
    "   \\tilde X_{[0]} =\n",
    "   \\left[ \\begin{matrix}\n",
    "      0 & 1 & 2 & 3 & 4 & 5 & 6 & 7\\\\\n",
    "      8 & 9 & 10 & 11 & 12 & 13 & 14 & 15\\\\\n",
    "      16 & 17 & 18 & 19 & 20 & 21 & 22 & 23\\\\\n",
    "   \\end{matrix} \\right]\n",
    "$$\n",
    "\n",
    "The 1-mode unfolding is given by:\n",
    "\n",
    "$$\n",
    "   \\tilde X_{[1]} =\n",
    "   \\left[ \\begin{matrix}\n",
    "      0 & 1 & 8 & 9 & 16 & 17\\\\\n",
    "      2 & 3 & 10 & 11 & 18 & 19\\\\\n",
    "      4 & 5 & 12 & 13 & 20 & 21\\\\\n",
    "      6 & 7 & 14 & 15 & 22 & 23\\\\\n",
    "   \\end{matrix} \\right]\n",
    "$$\n",
    "\n",
    "Finally, the 2-mode unfolding is the unfolding along the last axis:\n",
    "\n",
    "$$\n",
    "    \\tilde X_{[2]} =\n",
    "   \\left[ \\begin{matrix}\n",
    "      0 & 2 & 4 & 6 & 8 & 10 & 12 & 14 & 16 & 18 & 20 & 22\\\\\n",
    "      1 & 3 & 5 & 7 & 9 & 11 & 13 & 15 & 17 & 19 & 21 & 23\\\\\n",
    "   \\end{matrix} \\right]\n",
    "$$\n",
    "\n",
    "### In MXNet\n",
    "\n",
    "Let's define the unfolding function in MXNet. Given the mode $n$ along which to unfold, it will take a tensor, put the $n$-th dimension first, and matricize the result. Note that our definition of unfolding corresponds to a C-ordering of the elements. MXNet also has a C ordering of the elements, making that matricization a simple reshaping."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "def unfold(tensor, mode):\n",
    "    \"\"\"Returns the mode-`mode` unfolding of `tensor` with modes starting at `0`.\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    tensor : ndarray\n",
    "    mode : int, default is 0\n",
    "           indexing starts at 0, therefore mode is in ``range(0, tensor.ndim)``\n",
    "    \n",
    "    Returns\n",
    "    -------\n",
    "    ndarray\n",
    "        unfolded_tensor of shape ``(tensor.shape[mode], -1)``\n",
    "    \"\"\"\n",
    "    return nd.reshape(nd.moveaxis(tensor, mode, 0), (tensor.shape[mode], -1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[  0.   1.   2.   3.   4.   5.   6.   7.]\n",
       " [  8.   9.  10.  11.  12.  13.  14.  15.]\n",
       " [ 16.  17.  18.  19.  20.  21.  22.  23.]]\n",
       "<NDArray 3x8 @cpu(0)>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "unfold(X, mode=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[  0.   1.   8.   9.  16.  17.]\n",
       " [  2.   3.  10.  11.  18.  19.]\n",
       " [  4.   5.  12.  13.  20.  21.]\n",
       " [  6.   7.  14.  15.  22.  23.]]\n",
       "<NDArray 4x6 @cpu(0)>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "unfold(X, mode=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[  0.   2.   4.   6.   8.  10.  12.  14.  16.  18.  20.  22.]\n",
       " [  1.   3.   5.   7.   9.  11.  13.  15.  17.  19.  21.  23.]]\n",
       "<NDArray 2x12 @cpu(0)>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "unfold(X, mode=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3.2 Folding\n",
    "\n",
    "Folding is the inverse operation: we reshape the matrix into a tensor and move back the first dimension back to its original place."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def fold(unfolded_tensor, mode, shape):\n",
    "    \"\"\"Refolds the mode-`mode` unfolding into a tensor of shape `shape`\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    unfolded_tensor : ndarray\n",
    "        unfolded tensor of shape ``(shape[mode], -1)``\n",
    "    mode : int\n",
    "        the mode of the unfolding\n",
    "    shape : tuple\n",
    "        shape of the original tensor before unfolding\n",
    "    \n",
    "    Returns\n",
    "    -------\n",
    "    ndarray\n",
    "        folded_tensor of shape `shape`\n",
    "    \"\"\"\n",
    "    full_shape = list(shape)\n",
    "    mode_dim = full_shape.pop(mode)\n",
    "    full_shape.insert(0, mode_dim)\n",
    "    return nd.moveaxis(nd.reshape(unfolded_tensor, full_shape), 0, mode)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "[[[  0.   1.]\n",
       "  [  2.   3.]\n",
       "  [  4.   5.]\n",
       "  [  6.   7.]]\n",
       "\n",
       " [[  8.   9.]\n",
       "  [ 10.  11.]\n",
       "  [ 12.  13.]\n",
       "  [ 14.  15.]]\n",
       "\n",
       " [[ 16.  17.]\n",
       "  [ 18.  19.]\n",
       "  [ 20.  21.]\n",
       "  [ 22.  23.]]]\n",
       "<NDArray 3x4x2 @cpu(0)>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "unfolding = unfold(X, 1)\n",
    "original_shape = X.shape\n",
    "fold(unfolding, mode=1, shape=original_shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3.3 n-mode product\n",
    "\n",
    "Also known as **tensor contraction**. This is a natural generalization of matrix-vector and matrix-matrix product. When multiplying a tensor by a matrix or a vector, we now have to specify the **mode** $n$ along which to take the product.\n",
    "\n",
    "### Tensor times matrix\n",
    "\n",
    "In that case we are doing an operation analogous to a matrix multiplication on the $n$-th mode. Given a tensor $\\tilde X$ of size $(I_1, I_2, \\cdots, I_N)$, and a matrix $M$ of size $(D, I_n)$, the $n$-mode product of $\\tilde X$ by $M$ is written $\\tilde X \\times_n M$ and is of size $(D, I_1 \\times \\cdots \\times I_{n-1} \\times I_{n+1} \\cdots \\times I_N)$.\n",
    "\n",
    "\n",
    "One simple way to mathematically define the n-mode product is using the unfolding: if we write $\\tilde R = \\tilde X \\times_n M$, then we have:\n",
    "\n",
    "$$\n",
    "    \\tilde R_{[n]} = M \\times \\tilde X_{[n]}\n",
    "$$\n",
    "\n",
    "As a consequence, to get the n-mode product of $\\tilde X$ by $M$, we can simply take a matrix product between $M$ and the unfolding of $\\tilde X$ along the $n^{th}$ dimension, and refold the result into a tensor of shape $(I_1, \\cdots, I_{n-1}, D, I_{n+1}, \\cdots, I_N)$.\n",
    "\n",
    "### Tensor times vector\n",
    "\n",
    "In that case we are contracting over the $n$-th mode by multiplying it with a vector. Given a tensor $\\tilde X$ of size $(I_1, I_2, \\cdots, I_N)$, and a vector $v$ of size $(I_n)$, the $n$-mode product of $\\tilde X$ by $v$ is written $\\tilde X \\times_n v$ and is of size $(I_1 \\times \\cdots \\times I_{n-1} \\times I_{n+1} \\cdots \\times I_N)$ --we have essentially summed over (or contracted over) the $n$-th dimension--.\n",
    "\n",
    "![tensor_illustration](../img/tensor_contraction.png)\n",
    "\n",
    "\n",
    "### Example\n",
    "\n",
    "We will write a function `mode_dot` that works transparently for multiplying a tensor by a matrix or a vector, along a given mode."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def mode_dot(tensor, matrix_or_vector, mode):\n",
    "    \"\"\"n-mode product of a tensor by a matrix at the specified mode.\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    tensor : ndarray\n",
    "        tensor of shape ``(i_1, ..., i_k, ..., i_N)``\n",
    "    matrix_or_vector : ndarray\n",
    "        1D or 2D array of shape ``(J, i_k)`` or ``(i_k, )``\n",
    "        matrix or vectors to which to n-mode multiply the tensor\n",
    "    mode : int\n",
    "\n",
    "    Returns\n",
    "    -------\n",
    "    ndarray\n",
    "        `mode`-mode product of `tensor` by `matrix_or_vector`\n",
    "        * of shape :math:`(i_1, ..., i_{k-1}, J, i_{k+1}, ..., i_N)` if matrix_or_vector is a matrix\n",
    "        * of shape :math:`(i_1, ..., i_{k-1}, i_{k+1}, ..., i_N)` if matrix_or_vector is a vector\n",
    "    \"\"\"\n",
    "    # the mode along which to fold might decrease if we take product with a vector\n",
    "    fold_mode = mode\n",
    "    new_shape = list(tensor.shape)\n",
    "\n",
    "    # tensor times vector case: make sure the sizes are correct \n",
    "    # (we are contracting over one dimension which then disappearas)\n",
    "    if matrix_or_vector.ndim == 1: \n",
    "        if len(new_shape) > 1:\n",
    "            new_shape.pop(mode)\n",
    "            fold_mode -= 1\n",
    "        else:\n",
    "            new_shape = [1]\n",
    "\n",
    "    # This is the actual operation: we use the equivalent formulation of the n-mode-product using the unfolding\n",
    "    res = nd.dot(matrix_or_vector, unfold(tensor, mode))\n",
    "\n",
    "    # refold the result into a tensor and return it \n",
    "    return fold(res, fold_mode, new_shape)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Tensor times matrix\n",
    "\n",
    "With the tensor $\\tilde X$ of size (3, 4, 2) we defined previously, let's define a matrix M of size (5, 4) to multiply along the second mode:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(5, 4)\n"
     ]
    }
   ],
   "source": [
    "M = nd.arange(4*5).reshape((5, 4))\n",
    "print(M.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Keep in mind indexing starts at zero, so the second mode is represented by `mode=1`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "res = mode_dot(X, M, mode=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As expected the result is of shape (3, 5, 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(3, 5, 2)"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "res.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Tensor times vector\n",
    "\n",
    "Similarly, we can contract along mode 1 with a vector of size 4 (our tensor is of size (3, 4, 2).\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(4,)\n"
     ]
    }
   ],
   "source": [
    "v = nd.arange(4)\n",
    "print(v.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "res = mode_dot(X, v, mode=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since we have multiplied by a vector, we have effectively contracted out one mode of the tensor so the result is a matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(3, 2)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "res.shape"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}