Otimização é o problema de se encontrar uma solução numérica para a minimização ou a igualdade. O módulo scipy.optimize fornece algoritmos úteis para a função de minimização (escalar ou multi-dimensional), ajuste de curva e descoberta de raizes.

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

from scipy import optimize

[/pastacode]

Encontrando-se o mínimo de uma função escalar

Vamos definir a seguinte função:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

def f(x):
    return x**2 + 10*np.sin(x)

[/pastacode]

E graficá-la:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

x = np.arange(-10, 10, 0.1)
plt.plot(x, f(x)) 
plt.show()

[/pastacode]

scipy_optimize_example1Essa função tem um mínimo global em torno de -1,3 e um mínimo local em torno de 3,8. A forma geral e eficiente para encontrar um mínimo para esta função é a condução de um gradiente descendente a partir de um determinado ponto inicial. O algoritmo BFGS é uma boa maneira de fazer isso:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

optimize.fmin_bfgs(f, 0)

[/pastacode]

Um possível problema com esta abordagem é que, se a função tem mínimos locais o algoritmo pode encontrar estes mínimos locais, em vez de o mínimo global, dependendo do ponto inicial:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

optimize.fmin_bfgs(f, 3, disp=0)

[/pastacode]

Se não sabemos a vizinhança do mínimo global para escolher o ponto inicial, é preciso recorrer a uma otimização global mais cara. Para determinar o mínimo global, o algoritmo mais simples é o algoritmo de força bruta, em que a função é avaliada em cada ponto de uma determinada rede:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

grid = (-10, 10, 0.1)
xmin_global = optimize.brute(f, (grid,))
xmin_global

[/pastacode]

Para maiores tamanhos de grades, scipy.optimize.brute() torna-se bastante lento. scipy.optimize.anneal() fornece uma alternativa, usando arrefecimento simulado (simulated annealing). Algoritmos mais eficientes para diferentes classes de problemas de otimização global existem, mas estão fora do escopo do scipy. Alguns pacotes úteis para a otimização global são OpenOpt, IPOPT, PyGMO e PyEvolve.

Para encontrar o mínimo local, vamos restringir a variável no intervalo (0, 10) usando scipy.optimize.fminbound():

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

xmin_local = optimize.fminbound(f, 0, 10)    
xmin_local

[/pastacode]

Encontrar as raízes de uma função escalar

Para encontrar uma raiz, ou seja, um ponto em que f (x) = 0, da função f acima, podemos usar, por exemplo scipy.optimize.fsolve():

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

root = optimize.fsolve(f, 1)  # o chute inicial é 1
root

[/pastacode]

Note que apenas uma raiz é encontrada. Inspecionando o enredo, f revela que há uma segunda raiz em torno de -2,5. Nós encontramos o valor exato do mesmo, ajustando a nossa suposição inicial:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

root2 = optimize.fsolve(f, -2.5)
root2

[/pastacode]

Ajustes de curvas

Suponha que tenhamos dados amostrados de f com algum ruído:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

xdata = np.linspace(-10, 10, num=20)
ydata = f(xdata) + np.random.randn(xdata.size)

[/pastacode]

Agora, se nós sabemos a forma funcional da função a partir da qual as amostras foram retiradas (x ^ 2 + sin (x), nesse caso), mas não as amplitudes dos termos, podemos encontrar aqueles por mínimos quadrados do ajuste de curva. Primeiro temos que definir a função para ajustar:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

def f2(x, a, b):
    return a*x**2 + b*np.sin(x)

[/pastacode]

Então, podemos usar scipy.optimize.curve_fit() para encontrar a e b:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

guess = [2, 2]
params, params_covariance = optimize.curve_fit(f2, xdata, ydata, guess)
params

[/pastacode]

Agora, encontramos os mínimos e raízes de f e utilizamos o ajuste de curva nele, vamos colocar, então, todos essas resultados juntos em um único gráfico:

[pastacode lang=”python” message=”” highlight=”” provider=”manual”]

import numpy as np
from scipy import optimize
import matplotlib.pyplot as plt

x = np.arange(-10, 10, 0.1)
def f(x):
    return x**2 + 10*np.sin(x)

grid = (-10, 10, 0.1)
xmin_global = optimize.brute(f, (grid,))
xmin_local = optimize.fminbound(f, 0, 10)
root = optimize.fsolve(f, 1)  # chute inicial é 1
root2 = optimize.fsolve(f, -2.5)

xdata = np.linspace(-10, 10, num=20)
np.random.seed(1234)
ydata = f(xdata) + np.random.randn(xdata.size)

def f2(x, a, b):
    return a*x**2 + b*np.sin(x)

guess = [2, 2]
params, params_covariance = optimize.curve_fit(f2, xdata, ydata, guess)

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x, f(x), 'b-', label="f(x)")
ax.plot(x, f2(x, *params), 'r--', label="Resultado do Ajuste de Curva")
xmins = np.array([xmin_global[0], xmin_local])
ax.plot(xmins, f(xmins), 'go', label="Minimo")
roots = np.array([root, root2])
ax.plot(roots, f(roots), 'kv', label="Raizes")
ax.legend()
ax.set_xlabel('x')
ax.set_ylabel('f(x)')

[/pastacode]

Você pode encontrar algoritmos com as mesmas funcionalidades para problemas multi-dimensionais em scipy.optimize.