Algorytm Gaussa-Newtona - Gauss–Newton algorithm

Dopasowanie zaszumionej krzywej za pomocą asymetrycznego modelu szczytowego przy użyciu algorytmu Gaussa-Newtona ze zmiennym współczynnikiem tłumienia α.
U góry: surowe dane i model.
Dół: Ewolucja znormalizowanej sumy kwadratów błędów.

Algorytm Gaussa-Newtona służy do rozwiązania nieliniowych najmniejszych kwadratów problemy. Jest to modyfikacja metody Newtona dla znalezienia minimum o funkcji . W przeciwieństwie do metody Newtona, algorytm Gaussa-Newtona może być użyty tylko do zminimalizowania sumy kwadratów wartości funkcji, ale ma tę zaletę, że nie są wymagane drugie pochodne, które mogą być trudne do obliczenia.

Nieliniowe problemy najmniejszych kwadratów pojawiają się na przykład w regresji nieliniowej , gdzie parametry w modelu są poszukiwane tak, aby model był w dobrej zgodności z dostępnymi obserwacjami.

Metoda została nazwana na cześć matematyków Carla Friedricha Gaussa i Izaaka Newtona , a po raz pierwszy pojawiła się w pracy Gaussa z 1809 r. Theoria motus corporum coelestium in sectionibus conicis solem ambientum .

Opis

Podane funkcje (często nazywane resztami) zmiennych za pomocą algorytmu Gaussa-Newtona iteracyjnie znajdują wartość zmiennych, która minimalizuje sumę kwadratów

Zaczynając od wstępnego zgadywania minimum, metoda przebiega według iteracji

gdzie, jeśli r i βwektorami kolumnowymi , wpisy macierzy jakobianu to

a symbol oznacza transpozycję macierzy .

Jeśli m  =  n , iteracja upraszcza się do

co jest bezpośrednim uogólnieniem metody Newtona w jednym wymiarze.

W dopasowywaniu danych, gdzie celem jest znalezienie takich parametrów , aby dana funkcja modelu najlepiej pasowała do niektórych punktów danych , funkcje są resztami :

Następnie metodę Gaussa-Newtona można wyrazić w kategoriach jakobianu funkcji as

Należy pamiętać, że jest lewa pseudoinverse od .

Uwagi

Założenie m  ≥  n w oświadczeniu algorytmu jest konieczne, gdyż w przeciwnym razie macierz nie jest odwracalna, a równania normalne nie mogą być rozwiązane (przynajmniej jednoznacznie).

Algorytm Gaussa-Newtona można wyprowadzić przez liniową aproksymację wektora funkcji r i . Korzystając z twierdzenia Taylora możemy w każdej iteracji napisać:

z . Zadanie znalezienia minimalizacji sumy kwadratów po prawej stronie; tj,

jest liniowym problemem najmniejszych kwadratów , który można rozwiązać w sposób jawny, dając w algorytmie równania normalne.

Równania normalne są n równoczesnymi równaniami liniowymi w nieznanych przyrostach . Można je rozwiązać w jednym kroku, stosując dekompozycję Cholesky'ego lub, lepiej, faktoryzację QR dla . W przypadku dużych systemów metoda iteracyjna , taka jak metoda gradientu sprzężonego , może być bardziej wydajna. Czy występuje liniowa zależność między kolumnami J r iteracje nie uda, ponieważ staje się pojedynczą.

Gdy jest złożony, należy użyć formy sprzężonej: .

Przykład

Obliczona krzywa uzyskana za pomocą i (na niebiesko) w porównaniu z obserwowanymi danymi (na czerwono)

W tym przykładzie algorytm Gaussa-Newtona zostanie użyty do dopasowania modelu do niektórych danych poprzez zminimalizowanie sumy kwadratów błędów między danymi a przewidywaniami modelu.

W eksperymencie biologicznym badającym zależność między stężeniem substratu [ S ] a szybkością reakcji w reakcji enzymatycznej uzyskano dane przedstawione w poniższej tabeli.

i 1 2 3 4 5 6 7
[ S ] 0,038 0,194 0,425 0,626 1,253 2.500 3,740
Wskaźnik 0,050 0,127 0,094 0,2122 0,2729 0,2665 0,3317

Pożądane jest znalezienie krzywej (funkcji modelu) postaci

który najlepiej pasuje do danych w sensie najmniejszych kwadratów, z parametrami i do ustalenia.

Oznacz odpowiednio przez i wartości [ S ] i stawki za pomocą . Niech i . Znajdziemy i takie, że suma kwadratów reszt

jest zminimalizowany.

Jakobian wektora reszt względem niewiadomych jest macierzą z -tym wierszem zawierającym wpisy

