Optimización

La exposición en esta sección sigue de cerca el capítulo 4 del libro de Miranda y Fackrer.

En esta sección se describe cómo maximizar numéricamente una función respecto de un vector finito de variables. El problema a resolver es entonces:

$$\max_{x \in X} f(x)$$

con \(x \subseteq R^n\). En este problemas denominamos a \(f\) como la función objetivos, a \(X\) como el conjunto factible y a \(x^*\) como el máximo (si existe).

Los problemas de optimización finitos son muy comunes en economía.

Teorema de Wierstrass*

Si \(f\) es continua y \(X\) es un conjunto no vacío, cerrado y acotado, entonces \(f\) tiene un máximo en \(X\).

Condición de Primer Orden:

$$f'(x^*) = 0$$

Métodos Numéricos de Optimización: Ideas Básicas

El Método de Newton-Rapson

El método de Newton-Rapson usa una sucesión de aproximaciones cuadráticas de la función objetivo y el máximo de la aproximación debería converger al máximo de la función.

Este método está estrechamente relacionado con el método de Newton para buscar las raíces de una función.

El método es iterativo y se inicia con una conjetura \(x^{(0)}\) de la solución.

A continuación se actualiza la solución de \(x^{(k)}\) a \(x^{(k+1)}\) minimizando:

$$f(x) \approx f(x^{(k)}) + f'(x^{(k)})(x-x^{(k)}) + 0.5 (x-x^{(k)})^{T}f''(x^{(k)})(x-x^{(k)})$$

respecto de \(x\). La condición de primero orden es:

$$f'(x^{(k)}) + f''(x^{(k)})(x-x^{(k)})=0$$

y por tanto:

$$x^{(k+1)} \leftarrow x^{(k)} - [f''(x^{(k)})]^{-1} f'(x^{(k)})$$

El algoritmo termina cuando \(x^{(k+1)}\) es suficientemente cercano a \(x^{(k)}\).

Limitaciones: - Se requiere computar tanto la primera como la segunda derivadas de la función (el gradiente y la matriz hessiana). - No hay garantía que la función se incremente en la dirección de la actualización (esto se garantiza sólo si \(f''\) es definida negativa). - EL método de Newton-Rapson es usada sólo cuando la función es globalmente cóncava.

El Método de Quasi-Newton

Sigue la misma idea que el método de Newton-Rapson pero reemplaza la matriz hessiana con una aproximación que es definida negativa (garantizando incrementos en la función).

Por eficiencia, la aproximación se realiza directamente sobre la inversa de la matriz hessiana.

Dos métodos satisfacen está idea: - Davidson-Fletcher-Powell (DFP) - Broyden-Fletcher-Goldfarb-Shano (BFGS)

La diferencia entre ambos radica en cómo actualizan la aproximación de la inversa de la matriz hessiana.

El Toolbox de Optimización de Matlab

El toolbox de optimización de Matlab contiene funciones que permiten minimizar (o maximizar) funciones no lineales generales.

El toolbox incluye rutinas para realizar muchos tipos de procedimientos de optimización. Los que usaremos son: - Minimización no lineal no restringida fminunc - Minimización no lineal restringida fmincon

Las funciones fzero y fsolve ya vistas son también parte de este toolbox.

Nota: Matlab tiene implementados los métodos numéricos anteriores pero con el objetivo de minimizar una función. La maximización se realiza simplemente minimizando el negativo de la función objetivos.

Minimización No Restringida

Matlab minimiza (o maximiza) funciones no lineales sin restricciones usando la función fminunc.

La función busca iterativamente el mínimo de una función escalar con varias variables a partir de una conjetura inicial.

Esta función usa los métodos de Quasi-Newton si el gradiente no es provisto por el usuario.

La sintaxis general es:

[x,fval,exitflag,output,grad,hess] = fminunc('fun',x0,... opts)

Los inputs de la función son:

• 'fun': Función escalar objetivo. Puede estar escrita en un m-file o ser una función anónima.
• x0: Conjetura inicial.
• opts: estructura con las opciones de optimización (número de iteraciones, tolerancia, algoritmo, entre otros).

Los outputs de la función son:

• x: Solución (\(x^*\)).
• fval: Valor de la función objetivo en el mínimo.
• exitflag: Condición de salida de la función fminunc. La solución fue hallada si existflag>0.
• output: Estructura con información de la optimización (número de iteración, número de evaluaciones de la función, condición de primer orden, etc).
• grad: Gradiente evaluado en la solución.
• hess: Matriz hessiana evaluada en la solución.

Resolvamos el siguiente ejemplo:

$$f(x)=3x^2_1 +2x_1 x_2 +x^2_2 -4x_1 + 5x_2$$

Escribamos la función anónima.

clear all;
fun = @(x) 3*x(1)^2 + 2*x(1)*x(2) + x(2)^2 - 4*x(1) + 5*x(2)

fun =

  function_handle with value:

    @(x)3*x(1)^2+2*x(1)*x(2)+x(2)^2-4*x(1)+5*x(2)

