Gauß-Newton-Verfahren

Grundzüge des Verfahrens

• Die Tabelle der Messwerte hat $k$ Zeilen, es wurden also $k$ Messungen durchgeführt
• Bei jeder Messung wurden $n$ Stellgrößen $\, x_1, \dots , x_n$ vorgegeben und das Ergebnis $y$ gemessen
• Es gibt eine Modellfunktion $\, f(x_1, \dots , x_n) = y$ , welche den Zusammenhang zwischen $\, x_1, \dots , x_n$ und $y$ beschreibt. Diese Funktion hat $p$ verschiedene Parameter $\, a_1, \dots , a_p$, die nun so berechnet werden sollen, dass mit der Funktion y möglichst genau aus den $\, x_1, \dots , x_n$ berechnet werden kann. Das besondere an der Funktion ist, dass sie nichtlinear in den Parametern ist.

$\sum\limits_{i=1}^k (f(x_{i,1}, \dots , x_{i,n}) - y_i)^2 \rightarrow min!$

Algorithmus

Vorbereitung

• Gegeben ist die Funktion $\, f(x_1, \dots , x_n) = y$ mit n Variablen $x_1, \dots , x_n \,$ und p gesuchten Parametern $a_1, \dots , a_p \,$
• Aufstellen der Residuumsfunktion $\, r$ und des Residuenvektors $\mathbf{r}$

$\, r = f(x_1, \dots , x_n) - y$

$r_i = f(x_{i,1}, \dots , x_{i,n}) - y_i, i=1\, \, k \,$

• Berechnung der partiellen Ableitungen der Residuumsfunktion $\, r$ nach jedem Parameter $( a_1, \dots , a_p \, )$:

$r'_1 = \dfrac{\partial r}{\partial a_1}; \quad r'_2 = \dfrac{\partial r}{\partial a_2}; \quad \dots ; \quad r'_p = \dfrac{\partial r}{\partial a_p}$

$\, r'_{i,j} = r'_j(x_{i,1}, \dots , x_{i,n})$

• Aufbau der Matrizen für die Iteration:

$\mathbf{D} = \begin{pmatrix} r'_{1,1} & r'_{1,2} & \cdots & r'_{1,p} \\ r'_{2,1} & r'_{2,2} & \cdots & r'_{2,p} \\ \vdots & \vdots & \vdots & \vdots \\ r'_{k,1} & r'_{k,2} & \cdots & r'_{k,p} \end{pmatrix}$; $\mathbf{r} = \begin{pmatrix} r_1 \\ r_2 \\ \vdots \\ r_k \end{pmatrix}$; $\mathbf{a} = \begin{pmatrix} a_1 \\ a_2 \\ \vdots \\ a_p \end{pmatrix}$

Für den Aufbau der Matrizen ist folgendes zu beachten:

• Die Matrix $\bm{D}$ und der Spaltenvektor $\bm{r}$ haben $k$ Zeilen, also für jede Zeile der oben angegebenen Tabelle eine.
• Der Spaltenvektor $\bm{a}$ hat $p$ Zeilen, also für jeden Parameter $a_1, \dots , a_p \,$ eine
• Die Spalten in der Matrix $\bm{D}$ sind die partiellen Ableitungen nach den Parametern $a_1, \dots , a_p \,$. Die Reihenfolge der Spalten in $\bm{D}$ hängt mit der Reihenfolge der Parameter in $\bm{a}$ zusammen. Steht in Zeile 1 von $\bm{a}$ der Parameter $a_1 \,$, so muss in $\bm{D}$ die erste Spalte die Ableitungen nach $a_1 \,$ enthalten. Dementsprechend hat $\bm{D}$ $p$ Spalten, also für jeden Parameter $a_1, \dots , a_p \,$ eine.

• Die Anzahl der Variablen n hat keinen Einfluss auf den Aufbau der Matrix $\bm{D}$ und der beiden Vektoren $\bm{r},\, \bm{a}$.
• Zu Beginn der Iteration müssen Startwerte für die Parameter $a_1, a_2, \dots , a_p \,$ festgelegt werden

