SciPy, optimumigo

SciPy, optimumigo

SciPy (prononcita sai pie) estas matematika aplikaĵa pakaĵo bazita sur la etendaĵo Numpy Python. Kun SciPy, via interaga Python-sesio fariĝas la sama kompleta datumscienco kaj kompleksa sistema prototipa medio kiel MATLAB, IDL, Octave, R-Lab kaj SciLab. Hodiaŭ mi volas mallonge paroli pri kiel uzi iujn konatajn optimumigajn algoritmojn en la pako scipy.optimize. Pli detala kaj ĝisdatigita helpo pri uzado de funkcioj ĉiam povas esti akirita per la help() komando aŭ uzante Shift+Tab.

Enkonduko

Por savi vin kaj legantojn de serĉado kaj legado de ĉeffontoj, ligiloj al priskriboj de metodoj estos ĉefe en Vikipedio. Kiel regulo, ĉi tiu informo sufiĉas por kompreni la metodojn en ĝeneralaj terminoj kaj la kondiĉojn por ilia apliko. Por kompreni la esencon de matematikaj metodoj, sekvu la ligilojn al pli aŭtoritataj publikaĵoj, kiuj troviĝas ĉe la fino de ĉiu artikolo aŭ en via plej ŝatata serĉilo.

Do, la modulo scipy.optimize inkluzivas la efektivigon de la sekvaj proceduroj:

  1. Kondiĉa kaj senkondiĉa minimumigo de skalaj funkcioj de pluraj variabloj (minimumo) uzante diversajn algoritmojn (Nelder-Mead simplex, BFGS, Neŭtonaj konjugacigradientoj, COBYLA и SLSQP)
  2. Tutmonda optimumigo (ekzemple: basinhopping, diff_evoluo)
  3. Minimumigo de restaĵoj MNC (malplej_kvadratoj) kaj kurbaj konvenaj algoritmoj uzante neliniajn malplej kvadratojn (curve_fit)
  4. Minimumigante skalarajn funkciojn de unu variablo (minim_scalar) kaj serĉante radikojn (root_scalar)
  5. Plurdimensiaj solvantoj de sistemo de ekvacioj (radiko) uzante diversajn algoritmojn (hibrida Powell, Levenberg-Marquardt aŭ grandskalaj metodoj kiel ekz Neŭtono-Krylov).

En ĉi tiu artikolo ni konsideros nur la unuan eron el ĉi tiu tuta listo.

Senkondiĉa minimumigo de skalara funkcio de pluraj variabloj

La minimuma funkcio de la scipy.optimize-pakaĵo disponigas ĝeneralan interfacon por solvi kondiĉajn kaj senkondiĉajn minimumigproblemojn de skalaraj funkcioj de pluraj variabloj. Por pruvi kiel ĝi funkcias, ni bezonos taŭgan funkcion de pluraj variabloj, kiujn ni minimumigos en malsamaj manieroj.

Por tiuj celoj, la Rosenbrock-funkcio de N variabloj estas perfekta, kiu havas la formon:

SciPy, optimumigo

Malgraŭ tio, ke la Rosenbrock-funkcio kaj ĝiaj Jacobi kaj Hessianaj matricoj (la unua kaj dua derivaĵoj, respektive) estas jam difinitaj en la pako scipy.optimize, ni mem difinos ĝin.

import numpy as np

def rosen(x):
    """The Rosenbrock function"""
    return np.sum(100.0*(x[1:]-x[:-1]**2.0)**2.0 + (1-x[:-1])**2.0, axis=0)

Por klareco, ni desegnu en 3D la valorojn de la Rosenbrock-funkcio de du variabloj.

Desegna kodo

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib.ticker import LinearLocator, FormatStrFormatter

# Настраиваем 3D график
fig = plt.figure(figsize=[15, 10])
ax = fig.gca(projection='3d')

# Задаем угол обзора
ax.view_init(45, 30)

# Создаем данные для графика
X = np.arange(-2, 2, 0.1)
Y = np.arange(-1, 3, 0.1)
X, Y = np.meshgrid(X, Y)
Z = rosen(np.array([X,Y]))

# Рисуем поверхность
surf = ax.plot_surface(X, Y, Z, cmap=cm.coolwarm)
plt.show()

SciPy, optimumigo

Sciante anticipe, ke la minimumo estas 0 je SciPy, optimumigo, ni rigardu ekzemplojn de kiel determini la minimuman valoron de la Rosenbrock-funkcio uzante diversajn scipy.optimize procedurojn.

