Kurs:Numerik I/Lösung mit QS-Zerlegung

Diese Lernressoure zum Thema QS-Zerlegung kann als Wiki2Reveal Folien angezeigt werden. Einzelne Abschnitte werden als Folien betrachtet und Folien können in der Wiki2Reveal-Darstellung handschriftlich im Browser beschriftet werden.

Für eine QS-Zerlegung stellen wir zunächst fest, dass sich eine Matrix mit maximalen Spaltenrang in ein Produkt von

• einer orthogonalen quadratischen Matrix ${\displaystyle Q\in \mathbb {R} ^{n\times n}}$ und
• erweiterten oberen Dreiecksmatrix ${\displaystyle S\in \mathbb {R} ^{n\times k}}$

Sei ${\displaystyle A\in \mathbb {R} ^{n\times k}}$ mit ${\displaystyle n\geq k}$ und ${\displaystyle \operatorname {Rang} (A)=k}$ gegeben und ${\displaystyle A=QS}$ eine Zerlegung von ${\displaystyle A}$, wobei ${\displaystyle Q\in \mathbb {R} ^{n\times n}}$ eine orthogonale und ${\displaystyle S\in \mathbb {R} ^{n\times k}}$ eine Matrix der folgenden Gestalt ist mit ${\displaystyle \mathbf {0} :=(0)\in \mathbb {R} ^{(n-k)\times k}}$:

${\displaystyle S={\begin{pmatrix}R\\\mathbf {0} \end{pmatrix}}\in \mathbb {R} ^{n\times k},R:={\begin{pmatrix}*&\ldots &*\\\vdots &\ddots &\vdots \\0&\ldots &*\end{pmatrix}}\in \mathbb {R} ^{k\times k},}$

Dann gilt für die obere Dreiecksmatrix ${\displaystyle R}$:

(i) ${\displaystyle \det(R)\neq 0}$,

(ii) ${\displaystyle \operatorname {cond} _{2}(A^{T}A)={\big (}\operatorname {cond} _{2}(R){\big )}^{2}}$.

Der Beweis gliedert sich in folgende Schritte:

Da ${\displaystyle Q}$ eine orthogonale Matrix und ${\displaystyle A=QS}$ ist, gilt mit ${\displaystyle (QS)^{T}=S^{T}Q^{T}}$:

${\displaystyle A^{T}A=\underbrace {S^{T}Q^{T}} _{=A^{T}}\,\underbrace {QS} _{=A}=S^{T}S={\begin{pmatrix}R^{T}&\mathbf {0} \end{pmatrix}}{\begin{pmatrix}R\\\mathbf {0} \end{pmatrix}}=R^{T}R.}$

Dabei wurde u.a. verwendet, dass für orthogonale Matrizen ${\displaystyle Q}$ und ${\displaystyle s\in \mathbb {R} ^{n}}$ gilt:

${\displaystyle (Qs)^{T}Qs=\langle Qs,Qs\rangle =\langle s,s\rangle =s^{T}s}$

Wegen ${\displaystyle \operatorname {Rang} (A)=k}$ ist nach Lemma 4.1 insbesondere

${\displaystyle 0\neq \det(A^{T}A)=\det(R^{T}R)={\big (}\det(R){\big )}^{2}}$.

Damit gilt ${\displaystyle det(A^{T}A)>0}$ und ${\displaystyle det(R)\not =0}$ und (i) ist korrekt.

Ferner erhält man durch die Gleichheit von ${\displaystyle A^{T}A=R^{T}R}$ aus Beweisschritt 1 die Gleichheit der Kondition

${\displaystyle \operatorname {cond} _{2}(A^{T}A)=\operatorname {cond} _{2}(R^{T}R).}$

Damit kann man die Konditionszahl von ${\displaystyle A^{T}A}$ über die Konditionszahl der invertierbaren oberen Dreiecksmatrix ${\displaystyle R}$ berechnen.

Wendet man nun Lemma für die Konditionszahl auf die reguläre obere Dreiecksmatrix ${\displaystyle R}$ an, erhält man:

${\displaystyle \operatorname {cond} _{2}(R^{T}R)=(\operatorname {cond} _{2}(R))^{2}.}$

Nach dem Lemma über Eigenwerte positiv definiter Matrizen und mit der Invertierbarkeit von ${\displaystyle R}$ sind die Eigenwerte der Matrix ${\displaystyle R^{T}R}$ positiv. Insgesamt erhält man mit Schritt 4 die Eigenschaft (ii) aus dem Lemma.

${\displaystyle \operatorname {cond} _{2}(A^{T}A)=\operatorname {cond} _{2}(R^{T}R)={\big (}\operatorname {cond} _{2}(R){\big )}^{2}}$

q.e.d.

In dem obige Satz für die QS-Zerlegung wurde der maximale Spaltenrang ${\displaystyle k}$ der Matrix ${\displaystyle A\in \mathbb {R} ^{n\times k}}$ mit ${\displaystyle k\leq n}$ vorausgesetzt. Der folgende Satz gibt an, wie Fehlerquadratprobleme mittels einer QS-Zerlegung von ${\displaystyle A}$ gelöst werden können.

Voraussetzungen Es seien ${\displaystyle A\in \mathbb {R} ^{n\times k}}$ mit ${\displaystyle n\geq k}$, ${\displaystyle \operatorname {Rang} (A)=k}$, ${\displaystyle R\in \mathbb {R} ^{k\times k}}$ eine obere Dreiecksmatrix und ${\displaystyle Q\in \mathbb {R} ^{n\times n}}$eine orthogonale Matrix und ${\displaystyle S\in \mathbb {R} ^{n\times k}}$ gegeben wie in Lemma - QS-Zerlegung. Weiter sei der Vektor ${\displaystyle Q^{T}b\in \mathbb {R} ^{n}}$ dargestellt in der Form

${\displaystyle Q^{T}b={\begin{pmatrix}y^{(1)}\\y^{(2)}\end{pmatrix}}\in \mathbb {R} ^{n},y^{(1)}\in \mathbb {R} ^{k},y^{(2)}\in \mathbb {R} ^{n-k}.}$

Dann ist der Vektor ${\displaystyle x^{*}\in \mathbb {R} ^{k}}$ genau dann die eindeutige Lösung des linearen Ausgleichsproblems ${\displaystyle \min _{x\in \mathbb {R} ^{k}}\|Ax-b\|_{2},}$, wenn ${\displaystyle x^{*}}$ eindeutige die Lösung des Systems ${\displaystyle Rx=y^{(1)}}$ ist. Außerdem gilt für den Fehler ${\displaystyle \|Ax^{*}-b\|_{2}=\left\|y^{(2)}\right\|_{2}.}$

Ziel des Lösungsverfahrens für ein gesuchtest ${\displaystyle x}$ ist es, sich möglichst gut bzgl. der euklidschen Norm mit ${\displaystyle Ax}$ an den Vektor ${\displaystyle b}$ anzunähern. Die eindeutige Lösung des linearen Ausgleichsproblems liefert mit ${\displaystyle y^{(1)}}$ und ${\displaystyle y^{(2)}}$:

• ${\displaystyle Rx=y^{(1)}}$ ein Gleichungssystem mit einer einfach zu lösenden oberen Dreiecksmatrix ${\displaystyle R}$ und
• mit ${\displaystyle y^{(2)}}$ bzw. ${\displaystyle \|y^{(2)}\|_{2}}$ ein Maß für die minimale Abweichung von ${\displaystyle Ax}$ und ${\displaystyle b}$ über ${\displaystyle \|Ax^{*}-b\|_{2}=\|y^{(2)}\|_{2}}$ bzgl. der Lösung ${\displaystyle x^{*}}$.

Für einen beliebigen Vektor ${\displaystyle v\in \mathbb {R} ^{k}}$ erhalten wir für eine orthogonale Matrix ${\displaystyle Q}$ eine längetreue Abbildung über:

${\displaystyle \|Qv\|_{2}^{2}=\langle Qv,Qv\rangle =\langle v,v\rangle =\|v\|_{2}^{2}\quad (\ast )}$

