1.4. Cholesky-Verfahren#
In vielen Anwendungsgebieten der numerischen Mathematik tauchen symmetrische Matrizen auf, die häufig sogar positiv definit sind, wie zum Beispiel Korrelationsmatrizen bei statistischen Berechnungen, die Normalengleichungen bei linearen Approximationsproblemen oder Steifigkeitsmatrizen bei der Methode der finiten Elemente. Für diese spezielle Klasse von Matrizen lassen sich beim Lösen von linearen Gleichungssystemen der Form \(Ax = b\) ihre besonderen Eigenschaften für effizientere Algorithmen ausnutzen.
Zur Herleitung des Cholesky-Verfahrens im Folgenden benötigen wir zuerst einige grundlegende Definitionen.
Definition 1.7 (Hermitesche Matrizen)
Sei \(A\) eine \(n\times n\) Matrix. Wir nennen \(A\) hermitesch, falls die Matrix ihrer Adjungierten entspricht, d.h., es gilt \(A = A^*\).
Wir bezeichnen hierbei mit \(A^*\) diejenige Matrix, die durch Spiegelung der Einträge von \(A\) an der Hauptdiagonalen mit anschließender komplexen Konjugation entsteht, d.h., für alle \(i \neq j\) gilt \(a_{i,j} \ = \ \overline{a_{j,i}}\).
Lemma 1.4 (Eigenschaften hermitescher Matrizen)
Für hermitesche Matrizen lassen sich folgende Eigenschaften zeigen:
Jede hermitesche Matrix \(A\) mit rein reellen Einträgen \(a_{i,j} \in \R, 1\leq i,j\leq n\) ist symmetrisch.
Jede hermitesche Matrix \(A\) besitzt \(n\) reelle Eigenwerte und die zugehörigen Eigenvektoren bilden ein Orthogonalsystem des \(\mathbb{R}^n\).
Jede hermitesche Matrix \(A\) ist selbstadjungiert und es lässt sich mit dem komplexen Standardskalarprodukt zeigen, dass für alle \(x,y \in \mathbb{C}^n\) gilt:
\[\langle Ax, y \rangle \: = \: (Ax)^*y \: = \: x^*A^*y \: = \: x^*Ay \: = \: \langle x, Ay \rangle.\]
Nach der Einführung von hermiteschen bzw. symmetrischen Matrizen und ihren Eigenschaften können wir Aussagen über ihre Definitheit machen.
Definition 1.8 (Definitheit von Matrizen)
Sei \(A\) eine hermitesche Matrix.
Falls gilt \(\langle x, Ax \rangle > 0\) bzw. \(\langle x, Ax \rangle \geq 0\) für alle \(x \in \mathbb{C}^n, x\neq 0\), so nennen wir die Matrix \(A\) positiv definit bzw. positiv semi-definit.
Analog lassen sich die Begriffe negativ definit bzw. negativ semi-definit für den Fall \(\langle x, Ax \rangle < 0\) bzw. \(\langle x, Ax \rangle \leq 0\) definieren.
Ist A weder positiv noch negativ semi-definit so nennen wir \(A\) indefinit.
Das folgende Lemma sagt aus, dass die positive Definitheit einer Matrix eine starke Aussage ist, die sich auch auf ihre Untermatrizen überträgt.
Lemma 1.5 (Positive Definitheit von Untermatrizen)
Sei \(A\) eine positiv definite Matrix und sei \(m\in \mathbb{N}\) beliebig mit \(1 \leq m < n\). Dann sind alle Untermatrizen \(A_{i_1, \ldots, i_m} \in \mathbb{C}^{(n-m) \times (n-m)}\), die durch Streichung der Spalten und Zeilen von \(A\) mit Indizes \(\lbrace i_1, \ldots, i_m\rbrace\) auch positiv definit.
Proof. Wir definieren zunächst \(\overline{A} := A_{i_1, \ldots, i_m} \in \mathbb{C}^{(n-m) \times (n-m)}\) als diejenige Untermatrix, die durch Streichung der Spalten und Zeilen von \(A\) mit Indizes \(\lbrace i_1, \ldots, i_m\rbrace\) entsteht. Sei nun \(\overline{x} \in \mathbb{C}^{n-m}\) ein beliebiger Vektor mit \(\overline{x} \neq 0\). Wir müssen zeigen, dass dann schon gilt \(\langle \overline{x}, \overline{A} \overline{x} \rangle > 0\).
Betrachten wir nun die Einbettung \(\Phi \colon \mathbb{C}^{n-m} \rightarrow \mathbb{C}^n\), die jeden Vektor \(\overline{x} \in \mathbb{C}^{n-m}\) auf den eingebetteten Vektor \(x \in \mathbb{C}^n\) abbildet, d.h., auf denjenigen Vektor \(x\) dessen Einträge mit denen von \(\overline{x}\) übereinstimmen und zusätzlich Nullen als Einträge für die Indizes \(i_1, \ldots, i_m\) besitzt. Aus der positiven Definitheit von \(A\) folgt dann aber schon:
Da der Vektor \(\overline{x} \in \mathbb{C}^{n-m} / \lbrace 0 \rbrace\) und die Untermatrix \(\overline{A}\) beliebig gewählt waren, folgt hieraus schon die Aussage. ◻
Aus Lemma 1.5 lässt sich direkt folgendes Korollar ableiten.
Korollar 1.1 (Diagonaleinträge von positiv definiten Matrizen)
Sei \(A\) eine positiv definite Matrix. Dann sind alle Diagonalelemente von \(A\) positiv, d.h., es gilt \(a_{i,i} > 0\) für \(i=1,\ldots,n\).
Bemerkung 1.8
Man kann zeigen, dass jede positiv definite Matrix die schöne Eigenschaft hat nur reelle, positive Eigenwerte zu besitzen. Insbesondere sind also alle positiv definiten Matrizen folglich invertierbar.
Sei nun \(A \in \R^{n \times n}\) eine positiv definite Matrix, für welche eine LR-Zerlegung der Form \(A = LR\) existiert, wobei \(L\) eine normierte, linke, untere Dreiecksmatrix ist. Wir wissen nun auf Grund der Eigenschaften von hermiteschen Matrizen, dass folgender Zusammenhang gelten muss:
Man sieht also, dass \(R^*L^*\) ebenfalls eine Zerlegung von \(A\) in ein Produkt zweier Dreiecksmatrizen ist. Sei nun \(D\) die Matrix, die sich aus der Hauptdiagonalen von \(R\) ergibt. Dann sieht man ein, dass
ebenfalls eine LR-Zerlegung von \(A\) darstellt, die sich von \(A=R^*L^*\) nur durch die Normierung der Diagonalelemente der linken Seite unterscheidet. Da \(A\) als positiv definite Matrix regulär ist, folgt aber aus der Eindeutigkeit der LR-Zerlegung in Satz 1.4 schon, dass \(L = R^*({D^*})^{-1}\) und \(R = D^*L^*\) ist und wir somit eine eindeutige Darstellung als unitäre Diagonalisierung von \(A\) erhalten mit
erhalten. Da die Matrix \(A\) positiv definit ist, müssen die Elemente von \(D\) nach Bemerkung 1.8 reell und positiv sein. Wir können also \(D^* = D\) setzen. Definieren wir nun die Matrix \(\sqrt{D}\) als diejenige Matrix, die sich durch Wurzelziehen aus den Diagonalelementen von \(D\) ergibt, so können wir \(\tilde{L} \coloneqq L\sqrt{D}\) definieren und es ergibt sich die sogenannte Cholesky-Zerlegung der Matrix \(A\) als
Man beachte, dass die linke, untere Dreiecksmatrix \(\tilde{L}\) durch die Multiplikation mit \(\sqrt{D}\) im Allgemeinen nicht mehr normiert ist. Das Ergebnis lässt sich in folgendem Satz zusammenfassen.
Satz 1.5 (Existenz und Eindeutigkeit der Cholesky-Zerlegung)
Sei \(A\) eine positiv definite Matrix. Dann gibt es genau eine linke, untere Dreiecksmatrix \(\tilde{L}\) mit positiven Diagonalelementen, so dass eine eindeutige Cholesky-Zerlegung von \(A\) der folgenden Form existiert:
Der Beweis von Satz 1.5 ist konstruktiv und liefert zugleich den Algorithmus zur Berechnung der Cholesky-Zerlegung. Ohne Beschränkung der Allgemeinheit wollen wir diesen für den Fall von reellen Matrizen \(A \in \mathbb{R}^{n \times n}\) im Folgenden angeben.
Algorithmus 1.5 (Cholesky-Verfahren)
Zur Berechnung der Cholesky-Zerlegung einer symmetrischen Matrix \(A = LL^T\) betrachten wir die Gleichung zuerst elementweise:
Hierbei bezeichnet \(l^T_{k,j}\) ein Element der Matrix \(L^T\). Auf Grund der Symmetrie von \(A\) besteht das nichtlineare Gleichungssystem (1.14) nur aus \(n(n+1)/2\) Gleichungen. Da \(L\) eine linke, untere Dreiecksmatrix ist, gibt es genauso viele Unbekannte \(l_{i,j}\) für \(1 \leq j \leq i \leq n\).
Dieses Gleichungssystem lässt sich wie folgt rekursiv auflösen. Wir versuchen zunächst die erste Spalte von \(L\) zu berechnen indem wir \(j=1\) setzen. Da alle Summanden bis auf den ersten Null sind, ergeben sich die folgenden \(n\) Gleichungen:
Insbesondere ergibt sich für den Index \(i=1\) die Lösung
Hierbei haben wir ausgenutzt, dass \(l_{i,i} > 0\) für \(i=1,\ldots,n\) und die Diagonalelemente von \(A\) immer positiv sind nach Korollar 1.1. Die weiteren Elemente der ersten Spalte von \(L\) kann man durch Einsetzen von (1.15) lösen als
Um die zweite Spalte von \(L\) zu berechnen setzen wir wiederum \(j=2\) und erhalten die Gleichungen (1.14)
Für das Diagonalelement \(l_{2,2}\) erhalten wir im Fall \(i=2\) die Lösung
Hierbei haben wir angenommen, dass der Radikand positiv ist. Die weiteren Elemente der zweiten Spalte von \(L\) für \(3 \leq i \leq n\) kann man nun durch Einsetzen berechnen:
Wir wollen dieses Prinzip nun für den j-ten Schritt des Cholesky-Verfahrens verallgemeinern. Nehmen wir also an, dass die Spalten \(1\) bis \(j-1\) der Matrix \(L\) bereits bestimmt wurden. Dann lauten die Gleichungen (1.14) für die \(j\)-te Spalte von \(L\):
Für \(i=j\) ergibt sich für das Diagonalelement \(l_{i,i}\) der Matrix \(L\)
Auch hier nehmen wir an, dass der Radikand positiv ist. Die weiteren Elemente für \(i=j+1,\ldots,n\) lassen sich dann durch Einsetzen berechnen als
Wir sehen also, dass wir mithilfe des Cholesky-Verfahrens alle unbekannten Elemente von \(L\) berechnen können.
Wir haben in (1.16) von Algorithmus 1.5 angenommen, dass der Radikand immer positiv ist. Folgendes Lemma zeigt, dass diese Annahme gerechtfertigt ist.
Lemma 1.6 (Positive Radikanden im Cholesky-Verfahren)
Sei \(A\) eine positiv definite Matrix und \(j \in \lbrace 1,\ldots,n\rbrace\) beliebig. Dann gilt für den \(j\)-ten Schritt des Cholesky-Verfahrens in Algorithmus 1.5, dass der Radikand in (1.16) positiv ist, d.h., es gilt
Proof. Wir beweisen das Lemma durch vollständige Induktion. Der Induktionsanfang für \(j=1\) ist richtig, da \(a_{1,1} > 0\) nach Korollar 1.1 gilt.
Wir nehmen nun an, dass die Behauptung für die ersten \(j-1\) Schritte gilt. Dann lassen sich die ersten \(j-1\) Spalten der linken, unteren Dreiecksmatrix \(L\) mit dem Cholesky-Verfahren berechnen. Wir bezeichnen diese Spalten der Matrix \(L\) mit der \(n \times (j-1)\) Matrix \(L_j\) und es gilt
Es reicht die linke, untere Hälfte des Produkts \(L_jL_j^T\) zu betrachten, da sich die rechte Hälfte aus der Symmetrie von \(A\) ergibt. Außerdem ist klar, dass die Matrix \(L_j L_j^T\) nur maximal den Rang \(j-1\) haben kann.
Man kann ausrechnen, dass das Element \(x_{j,j}\) die folgende Gestalt hat
Wenn wir nun zeigen können, dass \(x_{j,j} < a_{j,j}\) gilt, so haben wir den Induktionsschritt vollzogen, denn dann ist der Radikand in (1.16) positiv.
Nehmen wir nun an, dass das Gegenteil gilt, d.h., \(x_{j,j} \geq a_{j,j}\). Da jede Untermatrix einer positiv definiten Matrix auch positiv definit ist nach Lemma 1.5, wäre die Matrix \(C \in \mathbb{R}^{j \times j}\) mit
auf Grund der Annahme \(x_{j,j} \geq a_{j,j}\) positiv definit. Dies sieht man durch folgendes Argument ein. Sei \(\delta \coloneqq x_{j,j} - a_{j,j}\), dann gilt \(\delta \geq 0\) nach Voraussetzung. Definieren wir nun eine Matrix \(B \in \mathbb{R}^{j\times j}\), die überall Nullen enthält außer dem Eintrag \(b_{j,j} = \delta \geq 0\) auf der Hauptdiagonalen. Es ist offensichtlich, dass die Matrix \(B\) positiv semi-definit ist und es gilt \(C = \overline{A} + B\), wobei \(\overline{A}\) die Untermatrix von \(A\) darstellt, die man erhält, wenn man sich auf die ersten \(j\) Zeilen und Spalten von \(A\) beschränkt. Dann gilt für alle Vektoren \(x \in \mathbb{R}^j / \lbrace \vec{0} \rbrace\):
Die Matrix \(C\) ist also auch positiv definit und somit hätte die Matrix auf der rechten Seite von (1.17) mindestens den Rang \(j\), während das Produkt \(L_jL_j^T\) aber maximal den Rang \(j-1\) haben kann. Dieser Widerspruch zeigt, dass \(x_{j,j} < a_{j,j}\) gelten muss. ◻
Abschließend wollen wir noch den Rechenaufwand zur Berechnung dieser besonderen Zerlegung angeben, der sich als effizienter als das Gaußsche Eliminationsverfahren herausstellt.
Bemerkung 1.9 (Rechenaufwand des Cholesky-Verfahrens)
Der numerische Rechenaufwand zur Berechnung der Cholesky-Zerlegung in Algorithmus 1.5 liegt in \(\frac{1}{6}n^3 + \mathcal{O}(n^2)\) und verhält sich dadurch asymptotisch doppelt so schnell wie das Gauss-Eliminationsverfahren in Algorithmus 1.2. Dies liegt vor allem daran, dass man beim Cholesky-Verfahren die Symmetrie der Matrix \(A\) geschickt ausnutzt.