Ausgleichsprobleme mit Nebenbedingungen

Lernziele

  • Sie verstehen das Konzept von Ausgleichsproblemen mit Nebenbedingungen.

  • Sie kennen eine approximative Lösungsmethode für diese Problemstellung und können diese implementieren.

  • Sie verstehen damit zusammenhängende, numerische Herausforderungen und berücksichtigen diese Aspekte bei der Lösung.

Theorie

Mitunter kann eine Fitaufgabenstellung nicht in der Form gestellt werden, welche auf eine direkte Minimierung von \(\Vert {\bf A} \, {\bf \xi} - {\bf b} \Vert^2 \to min\) führt. Eine Beispielsituation ist unten dargestetllt: von einem Prozess sei bekannt, dass die Daten in \(0 \le x < 5\) durch ein Polynom 3. Grades \(p_1(x)\) beschrieben werden, während im Interval \(5 \le x \le 10\) ein anderes Polynom 3. Grades \(p_2(x)\) vorliegt, dessen Koeffizienten verschieden sein können vom Polynom \(p_1(x)\). Nichtsdestotrotz verlangt man aber, dass die Funktion bestehend aus den beiden Polynomen an der Stelle \(x=5\) stetig ist und eine stetige 1. Ableitung aufweist!

datenLSmitNebenbed.png

Das Polynom \(p_1\) ist definiert durch die Koeffizienten \(\xi_i, \, i=1,\dots,4\), während das Polynom \(p_2\) die Koeffizienten \(\xi_i, \, i=5,\dots,8\) aufweist. Ein Moment der Überlegung zeigt, dass die Forderung nach der Stetigkeit der 0. und 1. Ableitung an der Stelle \(x=a\) zu linearen Bedingungsgleichungen für die Koeffizienten führt. Diese Bedingungen können in der Form \(\bf C \, \bf \xi = 0\) formuliert werden, wo \(\bf C\) eine zu bestimmende 2x8-Matrix darstellt. Die vorliegende Problemstellung lautet somit

\[\Vert {\bf A} \, {\bf \xi} - {\bf b} \Vert^2 \to min, \text{ unter der Bedingung } {\bf C} \, {\bf \xi} = 0\]

In dieser Situation können wir nicht direkt eine Normalengleichung formulieren und den gewohnten Lösungsweg einschlagen. Es gibt aber einen Trick, um zumindest approximativ wieder auf das gewohnte Minimierungsproblem zurückzukommen und dieses zu lösen. Wir betrachten dazu den Ausdruck

\[\Vert {\bf A} \, {\bf \xi} - {\bf b} \Vert^2 + \lambda \cdot \Vert {\bf C} \, {\bf \xi} \Vert^2 \to min, \text{ mit } \lambda > 0\]

Wenn \(\lambda\) genügend gross gewählt wird, führt die Minimierung dieses Ausdruckes dazu, dass die Nebenbedingung approximativ erfüllt wird, denn nur so kann der 2. Term, \(\Vert {\bf C} \, {\bf \xi} \Vert^2\), minimal werden. Durch diese (kontrollierbare) Approximation kann die obige Problemstellung auf eine lineare Ausgleichsrechnung zurückgeführt werden, denn bei konstantem \(\lambda\) sind alle Terme in dieser Gleichung quadratische Ausdrücke in den Vektorkomponenten von \(\bf \xi\).

Aufträge

  1. Daten einlesen

    Sie finden die entsprechenden Daten via den Link weiter unten. Die Daten lesen Sie in MATLAB ein mit:

    load('dataLSforMatlab')

    resp. unter PYTHON:

    import pickle as pkl
