{ "cells": [ { "cell_type": "markdown", "id": "4a485333", "metadata": {}, "source": [ "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
\n", "

MTH786 Machine Learning with Python

\n", " 		\n", " \n", "

Semester A, 2023/2024

\n", "
\n", "

Coursework 6

\n", " 		\n", " \n", "

Dr Nicola Perra

\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "id": "f1a1d188", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from numpy.testing import assert_array_almost_equal, assert_array_equal\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "id": "9d9f18fe", "metadata": {}, "source": [ "### Ridge regression\n", "By completing this exercise you will write a set of functions that are used for building a ridge regression for a given data samples. You will then learn how to perform K-cross validation for a hyperparameter selection. You will finish with applying the above methods to a real data, that needs to be standardised first.\n", "\n", "\n", "1. Implement function **ridge_regression_data** that computes (and outputs) the ridge regression data matrix $\\Phi\\left(\\mathbf{X}\\right)$ defined as\n", "- if a degree is given and it is larger than $1$ then the data matrix should coincide with the polynomial basis matrix\n", "$$\n", "\\Phi\\left(\\mathbf{X}\\right) = \n", "\\begin{pmatrix}\n", "1 & \\left(x^{(1)}\\right)^1 & \\left(x^{(1)}\\right)^2 & \\ldots & \\left(x^{(1)}\\right)^d \\\\\n", "1 & \\left(x^{(2)}\\right)^1 & \\left(x^{(2)}\\right)^2 & \\ldots & \\left(x^{(2)}\\right)^d \\\\\n", "\\vdots & \\vdots & \\vdots & \\ddots & \\vdots & \\\\\n", "1 & \\left(x^{(s)}\\right)^1 & \\left(x^{(s)}\\right)^2 & \\ldots & \\left(x^{(s)}\\right)^d \n", "\\end{pmatrix}\n", "$$\n", "- otherwise, if the degree is not provided or if it is equal to $1$, then the ridge regression data matrix should coincide with the linear regression data matrix, i.e.\n", "$$\n", "\\Phi\\left(\\mathbf{X}\\right) = \n", "\\begin{pmatrix}\n", "1 & x^{(1)}_1 & x^{(1)}_2 & \\ldots & x^{(1)}_d \\\\\n", "1 & x^{(2)}_1 & x^{(2)}_2 & \\ldots & x^{(2)}_d \\\\\n", "\\vdots & \\vdots & \\vdots & \\ddots & \\vdots & \\\\\n", "1 & x^{(s)}_1 & x^{(s)}_2 & \\ldots & x^{(s)}_d \\\\\n", "\\end{pmatrix}\n", "$$\n", "The function **ridge_regression_data**\n", "should take the $1$ compulsory argument *data_inputs* and $1$ optional argument *degree*\n", "- NumPy array *data_inputs* representing a list of inputs $x^{(1)}, x^{(2)}, \\ldots, x^{(s)}$. Each of $x^{(i)}$ is either a vector in the case of linear ridge regression ($d=1$) or a scalar in case of polynomial regression ($d>1$)\n", "- and integer *degree* representing the degree $d$ of a polynomial." ] }, { "cell_type": "code", "execution_count": 2, "id": "40ba602a", "metadata": {}, "outputs": [], "source": [ "def ridge_regression_data(data_inputs, degree=1):\n", "\n", " # this function does not produce a polynomial basis, if data_inputs.ndim>1 it assumes d=1\n", " if (degree > 1 and data_inputs.ndim > 1):\n", " raise NotImplementedError()\n", " \n", " # if the input is just a vector we use a polynomial regression\n", " if (data_inputs.ndim == 1):\n", " X_matrix = np.ones((len(data_inputs), 1))\n", " for i in range(degree):\n", " X_matrix = np.c_[X_matrix, np.power(data_inputs, i + 1)]\n", " else:\n", " first_column = np.ones((len(data_inputs), 1))\n", " X_matrix = np.c_[first_column, data_inputs]\n", " \n", " return X_matrix" ] }, { "cell_type": "markdown", "id": "d04619a1", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 3, "id": "72be9c11", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[ 1., 1., 1.],\n", " [ 1., 2., 4.],\n", " [ 1., 3., 9.],\n", " [ 1., 4., 16.]])" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "test_inputs = np.array([1, 2, 3, 4])\n", "test_degree = 2\n", "ridge_regression_data(test_inputs, test_degree)" ] }, { "cell_type": "code", "execution_count": 4, "id": "bf61232b", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[ 1., 10., 20.],\n", " [ 1., 12., 34.],\n", " [ 1., 43., 44.],\n", " [ 1., 44., 45.]])" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "test_inputs = np.array([[10, 20], [12, 34], [43, 44], [44, 45]])\n", "ridge_regression_data(test_inputs)" ] }, { "cell_type": "markdown", "id": "a98b7f3a", "metadata": {}, "source": [ "2. Write a function **ridge_regression** that takes three arguments *data_matrix*, *data_outputs* and *regularisation*, which computes and returns the solution $\\hat{\\mathbf{W}}$ of the normal equation\n", "$$\n", "\\left(\\Phi^{\\top}\\left(\\mathbf{X}\\right)\\Phi\\left(\\mathbf{X}\\right) +\\alpha I \\right)\\hat{\\mathbf{W}} = \\Phi^{\\top}\\left(\\mathbf{X}\\right)\\mathbf{Y}\n", "$$\n", "Here $\\Phi\\left(\\mathbf{X}\\right)$ is the mathematical representation of *data_matrix*,\n", "$\\mathbf{Y}$ is the mathematical representation of *data_outputs* and $\\alpha$ is the mathematical representation of *regularisation*, while $\\hat{\\mathbf{W}}$ is a mathematical representation for coefficients of the ridge regression." ] }, { "cell_type": "code", "execution_count": 5, "id": "1f617f3a", "metadata": {}, "outputs": [], "source": [ "def ridge_regression(data_matrix, data_outputs, regularisation=0):\n", " return np.linalg.solve(data_matrix.T @ data_matrix +regularisation * np.identity(len(data_matrix.T)),\n", " data_matrix.T @ data_outputs)" ] }, { "cell_type": "markdown", "id": "90653ea1", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 6, "id": "43665749", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhYAAAGdCAYAAABO2DpVAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAe30lEQVR4nO3df5DU9X348ddxJ3tKuW0xAke4UWRM9SSaCMGAWptUiMYhsTNJjdUMpuY7xdKIMkkqNQ1e2/RibR2bJpKGGGIGUSc/SGWiGNIW1KhBfnUk50QrRLDehSJ198Ry1rvP9w+Hqxfu4PbuvXssPB4z+8d97rPsy/egn6efz2d3a7IsywIAIIFRIz0AAHDsEBYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJBMXaVfsKenJ15++eUYO3Zs1NTUVPrlAYAhyLIsOjs7Y9KkSTFq1MDnJSoeFi+//HI0NTVV+mUBgAR2794dkydPHvD3FQ+LsWPHRsRbgzU0NFT65QGAISgWi9HU1NR7HB9IxcPi4OWPhoYGYQEAVeZItzG4eRMASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMhX/gCwAIL3uniw27twXezoPxPix9TFzyrioHVX57+QqKSzefPPNuPXWW+Pee++Njo6OaGxsjGuvvTa+8IUvHPYLSQCA8lm7vT1a1rRFe+FA77bGfH0sndccl05rrOgsJYXFbbfdFl//+tfjnnvuibPPPjs2bdoUn/rUpyKfz8eiRYvKNSMAMIC129vj+pVbIvu17R2FA3H9yi2x7JrzKhoXJYXFk08+GR/96Efj8ssvj4iI0047Le67777YtGlTWYYDAAbW3ZNFy5q2Q6IiIiKLiJqIaFnTFnOaJ1bsskhJ1y8uvPDC+Jd/+Zd47rnnIiLi3//93+Pxxx+PD3/4wwM+p6urK4rFYp8HADB8G3fu63P549dlEdFeOBAbd+6r2EwlnbH4sz/7sygUCnHmmWdGbW1tdHd3x5e+9KW46qqrBnxOa2trtLS0DHtQAKCvPZ0DR8VQ9kuhpDMWDzzwQKxcuTJWrVoVW7ZsiXvuuSf+7u/+Lu65554Bn7NkyZIoFAq9j927dw97aAAgYvzY+qT7pVDSGYvPfe5zcfPNN8cnPvGJiIh497vfHS+++GK0trbG/Pnz+31OLpeLXC43/EkBgD5mThkXjfn66Cgc6Pc+i5qImJh/662nlVLSGYvXX3/9kLeV1tbWRk9PT9KhAIAjqx1VE0vnNUfEWxHxdgd/XjqvuaKfZ1FSWMybNy++9KUvxY9+9KP45S9/GatXr4477rgjfv/3f79c8wEAh3HptMZYds15MTHf93LHxHx9xd9qGhFRk2VZf2dP+tXZ2Rl/8Rd/EatXr449e/bEpEmT4qqrroovfvGLMXr06EH9GcViMfL5fBQKhWhoaBjy4ADA/yn3J28O9vhdUlikICwAoPoM9vjtc7gBgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACCZksLitNNOi5qamkMeCxcuLNd8AEAVqStl56effjq6u7t7f96+fXvMmTMnPv7xjycfDACoPiWFxSmnnNLn5y9/+csxderUuPjii5MOBQBUp5LC4u3eeOONWLlyZSxevDhqamoG3K+rqyu6urp6fy4Wi0N9SQDgKDfkmzd/+MMfxquvvhrXXnvtYfdrbW2NfD7f+2hqahrqSwIAR7maLMuyoTzxQx/6UIwePTrWrFlz2P36O2PR1NQUhUIhGhoahvLSAECFFYvFyOfzRzx+D+lSyIsvvhg/+clP4gc/+MER983lcpHL5YbyMgBAlRnSpZAVK1bE+PHj4/LLL089DwBQxUoOi56enlixYkXMnz8/6uqGfO8nAHAMKjksfvKTn8SuXbvij/7oj8oxDwBQxUo+5TB37twY4v2eAMAxzneFAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMiWHxX/+53/GNddcEyeffHKcdNJJ8Z73vCc2b95cjtkAgCpTV8rO//3f/x0XXHBBfOADH4iHH344xo8fHy+88EL85m/+ZpnGAwCqSUlhcdttt0VTU1OsWLGid9tpp52WeiYAoEqVdCnkwQcfjBkzZsTHP/7xGD9+fLz3ve+N5cuXH/Y5XV1dUSwW+zwAgGNTSWGxY8eOWLZsWZxxxhnxyCOPxIIFC+KGG26I73znOwM+p7W1NfL5fO+jqalp2EMDAEenmizLssHuPHr06JgxY0Y88cQTvdtuuOGGePrpp+PJJ5/s9zldXV3R1dXV+3OxWIympqYoFArR0NAwjNEBgEopFouRz+ePePwu6YxFY2NjNDc399l21llnxa5duwZ8Ti6Xi4aGhj4PAODYVFJYXHDBBfGLX/yiz7bnnnsuTj311KRDAQDVqaSwuOmmm+Kpp56Kv/mbv4n/+I//iFWrVsU3vvGNWLhwYbnmAwCqSElh8b73vS9Wr14d9913X0ybNi3+6q/+Ku688864+uqryzUfAFBFSrp5M4XB3vwBABw9ynLzJgDA4QgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSqRvpAQAYed09WWzcuS/2dB6I8WPrY+aUcVE7qmakx6IKlRQWt956a7S0tPTZNmHChOjo6Eg6FACVs3Z7e7SsaYv2woHebY35+lg6rzkundY4gpNRjUq+FHL22WdHe3t77+OZZ54px1wAVMDa7e1x/cotfaIiIqKjcCCuX7kl1m5vH6HJqFYlXwqpq6uLiRMnlmMWACqouyeLljVtkfXzuywiaiKiZU1bzGme6LIIg1byGYvnn38+Jk2aFFOmTIlPfOITsWPHjsPu39XVFcVisc8DgJG3cee+Q85UvF0WEe2FA7Fx577KDUXVKykszj///PjOd74TjzzySCxfvjw6Ojpi9uzZ8corrwz4nNbW1sjn872PpqamYQ8NwPDt6Rw4KoayH0RE1GRZ1t9ZsEHZv39/TJ06NT7/+c/H4sWL+92nq6srurq6en8uFovR1NQUhUIhGhoahvrSAAzTky+8Elctf+qI+933/94fs6aeXIGJOJoVi8XI5/NHPH4P6+2mY8aMiXe/+93x/PPPD7hPLpeLXC43nJcBoAxmThkXjfn66Cgc6Pc+i5qImJh/662nMFjD+oCsrq6uePbZZ6Ox0duRAKpN7aiaWDqvOSLeioi3O/jz0nnNbtykJCWFxWc/+9nYsGFD7Ny5M372s5/Fxz72sSgWizF//vxyzQdAGV06rTGWXXNeTMzX99k+MV8fy645z+dYULKSLoW89NJLcdVVV8XevXvjlFNOife///3x1FNPxamnnlqu+QAos0unNcac5ok+eZMkhnXz5lAM9uYPAODoMdjjty8hAwCSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACRTN9IDAMeG7p4sNu7cF3s6D8T4sfUxc8q4qB1VM9JjARU2rDMWra2tUVNTEzfeeGOicYBqtHZ7e1x427/GVcufikX3b4urlj8VF972r7F2e/tIjwZU2JDD4umnn45vfOMbcc4556ScB6gya7e3x/Urt0R74UCf7R2FA3H9yi3iAo4zQwqL1157La6++upYvnx5/NZv/VbqmYAq0d2TRcuatsj6+d3BbS1r2qK7p789gGPRkMJi4cKFcfnll8cll1xyxH27urqiWCz2eQDHho079x1ypuLtsohoLxyIjTv3VW4oYESVfPPm/fffH1u2bImnn356UPu3trZGS0tLyYMBR789nQNHxVD2A6pfSWcsdu/eHYsWLYqVK1dGfX39oJ6zZMmSKBQKvY/du3cPaVDg6DN+7OD+OzDY/YDqV9IZi82bN8eePXti+vTpvdu6u7vj0Ucfja9+9avR1dUVtbW1fZ6Ty+Uil8ulmRY4qsycMi4a8/XRUTjQ730WNRExMf/WW0+B40NJZyx+7/d+L5555pnYtm1b72PGjBlx9dVXx7Zt2w6JCuDYVjuqJpbOa46ItyLi7Q7+vHRes8+zgONISWcsxo4dG9OmTeuzbcyYMXHyyScfsh04Plw6rTGWXXNetKxp63Mj58R8fSyd1xyXTmscwemASvPJm8CwXTqtMeY0T/TJm8Dww2L9+vUJxgCqXe2ompg19eSRHgMYYb6EDABIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkikpLJYtWxbnnHNONDQ0RENDQ8yaNSsefvjhcs0GAFSZksJi8uTJ8eUvfzk2bdoUmzZtig9+8IPx0Y9+NH7+85+Xaz4AoIrUZFmWDecPGDduXNx+++1x3XXXDWr/YrEY+Xw+CoVCNDQ0DOelAYAKGezxu26oL9Dd3R3f/e53Y//+/TFr1qwB9+vq6oqurq4+gwEAx6aSb9585pln4jd+4zcil8vFggULYvXq1dHc3Dzg/q2trZHP53sfTU1NwxoYADh6lXwp5I033ohdu3bFq6++Gt///vfjm9/8ZmzYsGHAuOjvjEVTU5NLIQBQRQZ7KWTY91hccsklMXXq1Pinf/qnpIMBAEePwR6/h/05FlmW9TkjAQAcv0q6efPP//zP47LLLoumpqbo7OyM+++/P9avXx9r164t13wAQBUpKSx+9atfxSc/+clob2+PfD4f55xzTqxduzbmzJlTrvkAgCpSUljcfffd5ZoDADgG+K4QACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIpqSwaG1tjfe9730xduzYGD9+fFxxxRXxi1/8olyzAQBVpqSw2LBhQyxcuDCeeuqpWLduXbz55psxd+7c2L9/f7nmAwCqSE2WZdlQn/xf//VfMX78+NiwYUP8zu/8zqCeUywWI5/PR6FQiIaGhqG+NABQQYM9ftcN50UKhUJERIwbN27Afbq6uqKrq6vPYADAsWnIN29mWRaLFy+OCy+8MKZNmzbgfq2trZHP53sfTU1NQ31JAOAoN+RLIQsXLowf/ehH8fjjj8fkyZMH3K+/MxZNTU0uhQBAFSnrpZDPfOYz8eCDD8ajjz562KiIiMjlcpHL5YbyMgBAlSkpLLIsi8985jOxevXqWL9+fUyZMqVccwEAVaiksFi4cGGsWrUq/vmf/znGjh0bHR0dERGRz+fjxBNPLMuAAED1KOkei5qamn63r1ixIq699tpB/RnebgoA1acs91gM4yMvAIDjgO8KAQCSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACRTN9IDQEREd08WG3fuiz2dB2L82PqYOWVc1I6qGemxAChRyWHx6KOPxu233x6bN2+O9vb2WL16dVxxxRVlGI3jxdrt7dGypi3aCwd6tzXm62PpvOa4dFrjCE4GQKlKvhSyf//+OPfcc+OrX/1qOebhOLN2e3tcv3JLn6iIiOgoHIjrV26JtdvbR2gyAIai5DMWl112WVx22WXlmIXjTHdPFi1r2iLr53dZRNRERMuatpjTPNFlEYAqUfabN7u6uqJYLPZ5QETExp37DjlT8XZZRLQXDsTGnfsqNxQAw1L2sGhtbY18Pt/7aGpqKvdLUiX2dA4cFUPZD4CRV/awWLJkSRQKhd7H7t27y/2SVInxY+uT7gfAyCv7201zuVzkcrlyvwxVaOaUcdGYr4+OwoF+77OoiYiJ+bfeegpAdfABWYyY2lE1sXRec0S8FRFvd/DnpfOa3bgJUEVKDovXXnsttm3bFtu2bYuIiJ07d8a2bdti165dqWfjOHDptMZYds15MTHf93LHxHx9LLvmPJ9jAVBlarIs6+8s9IDWr18fH/jABw7ZPn/+/Pj2t799xOcXi8XI5/NRKBSioaGhlJfmGOaTNwGOboM9fpd8j8Xv/u7vRoktAkdUO6omZk09eaTHAGCY3GMBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASEZYAADJCAsAIBlhAQAkIywAgGSEBQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAywgIASKZupAdIobsni40798WezgMxfmx9zJwyLmpH1Yz0WABw3BlSWNx1111x++23R3t7e5x99tlx5513xkUXXZR6tkFZu709Wta0RXvhQO+2xnx9LJ3XHJdOaxyRmQDgeFXypZAHHnggbrzxxrjlllti69atcdFFF8Vll10Wu3btKsd8h7V2e3tcv3JLn6iIiOgoHIjrV26JtdvbKz4TABzPSg6LO+64I6677rr49Kc/HWeddVbceeed0dTUFMuWLSvHfAPq7smiZU1bZP387uC2ljVt0d3T3x4AQDmUFBZvvPFGbN68OebOndtn+9y5c+OJJ57o9zldXV1RLBb7PFLYuHPfIWcq3i6LiPbCgdi4c1+S1wMAjqyksNi7d290d3fHhAkT+myfMGFCdHR09Puc1tbWyOfzvY+mpqahT/s2ezoHjoqh7AcADN+Q3m5aU9P3HRdZlh2y7aAlS5ZEoVDofezevXsoL3mI8WPrk+4HAAxfSe8Kecc73hG1tbWHnJ3Ys2fPIWcxDsrlcpHL5YY+4QBmThkXjfn66Cgc6Pc+i5qImJh/662nAEBllHTGYvTo0TF9+vRYt25dn+3r1q2L2bNnJx3sSGpH1cTSec0R8VZEvN3Bn5fOa/Z5FgBQQSVfClm8eHF885vfjG9961vx7LPPxk033RS7du2KBQsWlGO+w7p0WmMsu+a8mJjve7ljYr4+ll1zns+xAIAKK/kDsq688sp45ZVX4i//8i+jvb09pk2bFg899FCceuqp5ZjviC6d1hhzmif65E0AOArUZFlW0Q96KBaLkc/no1AoRENDQyVfGgAYosEev30JGQCQjLAAAJIRFgBAMsICAEhGWAAAyQgLACAZYQEAJCMsAIBkhAUAkEzJH+k9XAc/6LNYLFb6pQGAITp43D7SB3ZXPCw6OzsjIqKpqanSLw0ADFNnZ2fk8/kBf1/x7wrp6emJl19+OcaOHRs1Nem+KKxYLEZTU1Ps3r3bd5CUkXWuHGtdGda5MqxzZZRznbMsi87Ozpg0aVKMGjXwnRQVP2MxatSomDx5ctn+/IaGBn9pK8A6V461rgzrXBnWuTLKtc6HO1NxkJs3AYBkhAUAkMwxExa5XC6WLl0auVxupEc5plnnyrHWlWGdK8M6V8bRsM4Vv3kTADh2HTNnLACAkScsAIBkhAUAkIywAACSqaqwuOuuu2LKlClRX18f06dPj8cee+yw+2/YsCGmT58e9fX1cfrpp8fXv/71Ck1a3UpZ5x/84AcxZ86cOOWUU6KhoSFmzZoVjzzySAWnrV6l/n0+6Kc//WnU1dXFe97znvIOeAwpda27urrilltuiVNPPTVyuVxMnTo1vvWtb1Vo2upV6jrfe++9ce6558ZJJ50UjY2N8alPfSpeeeWVCk1bnR599NGYN29eTJo0KWpqauKHP/zhEZ9T8WNhViXuv//+7IQTTsiWL1+etbW1ZYsWLcrGjBmTvfjii/3uv2PHjuykk07KFi1alLW1tWXLly/PTjjhhOx73/tehSevLqWu86JFi7Lbbrst27hxY/bcc89lS5YsyU444YRsy5YtFZ68upS6zge9+uqr2emnn57NnTs3O/fccyszbJUbylp/5CMfyc4///xs3bp12c6dO7Of/exn2U9/+tMKTl19Sl3nxx57LBs1alT2D//wD9mOHTuyxx57LDv77LOzK664osKTV5eHHnoou+WWW7Lvf//7WURkq1evPuz+I3EsrJqwmDlzZrZgwYI+284888zs5ptv7nf/z3/+89mZZ57ZZ9sf//EfZ+9///vLNuOxoNR17k9zc3PW0tKSerRjylDX+corr8y+8IUvZEuXLhUWg1TqWj/88MNZPp/PXnnllUqMd8wodZ1vv/327PTTT++z7Stf+Uo2efLkss14rBlMWIzEsbAqLoW88cYbsXnz5pg7d26f7XPnzo0nnnii3+c8+eSTh+z/oQ99KDZt2hT/+7//W7ZZq9lQ1vnX9fT0RGdnZ4wbN64cIx4ThrrOK1asiBdeeCGWLl1a7hGPGUNZ6wcffDBmzJgRf/u3fxvvfOc7413veld89rOfjf/5n/+pxMhVaSjrPHv27HjppZfioYceiizL4le/+lV873vfi8svv7wSIx83RuJYWPEvIRuKvXv3Rnd3d0yYMKHP9gkTJkRHR0e/z+no6Oh3/zfffDP27t0bjY2NZZu3Wg1lnX/d3//938f+/fvjD/7gD8ox4jFhKOv8/PPPx8033xyPPfZY1NVVxb+2R4WhrPWOHTvi8ccfj/r6+li9enXs3bs3/uRP/iT27dvnPosBDGWdZ8+eHffee29ceeWVceDAgXjzzTfjIx/5SPzjP/5jJUY+bozEsbAqzlgc9Otfs55l2WG/er2//fvbTl+lrvNB9913X9x6663xwAMPxPjx48s13jFjsOvc3d0df/iHfxgtLS3xrne9q1LjHVNK+Tvd09MTNTU1ce+998bMmTPjwx/+cNxxxx3x7W9/21mLIyhlndva2uKGG26IL37xi7F58+ZYu3Zt7Ny5MxYsWFCJUY8rlT4WVsX/+rzjHe+I2traQ8p3z549h5TYQRMnTux3/7q6ujj55JPLNms1G8o6H/TAAw/EddddF9/97nfjkksuKeeYVa/Ude7s7IxNmzbF1q1b40//9E8j4q2DX5ZlUVdXFz/+8Y/jgx/8YEVmrzZD+Tvd2NgY73znO/t8PfRZZ50VWZbFSy+9FGeccUZZZ65GQ1nn1tbWuOCCC+Jzn/tcREScc845MWbMmLjooovir//6r51VTmQkjoVVccZi9OjRMX369Fi3bl2f7evWrYvZs2f3+5xZs2Ydsv+Pf/zjmDFjRpxwwgllm7WaDWWdI946U3HttdfGqlWrXB8dhFLXuaGhIZ555pnYtm1b72PBggXx27/927Ft27Y4//zzKzV61RnK3+kLLrggXn755Xjttdd6tz333HMxatSomDx5clnnrVZDWefXX389Ro3qewiqra2NiP/7P2qGb0SOhWW7LTSxg29luvvuu7O2trbsxhtvzMaMGZP98pe/zLIsy26++ebsk5/8ZO/+B99ic9NNN2VtbW3Z3Xff7e2mg1DqOq9atSqrq6vLvva1r2Xt7e29j1dffXWk/hGqQqnr/Ou8K2TwSl3rzs7ObPLkydnHPvax7Oc//3m2YcOG7Iwzzsg+/elPj9Q/QlUodZ1XrFiR1dXVZXfddVf2wgsvZI8//ng2Y8aMbObMmSP1j1AVOjs7s61bt2Zbt27NIiK74447sq1bt/a+rfdoOBZWTVhkWZZ97Wtfy0499dRs9OjR2XnnnZdt2LCh93fz58/PLr744j77r1+/Pnvve9+bjR49OjvttNOyZcuWVXji6lTKOl988cVZRBzymD9/fuUHrzKl/n1+O2FRmlLX+tlnn80uueSS7MQTT8wmT56cLV68OHv99dcrPHX1KXWdv/KVr2TNzc3ZiSeemDU2NmZXX3119tJLL1V46uryb//2b4f9b+7RcCz0tekAQDJVcY8FAFAdhAUAkIywAACSERYAQDLCAgBIRlgAAMkICwAgGWEBACQjLACAZIQFAJCMsAAAkhEWAEAy/x8UaK42XGMjCwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "test_inputs = np.array([0, 1 / 4, 1 / 2, 3 / 4, 1])\n", "test_outputs = np.array([0, 1, 4, 5, 8])\n", "plt.scatter(test_inputs,test_outputs)\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 7, "id": "54568bca", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhYAAAGdCAYAAABO2DpVAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAiCklEQVR4nO3df3BU1f3/8ddmAxul2W1BQxayX0BGqyGiAqKAqT8KCDKIZbBKwUGrn0+xVIOMWqhWTFsbLS0DVsWKijqAMmJiZVQU+ykYfyC/O2IctRIl4MYUKLsRm7VuzvePnUTWJJDdnN3NTZ6PmfvHnj2X++ZM4L5y7r3nuowxRgAAABZkZboAAADQdRAsAACANQQLAABgDcECAABYQ7AAAADWECwAAIA1BAsAAGANwQIAAFiTne4DNjY26rPPPlNubq5cLle6Dw8AAJJgjFF9fb369eunrKy25yXSHiw+++wzBQKBdB8WAABYUFNTo4KCgja/T3uwyM3NlRQrzOv1pvvwAAAgCeFwWIFAoPk83pa0B4umyx9er5dgAQCAwxzvNgZu3gQAANYQLAAAgDUECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgDcECAABYk/YFsgAAgH3RRqMt1YdUV9+gvNwcjRzUW+6s9L+TK6Fg8fXXX+vuu+/WqlWrVFtbK7/fr2uvvVZ33nnnMV9IAgAAUmf97qBK11UpGGpobvP7crRwcqEmFPnTWktCweK+++7Tww8/rCeffFJDhgzRtm3bdN1118nn86mkpCRVNQIAgDas3x3UjSt3yHyrvTbUoBtX7tCymcPSGi4SChZvv/22pkyZokmTJkmSBg4cqKefflrbtm1LSXEAAKBt0Uaj0nVVLUKFJBlJLkml66o0rjA/bZdFErp+ccEFF+hvf/ubPvzwQ0nSP/7xD73xxhu67LLL2twnEokoHA7HbQAAoOO2VB+Ku/zxbUZSMNSgLdWH0lZTQjMWv/zlLxUKhXT66afL7XYrGo3qnnvu0fTp09vcp6ysTKWlpR0uFAAAxKurbztUJNPPhoRmLNasWaOVK1dq9erV2rFjh5588kn98Y9/1JNPPtnmPgsWLFAoFGreampqOlw0AACQ8nJzrPazIaEZi9tuu03z58/X1VdfLUk688wz9emnn6qsrEyzZs1qdR+PxyOPx9PxSgEAQJyRg3rL78tRbaih1fssXJLyfbFHT9MloRmLL7/8ssVjpW63W42NjVaLAgAAx+fOcmnh5EJJsRBxtKbPCycXpnU9i4SCxeTJk3XPPffoxRdf1CeffKKKigotXrxYP/rRj1JVHwAAOIYJRX4tmzlM+b74yx35vpy0P2oqSS5jTGuzJ62qr6/Xr3/9a1VUVKiurk79+vXT9OnTddddd6lnz57t+jPC4bB8Pp9CoZC8Xm/ShQMAgG+keuXN9p6/EwoWNhAsAABwnvaev1mHGwAAWEOwAAAA1hAsAACANQQLAABgDcECAABYQ7AAAADWECwAAIA1BAsAAGANwQIAAFhDsAAAANYQLAAAgDUECwAAYE12pgsAAAAWRKNSZaUUDEp+v1RcLLndaS+DYAEAgNOVl0slJdK+fd+0FRRIS5dKU6emtRQuhQAA4GTl5dK0afGhQpL274+1l5entRyCBQAAThWNxmYqjGn5XVPb3LmxfmlCsAAAwKkqK1vOVBzNGKmmJtYvTQgWAAA4VTBot58FBAsAAJzK77fbzwKCBQAATlVcHHv6w+Vq/XuXSwoEYv3ShGABAIBTud2xR0qlluGi6fOSJWldz4JgAQCAk02dKq1dK/XvH99eUBBrT/M6FiyQBQCA002dKk2ZwsqbAADAErdbuuiiTFfBpRAAAGAPwQIAAFhDsAAAANYQLAAAgDUECwAAYA3BAgAAWEOwAAAA1iQULAYOHCiXy9VimzNnTqrqAwAADpLQAllbt25VNBpt/rx7926NGzdOV155pfXCAACA8yQULE4++eS4z/fee68GDx6sCy+80GpRAADAmZJe0vurr77SypUrNW/ePLnael2rpEgkokgk0vw5HA4ne0gAANDJJX3z5vPPP6/Dhw/r2muvPWa/srIy+Xy+5i0QCCR7SAAA0Mm5jDEmmR0vvfRS9ezZU+vWrTtmv9ZmLAKBgEKhkLxebzKHBgAAaRYOh+Xz+Y57/k7qUsinn36q1157TeXl5cft6/F45PF4kjkMAABwmKQuhaxYsUJ5eXmaNGmS7XoAAICDJRwsGhsbtWLFCs2aNUvZ2Unf+wkAALqghIPFa6+9pr179+qnP/1pKuoBAAAOlvCUw/jx45Xk/Z4AAKCL410hAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKxJOFjs379fM2fOVJ8+fXTiiSfq7LPP1vbt21NRGwAAcJjsRDr/+9//1pgxY3TxxRfr5ZdfVl5enj7++GN997vfTVF5AADASRIKFvfdd58CgYBWrFjR3DZw4EDbNQEAAIdK6FLICy+8oBEjRujKK69UXl6ezjnnHC1fvvyY+0QiEYXD4bgNAAB0TQkFiz179mjZsmU69dRT9corr2j27Nm6+eab9dRTT7W5T1lZmXw+X/MWCAQ6XDQAAOicXMYY097OPXv21IgRI/TWW281t918883aunWr3n777Vb3iUQiikQizZ/D4bACgYBCoZC8Xm8HSgcAAOkSDofl8/mOe/5OaMbC7/ersLAwru2MM87Q3r1729zH4/HI6/XGbQAAoGtKKFiMGTNGH3zwQVzbhx9+qAEDBlgtCgAAOFNCweKWW27R5s2b9fvf/17//Oc/tXr1aj3yyCOaM2dOquoDAAAOklCwOPfcc1VRUaGnn35aRUVF+u1vf6slS5ZoxowZqaoPAAA4SEI3b9rQ3ps/AABA55GSmzcBAACOhWABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsyc50AQCATiAalSorpWBQ8vul4mLJ7c50VXAgggUAdHfl5VJJibRv3zdtBQXS0qXS1KmZqwuOxKUQAOjOysuladPiQ4Uk7d8fay8vz0xdcCyCBQB0V9FobKbCmJbfNbXNnRvrB7QTwQIAuqvKypYzFUczRqqpifUD2imhYHH33XfL5XLFbfn5+amqDQCQSsGg3X6Akrh5c8iQIXrttdeaP7u5axgAnMnvt9sPUBLBIjs7m1kKAOgComMu0AHvSTo5fKDV6etGSXW+k3XymAvEr5Bor4Tvsfjoo4/Ur18/DRo0SFdffbX27NlzzP6RSEThcDhuAwBk3pa9Id11yf9KioWIozV9Xnjx/2jL3lBa64KzJRQszjvvPD311FN65ZVXtHz5ctXW1mr06NE6ePBgm/uUlZXJ5/M1b4FAoMNFAwA6rq6+Qa98f7RuvOJXqs09Ke672tyTdOMVv9Ir3x+tuvqGDFUIJ3IZ09pzRu1z5MgRDR48WLfffrvmzZvXap9IJKJIJNL8ORwOKxAIKBQKyev1JntoAEAHvf3xQU1fvlmSlNUY1ch97ynvi3+r7jvf05aCIWrMil0Aefp/zteowX0yWSo6gXA4LJ/Pd9zzd4dW3uzVq5fOPPNMffTRR2328Xg88ng8HTkMACAFRg7qLb8vR7WhBjVmubX5/w2N+94lKd+Xo5GDememQDhSh9axiEQiev/99+XnjmEAcBx3lksLJxdKioWIozV9Xji5UO6sb38LtC2hYHHrrbdq06ZNqq6u1jvvvKNp06YpHA5r1qxZqaoPAJBCE4r8WjZzmPJ9OXHt+b4cLZs5TBOK+MURiUnoUsi+ffs0ffp0HThwQCeffLLOP/98bd68WQMGDEhVfQCAFJtQ5Ne4wnxtqT6kuvoG5eXGLn8wU4FkdOjmzWS09+YPAADQebT3/M27QgAAgDUECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgDcECAABYQ7AAAADWECwAAIA1BAsAAGANwQIAAFhDsAAAANYQLAAAgDUECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgDcECAABYQ7AAAADWECwAAIA1BAsAAGANwQIAAFhDsAAAANYQLAAAgDUECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgTXamCwDQNUQbjbZUH1JdfYPycnM0clBvubNcmS4LQJp1aMairKxMLpdLc+fOtVQOACdavzuoC+77P01fvlklz+zS9OWbdcF9/6f1u4OZLg1AmiUdLLZu3apHHnlEQ4cOtVkPAIdZvzuoG1fuUDDUENdeG2rQjSt3EC6AbiapYPHFF19oxowZWr58ub73ve/ZrgmAQ0QbjUrXVcm08l1TW+m6KkUbW+sBoCtKKljMmTNHkyZN0tixY4/bNxKJKBwOx20AuoYt1YdazFQczUgKhhq0pfpQ+ooCkFEJ37z5zDPPaMeOHdq6dWu7+peVlam0tDThwgB0fnX1bYeKZPoBcL6EZixqampUUlKilStXKicnp137LFiwQKFQqHmrqalJqlAAnU9ebvv+H2hvPwDOl9CMxfbt21VXV6fhw4c3t0WjUb3++ut64IEHFIlE5Ha74/bxeDzyeDx2qgXQqYwc1Ft+X45qQw2t3mfhkpTviz16CqB7SGjG4oc//KHeffdd7dq1q3kbMWKEZsyYoV27drUIFQC6NneWSwsnF0qKhYijNX1eOLmQ9SyAbiShGYvc3FwVFRXFtfXq1Ut9+vRp0Q6ge5hQ5NeymcNUuq4q7kbOfF+OFk4u1IQifwarA5BurLwJoMMmFPk1rjCflTcBdDxYbNy40UIZAJzOneXSqMF9Ml0GgAzjJWQAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsyc50AQC6iGhUqqyUgkHJ75eKiyW3O9NVAUgzggWAjisvl0pKpH37vmkrKJCWLpWmTs1cXQDSjkshADqmvFyaNi0+VEjS/v2x9vLyzNQFICMIFgCSF43GZiqMafldU9vcubF+ALoFggWA5FVWtpypOJoxUk1NrB+AboFgASB5waDdfgAcj2ABIHl+v91+AByPYAEgecXFsac/XK7Wv3e5pEAg1g9At5BQsFi2bJmGDh0qr9crr9erUaNG6eWXX05VbQA6O7c79kip1DJcNH1esoT1LIBuJKFgUVBQoHvvvVfbtm3Ttm3bdMkll2jKlCl67733UlUfgM5u6lRp7Vqpf//49oKCWDvrWADdisuY1p4Ta7/evXtr0aJFuv7669vVPxwOy+fzKRQKyev1duTQADoTVt4EurT2nr+TXnkzGo3q2Wef1ZEjRzRq1Kg2+0UiEUUikbjCAHRBbrd00UWZrgJAhiV88+a7776r73znO/J4PJo9e7YqKipUWFjYZv+ysjL5fL7mLRAIdKhgAADQeSV8KeSrr77S3r17dfjwYT333HN69NFHtWnTpjbDRWszFoFAgEshAAA4SHsvhXT4HouxY8dq8ODB+stf/mK1MAAA0Hm09/zd4XUsjDFxMxIAAKD7SujmzV/96leaOHGiAoGA6uvr9cwzz2jjxo1av359quoDAAAOklCw+Pzzz3XNNdcoGAzK5/Np6NChWr9+vcaNG5eq+gAAgIMkFCwee+yxVNUBAAC6AN4VAgAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMCahIJFWVmZzj33XOXm5iovL09XXHGFPvjgg1TVBgAAHCahYLFp0ybNmTNHmzdv1oYNG/T1119r/PjxOnLkSKrqAwAADuIyxphkd/7Xv/6lvLw8bdq0ST/4wQ/atU84HJbP51MoFJLX60320AAAII3ae/7O7shBQqGQJKl3795t9olEIopEInGFAQCArinpmzeNMZo3b54uuOACFRUVtdmvrKxMPp+veQsEAskeEgAAdHJJXwqZM2eOXnzxRb3xxhsqKChos19rMxaBQIBLIQAAOEhKL4XcdNNNeuGFF/T6668fM1RIksfjkcfjSeYwAADAYRIKFsYY3XTTTaqoqNDGjRs1aNCgVNUFAAAcKKFgMWfOHK1evVp//etflZubq9raWkmSz+fTCSeckJICAQCAcyR0j4XL5Wq1fcWKFbr22mvb9WfwuCkAAM6TknssOrDkBQAA6AZ4VwgAALCGYAEAAKwhWAAAAGs6tKQ3YE00KlVWSsGg5PdLxcWS253pqgAACSJYIPPKy6WSEmnfvm/aCgqkpUulqVMzVxcAIGFcCkFmlZdL06bFhwpJ2r8/1l5enpm6AABJIVggc6LR2ExFa48xN7XNnRvrBwBwBIIFMqeysuVMxdGMkWpqYv0AAI5AsEDmBIN2+wEAMo5ggczx++32AwBkHMECmVNcHHv6o4130MjlkgKBWD8AgCMQLJA5bnfskVKpZbho+rxkCetZAICDECyQWVOnSmvXSv37x7cXFMTaWccCAByFBbKQeVOnSlOmsPImAHQBBAt0Dm63dNFFma4CANBBXAoBAADWECwAAIA1BAsAAGANwQIAAFhDsAAAANYQLAAAgDUECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgDcECAABYQ7AAAADWECwAAIA1BAsAAGBNdqYLACQp2mi0pfqQ6uoblJebo5GDesud5cp0WQCABCUcLF5//XUtWrRI27dvVzAYVEVFha644ooUlIbuYv3uoErXVSkYamhu8/tytHByoSYU+TNYGQAgUQlfCjly5IjOOussPfDAA6moB93M+t1B3bhyR1yokKTaUINuXLlD63cHM1QZACAZCc9YTJw4URMnTkxFLehmoo1GpeuqZFr5zkhySSpdV6VxhflcFgEAh0j5zZuRSEThcDhuAyRpS/WhFjMVRzOSgqEGbak+lL6iAAAdkvJgUVZWJp/P17wFAoFUHxIOUVffdqhIph8AIPNSHiwWLFigUCjUvNXU1KT6kHCIvNwcq/0AAJmX8sdNPR6PPB5Pqg8DBxo5qLf8vhzVhhpavc/CJSnfF3v0FADgDCyQhYxxZ7m0cHKhpFiIOFrT54WTC7lxEwAcJOFg8cUXX2jXrl3atWuXJKm6ulq7du3S3r17bdeGbmBCkV/LZg5Tvi/+cke+L0fLZg5jHQsAcBiXMaa1Weg2bdy4URdffHGL9lmzZumJJ5447v7hcFg+n0+hUEherzeRQ6MLY+VNAOjc2nv+Tvgei4suukgJZhHguNxZLo0a3CfTZQAAOoh7LAAAgDUECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgDcECAABYQ7AAAADWpPwlZGkRjUqVlVIwKPn9UnGx5HZnuioAALod5weL8nKppETat++btoICaelSaerUzNUFAEA35OxLIeXl0rRp8aFCkvbvj7WXl2emLgAAuinnBotoNDZT0dp7S5ra5s6N9QMAAGnh3GBRWdlypuJoxkg1NbF+AAAgLZwbLIJBu/0AAECHOTdY+P12+wEAgA5zbrAoLo49/eFytf69yyUFArF+AAAgLZwbLNzu2COlUstw0fR5yRLWswAAII2cGyyk2DoVa9dK/fvHtxcUxNpZxwIAgLRy/gJZU6dKU6aw8iYAAJ2A84OFFAsRF12U6SoAAOj2nH0pBAAAdCoECwAAYA3BAgAAWEOwAAAA1hAsAACANQQLAABgDcECAABYQ7AAAADWECwAAIA1BAsAAGBNl1jSO9potKX6kOrqG5SXm6ORg3rLndXG69QBAEDKJBUsHnroIS1atEjBYFBDhgzRkiVLVFxcbLu2dlm/O6jSdVUKhhqa2/y+HC2cXKgJRf6M1AQAQHeV8KWQNWvWaO7cubrjjju0c+dOFRcXa+LEidq7d28q6jum9buDunHljrhQIUm1oQbduHKH1u8Opr0mAAC6s4SDxeLFi3X99dfrhhtu0BlnnKElS5YoEAho2bJlqaivTdFGo9J1VTKtfNfUVrquStHG1noAAIBUSChYfPXVV9q+fbvGjx8f1z5+/Hi99dZbre4TiUQUDofjNhu2VB9qMVNxNCMpGGrQlupDVo4HAACOL6FgceDAAUWjUfXt2zeuvW/fvqqtrW11n7KyMvl8vuYtEAgkX+1R6urbDhXJ9AMAAB2X1OOmLlf8ExfGmBZtTRYsWKBQKNS81dTUJHPIFvJyc6z2AwAAHZfQUyEnnXSS3G53i9mJurq6FrMYTTwejzweT/IVtmHkoN7y+3JUG2po9T4Ll6R8X+zRUwAAkB4JzVj07NlTw4cP14YNG+LaN2zYoNGjR1st7HjcWS4tnFwoKRYijtb0eeHkQtazAAAgjRK+FDJv3jw9+uijevzxx/X+++/rlltu0d69ezV79uxU1HdME4r8WjZzmPJ98Zc78n05WjZzGOtYAACQZgkvkHXVVVfp4MGD+s1vfqNgMKiioiK99NJLGjBgQCrqO64JRX6NK8xn5U0AADoBlzEmrQs9hMNh+Xw+hUIheb3edB4aAAAkqb3nb15CBgAArCFYAAAAawgWAADAGoIFAACwhmABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKxJeEnvjmpa6DMcDqf70AAAIElN5+3jLdid9mBRX18vSQoEAuk+NAAA6KD6+nr5fL42v0/7u0IaGxv12WefKTc3Vy6XvReFhcNhBQIB1dTU8A6SFGKc04exTg/GOT0Y5/RI5TgbY1RfX69+/fopK6vtOynSPmORlZWlgoKClP35Xq+XH9o0YJzTh7FOD8Y5PRjn9EjVOB9rpqIJN28CAABrCBYAAMCaLhMsPB6PFi5cKI/Hk+lSujTGOX0Y6/RgnNODcU6PzjDOab95EwAAdF1dZsYCAABkHsECAABYQ7AAAADWECwAAIA1jgoWDz30kAYNGqScnBwNHz5clZWVx+y/adMmDR8+XDk5OTrllFP08MMPp6lSZ0tknMvLyzVu3DidfPLJ8nq9GjVqlF555ZU0Vutcif48N3nzzTeVnZ2ts88+O7UFdiGJjnUkEtEdd9yhAQMGyOPxaPDgwXr88cfTVK1zJTrOq1at0llnnaUTTzxRfr9f1113nQ4ePJimap3p9ddf1+TJk9WvXz+5XC49//zzx90n7edC4xDPPPOM6dGjh1m+fLmpqqoyJSUlplevXubTTz9ttf+ePXvMiSeeaEpKSkxVVZVZvny56dGjh1m7dm2aK3eWRMe5pKTE3HfffWbLli3mww8/NAsWLDA9evQwO3bsSHPlzpLoODc5fPiwOeWUU8z48ePNWWedlZ5iHS6Zsb788svNeeedZzZs2GCqq6vNO++8Y9588800Vu08iY5zZWWlycrKMkuXLjV79uwxlZWVZsiQIeaKK65Ic+XO8tJLL5k77rjDPPfcc0aSqaioOGb/TJwLHRMsRo4caWbPnh3Xdvrpp5v58+e32v/22283p59+elzbz372M3P++eenrMauINFxbk1hYaEpLS21XVqXkuw4X3XVVebOO+80CxcuJFi0U6Jj/fLLLxufz2cOHjyYjvK6jETHedGiReaUU06Ja7v//vtNQUFBymrsatoTLDJxLnTEpZCvvvpK27dv1/jx4+Pax48fr7feeqvVfd5+++0W/S+99FJt27ZN//3vf1NWq5MlM87f1tjYqPr6evXu3TsVJXYJyY7zihUr9PHHH2vhwoWpLrHLSGasX3jhBY0YMUJ/+MMf1L9/f5122mm69dZb9Z///CcdJTtSMuM8evRo7du3Ty+99JKMMfr888+1du1aTZo0KR0ldxuZOBem/SVkyThw4ICi0aj69u0b1963b1/V1ta2uk9tbW2r/b/++msdOHBAfr8/ZfU6VTLj/G1/+tOfdOTIEf34xz9ORYldQjLj/NFHH2n+/PmqrKxUdrYj/tl2CsmM9Z49e/TGG28oJydHFRUVOnDggH7+85/r0KFD3GfRhmTGefTo0Vq1apWuuuoqNTQ06Ouvv9bll1+uP//5z+koudvIxLnQETMWTb79mnVjzDFfvd5a/9baES/RcW7y9NNP6+6779aaNWuUl5eXqvK6jPaOczQa1U9+8hOVlpbqtNNOS1d5XUoiP9ONjY1yuVxatWqVRo4cqcsuu0yLFy/WE088wazFcSQyzlVVVbr55pt11113afv27Vq/fr2qq6s1e/bsdJTaraT7XOiIX31OOukkud3uFsm3rq6uRRJrkp+f32r/7Oxs9enTJ2W1Olky49xkzZo1uv766/Xss89q7NixqSzT8RId5/r6em3btk07d+7UL37xC0mxk58xRtnZ2Xr11Vd1ySWXpKV2p0nmZ9rv96t///5xr4c+44wzZIzRvn37dOqpp6a0ZidKZpzLyso0ZswY3XbbbZKkoUOHqlevXiouLtbvfvc7ZpUtycS50BEzFj179tTw4cO1YcOGuPYNGzZo9OjRre4zatSoFv1fffVVjRgxQj169EhZrU6WzDhLsZmKa6+9VqtXr+b6aDskOs5er1fvvvuudu3a1bzNnj1b3//+97Vr1y6dd9556SrdcZL5mR4zZow+++wzffHFF81tH374obKyslRQUJDSep0qmXH+8ssvlZUVfwpyu92SvvmNGh2XkXNhym4LtazpUabHHnvMVFVVmblz55pevXqZTz75xBhjzPz5880111zT3L/pEZtbbrnFVFVVmccee4zHTdsh0XFevXq1yc7ONg8++KAJBoPN2+HDhzP1V3CERMf523gqpP0SHev6+npTUFBgpk2bZt577z2zadMmc+qpp5obbrghU38FR0h0nFesWGGys7PNQw89ZD7++GPzxhtvmBEjRpiRI0dm6q/gCPX19Wbnzp1m586dRpJZvHix2blzZ/NjvZ3hXOiYYGGMMQ8++KAZMGCA6dmzpxk2bJjZtGlT83ezZs0yF154YVz/jRs3mnPOOcf07NnTDBw40CxbtizNFTtTIuN84YUXGkkttlmzZqW/cIdJ9Of5aASLxCQ61u+//74ZO3asOeGEE0xBQYGZN2+e+fLLL9NctfMkOs7333+/KSwsNCeccILx+/1mxowZZt++fWmu2ln+/ve/H/P/3M5wLuS16QAAwBpH3GMBAACcgWABAACsIVgAAABrCBYAAMAaggUAALCGYAEAAKwhWAAAAGsIFgAAwBqCBQAAsIZgAQAArCFYAAAAawgWAADAmv8PbTqXZ6dcIFYAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "test_degree = 2\n", "regularisation = 0.2\n", "test_data_matrix = ridge_regression_data(test_inputs, test_degree)\n", "weights=ridge_regression(test_data_matrix, test_outputs, regularisation)\n", "prediction=test_data_matrix@weights\n", "\n", "plt.scatter(test_inputs,test_outputs)\n", "plt.scatter(test_inputs,prediction,color='Red')\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "42a2cbe8", "metadata": {}, "source": [ "3. Write a function **prediction_function** that evaluates your predicted regression function at the given points $\\mathbf{X} = \n", "\\left\\{x^{(1)}, x^{(2)}, \\ldots, x^{(s)}\\right\\}$ for given coefficients $\\mathbf{W} = \\left(w^{(0)}, w^{(1)},\\ldots,w^{(d)}\\right)$. The function **prediction_function** takes the arguments *data_matrix* and *weights* as inputs and returns a value of the regression function evaluated for every $x \\in \\mathbf{X}$ via\n", "$$\n", "f_{\\mathbf{W}}\\left(\\mathbf{x}\\right)\n", "= w^{(0)}+w^{(1)}x_1+\\ldots+w^{(d)}x_d,\n", "$$\n", "or\n", "$$\n", "f_{\\mathbf{W}}\\left(x\\right)\n", "= w^{(0)}+w^{(1)}x+\\ldots+w^{(d)}x^d,\n", "$$\n", "where _data_matrix_ is a NumPy array representing data matrix $\\Phi\\left(\\mathbf{X}\\right)$, _weights_ is a NumPy representation of coefficients vector $\\mathbf{W}$. The function should return a vector of the regression function values $\\left(f_{\\mathbf{W}}\\left(x^{(1)}\\right), f_{\\mathbf{w}}\\left(x^{(2)}\\right), \\ldots, f_{\\mathbf{w}}\\left(x^{(s)}\\right)\\right)^{\\top}$.\n", "\n", "*Hint:* consider matrix $\\Phi\\left(\\mathbf{X}\\right)\\mathbf{W}$." ] }, { "cell_type": "code", "execution_count": 8, "id": "9a9ba731", "metadata": {}, "outputs": [], "source": [ "def prediction_function(data_matrix, weights):\n", " return data_matrix @ weights" ] }, { "cell_type": "markdown", "id": "117c1544", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 9, "id": "f4fa7a25", "metadata": {}, "outputs": [], "source": [ "test_inputs = np.array([[1, 2, 3]])\n", "test_weights = np.array([[-1, 1], [2, 2], [-3, 3], [4, 4]])\n", "test_data_matrix = ridge_regression_data(test_inputs)\n", "assert_array_almost_equal(prediction_function(test_data_matrix, test_weights),\n", " np.array([[7, 21]]))" ] }, { "cell_type": "markdown", "id": "4797723a", "metadata": {}, "source": [ "4. Write a function **prediction_error** that evaluates a mean-squared error over the set of data inputs and outputs. The function **prediction_error** takes the arguments _data_matrix_, _data_outputs_ and _weights_ as inputs and returns a mean squared error defined by\n", "$$\n", "\\mathrm{MSE} = \\frac{1}{2s} \\left\\|\\Phi\\left(\\mathbf{X}\\right)\\mathbf{W} - \\mathbf{Y} \\right\\|^2,\n", "$$\n", "where $\\Phi\\left(\\mathbf{X}\\right)$ is a mathematical representation of _data_matrix_, $\\mathbf{Y}$ is a mathematical representation of _data_outputs_ and $\\mathbf{W}$ is a mathematical representation of _weights_." ] }, { "cell_type": "code", "execution_count": 10, "id": "e78d17b4", "metadata": {}, "outputs": [], "source": [ "def prediction_error(data_matrix, data_outputs, weights):\n", " return 1 / (2 * len(data_matrix)) * (np.linalg.norm(data_matrix @ weights -data_outputs))**2" ] }, { "cell_type": "markdown", "id": "3b899f10", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 11, "id": "ff2177aa", "metadata": {}, "outputs": [], "source": [ "test_inputs = np.array([0, 1 / 4, 1 / 2, 3 / 4, 1])\n", "test_data_outputs = np.array([0, 1, 0, -1, 0])\n", "test_weights = np.array([2 / 5, -4 / 5, 6 / 5])\n", "test_data_matrix = ridge_regression_data(test_inputs, len(test_weights) - 1)\n", "\n", "assert_array_almost_equal(\n", " prediction_error(test_data_matrix, test_data_outputs, test_weights),\n", " 0.359125)" ] }, { "cell_type": "code", "execution_count": 12, "id": "dc7c6452", "metadata": {}, "outputs": [], "source": [ "test_inputs = np.array([[1, -1], [2, 2]])\n", "test_data_outputs = np.array([[-1, 2], [1, 3]])\n", "test_data_matrix = ridge_regression_data(test_inputs)\n", "test_weights = np.array([[0, 0], [1, 2], [3, 4]])\n", "assert_array_almost_equal(\n", " prediction_error(test_data_matrix, test_data_outputs, test_weights), 36.75)" ] }, { "cell_type": "markdown", "id": "c3243f38", "metadata": {}, "source": [ "### K-fold cross validation \n", "In this section you are asked to implement a number of functions that are used for a K-fold cross validation strategy as introduced in the lecture. This is a strategy for a calculation of a validation error using a smart data splitting and can be described as follows:\n", "1. split the data, joined inputs and outputs, into $K$ approximately equal chunks. Let us call them $D_1, D_2, \\ldots, D_K$;\n", "2. for every $i = 1,\\ldots,K$ evaluate optimal weights/coefficients $\\mathbf{W}_i$ for the ridge regression evaluated over the data set $D_1,D_2,\\ldots,D_{i-1},D_{i+1},\\ldots,D_K$ and a corresponding validation error $L_i$ evaluated over the set $D_i$;\n", "3. evaluate average of optimal weights \n", "$\\hat{\\mathbf{W}} = \\frac{1}{K}\\left(\\mathbf{W}_1+\\mathbf{W}_2+\\ldots+\\mathbf{W}_K\\right)$ and an average validation error \n", "$L = \\frac{1}{K}\\left(L_1+L_2+\\ldots+L_K\\right)$." ] }, { "cell_type": "markdown", "id": "b4afb1d3", "metadata": {}, "source": [ "1. Write a function **KFold_split** that takes two arguments *data_size* and *K* and outputs a random split of integer indexes $\\left[0,1,\\ldots,data\\_size\\right)$ into $K$ almost equal chunks. For example, for $K=2$ and $data\\_size = 5$\n", "it may return $[[3,0],[2,4,1]]$." ] }, { "cell_type": "code", "execution_count": 13, "id": "9c4e10e3", "metadata": {}, "outputs": [], "source": [ "def KFold_split(data_size, K):\n", " indexes = np.random.permutation(data_size)\n", " m, r = divmod(data_size, K) # m is the quotient and r the reminder of the division data_size/K\n", " indexes_split = [indexes[i * m + min(i, r):(i + 1) * m + min(i + 1, r)] for i in range(K)]\n", "\n", " return indexes_split" ] }, { "cell_type": "code", "execution_count": 14, "id": "82b318fe", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3 1\n", "---\n", "0 4\n", "4 7\n", "7 10\n" ] } ], "source": [ "# let's understand a bit better what is going on in the for loop\n", "K=3\n", "m, r = divmod(10, K)\n", "print (m,r)\n", "print ('---')\n", "for i in range(K):\n", " print (i * m + min(i, r),(i + 1) * m + min(i + 1, r))\n", "# hence this effectively tries to split the data into equally sized chuncks" ] }, { "cell_type": "markdown", "id": "79f8c85d", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 15, "id": "828a3fec", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "5 1\n", "[array([3, 2, 1, 4, 0])]\n" ] } ], "source": [ "test_data_size = np.random.randint(low=1, high=10)\n", "test_K = np.random.randint(low=1, high=5)\n", "print (test_data_size,test_K )\n", "indexes_split = KFold_split(test_data_size, test_K)\n", "print (indexes_split)\n", "data_indexes = np.array([])\n", "for i in range(test_K):\n", " data_indexes = np.append(data_indexes, indexes_split[i])\n", "data_indexes = np.sort(data_indexes)\n", "assert_array_almost_equal(data_indexes, np.array(range(test_data_size)))" ] }, { "cell_type": "markdown", "id": "d44452c9", "metadata": {}, "source": [ "2. Write a function **KFold_cross_validation** that takes $5$ arguments\n", "- *data_matrix* - a data matrix containing all inputs in an appropriate form (see **Ridge regression** section)\n", "- *data_outputs* - a data matrix containing all outputs in a matrix form\n", "- *K* - a positive integer number representing the number of chunks used for a validation algorithm\n", "- *model_evaluation* - a lambda-function that takes two parameters *data_matrix* and *data_outputs*, and evaluates optimal weights/coefficients/parameters of some ML model\n", "- *error_evaluation* - a lambda-function that takes three parameters *data_matrix*, *data_outputs*, and *weights* and evaluates a validation error.\n", "\n", "For an example of two last functions see the test below.\n", "\n", "The function **KFold_cross_validation** should return a matrix/vector of coefficients/weights and a validation error that are both evaluated as averages of optimal weights and validation error of every iteration step." ] }, { "cell_type": "code", "execution_count": 16, "id": "c93c0d0d", "metadata": {}, "outputs": [], "source": [ "def KFold_cross_validation(data_matrix, data_outputs, K, model_evaluation,error_evaluation):\n", "\n", " data_size = len(data_matrix)\n", " indexes_split = KFold_split(data_size, K)\n", " # each indexes_split[j] is an array with the ids of the data matrix for that split\n", " # all the others ids can be used for training and the j for validation\n", " for i in range(K):\n", " # hence this concatenates all the ids of the data matrix EXCEPT for those in the i split\n", " # those are used for validation!\n", " indexes = np.concatenate([indexes_split[j] for j in range(K) if (j != i)])\n", " \n", " # training the model in all the data EXCEPT for the i split\n", " weights = model_evaluation(data_matrix[indexes], data_outputs[indexes]) \n", " # the error is computed OUT OF SAMPLE in the i split\n", " error = error_evaluation(data_matrix[indexes_split[i]],data_outputs[indexes_split[i]], weights)\n", " # then the overall performance is measured as the average of the error for each split\n", " # and the weights are also averaged\n", " if (i == 0):\n", " optimal_weights = weights / K\n", " validation_error = error / K\n", " else:\n", " optimal_weights += weights / K\n", " validation_error += error / K\n", " \n", " return optimal_weights, validation_error" ] }, { "cell_type": "markdown", "id": "dee0adfb", "metadata": {}, "source": [ "Let's see how to use it" ] }, { "cell_type": "code", "execution_count": 17, "id": "239dcd42", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(array([ 0.37306987, -0.42083743, -0.47982548]), 0.48945014701218525)" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "model_evaluation = lambda data_matrix, data_outputs: ridge_regression(data_matrix, data_outputs, regularisation=0.1)\n", "error_evaluation = lambda data_matrix, data_outputs, weights: prediction_error(data_matrix, data_outputs, weights)\n", "\n", "inputs = np.array([0, 1 / 4, 1 / 2, 3 / 4, 1])\n", "data_outputs = np.array([0, 1, 0, -1, 0])\n", "\n", "data_matrix=ridge_regression_data(inputs, 2)\n", "\n", "KFold_cross_validation(data_matrix, data_outputs, 3, model_evaluation,error_evaluation)" ] }, { "cell_type": "markdown", "id": "5a702cf0", "metadata": {}, "source": [ "### Data standardisation \n", "In real-world problems we usually get a raw data in the form of $s$ samples each of which is described by numeric several values corresponding to different characteristics of the object. Such a data could be highly non-uniform. The goal of applying $\\textit{standardisation}$ is to make sure different features of objects are on almost on the same scale so that each feature is equally important and make it easier to process by most ML algorithms. The result of standardisation is that the features will be rescaled to ensure the mean and the standard deviation to be $0$ and $1$, respectively. This means that for a data given by $\\mathbf{X} = \\left(\n", "\\left(\\mathbf{x}^{(1)}\\right)^{\\top},\\left(\\mathbf{x}^{(2)}\\right)^{\\top},\\ldots,\\left(\\mathbf{x}^{(s)}\\right)^{\\top}\n", "\\right) \\in \\mathbb{R}^{s\\times d}$ we define a new, rescaled data as:\n", "$$\n", "\\hat{\\mathbf{x}}^{(i)}_k = \\frac{\\mathbf{x}^{(i)}_k - \\left\\langle \\mathbf{x}_k \\right\\rangle }{\\left(\\sigma_{\\mathbf{x}}\\right)_k},\n", "$$\n", "where $\\left\\langle \\mathbf{x}_k \\right\\rangle = \\frac{1}{s}\\sum\\limits_{j=1}^s \\mathbf{x}^{(j)}_k$, and\n", "$\\left(\\sigma_\\mathbf{x}\\right)_k = \\sqrt{\n", "\t\\frac{1}{s}\\sum\\limits_{j=1}^s \\left(\\mathbf{x}^{(j)}_k-\\left\\langle \\mathbf{x}_k \\right\\rangle\\right)^2}$\n", "are the mean and standard deviation of data vector $\\mathbf{x}$. \n", "\n", "Write two functions \n", "1. **standardise** to standardise the columns of a multi-dimensional array. The function **standardise**\ttakes the multi-dimensional array *data_matrix* as its input argument. It subtracts the means from each column and divides by the standard deviations. It returns the *standardised_matrix*, the *row_of_means* and the *row_of_standard_deviations*.\n", "2. **de_standardise** to de-standardise the columns of a multi-dimensional array. The function **de_standardise** reverses the above operation. It takes a *standardised_matrix*, the *row_of_means* and the *row_of_standard_deviations* as its arguments and returns a matrix for which the standardisation process is reversed." ] }, { "cell_type": "code", "execution_count": 18, "id": "a2271e9b", "metadata": {}, "outputs": [], "source": [ "def standardise(data_matrix):\n", " \n", " row_of_means = np.mean(data_matrix, axis=0)\n", " standardised_matrix = data_matrix - row_of_means\n", " row_of_stds = np.std(data_matrix, axis=0)\n", " \n", " return (standardised_matrix / row_of_stds), row_of_means, row_of_stds" ] }, { "cell_type": "markdown", "id": "1d224ad6", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 19, "id": "b4b23d2a", "metadata": {}, "outputs": [], "source": [ "test_data_matrix = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])\n", "test_standardise_data_matrix = np.array([[-1.224745, -1.224745, -1.224745],\n", " [0., 0., 0.],\n", " [1.224745, 1.224745, 1.224745]])\n", "test_row_of_means = np.array([4, 5, 6])\n", "test_row_of_stds = np.array(np.sqrt([6, 6, 6]))\n", "\n", "test_result_standardise_data_matrix, test_result_row_of_means, test_result_row_of_stds = standardise(\n", " test_data_matrix)\n", "assert_array_almost_equal(test_result_standardise_data_matrix,\n", " test_standardise_data_matrix)\n", "assert_array_almost_equal(test_result_row_of_means, test_row_of_means)\n", "assert_array_almost_equal(test_result_row_of_stds, test_row_of_stds)" ] }, { "cell_type": "code", "execution_count": 20, "id": "7522f7b3", "metadata": {}, "outputs": [], "source": [ "def de_standardise(standardised_matrix, row_of_means, row_of_stds):\n", " matrix = np.copy(standardised_matrix * row_of_stds)\n", " return matrix + row_of_means" ] }, { "cell_type": "markdown", "id": "a87eedee", "metadata": {}, "source": [ "### Boston housing price data\n", "Finally, we are going to apply all our tools for training a model predicting housing prices in Boston. In this exercise you will work with a real housing price data. The assignment folder contains $\\texttt{house_prices.csv}$ file which you will need to read the data from. This file contains an information about $N = 1200$ houses. The data columns are:\n", "- $\\texttt{StreetLength}$ - length of the street in front of the building\n", "- $\\texttt{Area}$ - total area of the lot\n", "- $\\texttt{Quality}$ - quality of building materials\n", "- $\\texttt{Condition}$ - condition of the building\n", "- $\\texttt{BasementArea}$ - area of the basement\n", "- $\\texttt{LivingArea}$ - total living area\n", "- $\\texttt{GarageArea}$ - a garage area\t\t\n", "- $\\texttt{SalePrice}$ - sale price\n", "\n", "Your task would be to build a ridge regression using $K$-fold cross validation strategy for validation and a grid search strategy for optimisation of the hyperparameter $\\alpha$ value.\n", "\n", "1. We start by loading the data. Please run the below cell to read the data from a .csv file. Make sure the .csv file is in the same folder with your Jupyter notebook." ] }, { "cell_type": "code", "execution_count": 21, "id": "2417729a", "metadata": {}, "outputs": [], "source": [ "housing_dataset_path = \"house_prices.csv\"\n", "housing_data = np.genfromtxt(housing_dataset_path,\n", " delimiter=\",\",\n", " skip_header=1,\n", " usecols=[0, 1, 2, 3, 4, 5, 6, 7])\n", "housing_data_input = housing_data[:, 0:7]\n", "housing_data_output = housing_data[:, 7]" ] }, { "cell_type": "markdown", "id": "3cf7d0bd", "metadata": {}, "source": [ "2. We then prepare the data for an analysis. You will need to standardise the inputs and outputs using the functions you developed before." ] }, { "cell_type": "code", "execution_count": 22, "id": "1e2033d2", "metadata": {}, "outputs": [], "source": [ "standardised_housing_data_input, means_housing_data_input,stds_housing_data_input = standardise(housing_data_input)\n", "standardised_housing_data_output, means_housing_data_output,stds_housing_data_output = standardise(housing_data_output)" ] }, { "cell_type": "markdown", "id": "708ff1b9", "metadata": {}, "source": [ "3. Write a function **grid_search** that performs a search for a minimum value of a given function on a given grid points. You function should have a signature *grid_search (objective, grid)*, where *objective* is a lambda-function to minimise and *grid* is a list of grid points. The function should return the grid point with the minimal value of objective function." ] }, { "cell_type": "code", "execution_count": 23, "id": "4795b97b", "metadata": {}, "outputs": [], "source": [ "def grid_search(objective, grid):\n", " values = np.array([])\n", " for point in grid:\n", " values = np.append(values, objective(point))\n", " return grid[np.argmin(values)]" ] }, { "cell_type": "markdown", "id": "2c7460dc", "metadata": {}, "source": [ "Test your function with the following unit tests" ] }, { "cell_type": "code", "execution_count": 24, "id": "2b73fab8", "metadata": {}, "outputs": [], "source": [ "test_function = lambda xy: xy[0]**2 + xy[1]**2 - 2 * xy[0] * xy[1]\n", "test_grid = [(0, 1), (1, 0), (2, 3), (5, 5)]\n", "assert_array_almost_equal(grid_search(test_function, test_grid), (5, 5))" ] }, { "cell_type": "markdown", "id": "c5d2894a", "metadata": {}, "source": [ "4. Implement a grid search algorithm to find an unknown hyperparameter $\\hat \\alpha$ such that\n", "$$\n", "\\hat{\\alpha} = \\arg\\min\\limits_{\\alpha\\geq 0} \\mathrm{Val}\\left( \\mathbf{W}_{\\alpha}\\right),\n", "$$\n", "where the validation error $\\mathrm{Val}\\left( \\mathbf{W}_{\\alpha}\\right)$ is evaluated using K-fold cross validation. Take $K=5$. Print out an optimal coefficients." ] }, { "cell_type": "code", "execution_count": 27, "id": "4f641d72", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "An optimal value of regularisation parameter is 11.0.\n", "For this value ofregularisation parameter one gets optimal weights of the form \n", "[-0.00875437 0.07058742 0.45124836 0.03148828 0.14961819 0.2659526\n", " 0.15911522]\n" ] } ], "source": [ "K = 5\n", "alpha_grid = np.append(np.array([i * 0.05 for i in range(20)]),\n", " np.array([i for i in range(1, 20)]))\n", "\n", "\n", "error_evaluation = lambda data_matrix, data_outputs, weights: prediction_error(data_matrix, data_outputs, weights)\n", "\n", "# the objective function is the validation error computed for each alpha value in the grid below\n", "validation_error = lambda alpha: KFold_cross_validation(standardised_housing_data_input,\n", " standardised_housing_data_output, \n", " K,lambda data_matrix,\n", " data_outputs: ridge_regression(data_matrix,\n", " data_outputs,\n", " regularisation=alpha),\n", " error_evaluation)[1]\n", "\n", "\n", "\n", "\n", "optimal_alpha = grid_search(validation_error, alpha_grid)\n", "\n", "optimal_weights = KFold_cross_validation(standardised_housing_data_input,\n", " standardised_housing_data_output, \n", " K,lambda data_matrix,\n", " data_outputs: ridge_regression(data_matrix, \n", " data_outputs,\n", " regularisation =optimal_alpha),\n", " error_evaluation)[0]\n", "\n", "print(\n", "\"An optimal value of regularisation parameter is {}.\\nFor this value ofregularisation parameter one gets optimal weights of the form \\n{}\".format(optimal_alpha, optimal_weights))" ] }, { "cell_type": "code", "execution_count": null, "id": "f7a0f5e1", "metadata": {}, "outputs": [], "source": [ "# let us get the names of the features from the data\n", "# they are the header of the file\n", "a=open(housing_dataset_path,'r')\n", "for i in a:\n", " features=i.strip().split(',')\n", " break\n", "a.close()\n", "# we need to delete the last since that is the target!\n", "features.pop()\n", "\n", "# now we can get the order of the weights in terms of their values\n", "c=0\n", "list_features=[]\n", "for i in optimal_weights:\n", " list_features.append([i,features[c]])\n", " c+=1\n", "\n", "list_features=sorted(list_features,reverse=True)\n", "for i in list_features:\n", " print (i)" ] }, { "cell_type": "markdown", "id": "3e36dbeb", "metadata": {}, "source": [ "## Bootstrapping" ] }, { "cell_type": "code", "execution_count": null, "id": "070ee6ee", "metadata": {}, "outputs": [], "source": [ "# to bootstrap we need to select N samples from the data with replacement\n", "# this means that we can just get N random numbers from 0 to len(data)-1 \n", "# indeed in this type of sampling the probability that one will be selected is constant 1/len(data)\n", "# we can pick N as 90% of the total sample size\n", "import random as rd\n", "data_size=len(standardised_housing_data_output)\n", "fraction=0.9\n", "samples_size=int(data_size*fraction)\n", "\n", "sample_input_list=[]\n", "sample_output_list=[]\n", "for i in range(samples_size):\n", " id_random=rd.randint(0,samples_size-1)\n", " sample_input_list.append(standardised_housing_data_input[id_random])\n", " sample_output_list.append(standardised_housing_data_output[id_random])\n", " \n", "sample_input=np.array(sample_input_list)\n", "sample_output=np.array(sample_output_list)\n", "\n", "# now can can create a function that does this plus for each sample compute the regression\n" ] }, { "cell_type": "code", "execution_count": null, "id": "063db0f5", "metadata": {}, "outputs": [], "source": [ "import random as rd\n", "def bootstrap_regression(standardised_data_input,standardised_data_output,fraction,M,alpha):\n", " # first we need to know what is N: the number of samples to extract\n", " data_size=len(standardised_data_output)\n", " samples_size=int(data_size*fraction)\n", "\n", " w_list=[]\n", " # then for each of the M extract\n", " for j in range(M):\n", " sample_input_list=[]\n", " sample_output_list=[]\n", " for i in range(samples_size):\n", " # we take N samples extract random numbers which are the id of the arrays that store the data\n", " id_random=rd.randint(0,samples_size-1)\n", " # note that we need to keep the X and Y correspondence hence the id_random is the same for each\n", " sample_input_list.append(standardised_data_input[id_random])\n", " sample_output_list.append(standardised_data_output[id_random])\n", "\n", " # convert the list to arrays\n", " sample_input=np.array(sample_input_list)\n", " sample_output=np.array(sample_output_list)\n", "\n", " # apply the regression, note that alpha is selected before with the model selection\n", " weights=ridge_regression(sample_input, sample_output, regularisation=alpha)\n", " # append the fitted values of Ws for the N samples in list\n", " w_list.append(weights)\n", " return w_list" ] }, { "cell_type": "code", "execution_count": null, "id": "fb17a546", "metadata": {}, "outputs": [], "source": [ "M=10000\n", "# note how the alpha is selected in the model selection before!\n", "w_list=bootstrap_regression(standardised_housing_data_input,standardised_housing_data_output,0.9,M,optimal_alpha)\n", "\n", "# now we have a list with all the values of w for each of the M samples\n", "# we can create a dictionary to get the values for each features\n", "dict_w={}\n", "for i in range(len(features)):\n", " dict_w[i]=[]\n", "\n", "# append the values to each feature as list\n", "for i in range(M):\n", " w=w_list[i]\n", " for j in range(len(w)):\n", " dict_w[j].append(w[j])\n", " \n", "for i in dict_w:\n", " dict_w[i]=np.array(dict_w[i]) # first convert the list into a nunpy array then we can manipulate them\n", " print (features[i],np.mean(dict_w[i])) # get for example the mean\n" ] }, { "cell_type": "code", "execution_count": null, "id": "4718691f", "metadata": {}, "outputs": [], "source": [ "# we can then use seaborn to plot, this requires pandas\n", "# let us first create the data for the df\n", "# this function takes the features and the list of w from the bootstrap, it write the data into a file for simplicity\n", "import pandas as pd\n", "\n", "def get_df(features,w_list):\n", " # let's first create the headers of the file\n", " a=open('data_for_plot.csv','w')\n", " for j in range(len(features)-1):\n", " a.write(features[j]+',')\n", " a.write(features[-1]+'\\n')\n", "\n", " # then we can write the data\n", " for j in w_list:\n", " for i in range(len(j)-1):\n", " a.write(str(j[i])+',')\n", " a.write(str(j[-1])+'\\n')\n", " a.close()\n", " df=pd.read_csv('data_for_plot.csv')\n", " \n", " return df" ] }, { "cell_type": "code", "execution_count": null, "id": "ee10efdb", "metadata": {}, "outputs": [], "source": [ "import seaborn as sns\n", "df=get_df(features,w_list) # get the dataframe\n", "\n", "\n", "sns.boxplot(data=df) # make the plot with seaborn\n", "plt.xticks(rotation=45)\n", "plt.ylabel('w')\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": null, "id": "a90cf58e", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.8.17" } }, "nbformat": 4, "nbformat_minor": 5 }