class04_Numpy.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "# Data Mining:<br>Statistical Modeling and Learning from Data\n",
    "\n",
    "## Dr. Ciro Cattuto<br>Dr. Laetitia Gauvin<br>Dr. André Panisson\n",
    "\n",
    "### Exercises - Numpy\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Indexing/Slicing Review"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* `s[i]` (indexing)\n",
    "* `s[i:j]` (slicing)\n",
    "* `s[i:j:k]` (step slicing)\n",
    "* meaning of negative indices\n",
    "* 0-base counting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## EXERCISE: Indexing Review"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = list(range(10)) #in Python3 range returns an iterator! Casting to list is needed\n",
    "m"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access the first position of the list\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access the last position of the list\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The triple [i:j:k] are in fact parameters of a **slice** object:\n",
    "\n",
    "### slice(start, stop[, step])\n",
    "> Return a slice object representing the set of indices specified by range(start, stop, step). The start and step arguments default to None. Slice objects have read-only data attributes start, stop and step which merely return the argument values (or their default). They have no other explicit functionality; however they are used by Numerical Python and other third party extensions. Slice objects are also generated when extended indexing syntax is used. For example: a[start:stop:step] or a[start:stop, i].\n",
    "\n",
    "http://docs.python.org/2/library/functions.html#slice\n",
    "\n",
    "For example, to return the first 3 elements in the even positions of a list, you can use:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m[slice(0,5,2)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "which is equivalent to"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m[0:5:2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## EXERCISE: Slicing Review"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access the first five elements of the list\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access the last five elements of the list\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access the list elements in reverse order\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## pylab mode\n",
    "\n",
    "These imports done for you in `pylab` mode.\n",
    "\n",
    "    import numpy as np\n",
    "    import matplotlib.pyplot as plt\n",
    "    \n",
    "The same can be done with the following command:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%pylab inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# NumPy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<http://www.numpy.org/>:\n",
    "\n",
    "NumPy is the **fundamental package for scientific computing with Python**. It contains among other things:\n",
    "\n",
    "* a powerful N-dimensional array object\n",
    "* sophisticated (**broadcasting**) functions [what is *broadcasting*?]\n",
    "* tools for integrating C/C++ and Fortran code\n",
    "* useful linear algebra, Fourier transform, and random number capabilities"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ndarray.ndim, ndarray.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# zero-dimensions\n",
    "\n",
    "a0 = array(5)\n",
    "a0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a0.ndim, a0.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 1-d array\n",
    "a1 = array([1,2])\n",
    "a1.ndim, a1.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 2-d array\n",
    "a2 = array(([1,2], [3,4]))\n",
    "a2.ndim, a2.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.dtype"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Array creation routines"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = array([1,2])\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = zeros((2,2))\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = ones((2,2))\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = empty((2,2))\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = eye(3,4,1)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = identity(3)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = diag(arange(4))\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = linspace(1,10,5) \n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = logspace(1,2,5) \n",
    "a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Type hierarchy: \n",
    "\n",
    "<img src=\"http://docs.scipy.org/doc/numpy/_images/dtype-hierarchy.png\">"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(10, dtype=float)\n",
    "a.dtype"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(10, dtype=byte)\n",
    "a.dtype"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a[0] = 128\n",
    "a[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a1 = a.astype(int16)\n",
    "a1[0] = 128\n",
    "a1[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### reshape, transpose"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(64)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# map a 0..63 1d array to a 8x8 2d array\n",
    "a1 = a.reshape(8,8)\n",
    "a1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.shape = (8,8)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.T"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### stacking & concatenation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = array([[1, 2], [3, 4]])\n",
    "b = array([[5, 6]])\n",
    "print (a.shape, b.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = concatenate((a, b), axis=0) # vertical stack\n",
    "print (x, x.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = concatenate((a, b.T), axis=1) # horizontal\n",
    "print (y, y.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (vstack((a,b)))\n",
    "print (hstack((a,b.T)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (append(a, b, axis=0))\n",
    "print (append(a, b.T, axis=1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (insert(a, 0, b, axis=0))\n",
    "print (insert(a, 0, b, axis=1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##Numpy operations\n",
    "\n",
    "example of [broadcasting](http://docs.scipy.org/doc/numpy/user/basics.broadcasting.html):\n",
    "\n",
    "> The term broadcasting describes how numpy treats arrays with different shapes during arithmetic operations. Subject to certain constraints, the smaller array is “broadcast” across the larger array so that they have compatible shapes. Broadcasting provides a means of vectorizing array operations so that looping occurs in C instead of Python. It does this without making needless copies of data and usually leads to efficient algorithm implementations. There are, however, cases where broadcasting is a bad idea because it leads to inefficient use of memory that slows computation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "2*a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a+2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "More broadcasting\n",
    "\n",
    "<img src=\"http://scipy-lectures.github.io/_images/numpy_broadcasting.png\" width=\"600\">"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = array([range(0,3)]*4)\n",
    "b = array([range(0,40,10)]*3).T\n",
    "a, b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a+b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b + arange(0,3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "arange(0,40,10).reshape(4,1) + arange(0,3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "w = arange(0,3)\n",
    "a = arange(0,12).reshape(4,3)\n",
    "print (w)\n",
    "print (a)\n",
    "print (w * a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Indexing/Slicing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a3 = arange(30) \n",
    "a3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (a3[0])\n",
    "print (a3[::-1])\n",
    "print (a3[2:5])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (a3[[2,3,4,6,5,2]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.mod(a3, 3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Select numbers divisible by 3."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# list comprehension\n",
    "[i for i in a3 if i % 3 == 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.mod(a3, 3) == 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "divisible_by_3 = np.mod(a3, 3) == 0\n",
    "a3[divisible_by_3]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2d, 3d slicing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(64).reshape(8,8)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a[0,:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a[:,0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a[:2,:2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a[::2,::2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = arange(27).reshape(3,3,3)\n",
    "b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b[0,:,:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b[:,0,:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b[:,:,0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise:  Calculate a series that holds all the squares less than 100"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Use arange, np.sqrt, astype\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# NumPy Functions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "http://docs.scipy.org/doc/numpy/reference/routines.math.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import random\n",
    "a = array([random.randint(0, 10) for i in range(10)])\n",
    "\n",
    "print (a)\n",
    "print (a.min())\n",
    "print (a.max())\n",
    "print (a.mean())\n",
    "print (a.std()) # standard deviation\n",
    "print (a.sum())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = arange(16).reshape(4, 4)\n",
    "print (b)\n",
    "print (b.T)\n",
    "print (b.trace())\n",
    "print (b.min())\n",
    "print (b.min(axis=0))\n",
    "print (b.min(axis=1))\n",
    "print (b.ravel())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(0,3)\n",
    "b = arange(1,4)\n",
    "print (a, b)\n",
    "print (np.dot(a,b))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Matrices"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = matrix([[3, 2, -1], [2, -2, 4], [-1, .5, -1]])\n",
    "print (a)\n",
    "print (a.trace())\n",
    "print (a.diagonal())\n",
    "print (a.T) # matrix transpose\n",
    "print (a.I) # matrix inverse\n",
    "print (a.H) # matrix conjugate transpose"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = matrix(arange(0,6).reshape(2,3))\n",
    "print (a * a) #it's ok to have an error, we need to use the transpose matrix!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (a * a.T)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: Implement the matrix product using python loops and compare the execution time with the Numpy implementation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = arange(50).reshape(5,10)\n",
    "\n",
    "def my_dot(a,b):\n",
    "    \n",
    "    for i in range(a.shape[0]):\n",
    "        \n",
    "    \n",
    "    # YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%timeit my_dot(a,a.T)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%timeit a*a.T"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise:  Reimplement the perceptron using NumPy (e.g., using the matrix product operation)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%file data.csv\n",
    "x1,x2,y\n",
    "0.4946,5.7661,0\n",
    "4.7206,5.7661,1\n",
    "1.2888,5.3433,0\n",
    "4.2898,5.3433,1\n",
    "1.4293,4.5592,0\n",
    "4.2286,4.5592,1\n",
    "1.1921,5.8563,0\n",
    "3.1454,5.8563,1\n",
    "1.063,5.7357,0\n",
    "5.1043,5.7357,1\n",
    "1.5079,5.8622,0\n",
    "3.9799,5.8622,1\n",
    "0.2678,6.9931,0\n",
    "4.5288,6.9931,1\n",
    "0.9726,3.6268,0\n",
    "4.106,3.6268,1\n",
    "2.5389,3.3884,0\n",
    "4.7555,3.3884,1\n",
    "2.473,5.6404,0\n",
    "4.7977,5.6404,1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = np.genfromtxt('data.csv', delimiter=',', skip_header=1)\n",
    "\n",
    "X = data[:,:-1]\n",
    "y = data[:,-1]\n",
    "\n",
    "print (X)\n",
    "print (y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c0 = y==0\n",
    "c1 = y==1\n",
    "\n",
    "plot(X[:,0][c0], X[:,1][c0], 'o', mec='r', mfc='none')\n",
    "plot(X[:,0][c1], X[:,1][c1], 'o', mec='g', mfc='none')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class Perceptron:\n",
    "        \n",
    "    def predict(self, x):\n",
    "        # YOUR CODE HERE\n",
    "\n",
    "    def fit(self, X, y):\n",
    "        self.w = zeros(len(x))\n",
    "        # YOUR CODE HERE, using dot product\n",
    "\n",
    "p = Perceptron()\n",
    "p.fit(X, y)\n",
    "\n",
    "for i,x in enumerate(X):\n",
    "    if p.predict(x) != y[i]:\n",
    "        print 'FAIL'\n",
    "        break\n",
    "else:\n",
    "    print 'SUCCESS!'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise:  Read a digits image and split it up in different training examples"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy import misc\n",
    "img = misc.imread('digits.png', flatten=1)\n",
    "\n",
    "print (img.shape)\n",
    "print (img.dtype)\n",
    "\n",
    "imshow(img, cmap=cm.Greys_r)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = []\n",
    "\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (data[20])\n",
    "imshow(data[20], cmap=cm.Greys_r)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Random Generators"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numpy import random\n",
    "data = random.normal(size=1000, loc=2)\n",
    "h=plt.hist(data, bins=50)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "http://docs.scipy.org/doc/numpy/reference/routines.random.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numpy import random\n",
    "data = random.exponential(size=1000)\n",
    "h = plt.hist(data, bins=50)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise:  Create a set of 100 points that follow the function $$f(x) = 0.5x + 1 $$ and add Gaussian white noise to the result. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = linspace(1,10,100)\n",
    "\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: Use polinomial fitting (numpy.polyfit) to fit the results in a line and verify that the error is Gaussian white noise."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c = polyfit(X, y, 1)\n",
    "\n",
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#Linear Regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##Exercise: Linear Regression with Ordinary Least Squares\n",
    "\n",
    "Find the weight values $\\mathbf{w}$ that minimize the error $E_{\\mathbf{in}}(\\mathbf{w}) = \\frac{1}{N} \\sum_{n=1}^n {(\\mathbf{w}^T \\mathbf{X}_n - \\mathbf{y}_n)^2}$.\n",
    "\n",
    "For this, implement Linear Regression and use the Ordinary Least Squares (OLS) closed-form expression to find the estimated values of $\\mathbf{w}$:\n",
    "\n",
    "$$\\mathbf{w} = (\\mathbf{X}^{\\rm T}\\mathbf{X})^{-1} \\mathbf{X}^{\\rm T}\\mathbf{y}$$\n",
    "\n",
    "You can use the `np.linalg.inv` function to invert a matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(np.linalg.inv)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following function generates some data to test the linear regression:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_linear_regression(c):\n",
    "    X = np.random.rand(100, len(c)-1)\n",
    "    noise = random.normal(size=len(X), loc=0, scale=0.1)\n",
    "    y = ((c[0] + X*c[1:]).sum(axis=1) + noise).reshape(100, 1)\n",
    "    \n",
    "    return X, y\n",
    "\n",
    "X, y = generate_linear_regression([2, 0.8])\n",
    "scatter(X, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To implement the linear regression, follow these steps:\n",
    "\n",
    "1. extend $\\mathbf{X}$ to $d + 1$ dimensions, setting the first column to 1\n",
    "2. calculate $\\mathbf{w}$ using the OLS closed form"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# YOUR CODE HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##Exercise: Linear Regression with nonlinear data\n",
    "\n",
    "Even if the data is nonlinear, linear regression can be used. \n",
    "The first step is to transform the data nonlinearly and append it as new features.\n",
    "\n",
    "The following function creates nonlinear data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_linear_regression_squared(c):\n",
    "    X = np.random.rand(100, len(c)-1)\n",
    "    noise = random.normal(size=len(X), loc=0, scale=0.1)\n",
    "    y = ((c[0] + c[1:]*X**2).sum(axis=1) + noise).reshape(100, 1)\n",
    "    \n",
    "    return X, y\n",
    "\n",
    "X, y = generate_linear_regression_squared([2, 2])\n",
    "scatter(X, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To apply the linear regression in nonlinear data, follow these steps:\n",
    "\n",
    "1. extend $\\mathbf{X}$ to add a nonlinear transformation of $\\mathbf{X}$ (e.g., $\\mathbf{Z} = [\\mathbf{X}, \\mathbf{X}^2$]).\n",
    "1. extend $\\mathbf{Z}$ to $d + 1$ dimensions, setting the first column to 1.\n",
    "2. calculate $\\mathbf{w}$ applying OLS to $\\mathbf{Z}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# YOUR CODE HERE"
   ]
  }
 ],
 "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}