Nelder-Mead simpla metodo

Estu komenca punkto x0 en 5-dimensia spaco. Ni trovu la minimuman punkton de la Rosenbrock-funkcio plej proksima al ĝi uzante la algoritmon Nelder-Mead simpla (la algoritmo estas specifita kiel la valoro de la metodoparametro):

from scipy.optimize import minimize
x0 = np.array([1.3, 0.7, 0.8, 1.9, 1.2])
res = minimize(rosen, x0, method='nelder-mead',
    options={'xtol': 1e-8, 'disp': True})
print(res.x)

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 339
         Function evaluations: 571
[1. 1. 1. 1. 1.]

La simpla metodo estas la plej simpla maniero por minimumigi eksplicite difinitan kaj sufiĉe glatan funkcion. Ĝi ne postulas kalkuli la derivaĵojn de funkcio; sufiĉas specifi nur ĝiajn valorojn. La metodo Nelder-Mead estas bona elekto por simplaj minimumigproblemoj. Tamen, ĉar ĝi ne uzas gradienttaksojn, povas daŭri pli longe por trovi la minimumon.

Powell-metodo

Alia optimumiga algoritmo en kiu nur la funkciovaloroj estas kalkulitaj estas La metodo de Powell. Por uzi ĝin, vi devas agordi metodo = 'powell' en la minimuma funkcio.

x0 = np.array([1.3, 0.7, 0.8, 1.9, 1.2])
res = minimize(rosen, x0, method='powell',
    options={'xtol': 1e-8, 'disp': True})
print(res.x)

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 19
         Function evaluations: 1622
[1. 1. 1. 1. 1.]

Broyden-Fletcher-Goldfarb-Shanno (BFGS) algoritmo

Por akiri pli rapidan konverĝon al solvo, la procedo BFGS uzas la gradienton de la objektiva funkcio. La gradiento povas esti specifita kiel funkcio aŭ kalkulita per unuaordaj diferencoj. Ĉiukaze, la BFGS-metodo tipe postulas malpli da funkciovokoj ol la simpla metodo.

Ni trovu la derivaĵon de la Rosenbrock-funkcio en analiza formo:

SciPy, optimumigo

SciPy, optimumigo

Ĉi tiu esprimo validas por la derivaĵoj de ĉiuj variabloj krom la unua kaj lasta, kiuj estas difinitaj kiel:

SciPy, optimumigo

SciPy, optimumigo

Ni rigardu la Python-funkcion, kiu kalkulas ĉi tiun gradienton:

def rosen_der (x):
    xm = x [1: -1]
    xm_m1 = x [: - 2]
    xm_p1 = x [2:]
    der = np.zeros_like (x)
    der [1: -1] = 200 * (xm-xm_m1 ** 2) - 400 * (xm_p1 - xm ** 2) * xm - 2 * (1-xm)
    der [0] = -400 * x [0] * (x [1] -x [0] ** 2) - 2 * (1-x [0])
    der [-1] = 200 * (x [-1] -x [-2] ** 2)
    return der

La gradienta kalkulfunkcio estas specifita kiel la valoro de la jac-parametro de la minimuma funkcio, kiel montrite sube.

res = minimize(rosen, x0, method='BFGS', jac=rosen_der, options={'disp': True})
print(res.x)

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 25
         Function evaluations: 30
         Gradient evaluations: 30
[1.00000004 1.0000001  1.00000021 1.00000044 1.00000092]

Konjuga gradientalgoritmo (Neŭtono)

Algoritmo Neŭtonaj konjugataj gradientoj estas modifita metodo de Neŭtono.
La metodo de Neŭtono baziĝas sur aproksimado de funkcio en loka ĉirkaŭaĵo per polinomo de la dua grado:

SciPy, optimumigo

kie SciPy, optimumigo estas la matrico de duaj derivaĵoj (hesia matrico, hesia).
Se la Hessian estas pozitiva difinita, tiam la loka minimumo de ĉi tiu funkcio povas esti trovita egaligante la nul gradienton de la kvadrata formo al nulo. La rezulto estos la esprimo:

SciPy, optimumigo

La inversa Hessian estas kalkulita uzante la konjugatan gradientmetodon. Ekzemplo de uzado de ĉi tiu metodo por minimumigi la Rosenbrock-funkcion estas donita malsupre. Por uzi la metodon Newton-CG, vi devas specifi funkcion, kiu kalkulas la Hessian.
La Hessian de la Rosenbrock-funkcio en analiza formo estas egala al:

