Redes Neuronales y Control difuso

Guia 3

Martin Noblía

Marcos Panizzo

Damián Presti


<span xmlns:dct="http://purl.org/dc/terms/" property="dct:title">Guia3 p3</span> por <a xmlns:cc="http://creativecommons.org/ns#" href="http://nbviewer.ipython.org/urls/raw.githubusercontent.com/elsuizo/Redes_neuronales_Fuzzy/master/guia1.ipynb?create=1" property="cc:attributionName" rel="cc:attributionURL">Martin Noblía Marcos Panizzo Damián Presti</a> se distribuye bajo una Licencia Creative Commons Atribución-CompartirIgual 4.0 Internacional.

3)

Minimizar la función:

$$F(x_{1},x_{2})=(x_{2} - 2 \, x_{1} + e^{-x_{1}})^{2} + (x_{1} - 2 \, x_{2} + e^{-x_{2}})^{2}$$

i) Utilizar el método gradiente descendente que parte de un $x_{0}$ al azar entre $[-1,1]$

ii) Utilizar el método steepest descent que utiliza información sobre el Gradiente y el Hessiano de la función $F$ para calcular el parámetro $\alpha$(learning rate)

iii) ¿Que algoritmo converge más rápido a la solución? ¿Cuales son los pros y contras de cada método?

i-ii)

Las definiciones básicas que necesitamos recordar son:

  • Minimo fuerte: El punto $\mathbf{x}^{*}$ es un minimo fuerte de $F(\mathbf{x})$ si $\epsilon$ $\delta > 0$ tal que $F(\mathbf{\mathbf{x}^{*}})< F(\mathbf{\mathbf{x}^{*}}+ \Delta \mathbf{x})$ tal que $\delta > \| \Delta \mathbf{x} \| > 0$

En otras palabras si nos movemos del punto en cualquier dirección una distancia pequeña la función tomará valores mayores

  • Mínimo Global: El punto $\mathbf{x}^{*}$ es un único mínimo global de $F(\mathbf{x})$ si $F(\mathbf{\mathbf{x}^{*}})< F(\mathbf{\mathbf{x}^{*}}+ \Delta \mathbf{x})$ para todo $\Delta \mathbf{x} \neq 0$

  • Mínimo débil: El punto $\mathbf{x}^{*}$ es un mínimo débil de $F(\mathbf{x})$ si no es un mínimo fuerte y existe un escalar $\delta > 0$ tal que $F(\mathbf{\mathbf{x}^{*}}) < F(\mathbf{\mathbf{x}^{*}}+ \Delta \mathbf{x})$ para todo $\Delta \mathbf{x}$ tal que $\delta > \| \Delta \mathbf{x} \| > 0$

Se demuestra a traves del desarrollo de series de Taylor que en un entorno reducido los terminos de orden superior no son importantes y por ello podemos resumir lo anterior en las dos condiciones siguientes:

Condición de primer orden:

El punto $\mathbf{x}^{*}$ es un candidato a mínimo si:

$$\nabla F (\mathbf{x}) |_{\mathbf{x}= \mathbf{x}^{*} } = 0$$

Condición de segundo orden:

En este caso suponemos que se cumple la primer condición y del análisis de la serie de Taylor alrededor del punto candidato queda:

$$\Delta \mathbf{x} \, \nabla^{2} F(\mathbf{x})|_{\mathbf{x} = \mathbf{x}^{*}} \Delta \mathbf{x} > 0 $$

Para que lo anterior se cumpla para un $\Delta \mathbf{x} \neq 0 $ se requiere que la matriz Hessiana sea definida positiva

Algoritmos de optimización

Los algoritmos de optimización que vamos a utilizar para hallar el/los minimo/s de $F(x)$ se basan en comenzar con algún iterante inicial $\mathbf{x}_{0}$ y actualizar de acuerdo a alguna ecuación de la forma:

$$\mathbf{x}_{k+1}=\mathbf{x}_{k} + \alpha_{k} \mathbf{p}_{k}$$

o $$\Delta \mathbf{x}_{k}=(\mathbf{x}_{k+1}-\mathbf{x}_{k})=\alpha_{k}\mathbf{p}_{k}$$

Donde el vector $\mathbf{p}_{k}$ representa una dirección de búsqueda y $\alpha_{k}$ es lo que se conoce como learning rate. Las diferentes estrategias para elegir a estos parámetros definen a cada algoritmo

Gradiente descendente

En este caso se elige a $\mathbf{p}_{k}$ con la dirección opuesta a la del gradiente o sea $\mathbf{p}_{k}= -\mathbf{g}_{k} = - \nabla F (\mathbf{x}) |_{\mathbf{x}= \mathbf{x}_{k} } $

Steepest descent

En este caso se minimiza con respecto a $\alpha_{k}$ en cada iteración, o sea elegir $\alpha_{k}$ para minimizar $F(\mathbf{x}_{k}+\alpha_{k} \mathbf{p}_{k})$ Para hacer esto en una función arbitraria se requiere:

$$\alpha_{k}=\frac{-\mathbf{g}_{k}\mathbf{p}_{k}}{\mathbf{p}_{k} \, A_{k} \, \mathbf{p}_{k}}$$

Donde la Matriz $A_{k} = \nabla^{2} F (\mathbf{x}) |_{\mathbf{x}= \mathbf{x}_{k} } $ (Que en el caso que $F$ sea una forma cuádratica es constante)

A continuación las implementaciones con el maravilloso y libre lenguaje de Programación Python :P


In [21]:
# Imports and setups
import numpy as np
import matplotlib.pyplot as plt
import sympy as sym
%matplotlib inline
sym.init_printing()  # Imprime la salida en LaTeX
plt.rcParams['figure.figsize'] = 8,6 #parámetros de tamaño

In [22]:
x_1, x_2 = sym.symbols('x_1 x_2')

In [23]:
# Funcion simbolica a minimizar
f = (x_2 - 2 * x_1 + sym.exp(-x_1))** 2 + (x_1 - 2 * x_2 + sym.exp(-x_2))** 2
f


Out[23]:
$$\left(- 2 x_{1} + x_{2} + e^{- x_{1}}\right)^{2} + \left(x_{1} - 2 x_{2} + e^{- x_{2}}\right)^{2}$$

In [24]:
# Funcion numerica a minimizar
f_numeric = sym.lambdify((x_1,x_2),f,'numpy')

In [25]:
# Definimos el gradiente simbolico
def grad(f):
    
    return sym.Matrix([f.diff(x_1),f.diff(x_2)])

In [26]:
#definimos el hessiano simbolico
def hess(f):
    
    return sym.Matrix([[f.diff(x_1,2),f.diff(x_1,x_2)],[f.diff(x_2,x_1),f.diff(x_2,2)]])

In [27]:
g = grad(f)
h = hess(f)
h


Out[27]:
$$\left[\begin{matrix}2 \left(\left(2 + e^{- x_{1}}\right)^{2} + \left(- 2 x_{1} + x_{2} + e^{- x_{1}}\right) e^{- x_{1}} + 1\right) & - 2 \left(4 + e^{- x_{2}} + e^{- x_{1}}\right)\\- 2 \left(4 + e^{- x_{2}} + e^{- x_{1}}\right) & 2 \left(\left(2 + e^{- x_{2}}\right)^{2} + \left(x_{1} - 2 x_{2} + e^{- x_{2}}\right) e^{- x_{2}} + 1\right)\end{matrix}\right]$$

In [28]:
g_numeric = sym.lambdify((x_1,x_2),g,'numpy')  # Gradiente numerico
h_numeric = sym.lambdify((x_1,x_2),h,'numpy')  # Hessiano numerico

Gradiente descendente


In [29]:
%%timeit
# gradiente descendente
ite1_max = 100
x_sol = np.empty((2,ite1_max))
x_0 = np.random.uniform(-1,1,size=2)

alpha = .03  # Learning rate
x_sol[:,0] = x_0  # Condicion inicial
tol1 = 1e-7 # tolerancia al error
delta1 = 1
i = 0
while (delta1 > tol1) and (i < (ite1_max - 1)):
    
    p = np.asarray(g_numeric(x_sol[0,i],x_sol[1,i])).reshape(2,)
    
    x_sol[:,i+1] = x_sol[:,i] - alpha * p
    
    delta1 = np.linalg.norm(x_sol[:,i+1] - x_sol[:,i])
    
    i += 1


100 loops, best of 3: 5.14 ms per loop

In [30]:
print 'cantidad de iteraciones: --> ', i


cantidad de iteraciones: -->  88

In [31]:
# Plots
X = np.arange(-2, 2, 0.1)
Y = np.arange(-2, 2, 0.1)
X, Y = np.meshgrid(X, Y)
Z = f_numeric(X, Y)
CS = plt.contour(X, Y, Z)
plt.clabel(CS, inline=1, fontsize=10)
plt.title('Curvas de Nivel', fontsize=20)
plt.plot(x_sol[0,0:i],x_sol[1,0:i],'--o')
plt.plot(x_0[0],x_0[1],'ro') # Punto inicial 
plt.plot(x_sol[0,i],x_sol[1,i],'go')  #Punto final


Out[31]:
[<matplotlib.lines.Line2D at 0x7f5dd868bdd0>]

In [32]:
# Verificamos la condición de primer orden
g = np.asarray(g_numeric(x_sol[0,i],x_sol[1,i])).reshape(2,)

In [33]:
# Es casi cero
g


Out[33]:
array([ -1.93748018e-06,  -1.93748018e-06])

In [34]:
# Verificamos la condicion de segundo orden(definida positiva)
H = np.asarray(h_numeric(x_sol[0,i],x_sol[1,i])).reshape(2,2)

In [35]:
H


Out[35]:
array([[ 15.18045235, -10.26857406],
       [-10.26857406,  15.18045235]])

In [36]:
# Verificamos que es un minimo fuerte ya que los autovalores del Hessiano en x* son positivos
np.linalg.eigvals(H)


