numpy_array_basics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Basics of NumPy arrays"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's quickly go over NumPy arrays, which are \"objects\" defined in the NumPy package that function extremely similar to Matlab's matrices. Any defined array comes with a series of \"attributes\" (i.e., values nested inside the object) and \"methods\" (i.e., functions nested inside the object) that are built into that object. All of the attributes and methods associated with a NumPy array for version 1.13 can be found [here](https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.ndarray.html).\n",
    "\n",
    "So, let's mess around with this for a bit."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preliminaries"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In order to use any NumPy functions, we first need to import the package (assuming you already have it installed using conda/pip/some package manager). It's traditional to import NumPy as \"np\" so we don't have to type as much later. Also, note that package imports are case-sensitive."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Creating and slicing NumPy arrays"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Unlike Matlab, arrays must be initialized before their values can be assigned. You can specify the datatype when you initialize (e.g., float, integer, string, complex, etc. -- default is float)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = np.empty((5, 2), dtype=float)  # without this line, the following lines would break!\n",
    "for i in range(5):\n",
    "    a[i, 0] = i ** 2  # assign square(i) to col 0\n",
    "    a[i, 1] = 2 * i  # assign double(i) to col 1\n",
    "a  # jupyter will print the last line if it's not assigned to anything"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Array slicing is extremely similar to Matlab, except for a few things:\n",
    "1. Slices are done with square brackets (curved brackets are only for function inputs)\n",
    "2. Indexing starts from 0\n",
    "3. The end index you specifiy is not included in the slice\n",
    "4. There is no `end`, but you can just omit the last index or use negative indexing\n",
    "5. Omitting the beginning index implies starting at 0; omitting the end implies ending at the end."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(a[0:3, 0])  # print elements 0, 1, and 2 in first column\n",
    "print(a[:3, 0])  # same as above, but we've omitted the zero in the indexing this time\n",
    "print(a[-2:, 1])  # print last two elements in second column\n",
    "print(a[:, 0])  # print all elements in first column\n",
    "print(a[[0, 1, 3], 1])  # print 1st, 2nd, and 4th elements in first column (uses broadcasting, more on that later)\n",
    "print(a[[0, 1, 3], [0, 1, 0]])  # print elements [0,0], [1,1], and [3,0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can create arrays with as many dimensions as you like"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = np.random.rand(10, 3, 20, 4)  # dim1 = 10, dim2 = 3, etc.\n",
    "print(b.shape)  # 'shape' is an attribute, it returns a tuple"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And you can use a lot of the methods \"built into\" the array. Remember, these are documented [here](https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.ndarray.html)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(a.mean())  # mean of entire matrix\n",
    "print(a.max(axis=0))  # max of each column\n",
    "print(a.min(axis=1))  # min of each row\n",
    "print(a.argmax())  # this is applied to the raveled matrix, so it returns an integer\n",
    "print(a.ravel()[a.argmax()])  # ...so we can either apply the index to the flattened matrix\n",
    "print(np.unravel_index(a.argmax(), a.shape))  # ...or convert the integer index to a tuple"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wait...what do we mean by a flattened array? And what does `ravel` do?\n",
    "\n",
    "In short, `flatten` will return a copy of the array, but `ravel` (and `reshape`) will not. This is important, because if you write something to the output of `flatten`, you will not update the original array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(a.ravel(), '...line 1')  # this is just a (contiguous) view of the array, it does not return a copy\n",
    "print(a.flatten(), '...line 2')  # this returns a copy of the array\n",
    "print(a.reshape((-1)), '...line 3')  # this will return a (non-contiguous) view of the array\n",
    "a.flatten()[8] = np.nan  # assign a nan to the output of flatten\n",
    "print(a.ravel(), '...line 4')  # ...but a is unchanged!\n",
    "a.ravel()[8] = np.nan  # assign a nan to the output of ravel\n",
    "print(a.ravel(), '...line 5')  # ...and now a is changed!\n",
    "a.ravel()[8] = 16.0  # assign the correct value back to a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Elementwise multiplication and array broadcasting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Matrix multiplication is one of the most common things we do in our analyses, but if you're not careful you can be tripped up by NumPy's matrix operations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "my_mat = np.array([[0, 1], [2, 3]])\n",
    "v1 = np.arange(2)  # 1d array of [0, 1]\n",
    "v2 = np.arange(2).reshape((-1, 1))  # v1 but as a column vector (note the '-1' indicates to infer the size along that axis)\n",
    "v3 = np.arange(2).reshape((1, -1))  # v1 but as a row vector"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Print their shapes just for clarity..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(my_mat.shape)\n",
    "print(v1.shape)  # note that the shape tuple for v1 has only one number, whereas v2 and v3 have two numbers\n",
    "print(v2.shape)\n",
    "print(v3.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What happens if we multiply things togeter? $\\mathrm{my\\_mat} \\cdot v_3$ should break, right?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(my_mat * v1)\n",
    "print(my_mat * v2)\n",
    "print(my_mat * v3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wait...this doesn't make sense at all! How are all these output matrices $2\\times2$?\n",
    "\n",
    "**Answer**: NumPy generally tries to do everything element-wise. So, `*` is element-wise multiplication.\n",
    "\n",
    "NumPy handles the mismatch in dimensions by \"broadcasting\" arrays into the other dimensions to try and make dimensions match so it can do things element-wise.\n",
    "\n",
    "**Broadcasting is done in two ways:**\n",
    "1. If a shape tuple has too few elements for an operation (e.g., `(2,)` vs. `(2,2)`), then NumPy will prepend \"1\" to the shape tuple until it has the right number of elements.\n",
    "2. If there are the right number of elements in the shape tuple, but 1) the shapes don't match and 2) one of the elements in the shape tuple is a 1, NumPy tries to repeat the array until it matches the size of the other array.\n",
    "\n",
    "If NumPy cannot get the dimensions to match using the above two steps, then it will throw an error."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's now explain what's going on in our example a little more clearly.\n",
    "\n",
    "#### Broadcasting of v1 array\n",
    "- v1 starts out as the following arrray: `[0, 1]` (shape `(2,)`)\n",
    "- A 1 is preprended to the shape tuple so it has the right number of elements: `[[0, 1]]` (shape `(1, 2)`)\n",
    "- The array is then tiled along the 0 axis (down rows) until its shape matches `my_mat`: `[[0, 1], [0, 1]]` (shape `(2, 2)`)\n",
    "- This broadcasted array is then multiplied by `my_mat` element-wise\n",
    "\n",
    "#### Broadcasting of v2 array\n",
    "- v2 starts out as the following arrray: `[[0], [1]]` (shape `(2, 1)`)\n",
    "- The array is then tiled along the 1 axis (along columns) until its shape matches `my_mat`: `[[0, 0], [1, 1]]` (shape `(2, 2)`)\n",
    "- This broadcasted array is then multiplied by `my_mat` element-wise\n",
    "\n",
    "#### Broadcasting of v3 array\n",
    "- v3 starts out as the following arrray: `[[0, 1]]` (shape `(1, 2)`)\n",
    "- The array is then tiled along the 0 axis (down rows) until its shape matches `my_mat`: `[[0, 1], [0, 1]]` (shape `(2, 2)`)\n",
    "- This broadcasted array is then multiplied by `my_mat` element-wise"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Matrix multiplication"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*I'm not sure how that's useful to know*, you are thinking. *I normally do matrix multiplication*. Well, let's figure that out, eh?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(my_mat @ v1)  # note this returns a 1D array\n",
    "print(my_mat @ v2)  # ...while this returns a 2D array\n",
    "print(np.dot(my_mat, v2))  # you can also use the numpy function\n",
    "# print(my_mat @ v3)  # this breaks with a 'dimension mismatch' error, just as we would expect"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Well. That was easy."
   ]
  }
 ],
 "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.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}