{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Introduction to Numpy\n", "## Adapted for M1 GGL Students (with BSc level in linear algebra)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**NumPy** is the fundamental package for scientific computing\n", "in Python. It is a Python library that provides a multidimensional array\n", "object. In this course, we will be using **NumPy** for linear algebra.\n", "\n", "If you are interested in learning more about **NumPy**, you can find the user\n", "guide and reference at the [documentation](https://docs.scipy.org/doc/numpy/index.html).\n", "\n", "Feel free to alter the code blocks of this notebook in case you want to try something on your own." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's first import the NumPy package" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np # we commonly use the np abbreviation when referring to numpy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Creating Numpy Arrays" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A Numpy array is similar to a Python list. However, it can only contain elements of the same type." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "New arrays can be made in several ways. We can take an existing list and convert it to a numpy array:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1 2 3]\n" ] } ], "source": [ "a = np.array([1,2,3])\n", "print(a)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1 2 3]\n", " [4 5 6]\n", " [7 8 9]]\n" ] } ], "source": [ "b = np.array([[1, 2, 3],\n", " [4, 5, 6],\n", " [7, 8, 9]])\n", "print(b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that for 2-dimensional Numpy arrays, all rows must have the same number of elements." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Shape of Numpy Arrays" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use [np.shape](https://numpy.org/doc/stable/reference/generated/numpy.shape.html) to get the shape of a Numpy array." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(3,)\n" ] } ], "source": [ "# This is a one dimensional array\n", "a = np.array([1, 2, 3])\n", "print(a.shape)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(3, 2)\n" ] } ], "source": [ "# This is a two dimensional array\n", "b = np.array([[1, 2],\n", " [3, 4],\n", " [5, 6]])\n", "print(b.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note the difference in shape between one dimensional and two dimensional arrays.\n", "\n", "For two dimensional arrays, you can consider the first element of the tuple to be the number of rows and the second element to be the number of columns." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Creating Numpy Arrays Based on Shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use the functions [np.zeros()](https://numpy.org/doc/stable/reference/generated/numpy.zeros.html) and [np.ones()](https://numpy.org/doc/stable/reference/generated/numpy.zeros.html) to create Numpy arrays of a specified shape. The input of the function is the shape of the Numpy array which you want to get." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[0. 0. 0. 0. 0.]\n" ] } ], "source": [ "a = np.zeros((5,))\n", "print(a)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0. 0. 0.]\n", " [0. 0. 0.]]\n" ] } ], "source": [ "b = np.zeros((2, 3))\n", "print(b)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1.]\n", " [1. 1.]\n", " [1. 1.]]\n" ] } ], "source": [ "b = np.ones((3, 2))\n", "print(b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Reshaping Arrays" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use the function [np.reshape()](https://numpy.org/doc/stable/reference/generated/numpy.reshape.html) to change the shape of a Numpy array." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One use of this function is to change one dimensional Numpy arrays to two dimensional row vectors or column vectors.\n", "\n", "You pass in the required dimension which the array needs to be reshaped to as input to the function. One of the values in the tuple can be -1. This value shall then be determined automatically. You shall get an error if you have more than one -1 in the shape.\n", "\n", "Reshape is often used in conjunction with [np.arange()](https://numpy.org/doc/stable/reference/generated/numpy.arange.html).\n", "\n", "`np.arange(x)` shall return a one dimensional array which contains the elements $[0, 1, \\dots, x - 1]$." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[0 1 2]\n", "(3,)\n" ] } ], "source": [ "a = np.arange(3)\n", "print(a)\n", "print(a.shape)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0 1 2]]\n", "(1, 3)\n" ] } ], "source": [ "# Reshaping it to a row vector\n", "# The first value of the tuple is 1 because a row vector only has 1 row\n", "a = a.reshape((1, -1))\n", "print(a)\n", "print(a.shape)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0]\n", " [1]\n", " [2]]\n", "(3, 1)\n" ] } ], "source": [ "# Reshaping it to a column vector\n", "# Why does this shape work?\n", "a = a.reshape((-1, 1))\n", "print(a)\n", "print(a.shape)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A: Original\n", "[0 1 2 3 4 5]\n", "(6,)\n", "A: Reshaped\n", "[[0 1]\n", " [2 3]\n", " [4 5]]\n", "(3, 2)\n" ] } ], "source": [ "# Reshaping an array we get from arange\n", "a = np.arange(6)\n", "print(\"A: Original\")\n", "print(a)\n", "print(a.shape)\n", "\n", "a = a.reshape((3, 2))\n", "print(\"A: Reshaped\")\n", "print(a)\n", "print(a.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Accessing Numpy Arrays" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use the common square bracket notation for accessing the elements of a Numpy array. \n", "\n", "If $A$ is a one dimensional array, then you can use `A[i]` to access the $i^{{th}}$ element of $A$.\n", "\n", "If $A$ is a two dimensional array, then you can use `A[i, j]` to access element $a_{{i, j}}$, where $i$ is the row and $j$ is the column.\n", "\n", "Note that Numpy arrays are 0 indexed." ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[0 1 2 3]\n", "The 3rd element of A\n", "2\n" ] } ], "source": [ "a = np.arange(4)\n", "\n", "print(a)\n", "\n", "print(\"The 3rd element of A\")\n", "print(a[2])" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0 1 2]\n", " [3 4 5]\n", " [6 7 8]]\n", "The 2nd row of the array:\n", "[3 4 5]\n", "Element at 0th row and 0th column:\n", "0\n", "Element at 3rd row and 2nd column:\n", "7\n" ] } ], "source": [ "a = np.arange(9).reshape(3, 3)\n", "\n", "print(a)\n", "\n", "print(\"The 2nd row of the array:\")\n", "print(a[1])\n", "\n", "print(\"Element at 0th row and 0th column:\")\n", "print(a[0, 0])\n", "\n", "print(\"Element at 3rd row and 2nd column:\")\n", "print(a[2, 1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Array Slicing" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can select more than one element along a certain dimenson by using the slicing notation. Suppose $A$ is a one dimensional array. If you want to get all elements from an array with $index >= start$ but $index < end$, you can use the following: `A[start:end]`\n", "\n", "To get all the elements, you can use: `A[:]`\n", "\n", "Read more about the slicing notation in the [documentation](https://numpy.org/doc/stable/reference/arrays.indexing.html)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can conveniently get the $i^{{th}}$ column of a two dimensional array $A$ by using: `A[:, i]`\n", "\n", "Think of this as getting all the rows of the $i^{{th}}$ column of $A$. Note that this returns the result as a one dimensional array. You need to reshape it if you want the answer as a column vector." ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 0 1 2 3]\n", " [ 4 5 6 7]\n", " [ 8 9 10 11]\n", " [12 13 14 15]]\n", "The 2th column of the array:\n", "[ 1 5 9 13]\n", "The 4th column of the array reshaped into a column vector:\n", "[[ 3]\n", " [ 7]\n", " [11]\n", " [15]]\n" ] } ], "source": [ "a = np.arange(16).reshape(4, 4)\n", "print(a)\n", "\n", "print(\"The 2th column of the array:\")\n", "print(a[:, 1])\n", "\n", "print(\"The 4th column of the array reshaped into a column vector:\")\n", "print(a[:, 3].reshape(-1, 1))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You shall now apply what you have learned so for in the following question." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Finding The Trace of a Matrix" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In the next function, you shall be given a square matrix $A$. You need to find the *trace* of $A$, which is defined to be the sum of all elements on the principal diagonal of $A$.\n", "\n", "**Hint**: Each element of the principal diagonal has an index of the form $(i, i)$." ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [], "source": [ "def trace(A):\n", " \"\"\"\n", " Find the trace of A\n", "\n", " Arguments:\n", " A: A two dimensional square matrix\n", "\n", " Returns:\n", " s: The trace of A\n", " \"\"\"\n", "\n", " s = 0\n", "\n", " # Get the number of rows in A\n", " n = A.shape[0] # Replace the 0 with the correct code to get the number of rows\n", "\n", " # Loop over all elements of the principal diagonal\n", " for i in range(n):\n", " s += A[i,i] # Replace the None with the required element of A\n", "\n", " return s" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In order to test your code, we shall be using the testing module of Numpy. More on this later." ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "17\n", "All tests passed!\n" ] } ], "source": [ "A = np.array([[3, 2, 7],\n", " [4, 9, 0],\n", " [1, 8, 5]])\n", "print(trace(A))\n", "\n", "# confront with the expected solution s_exp\n", "s_exp = 17\n", "np.testing.assert_allclose(trace(A), s_exp, rtol=1e-5)\n", "\n", "A = np.array([[12, -2, 31, 18],\n", " [32, -77, -24, 19],\n", " [87, 93, -53, 13],\n", " [49, 81, 63, 70]])\n", "s_exp = -48\n", "np.testing.assert_allclose(trace(A), s_exp, rtol=1e-5)\n", "\n", "print(\"All tests passed!\")" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-48\n", "-48\n" ] } ], "source": [ "print(trace(A))\n", "print(np.sum(np.diag(A))) # GW alternative way" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "- [GW] Few words on `teting.assert_allclose` from **NumPY**\n", " - Find out more at the end of the Notebook" ] }, { "cell_type": "code", "execution_count": 43, "metadata": {}, "outputs": [], "source": [ "# [GW] Test how np.teting.assert_allclose works (defaults are: atol=0, rtol=1e-7)\n", "# |actual - desired| <= atol + rtol * |desired|\n", "a = np.array([1.0, 2.0, 3.0])\n", "b = np.array([1.0000001, 2.0000001, 3.0000001])\n", "c = np.array([1.0, 2.0, 3.5])\n", "\n", "np.testing.assert_allclose(a, b) # shall pass (nothing is returned)\n", "#np.testing.assert_allclose(a, c) # shall raise an assertion error\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Operations on Numpy Arrays" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use the operations '`*`', '`+`', '`-`', '`**`' and '`/`' on numpy arrays and they operate element-wise.\n", "\n", "Generally, both the arrays must have the same shape, but see broadcasting later in the notebook." ] }, { "cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [], "source": [ "a = np.array([[4, 1, 2],\n", " [7, 2, 3]])\n", "\n", "b = np.array([[3, 6, 9],\n", " [7, 8, 2]])" ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 7 7 11]\n", " [14 10 5]]\n" ] } ], "source": [ "print(a + b)" ] }, { "cell_type": "code", "execution_count": 27, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 1 -5 -7]\n", " [ 0 -6 1]]\n" ] } ], "source": [ "print(a - b)" ] }, { "cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[12 6 18]\n", " [49 16 6]]\n" ] } ], "source": [ "print(a * b)" ] }, { "cell_type": "code", "execution_count": 29, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1.33333333 0.16666667 0.22222222]\n", " [1. 0.25 1.5 ]]\n" ] } ], "source": [ "print(a / b)" ] }, { "cell_type": "code", "execution_count": 30, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 64 1 512]\n", " [823543 256 9]]\n" ] } ], "source": [ "print(a ** b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Some other useful functions is Numpy are [np.sqrt()](https://numpy.org/doc/stable/reference/generated/numpy.sqrt.html), [np.exp()](https://numpy.org/doc/stable/reference/generated/numpy.exp.html), [np.log()](https://numpy.org/doc/stable/reference/generated/numpy.sqrt.html) which apply the corresponding operation to every element of the inputted Numpy array." ] }, { "cell_type": "code", "execution_count": 31, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1.73205081 2.44948974 3. ]\n", " [2.64575131 2.82842712 1.41421356]]\n" ] } ], "source": [ "print(np.sqrt(b))" ] }, { "cell_type": "code", "execution_count": 32, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 54.59815003 2.71828183 7.3890561 ]\n", " [1096.63315843 7.3890561 20.08553692]]\n" ] } ], "source": [ "print(np.exp(a))" ] }, { "cell_type": "code", "execution_count": 33, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1.38629436 0. 0.69314718]\n", " [1.94591015 0.69314718 1.09861229]]\n" ] } ], "source": [ "print(np.log(a))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Broadcasting" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Broadcasting allows you to perform element wise operations on Numpy arrays that are not of the same dimension but can be strected/duplicated so that they are of the same dimenson.\n", "\n", "The simplest example for this is when you have to multiply all elements of a Numpy array with a scaler or add a scaler to all elements of the array." ] }, { "cell_type": "code", "execution_count": 34, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 8 2 4]\n", " [14 4 6]]\n" ] } ], "source": [ "print(a * 2)" ] }, { "cell_type": "code", "execution_count": 35, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[6 3 4]\n", " [9 4 5]]\n" ] } ], "source": [ "print(a + 2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The way to think about this is to imagine that the number $2$ is being duplicated so that it fills all the elements of a Numpy array which has the same shape as that of $a$. Then, the regular element-wise operation is applicable.\n", "\n", "Broadcasting can also be applied during operation between two arrays of different shapes. Suppose that $a$ is an array of shape $(2, 3)$ whereas $c$ is an array of shape $(1, 3)$. Then, when an operation is performed between $a$ and $c$, you can think of $c$ being duplicated along the axis of the rows so that its shape is $(2, 3)$." ] }, { "cell_type": "code", "execution_count": 36, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "a:\n", " [[4 1 2]\n", " [7 2 3]]\n", "c:\n", " [[5 3 4]]\n", "a + c:\n", " [[ 9 4 6]\n", " [12 5 7]]\n" ] } ], "source": [ "c = np.array([5, 3, 4]).reshape((1, 3))\n", "\n", "print('a:\\n', a)\n", "print('c:\\n', c)\n", "print('a + c:\\n', a + c)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You shall get an error if you try to perform an element wise operation on two arrays which cannot be broadcasted together." ] }, { "cell_type": "code", "execution_count": 37, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0 1 2]\n", " [3 4 5]]\n", "[[0 1]\n", " [2 3]]\n" ] }, { "ename": "ValueError", "evalue": "operands could not be broadcast together with shapes (2,3) (2,2) ", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 5\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mb\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 6\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 7\u001b[0;31m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mb\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mValueError\u001b[0m: operands could not be broadcast together with shapes (2,3) (2,2) " ] } ], "source": [ "a = np.arange(6).reshape((2, 3))\n", "print(a)\n", "\n", "b = np.arange(4).reshape((2, 2))\n", "print(b)\n", "\n", "print(a + b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Broadcasting is a very powerful feauture of Numpy which often reduces the complexity of your code.\n", "\n", "To know more about broadcasting, you can read the [documentation](https://numpy.org/doc/stable/user/basics.broadcasting.html)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now you will write a function which makes use of broadcasting to perform operations on a Numpy array." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Sigmoid Function\n", "\n", "Sigmoid is a very important function in machine learning. Its formula is given by: $$\\operatorname{sigmoid}(x) = \\frac{1}{1+e^{-x}}$$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Write a function which computes the sigmoid of $x$. Note that $x$ might be a real number or a Numpy array." ] }, { "cell_type": "code", "execution_count": 38, "metadata": {}, "outputs": [], "source": [ "def sigmoid(x):\n", " \"\"\"\n", " Computes the sigmoid of x\n", "\n", " Arguments:\n", " x: A real number or a Numpy array\n", "\n", " Returns:\n", " s: The sigmoid of x\n", " \"\"\"\n", " s = 1 / (1 + np.exp(-x)) \n", " return s" ] }, { "cell_type": "code", "execution_count": 39, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "All tests passed!\n" ] } ], "source": [ "x = 2.1\n", "expected = 0.8909031\n", "np.testing.assert_allclose(sigmoid(x), expected, rtol=1e-5)\n", "\n", "x = np.array([[3.4, -7.1, 9.4],\n", " [-2.7, 8.882, -2.114]])\n", "expected = np.array([[9.67704535e-01, 8.24424686e-04, 9.99917283e-01],\n", " [6.29733561e-02, 9.99861153e-01, 1.07743524e-01]])\n", "np.testing.assert_allclose(sigmoid(x), expected, rtol=1e-5)\n", "\n", "print(\"All tests passed!\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Performing Operations Along an Axis" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Numpy also has other useful functions like [np.sum()](https://numpy.org/doc/stable/reference/generated/numpy.sum.html), [np.max()](https://numpy.org/doc/stable/reference/generated/numpy.amax.html) and [np.min()](https://numpy.org/doc/stable/reference/generated/numpy.amin.html).\n", "\n", "When you just pass the Numpy array to these functions, they simply return the answer of the operation over the entire array. Now suppose you want to find the maximum element of each column of an array. You can do this by passing an additional argument called `axis` to the function.\n", "\n", "If $A$ is a two dimensional array, if you want to find the maximum of each column of $A$, you can use: `np.max(A, axis=0)`\n", "\n", "If $A$ is a two dimensional array, if you want to find the maximum of each row of $A$, you can use: `np.max(A, axis=1)`\n", "\n", "A point to note is that these functions may return one dimensional arrays, which might cause errors with broadcasting. In order to make sure that two dimensional arrays are returned, you must pass the argument `keepdims=True` in the function: `np.max(A, axis=0, keepdims=True)`" ] }, { "cell_type": "code", "execution_count": 44, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[2 1 3]\n", " [4 5 6]] \n", "\n", "[1 4]\n", "[[1]\n", " [4]]\n" ] } ], "source": [ "# Minimum of all elements of array\n", "a = np.array([[2, 1, 3],\n", " [4, 5, 6]])\n", "\n", "print(a, '\\n')\n", "print(np.min(a, axis=1))\n", "print(np.min(a, axis=1, keepdims=True))" ] }, { "cell_type": "code", "execution_count": 47, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "21\n", "[6 6 9]\n" ] } ], "source": [ "# Sum of all columns of an array, returned as one dimensional array\n", "print(np.sum(a))\n", "print(np.sum(a, axis=0))" ] }, { "cell_type": "code", "execution_count": 48, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[6 6 9]]\n" ] } ], "source": [ "# Sum of all columns of an array, returned as a two dimensional array\n", "print(np.sum(a, axis=0, keepdims=True))" ] }, { "cell_type": "code", "execution_count": 49, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 6]\n", " [15]]\n" ] } ], "source": [ "# Sum of all rows of an array\n", "print(np.sum(a, axis=1, keepdims=True))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Normalizing All Columns of an Array" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Before a machine learning model is applied to data, it is very common to first normalize it. In the next exercise, you will need to normalize all the columns of an array. Suppose that $x_i$ is the $i^{{th}}$ column of the input array, and $c_i$ is the $i^{{th}}$ column of the output array, then: $$c_i = \\frac{x_i - \\operatorname{mean}(x_i)}{\\sigma(x_i)}$$ \n", "where $\\operatorname{mean}(x_i)$ is the algebraic mean of all the elements of $x_i$ and $\\sigma(x_i)$ is the standard deviation of all the elements of $x_i$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For this exercise, you may assume that $\\sigma(x)$ is never equal to $0$. \n", "\n", "Try to complete this exercise __without using any for__ loops. Numpy functions make use of a technique called __vectorization__ which make them __much faster than for loops__. Using vectorized implementations can often make magnitudes of difference in the training time of your models: where your model might initially take a couple of days to train, it would now be done in a couple of hours.\n", "\n", "If you are stuck, you may find [np.mean()](https://numpy.org/doc/stable/reference/generated/numpy.mean.html) and [np.std()](https://numpy.org/doc/stable/reference/generated/numpy.std.html) useful." ] }, { "cell_type": "code", "execution_count": 50, "metadata": {}, "outputs": [], "source": [ "def normalize(x):\n", " \"\"\"\n", " Normalize all the columns of x\n", "\n", " Arguments:\n", " x: A two dimensional Numpy array\n", "\n", " Returns:\n", " c: The normalized version of x\n", " \"\"\"\n", "\n", " # Calculate the mean of all columns\n", " col_means = np.mean(x, axis=0)\n", "\n", " # Calculate the standard deviation all columns\n", " col_std = np.std(x, axis=0)\n", "\n", " # Compute the final answer\n", " c = (x - col_means) / col_std\n", " return c\n" ] }, { "cell_type": "code", "execution_count": 51, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "All tests passed!\n" ] } ], "source": [ "# First test\n", "x = np.array([[1, 4],\n", " [3, 2]])\n", "expected = np.array([[-1, 1],\n", " [1, -1]])\n", "np.testing.assert_allclose(normalize(x), expected, rtol=1e-5)\n", "\n", "# Second test\n", "x = np.array([[324.33, 136.11, 239.38, 237.17],\n", " [123.43, 325.24, 243.52, 745.25],\n", " [856.36, 903.02, 430.25, 531.35]])\n", "expected = np.array([[-0.35694152, -0.97689929, -0.73023324, -1.28391946],\n", " [-1.00662188, -0.39712973, -0.68372539, 1.1554411 ],\n", " [ 1.3635634 , 1.37402902, 1.41395863, 0.12847837]])\n", "np.testing.assert_allclose(normalize(x), expected, rtol=1e-5)\n", "\n", "print(\"All tests passed!\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Linear Algebra" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this course, we use Numpy arrays for linear algebra. We shall use one dimensional arrays or two dimensional column vectors (two dimensional arrays of shape $(d, 1)$) to represent vectors and two dimensional arrays to represent matrices." ] }, { "cell_type": "code", "execution_count": 58, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A Vector:\n", " [0 1 2 3 4 5]\n", "A Matrix:\n", " [[1 2 3]\n", " [4 5 6]]\n" ] } ], "source": [ "a = np.arange(6)\n", "print('A Vector:\\n', a)\n", "\n", "a = np.array([[1, 2, 3],\n", " [4, 5, 6]]) # Can you think of a better way to create this array?\n", "a = np.arange(1,7).reshape((2,-1)) # Yes, I can...\n", "\n", "print('A Matrix:\\n', a)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Transposing a Matrix" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Suppose `A` is a matrix, then you can take the transpose of `A` with `A.T`" ] }, { "cell_type": "code", "execution_count": 59, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A\n", " [[5 2 9]\n", " [6 1 0]]\n", "A.T\n", " [[5 6]\n", " [2 1]\n", " [9 0]]\n", "B\n", " [[0 1 2]\n", " [3 4 5]\n", " [6 7 8]]\n", "B.T\n", " [[0 3 6]\n", " [1 4 7]\n", " [2 5 8]]\n" ] } ], "source": [ "A = np.array([[5, 2, 9],\n", " [6, 1, 0]])\n", "\n", "print('A\\n', A)\n", "print('A.T\\n', A.T)\n", "\n", "B = np.arange(9).reshape((3, 3))\n", "\n", "print('B\\n', B)\n", "print('B.T\\n', B.T)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that taking the transpose of a 1D array has **NO** effect." ] }, { "cell_type": "code", "execution_count": 60, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A:\n", " [1. 1. 1.]\n", "Shape of A:\n", " (3,)\n", "A.T:\n", " [1. 1. 1.]\n", "Shape of A.T:\n", " (3,)\n" ] } ], "source": [ "a = np.ones(3)\n", "print('A:\\n', a)\n", "print('Shape of A:\\n', a.shape)\n", "print('A.T:\\n', a.T)\n", "print('Shape of A.T:\\n', a.T.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "But it does work if you have a 2D array of shape $(d, 1)$\n" ] }, { "cell_type": "code", "execution_count": 61, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A:\n", " [[1.]\n", " [1.]\n", " [1.]]\n", "Shape of A:\n", " (3, 1)\n", "A.T:\n", " [[1. 1. 1.]]\n", "Shape of A.T:\n", " (1, 3)\n" ] } ], "source": [ "a = np.ones((3,1))\n", "print('A:\\n', a)\n", "print('Shape of A:\\n', a.shape)\n", "print('A.T:\\n', a.T)\n", "print('Shape of A.T:\\n', a.T.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Dot Product and Matrix Product" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can compute the dot product of two vectors or the matrix product between two matrices using [np.dot()](https://numpy.org/doc/stable/reference/generated/numpy.dot.html) or the `@` operator. Going forward, we shall only be using the `@` operator." ] }, { "cell_type": "code", "execution_count": 63, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "x:\n", " [1 2 3]\n", "y:\n", " [4 5 6]\n", "np.dot(x, y):\n", " 32\n", "x @ y:\n", " 32\n" ] } ], "source": [ "x = np.array([1, 2, 3])\n", "y = np.array([4, 5, 6])\n", "print('x:\\n', x)\n", "print('y:\\n', y)\n", "print('np.dot(x, y):\\n', np.dot(x, y))\n", "print('x @ y:\\n', x @ y)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use this operator to compute matrix-matrix and matrix-vector product as well." ] }, { "cell_type": "code", "execution_count": 64, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A:\n", " [[2 0 1]\n", " [1 3 4]\n", " [0 2 1]]\n", "B:\n", " [[0 1]\n", " [2 3]\n", " [4 5]]\n", "A @ B:\n", " [[ 4 7]\n", " [22 30]\n", " [ 8 11]]\n", "A:\n", " [[2 0 1]\n", " [1 3 4]\n", " [0 2 1]]\n", "C:\n", " [[3]\n", " [2]\n", " [8]]\n", "A @ C:\n", " [[14]\n", " [41]\n", " [12]]\n" ] } ], "source": [ "A = np.array([[2, 0, 1],\n", " [1, 3, 4],\n", " [0, 2, 1]])\n", "\n", "B = np.arange(6).reshape(3, 2)\n", "\n", "C = np.array([3, 2, 8]).reshape((-1, 1))\n", "\n", "print('A:\\n', A)\n", "print('B:\\n', B)\n", "print('A @ B:\\n', A @ B)\n", "\n", "print('A:\\n', A)\n", "print('C:\\n', C)\n", "print('A @ C:\\n', A @ C)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that if you try to multiply two matrices which are not compatible for multiplication, you get an error." ] }, { "cell_type": "code", "execution_count": 65, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(3, 2)\n", "(3, 3)\n" ] }, { "ename": "ValueError", "evalue": "shapes (3,2) and (3,3) not aligned: 2 (dim 1) != 3 (dim 0)", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mB\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mshape\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mA\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mshape\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mB\u001b[0m \u001b[0;34m@\u001b[0m \u001b[0mA\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mValueError\u001b[0m: shapes (3,2) and (3,3) not aligned: 2 (dim 1) != 3 (dim 0)" ] } ], "source": [ "print(B.shape)\n", "print(A.shape)\n", "print(B @ A)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### The Linalg Library of Numpy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The Numpy library ships with [np.linalg](https://numpy.org/doc/stable/reference/routines.linalg.html) package which lets us compute many properties of matrices.\n", "\n", "For example, we can compute the determinant of a square matrix by using [np.linalg.det()](https://numpy.org/doc/stable/reference/generated/numpy.linalg.det.html)." ] }, { "cell_type": "code", "execution_count": 78, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A:\n", " [[2 0 1]\n", " [1 3 4]\n", " [0 2 1]]\n", "Determinant of A : -7.999999999999998\n" ] } ], "source": [ "# This computes the determinant\n", "print('A:\\n', A)\n", "d = np.linalg.det(A)\n", "d_byHand = 2*(3-8)-1*(0-2)\n", "print('Determinant of A :', d)\n", "np.testing.assert_allclose(d, d_byHand, rtol=1e-5)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can compute the inverse of a matrix by using [np.linalg.inv()](https://numpy.org/doc/stable/reference/generated/numpy.linalg.inv.html)." ] }, { "cell_type": "code", "execution_count": 79, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 0.625 -0.25 0.375]\n", " [ 0.125 -0.25 0.875]\n", " [-0.25 0.5 -0.75 ]]\n" ] } ], "source": [ "# This computes the inverse\n", "print(np.linalg.inv(A))\n", "\n", "I = np.eye(3) # We can use this function to generate the identity matrix\n", "\n", "np.testing.assert_allclose(A @ np.linalg.inv(A), I, atol=1e-10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can compute the eigenvalues and eigenvectors of a matrix using [np.linalg.eig()](https://numpy.org/doc/stable/reference/generated/numpy.linalg.eig.html)." ] }, { "cell_type": "code", "execution_count": 81, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "The eigenvalues are:\n", " [ 5.10548262 1.77653793 -0.88202054]\n", "\n", "The eigenvectors are:\n", " [[-0.13964316 0.9724502 -0.24499029]\n", " [-0.8901906 -0.08437307 -0.66441635]\n", " [-0.43365942 -0.21730574 0.70606705]]\n" ] } ], "source": [ "# This computes the eigenvalues and eigenvectors\n", "eigenvalues, eigenvectors = np.linalg.eig(A)\n", "print(\"\\nThe eigenvalues are:\\n\", eigenvalues)\n", "print(\"\\nThe eigenvectors are:\\n\", eigenvectors)" ] }, { "cell_type": "code", "execution_count": 82, "metadata": {}, "outputs": [], "source": [ "np.linalg.eig?" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You shall now apply all the concepts that you have learned in the next sections." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Diagonalizing a Matrix" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this question, you shall be given a square matrix which you need to diagonalize. In particular, you will be given a diagonalizable matrix $A$ and you need to find matrices $S$ and $D$ such that: $$A = SDS^{{-1}}$$\n", "\n", "Recall that in order to do this, you must first find all the eigenvalues and eigenvectors of $A$. Then, $S$ is the matrix of all the eigenvectors arranged as columns, and $D$ is the matrix of the corresponding eigenvalues arranged along the diagonal.\n", "\n", "Suppose $A = \n", "\\begin{bmatrix}\n", "1 & 5 \\\\\n", "2 & 4 \\\\\n", "\\end{bmatrix} $\n", "\n", "Then, we can calculate $S = \n", "\\begin{bmatrix}\n", "-2.5 & 1 \\\\\n", "1 & 1 \\\\\n", "\\end{bmatrix} $\n", "\n", "And $D = \n", "\\begin{bmatrix}\n", "-1 & 0 \\\\\n", "0 & 6 \\\\\n", "\\end{bmatrix} $\n", "\n", "You might find [np.zeros()](https://numpy.org/devdocs/reference/generated/numpy.zeros.html), [np.linalg.eig()](https://numpy.org/devdocs/reference/generated/numpy.linalg.eig.html) and [np.linalg.inv()](https://numpy.org/devdocs/reference/generated/numpy.linalg.inv.html) useful. Note that for this exercise, you may assume that $A$ is always diagonalizable.\n", "\n", "For testing purposes, each eigenvector in $S$ must be of unit length. This shall always be the case if you use [np.linalg.eig()](https://numpy.org/devdocs/reference/generated/numpy.linalg.eig.html). However, if you do not use this function, then depending on your implementation, you might have to normalize the eigenvectors. Also, the eigenvalues must appear in non decreasing order." ] }, { "cell_type": "code", "execution_count": 83, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[-1.0000000e+00 0.0000000e+00]\n", " [ 4.4408921e-16 6.0000000e+00]]\n", "[[-1. 0.]\n", " [ 0. 6.]]\n" ] }, { "data": { "text/plain": [ "array([-1., 6.])" ] }, "execution_count": 83, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# GW check\n", "A = np.array([[1,5], [2,4]])\n", "eig_vals, S = np.linalg.eig(A)\n", "\n", "S_inv = np.linalg.inv(S)\n", "D = S_inv @ A @ S\n", "print(D)\n", "\n", "print(np.diag(eig_vals)) # directly from theory" ] }, { "cell_type": "code", "execution_count": 85, "metadata": {}, "outputs": [], "source": [ "def diagonalize(A):\n", " \"\"\"\n", " Diagonalizes the input matrix A\n", "\n", " Arguments:\n", " A: A two dimensional Numpy array which is guaranteed to be diagonalizable\n", "\n", " Returns:\n", " S, D, S_inv: As explained above\n", " \"\"\"\n", "\n", " # Retrieve the number of rows in A\n", " n = A.shape[0]\n", "\n", " # Get the eigenvalues and eigenvectors of A\n", " eig_vals, S = np.linalg.eig(A)\n", "\n", " # Start by initializing D to a matrix of zeros of the appropriate shape\n", " D = np.zeros((n,n))\n", " # Set the diagonal element of D to be the eigenvalues\n", " for i in range(n):\n", " D[i,i] = eig_vals[i]\n", "\n", " # Compute the inverse of S\n", " S_inv = np.linalg.inv(S)\n", "\n", " return S, D, S_inv" ] }, { "cell_type": "code", "execution_count": 86, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "All tests passed!\n" ] } ], "source": [ "A = np.array([[1, 5],\n", " [2, 4]])\n", "S_exp = np.array([[-0.92847669, -0.70710678],\n", " [ 0.37139068, -0.70710678]])\n", "D_exp = np.array([[-1, 0],\n", " [0, 6]])\n", "S_inv_exp = np.array([[-0.76930926, 0.76930926],\n", " [-0.40406102, -1.01015254]])\n", "\n", "\n", "S, D, S_inv = diagonalize(A)\n", "np.testing.assert_allclose(S_exp, S, rtol=1e-5, atol=1e-10)\n", "np.testing.assert_allclose(D_exp, D, rtol=1e-5, atol=1e-10)\n", "np.testing.assert_allclose(S_inv_exp, S_inv, rtol=1e-5, atol=1e-10)\n", "\n", "A = np.array([[4, -9, 6, 12],\n", " [0, -1, 4, 6],\n", " [2, -11, 8, 16],\n", " [-1, 3, 0, -1]])\n", "S_exp = np.array([[ -5.00000000e-01, 8.01783726e-01, 9.04534034e-01, 3.77964473e-01],\n", " [ -5.00000000e-01, 5.34522484e-01, 3.01511345e-01, 7.55928946e-01],\n", " [ 5.00000000e-01,-1.98636631e-14, 3.01511345e-01, 3.77964473e-01],\n", " [ -5.00000000e-01, 2.67261242e-01, -5.03145109e-15, 3.77964473e-01]])\n", "D_exp = np.array([[1, 0, 0, 0],\n", " [0, 2, 0, 0],\n", " [0, 0, 3, 0],\n", " [0, 0, 0, 4]])\n", "S_inv_exp = np.array([[ -2.00000000e+00, 1.00000000e+01, -4.00000000e+00, -1.40000000e+01],\n", " [ -3.74165739e+00, 2.24499443e+01, -1.12249722e+01, -2.99332591e+01],\n", " [ 3.31662479e+00, -1.32664992e+01, 6.63324958e+00, 1.65831240e+01],\n", " [ 2.74154909e-15, -2.64575131e+00, 2.64575131e+00, 5.29150262e+00]])\n", "\n", "S, D, S_inv = diagonalize(A)\n", "np.testing.assert_allclose(S_exp, S, rtol=1e-5, atol=1e-10)\n", "np.testing.assert_allclose(D_exp, D, rtol=1e-5, atol=1e-10)\n", "np.testing.assert_allclose(S_inv_exp, S_inv, rtol=1e-5, atol=1e-10)\n", "\n", "print(\"All tests passed!\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Introduction to Assert Statements and Testing" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "An [assert](w3schools.com/python/ref_keyword_assert.asp) statement lets you check if a condition in your program evaluates to true.\n", "\n", "If the condition does evaluate to true, then nothing happens and the program execution continues normally. However, if it evaluates to false, then the program execution is immediately terminated, an error is raised and an error message is printed." ] }, { "cell_type": "code", "execution_count": 90, "metadata": {}, "outputs": [], "source": [ "# An assert statement where condition evaluates to true\n", "x = 5\n", "assert x == 5" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Since the condition evaluates to true, nothing happens." ] }, { "cell_type": "code", "execution_count": 91, "metadata": {}, "outputs": [ { "ename": "AssertionError", "evalue": "", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mAssertionError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;31m# An assert statement where the condition evaluates to false\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mx\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m5\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0;32massert\u001b[0m \u001b[0mx\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m4\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mAssertionError\u001b[0m: " ] } ], "source": [ "# An assert statement where the condition evaluates to false\n", "x = 5\n", "assert x == 4" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This time, since the assert statement evaluates to false, an error is thrown and the line where the assert statement failed is printed to the standard output.\n", "\n", "You can also print a message which further explains what error took place." ] }, { "cell_type": "code", "execution_count": 92, "metadata": {}, "outputs": [ { "ename": "AssertionError", "evalue": "x does not store the intended value", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mAssertionError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;31m# An assert statement with an error message\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mx\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m5\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0;32massert\u001b[0m \u001b[0mx\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m4\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"x does not store the intended value\"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mAssertionError\u001b[0m: x does not store the intended value" ] } ], "source": [ "# An assert statement with an error message\n", "x = 5\n", "assert x == 4, \"x does not store the intended value\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Asserts are a very powerful way to find bugs in your code. Do not hesitate to use these statements in the functions you write to assist you in debugging them. Later, we shall also be using assert statements sometimes to test the functions that you shall write in the later labs.\n", "\n", "Let us see an example where using assert statements helps us spot bugs in our code." ] }, { "cell_type": "code", "execution_count": 109, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0 0\n", "1 1\n", "2 3\n", "3 6\n", "4 10\n", "4 6\n", "3 3\n", "2 1\n", "1 0\n" ] } ], "source": [ "def test_assert():\n", " \"\"\"\n", " This function demonstrates the use of assert statements in debugging\n", " \"\"\"\n", "\n", " A = np.arange(5)\n", " s = 0\n", "\n", " # We shall first add all elements of A to s\n", " for i in range(A.shape[0]):\n", " s += A[i]\n", " print(i, s)\n", " \n", " # We shall now subtract all the elements of A in the reverse order\n", " # Unfortunately, we have a bug\n", " for i in range(A.shape[0] - 1, 0, -1):\n", " s -= A[i]\n", " print(i, s)\n", " \n", " # Quite certainly, s must be equal to 0 at this point\n", " # Had our implementation been correct, this assert should pass\n", " assert s == 0\n", "\n", "test_assert()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Can you find the bug in the code and fix it? Have a look at the documentation of [range()](https://docs.python.org/3/library/stdtypes.html#typesseq-range)." ] }, { "cell_type": "code", "execution_count": 111, "metadata": {}, "outputs": [], "source": [ "# Yes, I can... The error was in the second argument, which should be changed to -1\n", "# Uncomment the print statements in the function for details." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The Numpy Testing Module" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Numpy has a very useful function called [np.testing.assert_allclose()](https://numpy.org/doc/stable/reference/generated/numpy.testing.assert_allclose.html) which allows us to test our functions.\n", "\n", "The function accepts two numbers or two Numpy arrays and checks them for equality. Note that you cannot compare floating point numbers using the `==` operator as **you have to account for a margin of error** which can be caused due to rounding. This function allows you to customize the error margin based on your needs.\n", "\n", "We highly advice you to read the documentation of this function." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The function takes two compulsory arguments, the first is the array which you want to test, and the second is the array which you want to test the array against. You can also configure the tolerance level if you wish.\n", "\n", "You can see examples of how this function is used earlier in the notebook." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In the next code block, we have written a function which tries to compute the inverse of a 2 by 2 matrix. However, the function is incorrect." ] }, { "cell_type": "code", "execution_count": 125, "metadata": {}, "outputs": [], "source": [ "def inverse(A):\n", " \"\"\"\n", " Computes (incorrectly) the inverse of A\n", "\n", " A must have shape (2, 2)\n", " \"\"\"\n", "\n", " return np.array([[A[1, 1], -A[0, 1]],\n", " [-A[1, 0], A[0, 0]]])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We have written a test for this function but it unfortunately passes. Can you write a test for this function which fails? Then, can you modify the function so that it is now correct and passes the test that you wrote?" ] }, { "cell_type": "code", "execution_count": 129, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 4 -5]\n", " [-2 1]]\n", "[[-0.66666667 0.83333333]\n", " [ 0.33333333 -0.16666667]]\n" ] }, { "ename": "AssertionError", "evalue": "\nNot equal to tolerance rtol=1e-05, atol=1e-10\n\n(mismatch 50.0%)\n x: array([-6, 0, 0, -6])\n y: array([1., 0., 0., 1.])", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mAssertionError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 12\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minverse\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mA\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 13\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mlinalg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0minv\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mA\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 14\u001b[0;31m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtesting\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0massert_allclose\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minverse\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mA\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m@\u001b[0m \u001b[0mA\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0meye\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrtol\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m1e-5\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0matol\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m1e-10\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 15\u001b[0m \u001b[0;31m#np.testing.assert_allclose(np.linalg.inv(A) @ A, np.eye(2), rtol=1e-5, atol=1e-10)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/opt/conda/lib/python3.6/site-packages/numpy/testing/_private/utils.py\u001b[0m in \u001b[0;36massert_allclose\u001b[0;34m(actual, desired, rtol, atol, equal_nan, err_msg, verbose)\u001b[0m\n\u001b[1;32m 1450\u001b[0m \u001b[0mheader\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m'Not equal to tolerance rtol=%g, atol=%g'\u001b[0m \u001b[0;34m%\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mrtol\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0matol\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1451\u001b[0m assert_array_compare(compare, actual, desired, err_msg=str(err_msg),\n\u001b[0;32m-> 1452\u001b[0;31m verbose=verbose, header=header, equal_nan=equal_nan)\n\u001b[0m\u001b[1;32m 1453\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1454\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/opt/conda/lib/python3.6/site-packages/numpy/testing/_private/utils.py\u001b[0m in \u001b[0;36massert_array_compare\u001b[0;34m(comparison, x, y, err_msg, verbose, header, precision, equal_nan, equal_inf)\u001b[0m\n\u001b[1;32m 787\u001b[0m \u001b[0mverbose\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mverbose\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mheader\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mheader\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 788\u001b[0m names=('x', 'y'), precision=precision)\n\u001b[0;32m--> 789\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mAssertionError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmsg\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 790\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mValueError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 791\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mtraceback\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;31mAssertionError\u001b[0m: \nNot equal to tolerance rtol=1e-05, atol=1e-10\n\n(mismatch 50.0%)\n x: array([-6, 0, 0, -6])\n y: array([1., 0., 0., 1.])" ] } ], "source": [ "A = np.array([[3, 5],\n", " [1, 2]])\n", "A_exp = np.array([[2, -5],\n", " [-1, 3]])\n", "\n", "np.testing.assert_allclose(inverse(A), A_exp, rtol=1e-5)\n", "\n", "# Add another test here\n", "A = np.array([[1, 5],\n", " [2, 4]])\n", "print(inverse(A))\n", "print(np.linalg.inv(A))\n", "\n", "np.testing.assert_allclose(inverse(A) @ A, np.eye(2), rtol=1e-5, atol=1e-10)\n", "#np.testing.assert_allclose(np.linalg.inv(A) @ A, np.eye(2), rtol=1e-5, atol=1e-10)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Debugging Your Code" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "While you are working through the rest of the labs of this course, you shall come across many situations where your code shall not work correctly. You are not alone if this happens to you. Debugging your code can be a difficult and daunting task, so in this last section, we shall give you **some practical guidelines to assist you in debugging** your code." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Dimension Errors in Matrix Multiplication" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This is one of the most frequent errors that you shall face. Let us first see the error message that Numpy prints if you try to multiply two matrices of incompatible dimensions." ] }, { "cell_type": "code", "execution_count": 134, "metadata": {}, "outputs": [ { "ename": "ValueError", "evalue": "shapes (3,2) and (4,5) not aligned: 2 (dim 1) != 4 (dim 0)", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0mA\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mones\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mB\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mones\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m4\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m5\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mA\u001b[0m \u001b[0;34m@\u001b[0m \u001b[0mB\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mValueError\u001b[0m: shapes (3,2) and (4,5) not aligned: 2 (dim 1) != 4 (dim 0)" ] } ], "source": [ "A = np.ones((3, 2))\n", "B = np.ones((4, 5))\n", "print(A @ B)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "That was a mouthful! However, in this course, we shall only be performing operations on arrays that have atmost a dimension of 2, which considerably simplifies things.\n", "\n", "The only important line in the error message is: `size 4 is different from 2`. This says that the 0th dimension of $B$ is 4 whereas the 1th dimension of $A$ is 2, and hence they are incompatible for matrix multiplication.\n", "\n", "One way to debug this is to print the dimensions of all the matrices before and after each matrix multiplication and track them, because a previous error which unfortunately passed dimension checks might be causing the problems here." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Errors can also take place when you try to multiply two one dimensional arrays together or try to multiply a one dimensional array with a two dimensional array.\n", "\n", "We would advice you to use [np.outer()](https://numpy.org/doc/stable/reference/generated/numpy.outer.html) and [np.inner()](https://numpy.org/doc/stable/reference/generated/numpy.outer.html) when computing the dot product of 1D arrays. If $X$ is a vector (represented as a 1D array in this course), then `np.inner(X, X)` calculates $X^T \\cdot X$ (the regular dot product) and `np.outer(X, X)` computes $X \\cdot X^T$.\n", "\n", "If you are performing matrix multiplication between a 2D array and a 1D array, we would advise you to first reshape the 1D array into a 2D array of shape $(d, 1)$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Practicing Code Debugging" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this question, we were trying to find the sum of the maximum element of each row in a 2D array, but we have unfortunately made a bug. Can you fix it so that our tests pass?" ] }, { "cell_type": "code", "execution_count": 136, "metadata": {}, "outputs": [], "source": [ "def sum_of_max(A):\n", " \"\"\"\n", " Computes the sum of the maximum element of each row of A\n", "\n", " A must be a 2D Numpy array\n", " \"\"\"\n", " #return np.sum(np.max(A, axis=0) # incorrect\n", " return np.sum(np.max(A, axis=1))" ] }, { "cell_type": "code", "execution_count": 137, "metadata": {}, "outputs": [], "source": [ "A = np.array([[1, 2],\n", " [3, 4]])\n", "np.testing.assert_allclose(sum_of_max(A), 6)\n", "\n", "A = np.array([[24, 69, 83],\n", " [74, 14, 27]])\n", "np.testing.assert_allclose(sum_of_max(A), 157)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Congratulations on making it to the end of this really long notebook! Now you are in a good place to proceed with the remaining labs of this course." ] } ], "metadata": { "interpreter": { "hash": "fa15421863e696996c586a7125639c86e606c31109b5cdd4c487ad3c1cbbc756" }, "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.4" } }, "nbformat": 4, "nbformat_minor": 2 }