Ferner gilt für orthogonale Matrizen ${\displaystyle Q^{T}Q=I}$, wobei ${\displaystyle I}$ die Einheitsmatrix als neutrales Element der Matrixmultiplikation ist.

unter Verwendung des QS-Zerlegungssatzes und der Ersetzung von ${\displaystyle A=QS}$ erhält man:

${\displaystyle {\begin{array}{rcl}\|b-Ax^{*}\|_{2}^{2}&=&\left\|(QQ^{T})b-QSx\right\|_{2}^{2}={\bigg \|}Q{\big (}\underbrace {Q^{T}b-Sx} _{=v}{\big )}{\bigg \|}_{2}^{2}\\&=&{\big \|}\underbrace {Q^{T}b-Sx} _{=v}{\big \|}_{2}^{2}{\mbox{ mit QS-Zerlegungssatz}}\\&=&\left\|{\begin{pmatrix}y^{1}\\y^{2}\end{pmatrix}}-{\begin{pmatrix}R\\\mathbf {0} \end{pmatrix}}x\right\|_{2}^{2}=\left\|y^{(1)}-Rx\right\|_{2}^{2}+\left\|y^{(2)}\right\|_{2}^{2}\\\end{array}}}$

Daraus folgt für einen beliebigen Vektor ${\displaystyle {\widehat {x}}\in \mathbb {R} ^{k}}$

${\displaystyle \|A{\widehat {x}}-b\|_{2}\geq \min _{x\in \mathbb {R} ^{k}}\|Ax-b\|_{2}=\left\|y^{(2)}\right\|_{2}}$

Ferner gilt

${\displaystyle \|b-Ax\|_{2}=\left\|y^{(2)}\right\|_{2}\Leftrightarrow \left\|y^{(1)}-Rx\right\|_{2}=0\Leftrightarrow Rx=y^{(1)}.}$

Damit ist alles gezeigt.

q.e.d.

Nach dem Satz und dem Lemma zur QS-Zerlegung besitzen das ${\displaystyle l_{2}}$-Approximationsproblem, wobei man versucht ${\displaystyle Ax}$ in der ${\displaystyle l_{2}}$-Norm möglichst gut anzunähern, eine eindeutige Lösung. Dies ist äquivalent zu der eindeutigen Lösbarkeit des Systems ${\displaystyle Rx=y^{1}}$. Abschließend geben wir zu dem letzten Satz noch ein Beispiel.

Wir betrachten das Fehlerquadratproblem mit den Daten

${\displaystyle A:={\begin{pmatrix}1&0\\1&3\\1&4\\1&7\end{pmatrix}},\quad b:={\begin{pmatrix}1\\2\\6\\4\end{pmatrix}}.}$

Wir suchen nun eine Vektor ${\displaystyle x^{*}\in \mathbb {R} ^{2}}$, der den Abstand ${\displaystyle \|Ax-b\|_{2}}$ zwischen ${\displaystyle Ax}$ und ${\displaystyle b}$ minimiert.

Der Spaltenrang von ${\displaystyle A}$ ist 2, da die beiden Spalten von ${\displaystyle A}$ nicht linear abhängig sind. Nun liefert die QS-Zerlegung von ${\displaystyle A}$ mit gleichzeitiger Berechnung von ${\displaystyle Q^{T}b}$ die folgende obere Dreiecksmatrix ${\displaystyle R}$ und die folgenden Vektoren ${\displaystyle y^{(1)}}$ und ${\displaystyle y^{(2)}}$:

${\displaystyle R:={\begin{pmatrix}-2&-7\\0&-5\end{pmatrix}},\quad y^{(1)}:={\begin{pmatrix}-{\frac {13}{2}}\\-{\frac {5}{2}}\end{pmatrix}},\quad y^{(2)}:={\begin{pmatrix}{\frac {99}{34}}\\-{\frac {5}{34}}\end{pmatrix}}.}$

Die Lösung von ${\displaystyle Rx=y^{(1)}}$ bzw.

${\displaystyle {\begin{pmatrix}-2&-7\\0&-5\end{pmatrix}}{\begin{pmatrix}x_{1}\\x_{2}\end{pmatrix}}={\begin{pmatrix}-{\frac {13}{2}}\\-{\frac {5}{2}}\end{pmatrix}}}$