data=pkl.load(open('dataLSforPython.pickle','rb'))

    Womit Ihnen ein 2D-Array \(data\) zur Verfügung steht, welcher spaltenweise die Werte der Variablen \(x_1\), \(x_2\), \(y_1\), \(y_2\) enthält.

  2. Fit mit einem einzelnen Polynom 3. Grades

    Zuerst versuchen Sie die Daten mit einem Fit an eine einzelnes Polynom 3. Grades. Bei einer graphische Visualisierung werden Sie feststellen, dass dieser Ansatz nicht befriedigend ist.

  3. Fit mit zwei unabhängigen Polynomen 3. Grades

    Nun testen Sie eine 2. Lösung, nämlich der Fit mit 2 unabhängigen Polynomen 3. Grades, je für die Daten \(x_1, \, y_1\) und \(x_2, \, y_2\). Stellen Sie diese Lösung im selben Graphen dar, wie das einzelne Polynom und notieren Sie, was Ihnen auffällt.

    Hinweis: Sie können diese beiden Fits mit einer einzigen, so genannten Blockdiagonalmatrix erstellen. Block oben links entspricht dem Fit der Daten zu \(x_1\), der Block unten rechts gehört zu den Daten für \(x_2\). Alle übrigen Matrixelemente sind gleich Null. Diese Matrixstruktur gibt dann entsprechend auch die Aufteilung des Vektors der Unbekannten \(\xi_i\), resp. der rechten Seite, bestehend aus den \(y\)-Datenwerten vor.

  4. Fit zwei Polynomen 3. Grades unter Berücksichtigung der Nebenbedingungen

    Weiter sollen Sie nun eine Prozedur für das Fitten mit linearen Nebenbedingungen mit der oben skizzierten, approximativen Methode erarbeiten. Die folgenden Teilfragen leiten Sie Schritt für Schritt durch den Prozess.

    • Matrix der Nebenbedingungen In einem ersten Schritte stellen Sie die Matrix \(\bf C\) auf, welche die Stetigkeit der 0. und 1. Ableitung mit der Gleichung \(\bf C \, \bf \xi = 0\) ausdrückt. Hinweis: Setzen Sie die Polynome mit den Koeffizienen aus \(\bf \xi\) an und formulieren Sie die entsprechenden Bedingungen an der Stelle \(x=5\).

    • Normalengleichung mit Nebenbedingungen Die Normalengleichung zum Finden des Minimums von \(\Vert {\bf A} \, {\bf \xi} - {\bf b} \Vert^2 + \lambda \cdot \Vert {\bf C} \, {\bf \xi} \Vert^2 = min\) kann analog zum 2. Schritt oben gefunden werden. Die zusätzlichen Terme, proportional zu \(\lambda\) im Funktional führen zu zusätzlichen Zeilen in der Matrix, welche mit Hilfe von \(\bf C\) gebildet werden können. Stellen Sie sich dazu vor, dass Sie zusätzliche Messpunkte vorliegen haben, deren Abstand zum vorgegebenen Wert (hier=0!) minimiert werden soll im Sinne der kleinsten Quadrate. Diese Überlegung führt Sie auch auf die richtige Spur zur Erweiterung der rechten Seite der modifzierten Normalengleichung. Mit der Lösung dieses modifizierten Systems kann der Fit nun unter Einhaltung der Nebenbedingungen ausgeführt werden. Der Parameter \(\lambda\) kann im Bereich \(0 < \lambda < 10^6\) gewählt werden. Experimentieren Sie mit dieser Grösse, indem Sie die Genauigkeit der Stetigkeitsbedingungen kontrollieren (Nebenfrage: Wie kann man diese Genauigkeit mit einer simplen Matrix-Vektor-Multiplikation testen, wenn die Lösung vorliegt?)

    • Numerische Aspekte Diese approximative Vorgehensweise ist algorithmisch simpel, da nur die Matrix, resp. die rechte Seite modifiziert und der Wert von \(\lambda\) korrekt gewählt werden muss. Aufgrund der obigen Überlegung ist ein möglichst grosser Wert von \(\lambda\) eine Garantie dafür, dass die Nebenbedingungen gut eingehalten werden. Aber leider hat die Sache auch einen Nachteil: grosse \(\lambda\)-Werte verschlechtern die Kondition der Systemmatrix. Untersuchen Sie diesen Zusammenhang, indem Sie \(\lambda\) systematisch variieren und die zugehörige Kondition der modifizierten Matrix \(\bf A\), resp. von \(\bf A^T\,A\) graphisch auftragen (logarithmisch) und Ihre Einsichten als Kommentare formulieren. Insbesondere zeigt diese Betrachtung, dass nicht die Normalengleichung gelöst werden sollte, sondern stattdessen z.B. die QR-Zerlegung genutzt wird.

Abgabe

Geben Sie Ihre Bearbeitung der Aufträge 2. - 4. bis in einer Woche ab.

Bemerkungen