SciPy, optimumigo

SciPy, optimumigo

kie SciPy, optimumigo и SciPy, optimumigo, difini la matricon SciPy, optimumigo.

La ceteraj ne-nulaj elementoj de la matrico estas egalaj al:

SciPy, optimumigo

SciPy, optimumigo

SciPy, optimumigo

SciPy, optimumigo

Ekzemple, en kvin-dimensia spaco N = 5, la Hessa matrico por la Rosenbrock-funkcio havas la formon de grupo:

SciPy, optimumigo

Kodo kiu kalkulas ĉi tiun Hessian kune kun kodo por minimumigi la Rosenbrock-funkcion uzante la konjugatan gradienton (Neŭtono) metodon:

def rosen_hess(x):
    x = np.asarray(x)
    H = np.diag(-400*x[:-1],1) - np.diag(400*x[:-1],-1)
    diagonal = np.zeros_like(x)
    diagonal[0] = 1200*x[0]**2-400*x[1]+2
    diagonal[-1] = 200
    diagonal[1:-1] = 202 + 1200*x[1:-1]**2 - 400*x[2:]
    H = H + np.diag(diagonal)
    return H

res = minimize(rosen, x0, method='Newton-CG', 
               jac=rosen_der, hess=rosen_hess,
               options={'xtol': 1e-8, 'disp': True})
print(res.x)

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 24
         Function evaluations: 33
         Gradient evaluations: 56
         Hessian evaluations: 24
[1.         1.         1.         0.99999999 0.99999999]

Ekzemplo kun la difino de la produktofunkcio de la Hessian kaj arbitra vektoro

En real-mondaj problemoj, komputado kaj stokado de la tuta hesia matrico povas postuli signifan tempon kaj memorresursojn. En ĉi tiu kazo, efektive ne necesas specifi la hesan matricon mem, ĉar la minimumiga proceduro postulas nur vektoron egalan al la produto de la Hessian kun alia arbitra vektoro. Tiel, de komputa vidpunkto, estas multe preferinde tuj difini funkcion kiu resendas la rezulton de la produkto de la Hessian kun arbitra vektoro.

Konsideru la hess-funkcion, kiu prenas la minimumigvektoron kiel la unuan argumenton, kaj arbitran vektoron kiel la duan argumenton (kune kun aliaj argumentoj de la minimumigota funkcio). En nia kazo, kalkuli la produkton de la Hessian de la Rosenbrock-funkcio kun arbitra vektoro ne estas tre malfacila. Se p estas arbitra vektoro, tiam la produkto SciPy, optimumigo aspektas kiel:

SciPy, optimumigo

La funkcio kiu kalkulas la produkton de la Hessian kaj arbitra vektoro estas pasita kiel la valoro de la hessp argumento al la minimumigi funkcio:

def rosen_hess_p(x, p):
    x = np.asarray(x)
    Hp = np.zeros_like(x)
    Hp[0] = (1200*x[0]**2 - 400*x[1] + 2)*p[0] - 400*x[0]*p[1]
    Hp[1:-1] = -400*x[:-2]*p[:-2]+(202+1200*x[1:-1]**2-400*x[2:])*p[1:-1] 
    -400*x[1:-1]*p[2:]
    Hp[-1] = -400*x[-2]*p[-2] + 200*p[-1]
    return Hp

res = minimize(rosen, x0, method='Newton-CG',
               jac=rosen_der, hessp=rosen_hess_p,
               options={'xtol': 1e-8, 'disp': True})

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 24
         Function evaluations: 33
         Gradient evaluations: 56
         Hessian evaluations: 66

Konjuga gradienta fidregionalgoritmo (Neŭtono)

Malbona kondiĉado de la hesia matrico kaj malĝustaj serĉdirektoj povas kaŭzi la konjugatan gradientalgoritmon de Neŭtono esti neefika. En tiaj kazoj oni preferas fidoregiona metodo (trust-region) konjugaci Neŭtonajn gradientojn.

Ekzemplo kun la difino de la Hessa matrico:

res = minimize(rosen, x0, method='trust-ncg',
               jac=rosen_der, hess=rosen_hess,
               options={'gtol': 1e-8, 'disp': True})
print(res.x)

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 20
         Function evaluations: 21
         Gradient evaluations: 20
         Hessian evaluations: 19
[1. 1. 1. 1. 1.]

Ekzemplo kun la produktofunkcio de la Hessian kaj arbitra vektoro:

res = minimize(rosen, x0, method='trust-ncg', 
                jac=rosen_der, hessp=rosen_hess_p, 
                options={'gtol': 1e-8, 'disp': True})
print(res.x)

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 20
         Function evaluations: 21
         Gradient evaluations: 20
         Hessian evaluations: 0
[1. 1. 1. 1. 1.]

Krylov-tipaj metodoj

Kiel la trust-ncg-metodo, Krylov-specaj metodoj estas bone konvenitaj por solvi grandskalajn problemojn ĉar ili uzas nur matric-vektorajn produktojn. Ilia esenco estas solvi problemon en konfida regiono limigita de detranĉita Krylov-subspaco. Por necertaj problemoj, estas pli bone uzi ĉi tiun metodon, ĉar ĝi uzas pli malgrandan nombron da neliniaj ripetoj pro la pli malgranda nombro da matrico-vektoraj produktoj per subproblemo, komparite kun la fido-ncg-metodo. Krome, la solvo al la kvadrata subproblemo estas trovita pli precize ol uzi la trust-ncg-metodon.
Ekzemplo kun la difino de la Hessa matrico:

res = minimize(rosen, x0, method='trust-krylov',
               jac=rosen_der, hess=rosen_hess,
               options={'gtol': 1e-8, 'disp': True})

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 19
         Function evaluations: 20
         Gradient evaluations: 20
         Hessian evaluations: 18

print(res.x)

    [1. 1. 1. 1. 1.]

Ekzemplo kun la produktofunkcio de la Hessian kaj arbitra vektoro:

res = minimize(rosen, x0, method='trust-krylov',
               jac=rosen_der, hessp=rosen_hess_p,
               options={'gtol': 1e-8, 'disp': True})

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 19
         Function evaluations: 20
         Gradient evaluations: 20
         Hessian evaluations: 0

print(res.x)

    [1. 1. 1. 1. 1.]

Algoritmo por proksimuma solvo en la konfida regiono

Ĉiuj metodoj (Neŭtono-CG, trust-ncg kaj trust-krylov) taŭgas por solvi grandskalajn problemojn (kun miloj da variabloj). Ĉi tio estas pro la fakto ke la subesta konjugacia gradientalgoritmo implicas proksimuman determinon de la inversa Hessa matrico. La solvo estas trovita ripete, sen eksplicita ekspansio de la hesa. Ĉar vi nur bezonas difini funkcion por la produkto de Hessian kaj arbitra vektoro, ĉi tiu algoritmo estas precipe bona por labori kun maldensaj (bendaj diagonalaj) matricoj. Ĉi tio provizas malaltajn memorkostojn kaj signifajn tempoŝparojn.

Por mezgrandaj problemoj, la kosto de stokado kaj fabrikado de la Hessian ne estas kritika. Ĉi tio signifas, ke eblas akiri solvon en malpli da ripetoj, solvante la subproblemojn de la fidregiono preskaŭ ekzakte. Por fari tion, kelkaj neliniaj ekvacioj estas solvitaj ripete por ĉiu kvadrata subproblemo. Tia solvo kutime postulas 3 aŭ 4 Cholesky-malkomponaĵojn de la hesia matrico. Kiel rezulto, la metodo konverĝas en pli malmultaj ripetoj kaj postulas pli malmultajn objektivan funkciokalkulojn ol aliaj efektivigitaj fidregionmetodoj. Ĉi tiu algoritmo nur implikas determini la kompletan Hessian matricon kaj ne subtenas la kapablon uzi la produktofunkcion de la Hessian kaj arbitra vektoro.

Ekzemplo kun minimumigo de la Rosenbrock-funkcio:

res = minimize(rosen, x0, method='trust-exact',
               jac=rosen_der, hess=rosen_hess,
               options={'gtol': 1e-8, 'disp': True})
res.x

Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 13
         Function evaluations: 14
         Gradient evaluations: 13
         Hessian evaluations: 14

array([1., 1., 1., 1., 1.])

Ni verŝajne haltos tie. En la sekva artikolo mi provos rakonti la plej interesajn aferojn pri kondiĉa minimumigo, la apliko de minimumigo en solvado de proksimumadproblemoj, minimumigo de funkcio de unu variablo, arbitraj minimumigiloj, kaj trovado de la radikoj de sistemo de ekvacioj uzante la scipy.optimize. pako.

fonto: https://docs.scipy.org/doc/scipy/reference/

fonto: www.habr.com

Aldoni komenton