# Import potrebných balíčkov
import pandas as pd
import numpy as np
import sympy as sp
from sympy.utilities.lambdify import lambdify
import matplotlib.pyplot as plt
from functools import reduce
from scipy.optimize import minimize
from scipy.optimize import curve_fit
# Uistíme sa, že máme všetky potrebné dáta
!mkdir -p data
!wget -nc -O data/sigmoid_regression_data.csv https://www.dropbox.com/s/p5q7gzupa2ndw55/sigmoid_regression_data.csv?dl=1
Ako ďalší príklad si uvedieme veľmi jednoduchú aplikáciu optimalizácie v rámci strojového učenia. Naším cieľom bude vykonať regresiu: dostaneme vstupné a výstupné dáta a budeme sa snažiť nájsť funkciu, ktorá produkuje takú závislosť.
Takúto úlohu je ľahké previesť na optimalizačný problém. Budeme predpokladať, že máme k dispozícii parametrickú funkciu $f_\theta(\mathbf{x})$, ktorej charakter je daný vektorom parametrov $\theta$. Naším cieľom je nájsť také parametre, pri ktorých bude $f_\theta(\mathbf{x})$ na našej dátovej množine robiť čo najmenšie chyby.
Dajme tomu, že dátová množina má vzorky v tvare $(\mathbf{x_i}, \mathbf{y_i})$, kde $\mathbf{x_i}$ je vstup a $\mathbf{y}_i$ je požadovaný výstup pre vzorku $i$. Potom sa dá náš cieľ opísať ako nasledujúci optimalizačný problém:
\begin{equation} \theta^* = \underset{\theta}{\arg\max} \sum_{(\mathbf{x}_i, \mathbf{y}_i)} \mathcal{L}(f_\theta(\mathbf{x}_i), \mathbf{y}_i) \end{equation}t.j. chceme nájsť taký vektor parametrov $\theta^*$, ktorý bude minimalizovať rozdiel medzi skutočnými a požadovanými výstupmi na celej dátovej množine v zmysle nejakej chybovej funkcie: napríklad kvadratickej chyby.
Definíciu nášho regresného problému začnime dátovou množinou, ktorú máme k dispozícii:
df = pd.read_csv("data/sigmoid_regression_data.csv")
X = df['x']
Y = df['y']
df.head()
Keďže dáta sú 2-rozmerné, môžeme si ich pre lepšiu predstavu vizualizovať.
plt.scatter(X, Y)
plt.xlabel('x')
plt.ylabel('y')
plt.grid(ls='--')
Zobrazená závislosť sa tvarom nápadne podobá na sigmoidnú (logistickú) kriku, ktorá je definovaná nasledovne: \begin{equation} \sigma(x) = \frac{1}{1 + e^{-x}}. \end{equation} Zdá sa však, že je trochu posunutá v smere x-ovej osi a možno nemá rovnakú strmosť. Zostavme preto regresný model tak, že vstup sigmoidnej funkcie prejde najprv jednoducho lineárnou transformáciou, ktorej parametre $a$ a $c$ sa naučíme z dát. Celý regresný model potom bude vyzerať nasledovne: \begin{align} u &= ax + c \\ \sigma(u) &= \frac{1}{1 + e^{-u}}. \end{align}
Alebo v rámci jednej funkcie: \begin{equation} \mathrm{f}(x, a, c) = \frac{1}{1 + e^{-ax - c}} \end{equation}
Pomocou balíčka sympy
symbolicky definujte náš regresný model ako funkciu $f(x, a, c)$. Výslednú funkciu vizualizujte (pre parametre $a=1$, $c=0$).
# Sem treba doplniť zdrojový kód.
f =
grad_f =
Vizualizácia regresnej funkciu f(x, a, c)
si kvôli kontrole vizualizujeme:
xx = np.linspace(-5, 5, 100)
a = 1; c = 0
yy = [f(x, a, c) for x in xx]
plt.plot(xx, yy)
plt.xlabel('x')
plt.ylabel('y')
plt.grid(ls='--')
Ako sme už povedali vyššie, naším cieľom je minimalizovať chybu na dátovej množine, t.j.
\begin{equation} \theta^* = \underset{\theta}{\arg\max} \sum_{(\mathbf{x}_i, \mathbf{y}_i)} \mathcal{L}(f_\theta(\mathbf{x}_i), \mathbf{y}_i) \end{equation}Vonkajšej sume sa pri určovaní gradientu nemusíme venovať. Z linearity derivácie vyplýva, že
\begin{equation} \nabla \sum_{(\mathbf{x}_i, \mathbf{y}_i)} \mathcal{L}(f_\theta(\mathbf{x}_i), \mathbf{y}_i) = \sum_{(\mathbf{x}_i, \mathbf{y}_i)} \nabla \mathcal{L}(f_\theta(\mathbf{x}_i), \mathbf{y}_i) \end{equation}stačí nám teda riešiť vnútornú časť sumy a potom gradienty pre jednotlivé vzorky z dátovej množiny sčítať dokopy.
Povedzme, že ako chybovú funkciu použijeme kvadratickú chybu, t.j.:
\begin{equation} \mathcal{L}(f_\theta(\mathbf{x}_i), \mathbf{y}_i) = \left( f_\theta(\mathbf{x}_i) - \mathbf{y}_i \right)^2. \end{equation}Zadefinujme si ju teraz symbolicky a určme gradient (len podľa parametrov $a$ a $c$ pretože tých optimálne hodnoty ideme hľadať):
symy = sp.symbols('y')
symL = (symf - symy)**2
L = lambdify((symx, symy, syma, symc), symL, "numpy")
sym_grad_L = sp.Matrix([symL]).jacobian([syma, symc])
grad_L = lambdify((symx, symy, syma, symc), sym_grad_L, "numpy")
sym_grad_L
Pri takejto definícii do funkcie grad_L
dosadíme za všetky potrebné argumenty ($x$, $y$, $a$, $c$), ale navracia sa len vektor s dvoma prvkami: parciálnou deriváciou podľa $a$ a podľa $c$. Napr.:
x = 0; y = 1; a = 1; c = 0
print("Hodnota chybovej funkcie: {}".format(L(x, y, a, c)))
print("Gradient: {}".format(grad_L(x, y, a, c)))
Ako už vieme, celkovú chybu dostaneme ako súčet chýb pre jednotlivé vzorky a rovnako aj celkový gradient bude súčtom všetkých jednotlivých gradientov. Zadefinujme si teda funkcie, ktoré nám pomôžu oba výsledky vyrátať:
def sumL(a, c):
L_sum = 0
for x, y in zip(X, Y):
L_sum += L(x, y, a, c)
return L_sum
def grad_sumL(a, c):
grad_sum = np.zeros(2)
for x, y in zip(X, Y):
grad_sum = grad_sum + grad_L(x, y, a, c)
return grad_sum
Že všetko správne funguje môžeme otestovať jednoducho:
print(sumL(a, c))
print(grad_sumL(a, c))
Ďalej už len použijeme funkciu minimize
, aby sme chybu minimalizovali a zobrazíme výslednú regresnú závislosť. Minimalizácia sa vykoná nasledovne:
res = minimize(fun=lambda X: sumL(*X),
x0=np.random.uniform(0, 1, 2),
method='L-BFGS-B',
jac=lambda X: grad_sumL(*X)
)
a, c = res.x
Zobrazíme si zvyškovú chybu na dátovej množine:
print("Zvyšková chyba: {}".format(res.fun))
xx = np.linspace(-20, 20, 100)
yy = [f(x, a, c) for x in xx]
plt.scatter(X, Y, c='r', label="pôvodné dáta")
plt.plot(xx, yy, label="regresná závislosť")
plt.xlabel('x')
plt.ylabel('y')
plt.grid(ls='--')
plt.legend()
V praxi sú samozrejme na regresiu k dispozícii aj jednoduchšie nástroje. Dobrým príkladom všeobecnej funkcie je scipy.optimize.curve_fit
. Existujú však aj špecializované funkcie – napríklad pre lineárnu regresiu, polynomickú regresiu a pod.
Pri použití týchto riešení stačí v predpísanej forme zadať regresnú funkciu (prípadne jej gradient) a nie je nutné sa zaoberať vecami ako je iterácia po dátovej množine a pod. – tu sme taký príklad uvádzali len kvôli názornosti.
Vyššie uvedený príklad by sme napr. pomocou scipy.optimize.curve_fit
vedeli zjednodušiť takto:
def sigmoid(x, a, c):
return 1 / (1 + np.exp(-a*x - c))
popt, _ = curve_fit(sigmoid, X, Y)
xx = np.arange(-20, 20, 0.05)
yy = sigmoid(xx, *popt)
plt.scatter(X, Y, c='r', label="pôvodné dáta")
plt.plot(xx, yy, label="regresná závislosť")
plt.xlabel('x')
plt.ylabel('y')
plt.grid(ls='--')
plt.legend()