Algorytm Levenberga-Marquardta

Algorytm Levenberga-Marquardta – algorytm optymalizacji nieliniowej. Jest to algorytm iteracyjny, łączący w sobie cechy metody największego spadku i metody Gaussa-Newtona.

Sformułowanie problemu

Mając daną serię danych $(t_i,y_i)\in {𝐑}^2,$ gdzie $i=1,2,\dots ,N,$ szukamy dopasowania $\overline{y}=f(t|𝐩),$ gdzie $𝐩\in {𝐑}^n$ – wektor parametrów. Zakładamy, że najlepszym dopasowaniem jest to minimalizujące funkcjonał:

${\chi }^2(f)={\chi }^2(𝐩)=\sum _{i=1}^N[y_i-f(t_i|𝐩){]}^2.$

Algorytm Levenberga-Marquardta w ogólności znajduje rozwiązanie zadania optymalizacji nieliniowej funkcji dającej się zapisać w postaci:

$\mathrm{\Phi }(𝐱)=\frac{1}{2}\sum _{i=1}^N{r}_i^2(𝐱),$

gdzie $𝐱\in {𝐑}^n$ i zakładamy, że $N⩾n.$ Jak łatwo zauważyć, funkcjonał ${\chi }^2$ daje się zapisać w taki sposób. Dla uproszczenia, przedstawmy funkcje $r_i$ jako wektor $𝐫(𝐱)=(r_1(𝐱),\dots ,r_N(𝐱))$ (zwany wektorem rezydualnym). Wtedy $\mathrm{\Phi }(𝐱)=\Vert 𝐫(𝐱){\Vert }^2.$ Pochodne funkcji $\mathrm{\Phi }$ można zapisać przy użyciu Macierzy Jacobiego funkcji $𝐫,$ zdefiniowanego jako $\{𝖩(𝐱){\}}_{ij}=\frac{\partial r_i}{\partial x_j}(𝐱).$ W ogólnym przypadku gradient funkcji $\mathrm{\Phi }$ można zapisać:

$\mathrm{\nabla }\mathrm{\Phi }(𝐱)=\sum _{i=1}^N{r}_i(𝐱)\mathrm{\nabla }{r}_i(𝐱)=𝖩(𝐱{)}^{𝖳}𝐫(𝐱),$

a jej Macierz Hessego:

${\mathrm{\nabla }}^2\mathrm{\Phi }(𝐱)=𝖩(𝐱{)}^{𝖳}𝖩(𝐱)+\sum _{i=1}^N{r}_j(𝐱){\mathrm{\nabla }}^2{r}_j(𝐱).$

W przypadku, gdy funkcje $r_j$ można aproksymować funkcjami liniowymi w otoczeniu interesującego nas punktu (wtedy ${\mathrm{\nabla }}^2{r}_j(𝐱)$ jest bliskie zeru), lub gdy $r_j(𝐱)$ jest małe, hesjan funkcji $\mathrm{\Phi }$ przyjmuje prostszą postać:

${\mathrm{\nabla }}^2\mathrm{\Phi }(𝐱)=𝖩(𝐱{)}^{𝖳}𝖩(𝐱),$

a więc hesjan można otrzymać wprost mając dany jakobian wektora rezydualnego $𝐫(𝐱),$ co jest charakterystyczne dla zadania najmniejszych kwadratów.

Opis metody

Najprostszym podejściem do problemu minimalizacji funkcji $\mathrm{\Phi }$ jest metoda największego spadku, opisana schematem:

${𝐱}_{i+1}={𝐱}_i-\lambda \mathrm{\nabla }\mathrm{\Phi }({𝐱}_i),$

która jest, w ogólnym przypadku, wolno zbieżna. Aby poprawić jej zbieżność, można skorzystać z wiedzy o drugiej pochodnej minimalizowanej funkcji w badanym punkcie. Jednym z możliwych podejść jest rozwinięcie gradientu minimalizowanej funkcji w szereg Taylora:

$\mathrm{\nabla }\mathrm{\Phi }(𝐱)=\mathrm{\nabla }\mathrm{\Phi }({𝐱}_0)+(𝐱-{𝐱}_0{)}^{𝖳}{\mathrm{\nabla }}^2\mathrm{\Phi }({𝐱}_0)+\dots$