lautet

${\displaystyle x^{*}=(x_{1}^{*},x_{2}^{*})^{T}:=\left({\frac {3}{2}},{\frac {1}{2}}\right)^{T}.}$

Der zugehörige Approximationsfehler in der Euklidischen Norm ist gegeben durch

${\displaystyle \|b-Ax^{*}\|_{2}=\left\|y^{2}\right\|_{2}=\left\|{\begin{pmatrix}{\frac {99}{34}}\\-{\frac {5}{34}}\end{pmatrix}}\right\|_{2}={\frac {\sqrt {9826}}{34}}={\frac {17{\sqrt {34}}}{34}}={\frac {1}{2}}{\sqrt {34}}\approx 2.915.}$

Wir wollen nun die beiden beschriebenen Lösungswege für das ${\displaystyle l_{2}}$-Problem, die Lösung der Normalgleichungen mittels Cholesky-Zerlegung und die Lösung über den in QS-Zerlegungssatz dargestellten Weg, vergleichen.

In beiden Fällen muss die rechte Seite ${\displaystyle b}$ des Systems ${\displaystyle Ax=b}$ von links mit einer Matrix multipliziert werden. Weiter sind bei der Cholesky-Zerlegung zwei und bei dem Weg über die QS-Zerlegung ein gestaffeltes lineares Gleichungssystem zu lösen. Da die Lösung eines solchen Systems nur etwa ${\displaystyle {\frac {k^{2}}{2}}}$ Multiplikationen und Divisionen erfordert, vernachlässigen wir hier diesen zusätzlichen Aufwand bei der Cholesky-Zerlegung.

Der Lösungsvektor ${\displaystyle x\in \mathbb {R} ^{k}}$ hat in der Regel eine feste Dimension ${\displaystyle k\leq n}$ während ${\displaystyle n}$ die Anzahl der Gleichungen für Daten angibt, die bzgl. ${\displaystyle \|Ax-b\|_{2}}$ in der euklidischen Norm minimiert werden soll.

Im Fall der Lösung der Normalgleichungen benötigt man ferner zur Berechnung der symmetrischen Matrix ${\displaystyle A^{T}A}$ etwa ${\displaystyle {\frac {1}{2}}nk^{2}}$ und für die Cholesky-Zerlegung etwa ${\displaystyle {\frac {1}{6}}k^{3}}$ wesentliche Rechenoperationen. Daneben erfordert die Berechnung des Residuenvektors ${\displaystyle Ax-b}$ weitere ${\displaystyle n\cdot k}$ Multiplikationen, so dass die zuletzt genannten Teilaufgaben zusammen ungefähr

${\displaystyle {\frac {1}{2}}n\cdot k^{2}+{\frac {1}{6}}\cdot k^{3}+n\cdot k}$

wesentliche Rechenoperationen erfordern, das sind

${\displaystyle \sim {\frac {1}{2}}nk^{2}}$ für ${\displaystyle n\gg k}$ und ${\displaystyle \sim {\frac {2}{3}}n^{3}}$ für ${\displaystyle n\approx k}$.

Für das sich aus dem QS-Zerlegungssatz ergebende Vorgehen benötigt man über die Matrix-Vektor-Multiplikation und die Lösung eines Dreieckssystems hinaus nur die Berechnung der QS-Zerlegung von ${\displaystyle A}$.

Wenn ${\displaystyle n}$ nahe bei ${\displaystyle k}$ liegt (${\displaystyle n\approx k}$) benötigen im Vergleich der beiden Verfahren demnach für beide Wege etwa gleich viele wesentliche Rechenoperationen.

Ist ${\displaystyle n}$ sehr groß im Vergleich zu ${\displaystyle k}$ müssen über den Weg der QS-Zerlegung etwa doppelt so viele wesentliche Rechenoperationen durchgeführt werden, wie der Lösung über die Normalgleichungen (Cholesky).