Out[36]:
array([ 25.4490264 ,   4.91187829])

In [37]:
%%timeit
# Steepest descent
ite2_max = 100
x_sol2 = np.empty((2,ite2_max))
#x_02 = np.random.uniform(-1,1,size=2)

k = 0
x_sol2[:,0] = x_0
tol = 1e-7
delta2 = 1

while (delta2 > tol) and (k < (ite2_max-1)):
    
    grad_k = np.asarray(g_numeric(x_sol2[0,k],x_sol2[1,k])).reshape(2,)
   
    p_k = -grad_k
   
    hess_k = np.asarray(h_numeric(x_sol2[0,k],x_sol2[1,k])).reshape(2,2)
   
    alpha_k = -(np.dot(grad_k,p_k) ) / (np.dot(p_k,hess_k).dot(p_k))   
    
    x_sol2[:,k+1] = x_sol2[:,k] + alpha_k * p_k
    
    delta2 = np.linalg.norm(x_sol2[:,k] - x_sol2[:,k-1])
    
    k += 1


100 loops, best of 3: 2.44 ms per loop

In [38]:
print 'Cantidad de iteraciones: --> ', k


Cantidad de iteraciones: -->  42

In [39]:
X = np.arange(-2, 2, 0.1)
Y = np.arange(-2, 2, 0.1)
X, Y = np.meshgrid(X, Y)
Z = f_numeric(X, Y)
CS = plt.contour(X, Y, Z)
plt.clabel(CS, inline=1, fontsize=10)
plt.title('Curvas de Nivel', fontsize=20)
plt.plot(x_sol2[0,0:k],x_sol2[1,0:k],'--o')
plt.plot(x_0[0],x_0[1],'ro') # Punto inicial 
plt.plot(x_sol2[0,k],x_sol2[1,k],'go')  #Punto final


Out[39]:
[<matplotlib.lines.Line2D at 0x7f5dd8dd6590>]

Graficamos la función en el espacio para observar ese hermoso valle que posee alrededor del/los minimos


In [40]:
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
from matplotlib.ticker import LinearLocator, FormatStrFormatter
import matplotlib.pyplot as plt
import numpy as np

fig = plt.figure()
ax = fig.gca(projection='3d')
X = np.arange(-2, 2, 0.1)
Y = np.arange(-2, 2, 0.1)
X, Y = np.meshgrid(X, Y)
Z = f_numeric(X, Y)
#ax.scatter(x_sol2[0],x_sol2[1],'ro')
surf = ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=cm.YlGnBu_r,
        linewidth=0, antialiased=False)
#ax.set_zlim(-1.01, 1.01)

ax.zaxis.set_major_locator(LinearLocator(10))
ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))
ax.view_init(elev=30., azim=30.)
fig.colorbar(surf, shrink=0.5, aspect=5)


Out[40]:
<matplotlib.colorbar.Colorbar instance at 0x7f5dd8e4dc20>

iii)

Como vemos para este caso el algoritmo que converge más rápido(por muy poco en cuanto a tiempo) es el de Steepest descent y además realiza la tarea en menos iteraciones(con la misma tolerancia). Los pros de cada método son:

Gradiente descendente:

  • Si el iterante inicial esta cerca de un mínimo su desempeño es bueno.
  • El $\alpha$ es una constante y no requiere ser calculada iteración a iteración(mejora la velocidad de procesamiento)
  • Solo se necesita la información del gradiente de la función ( $\nabla F$)
  • Si la función $F$ es una forma cuadrática sabemos de antemano cual es la cota para $\alpha$, $\alpha_{max}$ dado por:
$$\alpha_{max} < \frac{2}{\lambda_{max}}$$

Donde $\lambda_{max}$ es el autovalor más grande de la representación matricial de $F$ y esta relacionado con la curvatura de la función

Steepest descent:

  • Si la Función $F$ es una forma cuadrática mejora el rendimiento con respecto al anterior
  • El calculo de $\alpha$ no depende de prueba y error

Desventajas:

Gradiente descendente:

  • Si la función $F$ no es una forma cuadrática no tenemos cotas para el alpha y nuestra elección es por prueba y error
  • Si la función $F$ no es una forma cuadrática, tiene curvatura elevada y no elegimos correctamente $\alpha$ diverge
  • Si el minimo de $F$ esta relativamente "lejos" del iterante inicial el camino no es optimo

Steepest descent:

  • Si la función $F$ no es una forma cuadrática el cálculo de $\alpha_{k}$ se vuelve más costoso computacionalmente ya que por cada iteración debemos evaluar el Hessiano de $F$ ($\nabla^{2}F$)

  • Si la función $F$ posee valles como por ejemplo el famoso caso de la función de Rosenbrock( $F(x_{1},x_{2})=(1-x_{2})^{2}+100(x_{1}-x_{2}^{2})^{2}$ ) el método tardará mucho en converger

Referencias

[1] Martin T. Hagan, Howard B. Demuth, Mark Beale, Neural Network design

[2] Alfio Quateroni, Ricardo Sacco, Fausto Saleri, Numerical Mathematics