Ahora usamos la función fminunc para hallar el mínimo partiendo de la conjetura $[1,1]$.

warning('off');
x0 = [1,1];
[x,fval, eflag, out] = fminunc(@(x) fun(x),x0)

Local minimum found.

Optimization completed because the size of the gradient is less than
the default value of the optimality tolerance.

x =

    2.2500   -4.7500


fval =

  -16.3750


eflag =

     1


out =

  struct with fields:

       iterations: 8
        funcCount: 27
         stepsize: 7.3113e-05
     lssteplength: 1
    firstorderopt: 3.8147e-06
        algorithm: 'quasi-newton'
          message: 'Local minimum found....'

x0 = [1,1];
[x,fval, eflag] = fminsearch(@(x) fun(x),x0)

x =

    2.2500   -4.7500


fval =

  -16.3750


eflag =

     1

Minimización Restringida

Ahora el objetivo es encontrar $x$ que minimiza la función no lineal objetivo $f(x)$ sujeto a restricciones, que pueden ser lineales o no lineales y pueden de igualdad o de desigualdad.

El problema general es:

\begin{eqnarray} \min_{x} && f(x) \ s.a & & \ && Ax = b \ && C(x) = 0 \ && AIx \leq bI \ && CI(x) = 0 \ && lb \leq x \leq ub \end{eqnarray}

donde $A$ y $AI$ son matrices, $C(x)$ y $CI(x)$ son funciones no lineales, $b$ y $bI$ son vectores columna. La dimensiones de estos objetos dependen del número de restricciones. Finalmente, $lb$ y $ub$ son vectores columnas de la misma dimensión que $x$.

La función para resolver problemas de optimización con restricciones en Matlab es fmincon. Funciones de la misma forma que fminunc pero tiene una mayor cantidad de inputs.

Las sintaxis general es:

[x,fval,existflag,output,lambda,grad,hess] = fmincon('fun',x0,AI,bI,A,b,lb,ub,'NonConFunc',Options)

• 'fun': Función escalar objetivo. Puede estar escrita en un m-file o ser una función anónima.
• x0: Conjetura inicial.
• AI y bI caracterizan las restricciones de desigualdad (si no existen éstas usar AI=[] y bI=[]).
• A y b caracterizan las restricciones de igualdad (si no existen éstas usar A=[] y b=[]).
• lb y ub caracterizan las cotas inferiores y superiores de $x$ (si no existen éstas usar lb=[] y ub=[]).
• 'NonConFunc': es una función que tiene x como input y retorna vectores CI(x) y C(x) asociados a las restricciones de desigualdad e igualdad no lineales.
• opts: estructura con las opciones de optimización (número de iteraciones, tolerancia, algoritmo, entre otros).

Los outputs de la función son:

• x: Solución ($x^*$).
• fval: Valor de la función objetivo en el mínimo.
• exitflag: Condición de salida de la función fminunc. La solución fue hallada si existflag>0.
• output: Estructura con información de la optimización (número de iteración, número de evaluaciones de la función, condición de primer orden, etc).
• lambda: Estructura con los multiplicadores de lagrange asociados a las restricciones.
• grad: Gradiente evaluado en la solución.
• hess: Matriz hessiana evaluada en la solución.

A modo de ejemplo, resolvamos el siguiente problema:

\begin{eqnarray} \min_{x_1,x_2} && f(x_1,x_2) = 100(x_2-x^2_1)^2 + (1-x_1)^2 \ s.a & & \ && x_1 + 2x_2 \leq 1\ && 2x_1 + x_2= 1
\end{eqnarray}

Definimos la función objetivo:

function f = rosenbrock(x)

x1 = x(1);
x2 = x(2);

f = 100*(x2-x1^2)^2 + (1-x1)^2;

end

Definimos las restricciones:

clear all;
x0 = [0.5,0];
A = [1,2];
b = 1;
Aeq = [2,1];
beq = 1;

Optimizamos:

[x, fval] = fmincon(@rosenbrock,x0,A,b,Aeq,beq)

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the default value of the optimality tolerance,
and constraints are satisfied to within the default value of the constraint tolerance.

x =

    0.4149    0.1701


fval =

    0.3427

Ahora resolvamos el siguiente problema:

\begin{eqnarray} \min_{x_1,x_2} && f(x_1,x_2) = 1+\frac{x_1}{(1+x_2)} - 3x_1x_2 + x_2(1+x_1) \ s.a & & \ && 0 \leq x_1 \leq 1\ && 0 \leq x_2 \leq 2\ \end{eqnarray}

FunObj = @(x) 1+x(1)/(1+x(2)) - 3*x(1)*x(2) + x(2)*(1+x(1));

Restricciones:

lb = [0,0];
ub = [10,10];

Optimizamos:

x0 = [0.5,1];
x = fmincon(FunObj,x0,[],[],[],[],lb,ub)

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the default value of the optimality tolerance,
and constraints are satisfied to within the default value of the constraint tolerance.

x =

   10.0000   10.0000