Letzterem Nachteil der QS-Zerlegung steht jedoch entgegen, dass bei ihr nach Lemma 4.6 die Konditionszahl ${\displaystyle \kappa :=\left[\operatorname {cond} _{2}\left(A^{T}A\right)\right]^{1/2}}$ bzw. für quadratisches ${\displaystyle A}$ (vgl. Satz 4.5) die Zahl ${\displaystyle \kappa :=\operatorname {cond} _{2}(A)}$ den Rundungsfehlereinfluss bei der Lösung des Dreieckssystems bestimmt, während dies bei dem Weg über die Normalgleichungen die Zahl ${\displaystyle \kappa ^{2}}$ ist. Wir geben dazu ein (allerdings etwas extremes) Beispiel:

Es sei ${\displaystyle n:=6,k:=5}$ und ${\displaystyle A\in \mathbb {R} ^{6\times 5}}$ mit einem ${\displaystyle \varepsilon >0}$ gegeben durch

${\displaystyle A:={\begin{pmatrix}1&1&1&1&1\\\varepsilon &&&&\\&\varepsilon &&&\\&&\varepsilon &&\\&&&\varepsilon &\\&&&&\varepsilon \end{pmatrix}}.}$

Damit ergibt sich für die Berechnung von ${\displaystyle A^{T}A}$:

${\displaystyle A^{T}A={\begin{pmatrix}1+\varepsilon ^{2}&1&1&1&1\\1&1+\varepsilon ^{2}&1&1&1\\1&1&1+\varepsilon ^{2}&1&1\\1&1&1&1+\varepsilon ^{2}&1\\1&1&1&1&1+\varepsilon ^{2}\end{pmatrix}}.}$

Es seien nun ${\displaystyle b:=10}$ und ${\displaystyle r:=10}$ die Basis und Mantissenlänge des verwendeten Computers, so dass ${\displaystyle eps=5\cdot 10^{-10}}$ die zugehörige Maschinengenauigkeit ist. Weiter sei ${\displaystyle \varepsilon :=0.5_{10}-5}$ und damit ${\displaystyle \varepsilon ^{2}=.25_{10}-10}$. Damit liegt ${\displaystyle \varepsilon ^{2}}$ unterhalb der zugehörigen Maschinengenauigkeit und wird als 0 gespeichert.

Dann erhält man Gleitkomma-Darstellung für ${\displaystyle 1+\varepsilon ^{2}}$

${\displaystyle gl(1+\varepsilon ^{2})=1,\quad gl(A^{T}A)=(a_{ij})}$ mit ${\displaystyle a_{ij}:=1}$ ${\displaystyle (i,j=1,\ldots ,5)}$.

Die Matrix ${\displaystyle gl(A^{T}A)}$ hat offenbar den Rang 1 und ist singulär.

Man errechnet hier für die Konditionszahl:

${\displaystyle \operatorname {cond} _{2}(A^{T}A)={\frac {5+\varepsilon ^{2}}{\varepsilon ^{2}}}\approx 2_{10}+11,}$

${\displaystyle \left[\operatorname {cond} _{2}(A^{T}A)\right]^{1/2}=\left[{\frac {5+\varepsilon ^{2}}{\varepsilon ^{2}}}\right]^{1/2}\approx 4.5_{10}+5.}$

• Dreiecksmatrix
• orthogonale Matrix
• Konditionszahl
• Lemma - Eigenwerte positiv definiter Matrizen
• Lemma - Konditionszahl Normalengleichung

Diese Lernresource können Sie als Wiki2Reveal-Foliensatz darstellen.

Dieser Wiki2Reveal Foliensatz wurde für den Lerneinheit Kurs:Numerik I' erstellt der Link für die Wiki2Reveal-Folien wurde mit dem Wiki2Reveal-Linkgenerator erstellt.

• Die Seite wurde als Dokumententyp PanDocElectron-SLIDE erstellt.
• Link zur Quelle in Wikiversity: https://de.wikiversity.org/wiki/Kurs:Numerik%20I/L%C3%B6sung%20mit%20QS-Zerlegung
• siehe auch weitere Informationen zu Wiki2Reveal und unter Wiki2Reveal-Linkgenerator.