Iteration

• Die Iteration wird mit folgender Matrixgleichung durchgeführt:

$\mathbf{a}_{i+1} = \mathbf{a}_i - (\mathbf{D}^T \cdot \mathbf{D})^{-1} \cdot \mathbf{D}^T \cdot \mathbf{r}$

Dabei wird in jedem Schritt der Vektor $\mathbf{a}_i$, der die Parameter $a_1, a_2, \dots , a_p \,$ enthält, verbessert.

Die Matrix $\bm{D}$ wird berechnet, indem man zunächst alle Werte in den $k$ Zeilen der Tabelle in die Funktion $r'_1 \,$ einsetzt. Das Ergebnis schreibt man untereinander in die Spalte 1 von $\bm{D}$. Danach setzt man alle Werte in den $k$ Zeilen der Tabelle in die Funktion $r'_2 \,$ ein und schreibt sie in Spalte 2 der Matrix $\bm{D}$ usw.

Um den Vektor $\bm{r}$ zu berechnen, setzt man alle Werte in den $k$ Zeilen der Tabelle in die Funktion $r_i \,$ und schreibt die Ergebnisse jeweils untereinander als Vektor $\bm{r}$ auf.

Für die numerische Berechnung empfiehlt sich eine Aufspaltung der Berechnung, damit die Matrixinversion durch die Lösung eines linearen Gleichungssystems für den unbekannten Lösungsvektor $\bm{s}$ ersetzt werden kann:

$(\mathbf{D}^T \cdot \mathbf{D}) \cdot \mathbf{s} = \mathbf{D}^T \cdot \mathbf{r}; \qquad \mathbf{a}_{i+1} = \mathbf{a}_i - \mathbf{s}$

Die Vorteile liegen in einer schnelleren Berechnung bei höherer Genauigkeit.

• Sobald $\mathbf{a}_{i+1}$ berechnet wurde, müssen auch Matrizen neu berechnet werden, um den nächsten Iterationsschritt vorzubereiten. Um den Rechenaufwand zu verringern, kann auch mehrfach ohne Neuberechnung von $\bm{s}$ iteriert werden. Dieses Vorgehen wird beim Newtonschen Verfahren häufig empfohlen, reduziert aber die Konvergenzgeschwindigkeit und sollte erst angewendet werden, wenn sich $\bm{a}$ nur noch wenig ändert.
• Die Iteration wird abgebrochen, falls $\mathbf{a}_{i+1} = \mathbf{a}_i$, also bei den $a_1, a_2, \dots , a_p \,$ keine Änderung mehr eintritt.

Anmerkungen

• Da $\mathbf{D}^T \cdot \mathbf{D}$ symmetrisch und positiv definit ist, eignet sich die Cholesky-Zerlegung besonders gut für die Auflösung des Gleichungssystems.
• Die Anzahl der Zeilen in der Tabelle ($k$) muss stets größer, als die Anzahl der Parameter $a_1, a_2, \dots , a_p \,$ sein. Falls die Anzahl der Paramater gleich $k$ ist, bestimmt das Verfahren die Parameter exakt (im Rahmen der Genauigkeit der Iteration), es ist also nicht nur die optimale Lösung im Sinne der Fehlerquadrate. Das System ist unterbestimmt, wenn die Anzahl der Parameter größer als $k$ ist.
• Durch Einführung eines Schrittweiteparameters $\gamma$ lässt sich ein Abstieg, d.h. die Verringerung der Fehlerquadratsumme erreichen.

Literatur

• Ralf Pfeifer: Effektive Messauswertung mit der Gauß'schen Fehlerquadratmethode, ISBN 3-89001-251-5

Ein Mathematiker ist eine Maschine, die Kaffee in Theoreme verwandelt.

Paul Erdös

Optimierung

• Lineare Optimierung
• Konvexe Optimierung
• Ausgleichungsrechnung
  • Methode der kleinsten Quadrate
  • Gauß-Newton-Verfahren
• Optimale Steuerung