{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## Uczenie maszynowe UMZ 2020/2021\n",
"### 10 i 17 marca 2021\n",
"# 2. Regresja liniowa"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## 2.1. Funkcja kosztu"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"### Zadanie\n",
"Znając $x$ – ludność miasta (w dziesiątkach tysięcy mieszkańców),\n",
"należy przewidzieć $y$ – dochód firmy transportowej (w dziesiątkach tysięcy dolarów).\n",
"\n",
"(Dane pochodzą z kursu „Machine Learning”, Andrew Ng, Coursera)."
]
},
{
"cell_type": "code",
"execution_count": 48,
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"outputs": [],
"source": [
"import numpy as np\n",
"import matplotlib\n",
"import matplotlib.pyplot as pl\n",
"import ipywidgets as widgets\n",
"\n",
"%matplotlib inline\n",
"%config InlineBackend.figure_format = 'svg'\n",
"\n",
"from IPython.display import display, Math, Latex"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"### Dane"
]
},
{
"cell_type": "code",
"execution_count": 49,
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"6.1101,17.592\n",
"\n",
"5.5277,9.1302\n",
"\n",
"8.5186,13.662\n",
"\n",
"7.0032,11.854\n",
"\n",
"5.8598,6.8233\n",
"\n",
"8.3829,11.886\n",
"\n",
"7.4764,4.3483\n",
"\n",
"8.5781,12\n",
"\n",
"6.4862,6.5987\n",
"\n",
"5.0546,3.8166\n",
"\n"
]
}
],
"source": [
"with open('data01_train.csv') as data:\n",
" for line in data.readlines()[:10]:\n",
" print(line)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"### Wczytanie danych"
]
},
{
"cell_type": "code",
"execution_count": 50,
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"x = [6.1101, 5.5277, 8.5186, 7.0032, 5.8598, 8.3829, 7.4764, 8.5781, 6.4862, 5.0546]\n",
"y = [17.592, 9.1302, 13.662, 11.854, 6.8233, 11.886, 4.3483, 12.0, 6.5987, 3.8166]\n"
]
}
],
"source": [
"import csv\n",
"\n",
"reader = csv.reader(open('data01_train.csv'), delimiter=',')\n",
"\n",
"x = list()\n",
"y = list()\n",
"for xi, yi in reader:\n",
" x.append(float(xi))\n",
" y.append(float(yi)) \n",
" \n",
"print('x = {}'.format(x[:10])) \n",
"print('y = {}'.format(y[:10]))"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"### Hipoteza i parametry modelu"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Jak przewidzieć $y$ na podstawie danego $x$? W celu odpowiedzi na to pytanie będziemy starać się znaleźć taką funkcję $h(x)$, która będzie najlepiej obrazować zależność między $x$ a $y$, tj. $y \\sim h(x)$.\n",
"\n",
"Zacznijmy od najprostszego przypadku, kiedy $h(x)$ jest po prostu funkcją liniową."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Ogólny wzór funkcji liniowej to\n",
"$$ h(x) = a \\, x + b $$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Pamiętajmy jednak, że współczynniki $a$ i $b$ nie są w tej chwili dane z góry – naszym zadaniem właśnie będzie znalezienie takich ich wartości, żeby $h(x)$ było „możliwie jak najbliżej” $y$ (co właściwie oznacza to sformułowanie, wyjaśnię potem).\n",
"\n",
"Poszukiwaną funkcję $h$ będziemy nazywać **funkcją hipotezy**, a jej współczynniki – **parametrami modelu**.\n",
"\n",
"W teorii uczenia maszynowego parametry modelu oznacza się na ogół grecką literą $\\theta$ z odpowiednimi indeksami, dlatego powyższy wzór opisujący liniową funkcję hipotezy zapiszemy jako\n",
"$$ h(x) = \\theta_0 + \\theta_1 x $$\n",
"\n",
"**Parametry modelu** tworzą wektor, który oznaczymy po prostu przez $\\theta$:"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"$$ \\theta = \\left[\\begin{array}{c}\\theta_0\\\\ \\theta_1\\end{array}\\right] $$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Żeby podkreślić fakt, że funkcja hipotezy zależy od parametrów modelu, będziemy pisać $h_\\theta$ zamiast $h$:"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"$$ h_{\\theta}(x) = \\theta_0 + \\theta_1 x $$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Przyjrzyjmy się teraz, jak wyglądają dane, które mamy modelować:"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"Na poniższym wykresie możesz spróbować ręcznie dopasować parametry modelu $\\theta_0$ i $\\theta_1$ tak, aby jak najlepiej modelowały zależność między $x$ a $y$:"
]
},
{
"cell_type": "code",
"execution_count": 51,
"metadata": {
"slideshow": {
"slide_type": "skip"
}
},
"outputs": [
{
"data": {
"image/svg+xml": [
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"# Funkcje rysujące wykres kropkowy oraz prostą regresyjną\n",
"\n",
"def regdots(x, y): \n",
" fig = pl.figure(figsize=(16*.6, 9*.6))\n",
" ax = fig.add_subplot(111)\n",
" fig.subplots_adjust(left=0.1, right=0.9, bottom=0.1, top=0.9)\n",
" ax.scatter(x, y, c='r', s=50, label='Dane')\n",
" \n",
" ax.set_xlabel(u'Wielkość miejscowości [dzies. tys. mieszk.]')\n",
" ax.set_ylabel(u'Dochód firmy [dzies. tys. dolarów]')\n",
" ax.margins(.05, .05)\n",
" pl.ylim(min(y) - 1, max(y) + 1)\n",
" pl.xlim(min(x) - 1, max(x) + 1)\n",
" return fig\n",
"\n",
"def regline(fig, fun, theta, x):\n",
" ax = fig.axes[0]\n",
" x0, x1 = min(x), max(x)\n",
" X = [x0, x1]\n",
" Y = [fun(theta, x) for x in X]\n",
" ax.plot(X, Y, linewidth='2',\n",
" label=(r'$y={theta0}{op}{theta1}x$'.format(\n",
" theta0=theta[0],\n",
" theta1=(theta[1] if theta[1] >= 0 else -theta[1]),\n",
" op='+' if theta[1] >= 0 else '-')))\n",
"\n",
"def legend(fig):\n",
" ax = fig.axes[0]\n",
" handles, labels = ax.get_legend_handles_labels()\n",
" # try-except block is a fix for a bug in Poly3DCollection\n",
" try:\n",
" fig.legend(handles, labels, fontsize='15', loc='lower right')\n",
" except AttributeError:\n",
" pass\n",
"\n",
"fig = regdots(x,y)\n",
"legend(fig)"
]
},
{
"cell_type": "code",
"execution_count": 52,
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"outputs": [],
"source": [
"# Hipoteza: funkcja liniowa jednej zmiennej\n",
"\n",
"def h(theta, x):\n",
" return theta[0] + theta[1] * x"
]
},
{
"cell_type": "code",
"execution_count": 53,
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"outputs": [],
"source": [
"# Przygotowanie interaktywnego wykresu\n",
"\n",
"sliderTheta01 = widgets.FloatSlider(min=-10, max=10, step=0.1, value=0, description=r'$\\theta_0$', width=300)\n",
"sliderTheta11 = widgets.FloatSlider(min=-5, max=5, step=0.1, value=0, description=r'$\\theta_1$', width=300)\n",
"\n",
"def slide1(theta0, theta1):\n",
" fig = regdots(x, y)\n",
" regline(fig, h, [theta0, theta1], x)\n",
" legend(fig)"
]
},
{
"cell_type": "code",
"execution_count": 54,
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "a351866ebbb34b099eade2549da1f5d0",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"interactive(children=(FloatSlider(value=0.0, description='$\\\\theta_0$', max=10.0, min=-10.0), FloatSlider(valu…"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
""
]
},
"execution_count": 54,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"widgets.interact_manual(slide1, theta0=sliderTheta01, theta1=sliderTheta11)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Skąd wiadomo, że przewidywania modelu (wartości funkcji $h(x)$) zgadzaja się z obserwacjami (wartości $y$)?\n",
"\n",
"Aby to zmierzyć wprowadzimy pojęcie funkcji kosztu."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"### Funkcja kosztu"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Funkcję kosztu zdefiniujemy w taki sposób, żeby odzwierciedlała ona różnicę między przewidywaniami modelu a obserwacjami.\n",
"\n",
"Jedną z możliwosci jest zdefiniowanie funkcji kosztu jako wartość **błędu średniokwadratowego** (metoda najmniejszych kwadratów, *root-mean-square error, RMSE*):"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"$$ J(\\theta) \\, = \\, \\frac{1}{2m} \\sum_{i = 1}^{m} \\left( h_{\\theta} \\left( x^{(i)} \\right) - y^{(i)} \\right) ^2 $$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"gdzie $m$ jest liczbą wszystkich przykładów (obserwacji), czyli wielkością zbioru danych uczących.\n",
"\n",
"W powyższym wzorze sumujemy kwadraty różnic między przewidywaniami modelu ($h_\\theta \\left( x^{(i)} \\right)$) a obserwacjami ($y^{(i)}$) po wszystkich przykładach $i$."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"Teraz nasze zadanie sprowadza się do tego, że będziemy szukać takich parametrów $\\theta = \\left[\\begin{array}{c}\\theta_0\\\\ \\theta_1\\end{array}\\right]$, które minimalizują fukcję kosztu $J(\\theta)$:"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"$$ \\hat\\theta = \\mathop{\\arg\\min}_{\\theta} J(\\theta) $$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"$$ \\theta \\in \\mathbb{R}^2, \\quad J \\colon \\mathbb{R}^2 \\to \\mathbb{R} $$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Proszę zwrócić uwagę, że dziedziną funkcji kosztu jest zbiór wszystkich możliwych wartości parametrów $\\theta$."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"$$ J(\\theta_0, \\theta_1) \\, = \\, \\frac{1}{2m} \\sum_{i = 1}^{m} \\left( \\theta_0 + \\theta_1 x^{(i)} - y^{(i)} \\right) ^2 $$"
]
},
{
"cell_type": "code",
"execution_count": 55,
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [],
"source": [
"def J(h, theta, x, y):\n",
" \"\"\"Funkcja kosztu\"\"\"\n",
" m = len(y)\n",
" return 1.0 / (2 * m) * sum((h(theta, x[i]) - y[i])**2 for i in range(m))"
]
},
{
"cell_type": "code",
"execution_count": 56,
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"outputs": [],
"source": [
"# Oblicz wartość funkcji kosztu i pokaż na wykresie\n",
"\n",
"def regline2(fig, fun, theta, xx, yy):\n",
" \"\"\"Rysuj regresję liniową\"\"\"\n",
" ax = fig.axes[0]\n",
" x0, x1 = min(xx), max(xx)\n",
" X = [x0, x1]\n",
" Y = [fun(theta, x) for x in X]\n",
" cost = J(fun, theta, xx, yy)\n",
" ax.plot(X, Y, linewidth='2',\n",
" label=(r'$y={theta0}{op}{theta1}x, \\; J(\\theta)={cost:.3}$'.format(\n",
" theta0=theta[0],\n",
" theta1=(theta[1] if theta[1] >= 0 else -theta[1]),\n",
" op='+' if theta[1] >= 0 else '-',\n",
" cost=cost)))\n",
"\n",
"sliderTheta02 = widgets.FloatSlider(min=-10, max=10, step=0.1, value=0, description=r'$\\theta_0$', width=300)\n",
"sliderTheta12 = widgets.FloatSlider(min=-5, max=5, step=0.1, value=0, description=r'$\\theta_1$', width=300)\n",
"\n",
"def slide2(theta0, theta1):\n",
" fig = regdots(x, y)\n",
" regline2(fig, h, [theta0, theta1], x, y)\n",
" legend(fig)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Poniższy interaktywny wykres pokazuje wartość funkcji kosztu $J(\\theta)$. Czy teraz łatwiej jest dobrać parametry modelu?"
]
},
{
"cell_type": "code",
"execution_count": 57,
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "55119e223ecc4a22b89b532de485a76c",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"interactive(children=(FloatSlider(value=0.0, description='$\\\\theta_0$', max=10.0, min=-10.0), FloatSlider(valu…"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
""
]
},
"execution_count": 57,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"widgets.interact_manual(slide2, theta0=sliderTheta02, theta1=sliderTheta12)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"### Funkcja kosztu jako funkcja zmiennej $\\theta$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"source": [
"Funkcja kosztu zdefiniowana jako RMSE jest funkcją zmiennej wektorowej $\\theta$, czyli funkcją dwóch zmiennych rzeczywistych: $\\theta_0$ i $\\theta_1$.\n",
" \n",
"Zobaczmy, jak wygląda jej wykres."
]
},
{
"cell_type": "code",
"execution_count": 58,
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"outputs": [],
"source": [
"# Wykres funkcji kosztu dla ustalonego theta_1=1.0\n",
"\n",
"def costfun(fun, x, y):\n",
" return lambda theta: J(fun, theta, x, y)\n",
"\n",
"def costplot(hypothesis, x, y, theta1=1.0):\n",
" fig = pl.figure(figsize=(16*.6, 9*.6))\n",
" ax = fig.add_subplot(111)\n",
" fig.subplots_adjust(left=0.1, right=0.9, bottom=0.1, top=0.9)\n",
" ax.set_xlabel(r'$\\theta_0$')\n",
" ax.set_ylabel(r'$J(\\theta)$')\n",
" j = costfun(hypothesis, x, y)\n",
" fun = lambda theta0: j([theta0, theta1])\n",
" X = np.arange(-10, 10, 0.1)\n",
" Y = [fun(x) for x in X]\n",
" ax.plot(X, Y, linewidth='2', label=(r'$J(\\theta_0, {theta1})$'.format(theta1=theta1)))\n",
" return fig\n",
"\n",
"def slide3(theta1):\n",
" fig = costplot(h, x, y, theta1)\n",
" legend(fig)\n",
"\n",
"sliderTheta13 = widgets.FloatSlider(min=-5, max=5, step=0.1, value=1.0, description=r'$\\theta_1$', width=300)"
]
},
{
"cell_type": "code",
"execution_count": 59,
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "4b65a40539ad489796b42d92ce6560c9",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"interactive(children=(FloatSlider(value=1.0, description='$\\\\theta_1$', max=5.0, min=-5.0), Button(description…"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
""
]
},
"execution_count": 59,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"widgets.interact_manual(slide3, theta1=sliderTheta13)"
]
},
{
"cell_type": "code",
"execution_count": 60,
"metadata": {
"slideshow": {
"slide_type": "notes"
}
},
"outputs": [],
"source": [
"# Wykres funkcji kosztu względem theta_0 i theta_1\n",
"\n",
"from mpl_toolkits.mplot3d import Axes3D\n",
"import pylab\n",
"\n",
"%matplotlib inline\n",
"\n",
"def costplot3d(hypothesis, x, y, show_gradient=False):\n",
" fig = pl.figure(figsize=(16*.6, 9*.6))\n",
" ax = fig.add_subplot(111, projection='3d')\n",
" fig.subplots_adjust(left=0.0, right=1.0, bottom=0.0, top=1.0)\n",
" ax.set_xlabel(r'$\\theta_0$')\n",
" ax.set_ylabel(r'$\\theta_1$')\n",
" ax.set_zlabel(r'$J(\\theta)$')\n",
" \n",
" j = lambda theta0, theta1: costfun(hypothesis, x, y)([theta0, theta1])\n",
" X = np.arange(-10, 10.1, 0.1)\n",
" Y = np.arange(-1, 4.1, 0.1)\n",
" X, Y = np.meshgrid(X, Y)\n",
" Z = np.array([[J(hypothesis, [theta0, theta1], x, y) \n",
" for theta0, theta1 in zip(xRow, yRow)] \n",
" for xRow, yRow in zip(X, Y)])\n",
" \n",
" ax.plot_surface(X, Y, Z, rstride=2, cstride=8, linewidth=0.5,\n",
" alpha=0.5, cmap='jet', zorder=0,\n",
" label=r\"$J(\\theta)$\")\n",
" ax.view_init(elev=20., azim=-150)\n",
"\n",
" ax.set_xlim3d(-10, 10);\n",
" ax.set_ylim3d(-1, 4);\n",
" ax.set_zlim3d(-100, 800);\n",
"\n",
" N = range(0, 800, 20)\n",
" pl.contour(X, Y, Z, N, zdir='z', offset=-100, cmap='coolwarm', alpha=1)\n",
" \n",
" ax.plot([-3.89578088] * 2,\n",
" [ 1.19303364] * 2,\n",
" [-100, 4.47697137598], \n",
" color='red', alpha=1, linewidth=1.3, zorder=100, linestyle='dashed',\n",
" label=r'minimum: $J(-3.90, 1.19) = 4.48$')\n",
" ax.scatter([-3.89578088] * 2,\n",
" [ 1.19303364] * 2,\n",
" [-100, 4.47697137598], \n",
" c='r', s=80, marker='x', alpha=1, linewidth=1.3, zorder=100, \n",
" label=r'minimum: $J(-3.90, 1.19) = 4.48$')\n",
" \n",
" if show_gradient:\n",
" ax.plot([3.0, 1.1],\n",
" [3.0, 2.4],\n",
" [263.0, 125.0], \n",
" color='green', alpha=1, linewidth=1.3, zorder=100)\n",
" ax.scatter([3.0],\n",
" [3.0],\n",
" [263.0], \n",
" c='g', s=30, marker='D', alpha=1, linewidth=1.3, zorder=100)\n",
"\n",
" ax.margins(0,0,0)\n",
" fig.tight_layout()"
]
},
{
"cell_type": "code",
"execution_count": 61,
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"outputs": [
{
"data": {
"image/svg+xml": [
"\r\n",
"\r\n",
"\r\n",
"