Idea: a partir de un punto obtenido, escoger una dirección para dar el próximo paso
Definición (dirección de descenso): sean $f:\mathbb{R}^n \rightarrow \mathbb{R}$, $\bar{x}\in\mathbb{R}^n$ y $d\in\mathbb{R}^n-\{0\}$, diremos que $d$ es una dirección de descenso para $f$ a partir de $\bar{x}$ si existe $\delta>0$ tal que $f(\bar{x}+td)<f(\bar{x}) \quad \forall t\in(0,\delta)$
Teorema: Si $\nabla f(\bar{x})d < 0$, entonces $d$ es dirección de descenso para $f$ a partir de $\bar{x}$
Dados: $f,\; x^0 \in \mathbb{R}^n,\; \varepsilon>0,\; k_{MAX}>0$
k = 0
REPETIR mientras $\nabla f(x^k) > \varepsilon$ y $k<k_{MAX}$ :
Calcular $d^k$ tal que $\nabla f(x^k)^Td^k < 0$
Escoger $t_k>0$ tal que $f(x^k+t_kd^k)<f(x^k)$
Hacer $x^{k+1}=x^k + t_kd^k$
$k = k+1$
Lo que nosotros querríamos hacer es resolver el siguiente problema:
De aquí en adelante definimos $\varphi\colon [0,+\infty)\rightarrow \mathbb{R}$ de la siguiente manera: $$\varphi(t)=f(x+td)$$ Definición (Función unimodal): una función continua $\varphi\colon [0,+\infty)\rightarrow \mathbb{R}$ se dice unimodal si admite un conjunto de minimizadores $[t_1, t_2]$ y es estrictamente decreciente en $[0,t_1]$ y estrictamente creciente en $[t_2, +\infty)$
Cómo encontrar un minimizador de una función unimodal
Supongamos que sabemos que un minimizador se encuentra en $[a,b]$.
El método de la sección áurea consiste en tomar: $$u = a + \frac{3-\sqrt{5}}{2}(b-a) \qquad v = a + \dfrac{\sqrt{5}-1}{2}(b-a)$$
Lema: en el método de la sección áurea, si $[v,b]$ es descartado, entonces $v^+=u$. Análogamente, si $[a,u]$ es descartado, entonces $u^+=v$.
El algoritmo de la sección áurea consiste en 2 fases: en la primera, se busca un intervalo $[a,b]$ que contenga un minimizador de $\varphi$ ; en la segunda, se reduce el intervalo $[a,b]$ hasta que la precisión deseada $\varepsilon$ es alcanzada.
Algoritmo de sección áurea
Dados: $f, \; x\in\mathbb{R}^n,\; d\in\mathbb{R}^n,\; \varepsilon > 0,\; \rho>0$
Definir $\theta_1 = \frac{3-\sqrt{5}}{2}$, $\theta_2 = 1 - \theta_1$
Definir $a = 0$, $s = \rho$, $b=2\rho$
Definir $\varphi_b = f(x+b*d), \; \varphi_s = f(x+s*d)$
REPETIR mientras $\varphi_b<\varphi_s$:
$a=s$, $s=b$, $b=2b$
$\varphi_s = \varphi_b, \; \varphi_b = f(x+b*d)$
Definir $u = a+\theta_1(b-a), \; v=a+\theta_2(b-a)$
Definir $\varphi_u = f(x+u*d), \; \varphi_v = f(x+v*d)$
REPETIR mientras $(b-a)>\varepsilon$:
Si $\varphi_u<\varphi_v$:
$b=v, \; v=u,\; u=a+\theta_1(b-a)$
$\varphi_v = \varphi_u,\; \varphi_u = f(x+u*d)$
Si no:
$a=u,\; u=v,\; v= a+\theta_2(b-a)$
$\varphi_u = \varphi_v,\; \varphi_v = f(x+v*d)$
Definir $\bar{t} = \dfrac{u+v}{2}$
Valores usuales para los parámetros: $\varepsilon = 10^{-5},\; \rho=1$
Cuando $\varphi$ es unimodal, la sección áurea funciona muy bien para encontrar una aproximación al minimizador. Si $\varphi$ no es unimodal, este algoritmo puede no ser eficiente.
A diferencia de la sección áurea, este algoritmo no busca minimizar $\varphi(t)$ sino encontrar $t$ tal que haya una buena reducción en la función objetivo. Más específicamente, dados $\bar{x}\in\mathbb{R}^n$, $d$ dirección de descenso y $\eta\in(0,1)$ busca $\bar{t}$ tal que: $$f(\bar{x}+\bar{t}d) \leq f(\bar{x})+\eta\bar{t}\nabla f(\bar{x})^Td$$ es decir, la reducción debe ser proporcional al tamaño del paso.
Algoritmo de Búsqueda de Armijo
Dados: $f,\;\bar{x}\in\mathbb{R}^n$, $d\in \mathbb{R^n}$ dirección de descenso, $\gamma,\;\eta\in(0,1)$
t=1
REPETIR mientras $f(\bar{x}+td) > f(\bar{x})+\eta t\nabla f(\bar{x})^Td$
$t=\gamma t$
Valores usuales para los parámetros: $\gamma = 0.7,\; \eta=0.45$
La condicion de Armijo impone cotas a cuan grande puede ser el paso en la direccion $d$. Sin embargo, puede ocurrir que el paso sea tan pequeño que $x^k$ no converja a un minimizador local.
Por ejemplo, si $f(x)=x^2$, si comenzamos en $x^0=2$, $d=-1$ es dirección de descenso en $x^k = 1+2^{-k}$ y $t_k = 2^{-k-1}$ cumple con la condición y efectivamente se logra un descenso en $f$. Sin embargo, $x_k\rightarrow 1$ que no es el minimizador de $f$. No se converge al mínimo pues los pasos son muy pequeños.
Entonces, necesitamos otra condición que acote inferiormente a $t$.
Wolfe agrega otra condición sobre $t$, la condición de curvatura: $$\nabla f (\bar{x} + \bar{t}d)^Td \geq \zeta \nabla f(\bar{x})^T d$$ donde $\zeta\in(\eta,1)$ con $\eta$ siendo la constante de la condición de Armijo. El lado izquierdo es $\varphi'(\bar{t})$, por lo que la condición impone que la pendiente de $\varphi$ en $\bar{t}$ sea mayor que $\zeta$ veces la pendiente inicial.
Si $\varphi'$ es muy negativa $\Rightarrow$ se puede decrecer mucho en esta dirección
Si $\varphi'$ no es muy negativa o es positiva $\Rightarrow$ terminar la búsqueda lineal, no se pueden lograr (muchas) mejoras
Condiciones de Wolfe:
$$\begin{array}{rcl} f(\bar{x}+\bar{t}d) &\leq& f(\bar{x})+c_1\bar{t}\nabla f(\bar{x})^Td \\ \nabla f (\bar{x} + \bar{t}d)^Td &\geq& c_2 \nabla f(\bar{x})^T d \end{array}$$Con $0<c_1<c_2<1$.
Algoritmo de Búsqueda de Wolfe
Dados: $f,\; \bar{x}\in\mathbb{R}^n,\; d$ dirección de descenso,$\; 0<c_1<c_2<1$
Definir $\alpha = 0,\; t=1,\; \beta = +\infty$ (beta=np.inf
)
REPETIR
Si $f(\bar{x}+td) > f(\bar{x})+c_1t\nabla f(\bar{x})^Td$ :
Definir $\beta=t,\; t=\frac{1}{2}(\alpha+\beta)$
De lo contrario, si $\nabla f (\bar{x} + td)^Td < c_2 \nabla f(\bar{x})^T d$:
Definir $\alpha=t ,\; t=\begin{cases} 2\alpha \quad \text{si } \beta=+\infty \\ \frac{1}{2}(\alpha+\beta) \quad \text{c.c.} \end{cases} $
Si no: PARAR
Valores usuales para parámetros: $c_1= 0.5,\; c_2=0.75$
El comando break
de Python permite salir de un ciclo while
o for
Lema: sean $f:\mathbb{R}^n\rightarrow \mathbb{R}$, $f\in C^1$, $\bar{x}\in\mathbb{R}^n$ y $d$ una dirección de descenso, entonces una de las siguientes dos situaciones pueden ocurrir con el método antes expuesto para las condiciones de Wolfe:
i) El procedimiento temrina en una cantidad finita de pasos, devolviendo un valor $\bar{t}$ que satisface las condiciones de Wolfe
ii) El procedimiento no termina: el parámetro $\beta$ nunca toma un valor finito, $\alpha$ se vuelve positivo en la primera iteración y se duplica con las iteraciones siguientes, y $f(x+td)\rightarrow -\infty$
En este método, se toma como dirección de descenso: $$d = -\nabla f(x)$$ Observar que, si $\nabla f(x) \neq 0$, efectivamente $d$ es dirección de descenso: $$\nabla f(x)^T d = -\left\Vert \nabla f(x) \right\Vert ^2 < 0$$
Algoritmo del Método del Gradiente
Dados: $f, x^0 \in \mathbb{R}^n,\; \varepsilon>0,\; k_{MAX}>0 $
$k = 0$
REPETIR mientras $\lVert\nabla f(x^{k})\rVert > \varepsilon$ y $k<k_{MAX}$ :
Definir $d^k = -\nabla f(x^k)$
Obtener $t_k>0$ tal que $f(x^k+t_kd^k)<f(x^k)$ (con sección áurea, Armijo o Wolfe)
Hacer $x^{k+1} = x^k + t_kd^k$
$k = k+1$
Lema: en el algoritmo anterior, sea $t_k$ el obtenido por una minimiación local de $f(x^k+td^k)$, entonces $<d^{k+1}, d^k> = 0$
Demo: como $t_k$ es minimizador local de $\varphi$ tenemos que: $$0 =\varphi'(t_k) = \nabla f(x^k +t_kd^k)^Td^k=\nabla f(x^{k+1})^Td^k= (-d^{k+1})^Td^k$$
En este método, se elige como dirección de descenso: $$d = -Hf(x)^{-1}\nabla f(x)$$ Notar que, si $Hf(x)$ es definida positiva, $d$ es dirección descenso: $$\nabla f(x)^T d= -\nabla f(x)^THf(x)^{-1}\nabla f(x) < 0$$
Observación: $d$ puede no estar bien definido si $Hf(x)$ es singular o puede no ser una dirección de descenso si $Hf(x)$ no es definida positiva.
Teorema: si $Hf(x)$ es definida positiva para todo $x\in\mathbb{R}^n$, el Método de Newton es globalmente convergente.
Algoritmo del Método de Newton
Dados: $f, x^0 \in \mathbb{R}^n,\; \varepsilon>0,\; k_{MAX}>0$
$k=0$
REPETIR mientras $\lVert\nabla f(x^{(k)})\rVert > \varepsilon$ y $k<k_{MAX}$:
Definir $d^{(k)}$ como la solución de $Hf(x^{(k)})d^{(k)} = -\nabla f(x^{(k)})$
Determinar el tamaño del paso $t_k > 0$
Hacer $x^{(k+1)} = x^{(k)} + t_kd^{(k)}$
k = k + 1
Considerar $$x^{(k+1)} = x^{(k)} - t_k(Hf(x^{(k)})+\mu I)^{-1}\nabla f(x^{(k)})$$ donde $\mu>0$ es lo suficientemente grande como para que $B^{(k)}= Hf(x^{(k)})+\mu I$ sea definida positiva.
Obs: si $\mu$ es demasiado grande, $d^{(k)}$ tiende a estar en la dirección de $-\nabla f(x^{(k)})$, pues $(Hf(x^{(k)})+\mu I)^{-1}\nabla f(x^{(k)}) \approx \frac{1}{\mu}\nabla f(x^{(k)})$ . En cambio, si $\mu$ es pequeño, $d^{(k)}$ se parece más a la dirección de descenso del Método de Newton
Con la modificación de L-M, se obtiene un método que utiliza las ventajas del Método del Gradiente y del Método de Newton. Pensándolo con ideas generales:
Obs: Como pedíamos que $f\in C^3$, en particular resulta que $Hf(x)$ es simétrica $\forall x \in \mathbb{R}^n$, entonces los autovalores de $Hf(x)$ son reales $\forall x \in \mathbb{R}^n$.
Obs: Si $\lambda_1, \dots, \lambda_n$ son los autovalores de $Hf(x)$, entonces $\lambda_1+\mu, \dots, \lambda_n+\mu$ son los autovalores de $B=Hf(x)+\mu I$, y $B$ tiene los mismos autovectores que $Hf(x)$. Sea $v_i$ autovector de $Hf(x)$: $$\begin{array}{rcl} Bv_i &=& (Hf(x)+\mu I)v_i \\ &=& Hf(x)v_i+\mu Iv_i \\ &=& \lambda_iv_i+\mu v_i \\ &=& (\lambda_i+\mu)v_i \end{array}$$ Entonces, en caso de que $Hf(x)$ no sea definida positiva, se puede utilizar su mínimo autovalor como una aproximación a $-\mu$
Algoritmo del Método de Newton con modificación de Levenberg-Marquardt
Dados: $f, x^0 \in \mathbb{R}^n,\; k_{MAX}>0,\; \varepsilon>0,\; \gamma>0$
$k=0$
REPETIR mientras $\lVert\nabla f(x^{(k)})\rVert > \varepsilon$ y $k<k_{MAX}$:
$B = Hf(x^{(k)})$
$\mu = \min \{ \lambda \colon \lambda \text{ es autovalor de } B\}$
Si $\mu \leq 0$:
$B = B+(-\mu+\gamma)I$
Definir $d^{(k)}$ como la solución de $Bd^{(k)} = -\nabla f(x^{(k)})$
Determinar el tamaño del paso $t_k > 0$
Hacer $x^{(k+1)} = x^{(k)} + t_kd^{(k)}$
k = k + 1
Utilizar la función eigvals
de numpy.linalg
para calcular los autovalores.
Comparación de los métodos, aplicados a: $$f(x, y) = 2(e^{-x^2-y^2}- e^{-(x-1)^2 - (y-1)^2})$$
Implementar una función para cada una de las siguientes tareas:
Las últimas tres deben tener metodo
como un input opcional que permita elegir qué metodo utilizar (sección áurea, Armijo o Wolfe) para encontrar la longitud del paso. Deben imprimir la cantidad de iteraciones que realizaron, el minimizador y el valor de la función en el mismo.
Aplicar los tres métodos sobre la función de Rosenbrock en 2 dimensiones y comparar los resultados. El minimizador de esta función es $x^*=(1,1)$ y $f(x^*)=0$
Se pueden aplicar a otras funciones del archivo fx_test.py
. Se pueden importar todas las funciones de ese archivo con el siguiente comando
from fx_test import *
Para utilizar la librería plotly
hay que instalar también jupyter
. Para instalar ambos paquetes, desde la terminal se ejecuta, para Python 2.7 : pip install jupyter
y luego pip install plotly
.
Para Python 3+: pip3 install jupyter
y luego pip3 install plotly
Dependiendo de su configuración, podría ser necesario anteponer el comando sudo
Para graficar las funciones del archivo fx_test
pueden usar la función plot_function
que se encuentra en ese archivo. Por ejemplo, si se desea graficar la función de ackley
:
plot_function(ackley)
Pueden modificar los límites del gráfico con el input opcional limits
e ingresando una tupla $(x_1, x_2, y_1, y_2)$ para que se grafique la función en el dominio $[x_1, x_2]\times [y_1, y_2]$. Por ejemplo:
plot_function(ackley, limits=(-1,1,-3,2))
grafica la función de Ackley en $[-1,1]\times [-3,2]$. El argumento opcional points
toma una lista de puntos $(x,y)$ y los grafica sobre la superficie de la función (se ve mejor utilizando plotly