i przyjęcie przybliżenia kwadratowego funkcji $\mathrm{\Phi }$ w otoczeniu ${𝐱}_0$ do rozwiązania równania $\mathrm{\nabla }\mathrm{\Phi }(\overline{𝐱})=0.$ W ten sposób otrzymujemy metodę Gaussa-Newtona opisaną schematem:

${𝐱}_{i+1}={𝐱}_i-({\mathrm{\nabla }}^2\mathrm{\Phi }({𝐱}_i){)}^{-1}\mathrm{\nabla }\mathrm{\Phi }({𝐱}_i),$

gdzie hesjan funkcji $\mathrm{\Phi }$ nie musi być znany dokładnie i często wystarczy podane wcześniej przybliżenie. Niestety, szybkość zbieżności tej metody zależy od wyboru punktu początkowego, a konkretnie od liniowości minimalizowanej funkcji w otoczeniu punktu startowego. Kenneth Levenberg zauważył, że opisane metody (największego spadku i Gaussa-Newtona) nawzajem się uzupełniają i zaproponował następującą modyfikację kroku metody:

${𝐱}_{i+1}={𝐱}_i-(𝖧({𝐱}_i)+\lambda 𝖨{)}^{-1}\mathrm{\nabla }\mathrm{\Phi }({𝐱}_i),$(*)

wraz z następującym algorytmem:

1. oblicz wartość ${𝐱}_{i+1}$ na podstawie ${𝐱}_i$ i równania (*),
2. oblicz wartość błędu w punkcie ${𝐱}_{i+1},$
3. jeśli błąd wzrósł, wróć do wartości ${𝐱}_i,$ zwiększ wartość $\lambda$ $k$-krotnie i wróć do kroku 1 (przybliżenie liniowe minimalizowanej funkcji w otoczeniu ${𝐱}_i$ okazało się nie dość ścisłe, więc zwiększamy „wpływ” metody największego spadku),
4. jeśli błąd zmalał, zaakceptuj ten krok i zmniejsz wartość $\lambda$ $k$-krotnie (założenie o liniowości minimalizowanej funkcji w otoczeniu ${𝐱}_i$ okazało się wystarczająco ścisłe, więc zwiększamy „wpływ” metody Gaussa-Newtona).

W typowych zastosowaniach $k=10.$ W przypadku, gdy $\lambda$ jest duże, hesjan w zasadzie nie jest wykorzystywany. Donald Marquardt zauważył, że nawet w takiej sytuacji można wykorzystać informację zawartą w drugiej pochodnej minimalizowanej funkcji, poprzez skalowanie każdego komponentu wektora gradientu w zależności od krzywizny w danym kierunku (co pomaga w źle uwarunkowanych zadaniach minimalizacji typu error valley). Po uwzględnieniu poprawki Marquardta otrzymujemy następującą postać kroku metody:

${𝐱}_{i+1}={𝐱}_i-(𝖧({𝐱}_i)+\lambda \text{diag}[𝖧]{)}^{-1}\mathrm{\nabla }\mathrm{\Phi }({𝐱}_i),$

gdzie:

$\text{diag}[𝖧]=\left[\begin{array}{cccc}{h}_{11}& 0& \dots & 0\\ 0& h_{22}& \dots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \dots & h_{nn}\end{array}\right].$

Największą zaletą algorytmu Levenberga-Marquardta jest jego szybka zbieżność, w porównaniu z konkurencyjnymi metodami. Najkosztowniejszą operacją jest natomiast wyznaczenie macierzy odwrotnej, które w praktyce jest przeprowadzane w sposób przybliżony, na przykład przy użyciu metody SVD. Tym niemniej, nawet w najoszczędniejszych przypadkach koszt jednego kroku rośnie niedopuszczalnie wraz ze wzrostem rozmiaru zadania powyżej tysiąca parametrów. Z drugiej dla zadań o umiarkowanej ilości parametrów (rzędu kilkuset), metoda Levenberga-Marquardta jest dużo szybsza od metody największego spadku.

Zobacz też

• optymalizacja (matematyka)

Linki zewnętrzne

• Szablon:Cytuj stronę