• y(x;a1, . . . , aM)hängt nicht linear von seinen Parameternai ab. Die funktionale Form vony ist komplett frei (dem Nutzer überlassen).
• Vorsicht bei Scheinparametern! Das Modell
y(x) =aexp(−bx+d) (11.46)
mit (scheinbar) drei Parametern (a, b, d) enthält tatsächlich nur zwei (b, c):
y(x) =aexp(−bx+d) = aexp(d) exp(−bx) =cexp(−bx). (11.47) Versucht man trotzdem, das Modell in der ursprünglichen Form mit drei Parametern zu fitten, wird das resultierende Gleichungssystem exakt singulär und einfache Verfahren versagen (nicht jedoch SVD).
• Wiederum wird eine Fehlerquadratsummeχ2minimiert, aber Ausgangspunkt und Strategie sind i. A. etwas anders als bei linearer Regression:
– Das Modell,y(x, a1, . . . , aM), kann beliebig allgemein sein und nichtlinear von seinen Parameternaiabhängen.
– Im Gegensatz zur einzigen, eindeutigen Lösung bei der linearen Regression (M Glei-chungen fürM Unbekannte) gibt es hier einen wesentlich komplizierteren funktionalen Verlauf vonχ2und eine Vielzahl möglicher Lösungen. Daher
∗ braucht man einiteratives Verfahren(basierend auf einem einfachen, approximativen Modell an den kompliziertenχ2-Verlauf),
∗ ist es nicht möglich, in einem Schritt die „richtige“ (welche?) Lösung zu finden,
∗ ist eine gefundene Lösung evtl. weder die beste noch die eigentlich gesuchte; ver-schiedene Startpunkteder Iteration liefern evtl.verschiedene Lösungen.
• Daher: Vorgang analog zu Abschnitt 6.4: Die Approximation einer beliebigen mehrdimensio-nalen Funktionf(a)an der Stellea= (a1, a2, . . . , an)T geschieht mit einer Taylorreihe um einen Entwicklungspunkta0 (hier zur Vereinfachunga0 =0):
f(a) =f(a0) +X
• Analog zur Prozedur in Abschnitt 6.4 ergibt Abbruch der Reihenentwicklung nach 2. Ordnung die iterative Formel
ai−am =D−1 ∇f|a
i. (11.50)
• Gegeben sei nun eine beliebige Modellfunktiony=y(x,a)mit linearen und/oder nichtlinearen Parameternak,k = 1,2, . . . , M.
• Wie bei der linearen Regression versuchen wir, diese Parameter dadurch zu bestimmen, dass wir das Minimum des Fit-Maßesχ2 suchen, gegenüberN gegebenen Datenpunkten(xi, yi), i= 1,2, . . . , N:
• Diese im Praxisfall vermutlich komplizierte Funktion approximieren wir wie oben mit einer nach dem quadratischen Term abgebrochenen Taylorreihe:
χ2(a)≈γ−dTa+ 1
2aTDa. (11.52)
• Wichtig: Im Gegensatz zu dem Normalfall einer mehrdimensionalen Minimierung ist hier die Hesse-Matrixbekannt.
• Wenn die quadratische Approximation also vernünftig ist, ergibt sich, ausgehend von einem Startvektora0, iterativ die bessere Lösung
aj =ai−D−1 ∇χ2
ai. (11.53)
• Diese Näherung muss aber nicht notwendigerweise gut sein! In diesem Fall ist ein kleiner steepest-descent-artiger Schritt besser:
aj =ai−const.· ∇χ2
ai. (11.54)
• Nur Anwendung von Gl.(11.53)(Newton-Verfahren) sollte in der Praxis sonicht angewendet werden, da es nur in der Nähe von Minima wie gewünscht funktioniert. An anderen Stellen ergibt sich i. A.keinerleiKonvergenz zu einem Minimum.
• Auch wenn wir das Modelly(x,a)überhaupt nicht näher spezifizieren, können wir dennoch aufgrund der allgemeinen Form vonχ2nach Gl.(11.51)den Gradienten und die Hesse-Matrix angeben, vorausgesetzt, wir kennen auch die Ableitungen vony(x,a)nacha.
• Für den Gradienten ergibt sich
∂χ2 und für die Elemente der Hesseschen Matrix (Produktregel):
∂2χ2
• Da wir ohnehin nur approximativ-iterativ rechnen, können wir es uns hier leisten, den zweiten Term gegenüber dem ersten zu vernachlässigen, was auch deshalb berechtigt ist, weil er den (vermutlich kleinen) Faktor(yi−y(xi,a))enthält, der nicht mit dem Modell korreliert und bei der Summenbildung zur Kompensation tendiert. Einbeziehung des Terms führt daher oftmals zuschlechterenKonvergenz.
• Wichtig: Da der Gradient nicht geändert wird, führen Variationen an der Hesse-Matrix dennoch zum (lokalen) Minimum.
157
• Mit dieser „Näherung“ und mit den Abkürzungen
wird aus Gl.(11.53)für den Newton-artigen Schritt:
M
X
`=1
αk`δa` =βk , k = 1,2, . . . , M (11.60)
• Dies ist eine lineare Matrix-Vektor-Gleichung α
αα δa=βββ. (11.61)
• Gleichzeitig wird aus Gl.(11.54)für den steepest-descent-Schritt
δa`=const.·β`. (11.62)
• Für die Regression wird also sowohl die Modellfunktion nebst Gradienten inaals auch die Jacobi-Matrix der Modellfunktion benötigt.
11.3.1 Levenberg-Marquardt-Methode
• Wir wissen qualitativ, wann wir Newton-Schritte und wann wir steepest-descent-Schritte machen sollten.
• Dies würden wir jedoch gerne automatisiert entscheiden lassen und dabei eine Art „fließenden Übergang“ zwischen beiden Schritt-Typen konstruieren.
• Es zeigt sich, daß wir zu diesem Zweck die noch freie Konstante in Gl. (11.62)ausnutzen können.
• Setzen wir für diese Konstante (diagonale Approximation der Hesse-Matrix) const.= 1
λ α`` (11.63)
mit einem adjustierbaren Parameterλan, wird dadurch aus Gl.(11.62)
λ α``δa` =β`. (11.64)
• Dadurch können wir diese modifizierte Gleichung für den steepest-descent-Schritt „auf die Diagonale“ der Matrixgleichung Gl.(11.61)für den Newton-Schritt bringen, indem wir eine neue Matrixααα0definieren, mit den Elementen
α0jk = αjk, für j 6=k, (11.65)
α0jj = αjj(1 +λ). (11.66)
• Dadurch können wir beide Schritt-Typen in einem einzigen linearen Gleichungssystem vereini-gen:
ααα0δa=βββ. (11.67)
• Wenn dabei der adjustierbare Parameterλgroß wird, wird die Matrixααα0 diagonal dominant und wir führen steepest-descent-Schritte aus.
• Fürλ →0, gehtααα0gegen die unmodifizierte Matrixαααund wir machen Newton-Schritte.
Algorithmus:
• Die grobe Programmstruktur für den Levenberg-Marquardt-Algorithmus sieht folgendermaßen aus:
– Rate einen Satz von Parameterna(Startpunkt).
– Berechne das zugehörigeχ2(a).
– Rate einen Startwert fürλ(z.B.λ= 0.001).
∗ Stelle die Matrix und den Rechte-Seiten-Vektor für das Gleichungssystem Gl.(11.67)auf.
– Löse dieses Gleichungssystem;
die Lösung ist der Schrittvektorδa. – Berechneχ2(a+δa).
– Wennχ2(a+δa)≥χ2(a), dann war der Schritt schlecht und wird verworfen;
vergrößereλ(z.B. Multiplikation mit 10) und gehe zurück zu∗.
– Anderenfalls war der Schritt in Ordnung und wird akzeptiert (d.h. wir setzena=a+δa);
für den nächsten Schritt probieren wir eine Verringerung vonλ(z.B. Division durch 10) und gehen zu∗.
• Die Konstruktion der Matrix und des Rechte-Seiten-Vektors ist sehr ähnlich zum Fall der linearen Regression.
• Das Verhalten des Algorithmus in dieser Form ist allerdings in der Praxis nicht unkritisch:
Großeλ-Werte garantieren, dassααα0positiv definit ist (was im Minimum der Fall sein sollte) und nicht singulär wird, bedeuten aber auch kleine steepest-descent-Schrittchen – und umgekehrt.
• Daher gilt auch hier, dass gegebenenfalls Singulärwertzerlegung verwendet werden sollte.
• Professionelle Implementationen haben noch raffiniertere Stabilisierungsmethoden.
• Die Iteration sollte beendet werden, wenn folgendes eintritt:
δχ2 =|χ2(a+δa)−χ2(a)| < 0.1 (absolut), (11.68) oder δχ2
χ2 < 10−3 (relativ). (11.69)
• Ein lediglich steigendesχ2ist kein Grund zum Abbruch, da diese Situation in den folgenden Iterationsschritten in der Regel vonλ-Änderungen ( = Schritt-Typ-Änderungen) abgefangen wird.
159
• Nach Konvergenz ergibt sich die Covarianz-Matrix wie bei der linearen Regression durch (λ= 0!)
C=ααα−1. (11.70)
• Die Optimierung sollteimmer mitverschiedenenStartparametern adurchgeführt werden.
Je nach Modell und Parametersatz ist die Regression sehr schwer und es gibt viele mögli-che Lösungen. Durch Variation der Startparameter lässt sich auch die Empfindlichkeit der Fitparameteraiabschätzen.
• Wichtig: Das Verfahren kann auch scheitern! (Keine Konvergenz.)