In [0]:
# Import potrebných balíčkov
import numpy as np
import matplotlib.pyplot as plt
from sympy.utilities.lambdify import lambdify
import sympy as sp
In [0]:
# Funkcie z predchádzajúceho notebook-u
def plot_func(xx, yy, f, X=None):
    if not X is None:
        Xmin, Xmax = X[:, 0].min(), X[:, 0].max() 
        Ymin, Ymax = X[:, 1].min(), X[:, 1].max()
        
        if (Xmin < xx.min() or Xmax > xx.max() or
                Ymin < yy.min() or Ymax > yy.max()):            
            xx = np.linspace(Xmin, Xmax, 100)
            yy = np.linspace(Ymin, Ymax, 100)
            xx, yy = np.meshgrid(xx, yy)
            
        plt.scatter(X[:, 0], X[:, 1], zorder=10)
        
    zz = f(xx, yy)
    plt.contour(xx, yy, zz, cmap='Spectral')
    # obe osi v rovnakej škále + legenda
    plt.gca().set_aspect('equal')
    plt.xlabel('x'); plt.ylabel('y')
    plt.colorbar(label='z')
        
def grad_desc(grad_f, init_point,
              num_steps, learning_rate):
    X = np.zeros((num_steps + 1, 2))
    X[0] = init_point
    
    for i in range(num_steps):
        X[i+1] = X[i] - learning_rate * grad_f(*X[i])
    
    return X

Metóda klesajúceho gradientu na roztiahnutom paraboloide

Skúsme teraz ten istý postup aplikovať na o trochu náročnejší problém: na roztiahnutý paraboloid. Problém bude v tom, že nebude možné nájsť rýchlosť učenia, ktorá by dobre fungovala v jednom aj v druhom smere.

Povedzme, že funkcia, ktorú chceme minimalizovať, bude mať nasledujúci tvar:

\begin{equation} z = f(x, y) = (5x)^2 + y^2 \end{equation}

Vizualizácia funkcie

Na začiatok si funkciu znovu zadefinujme a vizualizujme. Funkciu zadefinujeme najprv symbolicky, aby sme neskôr mohli využiť automatický spôsob na symbolický výpočet gradientu:

In [3]:
symx, symy = sp.symbols('x y')
symf = (5*symx)**2 + symy ** 2
symf
Out[3]:
25*x**2 + y**2
In [0]:
f = lambdify((symx, symy), symf, "numpy")

Vizualizáciu funkcie realizujeme rovnako, ako predtým:

In [5]:
xx, yy = np.mgrid[-10:10.2:0.2, -10:10.2:0.2]
zz = f(xx, yy)
plot_func(xx, yy, f)

Úloha 1: Automatické určenie gradientu

Pomocou balíčka sympy znovu automaticky určte gradient a jeho symbolické vyjadrenie preveďte na klasickú numerickú funkciu grad_f:

In [0]:

Minimalizácia funkcie metódou klesajúceho gradientu

Ak sa teraz pokúsime danú účelovú funkciu minimalizovať pomocou metódy klesajúceho gradientu, narazíme na problém: bude ťažké nastaviť rýchlosť učenia, ktorá by fungovala v oboch smeroch. Buď nastavíme rýchlosť učenia vysokú a dôjde ku oscilácii v smere, kde funkcia klesá prudšie alebo nastavíme rýchlosť učenia nízku a minimalizácia v druhom smere bude postupovať extrémne pomaly.

Môžeme si to overiť aj prakticky. S $\gamma = 0.1$ dochádza ku oscilácii:

In [7]:
X = grad_desc(grad_f, init_point=[-9, -8],
              num_steps=20, learning_rate=0.1)
plot_func(xx, yy, f, X)

Ak zvolíme $\gamma=0.01$, minimalizácia v smere $y$ postupuje pomaly:

In [8]:
X = grad_desc(grad_f, init_point=[-9, -8],
              num_steps=20, learning_rate=0.01)
plot_func(xx, yy, f, X)

Klesajúci gradient s hybnosťou

Jeden spôsob, ako sa brániť voči problému, na ktorý sme narazili pri použití metódy klesajúceho gradientu, je zaviesť do minimalizačného pravidla ďalší člen – hybnosť. V tom prípade berieme pri výpočte nového bodu do úvahy aj veľkosť zmeny z predchádzajúcej iterácie:

\begin{align} \Delta \mathbf{x}_{i+1} &= \alpha \Delta \mathbf{x}_i - \gamma \nabla f(\mathbf{x}_i) \\[0.75em] \mathbf{x}_{i+1} &= \mathbf{x}_i + \Delta \mathbf{x}_{i+1}, \end{align}

kde $\alpha$ je koeficien hybnosti (určuje aká váha sa pri výpočte prikladá zmene z predchádzajúcej iterácie $\Delta \mathbf{x}_i$).

Výhodou je, že ak minimalizácia dlho postupuje jedným smerom, hybnosť sa bude kumulovať a kroky v tom smere sa budú postupne zväčšovať. Naopak v smere, kde dochádza ku oscilácii sa bude neustále meniť znamienko zmeny, čo bude mať tendenciu osciláciu tlmiť.

Úloha 2: Doplnenie hybnosti

Doplňte do nasledujúcej bunky kód metódy klesajúceho gradientu s hybnosťou (podľa vzorcov uvedených vyššie):

In [0]:
def grad_desc_momentum(grad_f, init_point,
              num_steps, learning_rate, alpha):
    X = np.zeros((num_steps + 1, 2))
    X[0] = init_point
    deltaX = 0
    
    for i in range(num_steps):
        
        
        # sem treba doplniť kód
        
    
    return X

Pozrime sa teraz, ako sa bude nášmu novému algoritmu dariť na probléme roztiahnutého paraboloidu:

In [10]:
X = grad_desc_momentum(grad_f, init_point=[-9, -8],
              num_steps=20, learning_rate=0.01, alpha=0.8)
plot_func(xx, yy, f, X)

Ako vidno, upravený algoritmus nemá problém nájsť minimum, pretože hybnosť akceleruje jeho postup v smere, kde sa znamienko gradientu dlhodobo nemení a naopak tlmí osciláciu v smere, kde minimum preskakujeme.