Rozpoczynając od wstępnych oszacowań i , po pięciu iteracjach algorytmu Gaussa-Newtona uzyskuje się wartości optymalne i . Suma kwadratów reszt zmniejszyła się z wartości początkowej 1,445 do 0,00784 po piątej iteracji. Wykres na rysunku po prawej przedstawia krzywą wyznaczoną przez model dla optymalnych parametrów z obserwowanymi danymi.

Właściwości konwergencji

Można wykazać, że przyrost Δ jest kierunek zejście do S , a gdy zbieżny algorithm to granica jest stacjonarny punkt z S . Jednak zbieżność nie jest gwarantowana, nawet lokalna zbieżność, jak w metodzie Newtona , lub zbieżność w zwykłych warunkach Wolfe'a.

Szybkość zbieżności algorytmu Gaussa-Newtona może zbliżyć się do kwadratu . Algorytm może być zbieżny powoli lub wcale, jeśli początkowe domysły są dalekie od minimum lub macierz jest źle uwarunkowana . Rozważmy na przykład problem z równaniami i zmienną, podaną przez

Optimum jest przy . (Właściwie optimum to dla , ponieważ , ale .) Jeśli , to problem jest w rzeczywistości liniowy i metoda znajduje optimum w jednej iteracji. Jeżeli |λ| < 1, to metoda zbiega się liniowo i błąd maleje asymptotycznie ze współczynnikiem |λ| w każdej iteracji. Jednakże, jeśli |λ| > 1, to metoda nie jest nawet lokalnie zbieżna.

Wyprowadzenie z metody Newtona

W dalszej części algorytm Gaussa-Newtona zostanie wyprowadzony z metody Newtona do optymalizacji funkcji przez aproksymację. W konsekwencji szybkość zbieżności algorytmu Gaussa-Newtona może być kwadratowa w pewnych warunkach regularności. Generalnie (w słabszych warunkach) współczynnik konwergencji jest liniowy.

Relacja rekurencyjna dla metody Newtona minimalizacji funkcji S parametrów wynosi

w którym g oznacza wektor gradientu z S , a H oznacza Hessian matrycę z S .

Ponieważ gradient jest podany przez

Elementy hesu są obliczane przez zróżnicowanie elementów gradientu , w odniesieniu do :

Metodę Gaussa-Newtona uzyskuje się przez zignorowanie terminów pochodnych drugiego rzędu (drugi termin w tym wyrażeniu). Oznacza to, że hes jest przybliżony przez

gdzie są wpisy jakobian J r . Gradient i przybliżony hesjan można zapisać w notacji macierzowej jako

Wyrażenia te są podstawiane do powyższej relacji rekurencyjnej w celu uzyskania równań operacyjnych

Zbieżność metody Gaussa-Newtona nie jest gwarantowana we wszystkich przypadkach. Przybliżenie

które muszą być spełnione, aby móc zignorować terminy pochodne drugiego rzędu, mogą być ważne w dwóch przypadkach, dla których należy oczekiwać zbieżności:

  1. Wartości funkcji są małe, przynajmniej w okolicach minimum.
  2. Funkcje są tylko „łagodnie” nieliniowe, więc ich wielkość jest stosunkowo niewielka.

Ulepszone wersje

Z metodą Gaussa-Newtona suma kwadratów reszt S może nie zmniejszać się w każdej iteracji. Jednakże, ponieważ Δ jest kierunkiem opadania, chyba że jest punktem stacjonarnym, to dla wszystkich wystarczająco małych . Zatem w przypadku wystąpienia dywergencji jednym z rozwiązań jest zastosowanie ułamka wektora przyrostu Δ we wzorze aktualizacyjnym:

.

Innymi słowy, wektor przyrostu jest zbyt długi, ale nadal wskazuje „w dół”, więc przejście tylko części drogi zmniejszy funkcję celu S . Optymalną wartość dla można znaleźć za pomocą algorytmu wyszukiwania linii , to znaczy, że wielkość jest określana przez znalezienie wartości, która minimalizuje S , zwykle przy użyciu metody wyszukiwania bezpośredniego w przedziale lub wyszukiwania linii z cofaniem, takiego jak wyszukiwanie linii Armijo . Zazwyczaj powinien być tak dobrany, aby spełniał warunki Wolfe'a lub Goldsteina .

W przypadkach, gdy kierunek wektora przesunięcia jest taki, że optymalny ułamek α jest bliski zeru, alternatywną metodą obsługi dywergencji jest użycie algorytmu Levenberga-Marquardta , metody regionu zaufania . Równania normalne są modyfikowane w taki sposób, że wektor przyrostu jest obracany w kierunku najbardziej stromego opadania ,

gdzie D jest dodatnią macierzą diagonalną. Należy zauważyć, że jeśli D oznacza macierz identyczności I i , następnie , w związku z czym kierunek hemibursztynianu zbliża się do kierunku gradientu ujemnym .

Tak zwany parametr Marquardt może być również zoptymalizowany przez przeszukiwanie linii, ale jest to nieefektywne, ponieważ wektor przesunięcia musi być przeliczany przy każdej zmianie. Bardziej wydajna strategia jest następująca: Gdy wystąpi dywergencja, zwiększaj parametr Marquardt, aż nastąpi spadek S . Następnie zachowaj wartość z jednej iteracji do następnej, ale zmniejszaj ją, jeśli to możliwe, aż do osiągnięcia wartości odcięcia, kiedy parametr Marquardt można ustawić na zero; minimalizacja S staje się wtedy standardową minimalizacją Gaussa-Newtona.

Optymalizacja na dużą skalę

W przypadku optymalizacji na dużą skalę metoda Gaussa-Newtona jest szczególnie interesująca, ponieważ często (choć na pewno nie zawsze) prawdą jest, że macierz jest bardziej rzadka niż przybliżona Hessian . W takich przypadkach samo obliczenie kroku będzie zazwyczaj musiało być wykonane przy użyciu przybliżonej metody iteracyjnej odpowiedniej dla dużych i rzadkich problemów, takiej jak metoda gradientu sprzężonego .

Aby takie podejście zadziałało, potrzebna jest przynajmniej wydajna metoda obliczania produktu

dla niektórych wektorów p . W przypadku przechowywania macierzy rzadkiej , ogólnie praktyczne jest przechowywanie wierszy w postaci skompresowanej (np. bez wpisów zerowych), co utrudnia bezpośrednie obliczenie powyższego produktu ze względu na transpozycję. Jeśli jednak zdefiniujemy c i jako wiersz i macierzy , zachodzi następująca prosta zależność:

aby każdy rząd wnosił addytywny i niezależny wkład do produktu. Poza przestrzeganiem praktycznej, rzadkiej struktury pamięci, wyrażenie to dobrze nadaje się do obliczeń równoległych . Zauważ, że każdy wiersz c i jest gradientem odpowiadającej mu reszty r i ; mając to na uwadze, powyższy wzór podkreśla fakt, że reszty przyczyniają się do problemu niezależnie od siebie.

Powiązane algorytmy

W metodzie quasi-Newtona , takiej jak metoda Davidona, Fletchera i Powella lub Broydena-Fletchera-Goldfarba-Shanno ( metoda BFGS ), oszacowanie pełnego Hessu jest budowane numerycznie przy użyciu tylko pierwszych pochodnych, tak że po n cyklach udokładniania Metoda jest bardzo zbliżona do metody Newtona pod względem wydajności. Zauważ, że metody quasi-Newtona mogą minimalizować ogólne funkcje o wartościach rzeczywistych, podczas gdy Gauss-Newton, Levenberg-Marquardt itp. pasują tylko do nieliniowych problemów najmniejszych kwadratów.

Inną metodą rozwiązywania problemów minimalizacji przy użyciu tylko pierwszych pochodnych jest metoda gradientu . Jednak ta metoda nie uwzględnia drugich pochodnych nawet w przybliżeniu. W konsekwencji jest wysoce nieefektywny dla wielu funkcji, zwłaszcza jeśli parametry mają silne interakcje.

Uwagi

  1. ^ Mittelhammer, Ron C.; Miller, Douglas J.; Sędzia George G. (2000). Podstawy ekonometryczne . Cambridge: Wydawnictwo Uniwersytetu Cambridge. s. 197–198. Numer ISBN 0-521-62394-4.
  2. ^ Floudas, Christodoulos A. ; Pardalos, Panos M. (2008). Encyklopedia Optymalizacji . Skoczek. P. 1130. Numer ISBN 9780387747583.
  3. ^ B Björck (1996)
  4. ^ Björck (1996), s. 260.
  5. ^ Mascarenhas (2013), "Rozbieżność metod BFGS i Gaussa Newtona", Programowanie matematyczne , 147 (1): 253-276, arXiv : 1309.7922 , doi : 10.1007/s10107-013-0720-6
  6. ^ Björck (1996), s. 341, 342.
  7. ^ Fletcher (1987), s. 113.
  8. ^ „Zarchiwizowana kopia” (PDF) . Zarchiwizowane z oryginału (PDF) w dniu 2016-08-04 . Pobrano 2014-04-25 .CS1 maint: zarchiwizowana kopia jako tytuł ( link )
  9. ^ Nocedal (1999), s. 259.
  10. ^ Nocedal, Jorge. (1999). Optymalizacja numeryczna . Wright, Stephen J., 1960-. Nowy Jork: Springer. Numer ISBN 0387227423. OCLC  54849297 .

Bibliografia

Linki zewnętrzne

Realizacje

  • Artelys Knitro to nieliniowy solwer z implementacją metody Gaussa-Newtona. Jest napisany w C i posiada interfejsy do C++/C#/Java/Python/MATLAB/R.