Pseudoinverse

Die Pseudoinverse einer Matrix ist ein Begriff aus dem mathematischen Teilgebiet der linearen Algebra, der auch in der numerischen Mathematik eine wichtige Rolle spielt. Sie ist eine Verallgemeinerung der inversen Matrix auf singuläre und nichtquadratische Matrizen, weshalb sie häufig auch als verallgemeinerte Inverse bezeichnet wird. Der häufigste Anwendungsfall für Pseudoinversen ist die Lösung linearer Gleichungssysteme und linearer Ausgleichsprobleme.

Eine erste Form wurde von E. H. Moore (1920)[1] und Roger Penrose (1955)[2] beschrieben. Die nach ihnen benannte Moore-Penrose-Inverse ist nicht die einzige Möglichkeit, eine Pseudoinverse zu definieren, häufig wird aber Pseudoinverse synonym mit Moore-Penrose-Inverse benutzt.[3] Die Moore-Penrose-Inverse ist für alle Matrizen mit Einträgen aus den reellen oder komplexen Zahlen definiert und eindeutig. Mit ihr kann man bei linearen Ausgleichsproblemen die optimale Lösung hinsichtlich der kleinsten Summe quadrierter Abweichungen der euklidischen Normen berechnen.

Eine numerisch robuste Methode zur Bestimmung der Moore-Penrose-Inversen baut auf der Singulärwertzerlegung auf.

Allgemeine Pseudoinversen

Die Verallgemeinerung der Bildung der Inversen einer Matrix auf singuläre Matrizen wird in der Literatur nicht einheitlich gehandhabt und orientiert sich oftmals an der zu lösenden Aufgabenstellung (einige Beispiele solcher Verallgemeinerungen sind weiter unten aufgeführt).

Nach Adi Ben-Israel[4] sollte eine Definition von verallgemeinerten Inversen zumindest die folgenden drei Forderungen erfüllen:

  1. Für reguläre Matrizen sollte sich eindeutig die gewöhnliche Inverse ergeben.
  2. Im verallgemeinerten Sinne sollten auch singuläre Matrizen invertierbar sein (wenigstens einige, nicht notwendigerweise alle).
  3. Für singuläre Matrizen sollten die verallgemeinerten Inversen ähnliche Eigenschaften haben wie gewöhnliche Inverse regulärer Matrizen.

Als Ausgangspunkt für die Konstruktion von verschiedenen Pseudoinversen schwächt Adi Ben-Israel[4] dann die vier definierenden Aussagen für die im nächsten Abschnitt beschriebene Moore-Penrose-Inverse in verschiedene Richtungen ab und ergänzt sie durch andere Bedingungen. Die Mindestforderung an eine Pseudoinverse ist die folgende: Eine Matrix B {\displaystyle B} ist genau dann Pseudoinverse von A {\displaystyle A} , wenn gilt:

A B A = A {\displaystyle ABA=A}

Dagegen bezeichnet Max Koecher[5] eine Matrix B {\displaystyle B} genau dann als Pseudoinverse von A {\displaystyle A} , wenn für sie die folgenden beiden Aussagen

A B A = A {\displaystyle ABA=A} und
B A B = B {\displaystyle BAB=B}

zutreffen.

Die erste Bedingung sichert dabei, dass die Spalten y {\displaystyle y} von A {\displaystyle A} durch B {\displaystyle B} auf Lösungen x {\displaystyle x} des Gleichungssystems y = A x {\displaystyle y=Ax} abgebildet werden. Durch die zweite Aussage können keine vom Nullvektor verschiedene Spalten von B {\displaystyle B} im Kern von A {\displaystyle A} liegen.

Die Moore-Penrose-Inverse

Die Moore-Penrose-Inverse (auch einfach Pseudoinverse) einer Matrix A C m × n {\displaystyle A\in \mathbb {C} ^{m\times n}} ist die eindeutig bestimmte Matrix A + C n × m {\displaystyle A^{+}\in \mathbb {C} ^{n\times m}} , welche die folgenden vier Eigenschaften („Moore-Penrose-Bedingungen“) erfüllt:

  • A A + A = A {\displaystyle AA^{+}A\;\,\,=A}
( A + {\displaystyle A^{+}} ist eine verallgemeinerte Inverse.)
  • A + A A + = A + {\displaystyle A^{+}AA^{+}=A^{+}}
( A + {\displaystyle A^{+}} verhält sich wie eine schwache Inverse.)
  • ( A A + ) = A A + {\displaystyle (AA^{+})^{*}\,=AA^{+}}
(Die Matrix A A + {\displaystyle AA^{+}} ist hermitesch.)
  • ( A + A ) = A + A {\displaystyle (A^{+}A)^{*}\,=A^{+}A}
(Die Matrix A + A {\displaystyle A^{+}A} ist ebenfalls hermitesch.)

Dabei bezeichnet M {\displaystyle M^{*}} die adjungierte Matrix zu einer Matrix M {\displaystyle M} . Bei Matrizen mit Einträgen aus den reellen Zahlen ist diese identisch mit der zu M {\displaystyle M} transponierten Matrix M T {\displaystyle M^{T}} .

Die Moore-Penrose-Inverse kann auch durch einen Grenzwert definiert werden:

A + = lim δ 0 ( A A + δ E n ) 1 A = lim δ 0 A ( A A + δ E m ) 1 {\displaystyle {\begin{aligned}A^{+}&=\lim _{\delta \searrow 0}\,\;(A^{*}A+\delta E_{n})^{-1}A^{*}\\&=\lim _{\delta \searrow 0}\;A^{*}(AA^{*}+\delta E_{m})^{-1}\end{aligned}}}

mit E k {\displaystyle E_{k}} als der Einheitsmatrix in C k × k {\displaystyle \mathbb {C} ^{k\times k}} . Dieser Grenzwert existiert auch dann, wenn ( A A ) 1 {\displaystyle (AA^{*})^{-1}} und ( A A ) 1 {\displaystyle (A^{*}A)^{-1}} nicht existieren.

Rechenregeln

( A + ) + = A {\displaystyle (A^{+})^{+}=A}
A ¯ + = A + ¯ {\displaystyle {\overline {A}}^{+}={\overline {A^{+}}}}
( A ) + = ( A + ) {\displaystyle (A^{*})^{+}=(A^{+})^{*}}
( λ A ) + = λ 1 A + {\displaystyle (\lambda A)^{+}=\lambda ^{-1}A^{+}} für λ 0 {\displaystyle \lambda \neq 0}

Spezialfälle

Sind die Spalten der Matrix A {\displaystyle A} linear unabhängig, dann ist A A {\displaystyle A^{*}A} invertierbar. In diesem Fall gilt die folgende Gleichung[4]

A + = ( A A ) 1 A . {\displaystyle A^{+}=(A^{*}A)^{-1}A^{*}.}

Nimmt man die erste Grenzwertdefinition für die Moore-Penrose-Inverse, so verschwindet der Summand δ E n {\displaystyle \delta E_{n}} . Daraus folgt, dass A + {\displaystyle A^{+}} eine Linksinverse zu A {\displaystyle A} ist.

A + A = E n {\displaystyle A^{+}A=E_{n}}

Sind die Zeilen der Matrix A {\displaystyle A} linear unabhängig, dann ist A A {\displaystyle AA^{*}} invertierbar. In diesem Fall gilt die folgende Gleichung

A + = A ( A A ) 1 . {\displaystyle A^{+}=A^{*}(AA^{*})^{-1}.}

Nimmt man die zweite Grenzwertdefinition für die Moore-Penrose-Inverse, so verschwindet der Summand δ E m {\displaystyle \delta E_{m}} . Daraus folgt, dass A + {\displaystyle A^{+}} eine Rechtsinverse zu A {\displaystyle A} ist.

A A + = E m {\displaystyle AA^{+}=E_{m}}

Sind sowohl Spalten als auch die Zeilen einer Matrix unabhängig, dann ist die Matrix invertierbar, und die Pseudoinverse stimmt mit der Inversen überein.

Ist das Produkt A B {\displaystyle AB} zweier Matrizen definiert und eine der beiden eine unitäre Matrix, dann gilt

( A B ) + = B + A + . {\displaystyle (AB)^{+}=B^{+}A^{+}.}

Man kann die Pseudoinverse auch für Skalare und Vektoren definieren, indem man diese als Matrizen betrachtet. Bei Skalaren ist die Pseudoinverse von null wieder null, und für alle anderen Werte x {\displaystyle x} ist sie x 1 {\displaystyle x^{-1}} . Für Vektoren gilt

x + = { 0 T , wenn  x = 0 , x x x , sonst . {\displaystyle x^{+}\!={\begin{cases}0^{T},&{\text{wenn }}x=0,\\{\frac {x^{*}}{x^{*}x}},&{\text{sonst}}.\end{cases}}}

Diese Behauptungen lassen sich überprüfen, indem man die Kriterien für die Moore-Penrose-Inverse nachprüft.

Ist die Matrix A {\displaystyle A} hermitesch (oder symmetrisch im reellen Fall), dann ist A + {\displaystyle A^{+}} ebenfalls hermitesch (symmetrisch). Aus dem Spektralsatz folgt in diesem Fall die Zerlegung

U A U = D {\displaystyle UAU^{*}=D}

und damit

A + = U D + U {\displaystyle A^{+}=U^{*}D^{+}U} ,

wobei die Pseudoinverse D + {\displaystyle D^{+}} der Diagonalmatrix D {\displaystyle D} durch

d i i + = { 0 , wenn  d i i = 0 , 1 d i i , sonst {\displaystyle d_{ii}^{+}={\begin{cases}0,&{\text{wenn }}d_{ii}=0,\\{\frac {1}{d_{ii}}},&{\text{sonst}}\end{cases}}}

für alle Diagonaleinträge gegeben ist.

Berechnung

  • Ist k {\displaystyle k} der Rang der m × n {\displaystyle m\times n} -Matrix A {\displaystyle A} , dann kann A {\displaystyle A} in das Produkt A = B C {\displaystyle A=BC} einer m × k {\displaystyle m\times k} -Matrix B {\displaystyle B} und einer k × n {\displaystyle k\times n} -Matrix C {\displaystyle C} zerlegt werden. Es gilt
A + = C ( C C ) 1 ( B B ) 1 B . {\displaystyle A^{+}=C^{*}(CC^{*})^{-1}(B^{*}B)^{-1}B^{*}.}
Hat A {\displaystyle A} vollen Zeilenrang, das heißt, es gilt k = m {\displaystyle k=m} , dann kann für B {\displaystyle B} die Einheitsmatrix gewählt werden und obige Formel reduziert sich zu
A + = A ( A A ) 1 . {\displaystyle A^{+}=A^{*}(AA^{*})^{-1}.}
In ähnlicher Weise gilt für eine Matrix A {\displaystyle A} mit vollem Spaltenrang, das heißt, es gilt k = n {\displaystyle k=n} , die Gleichung
A + = ( A A ) 1 A . {\displaystyle A^{+}=(A^{*}A)^{-1}A^{*}.}
  • Mit der Singulärwertzerlegung existiert ein anderes Verfahren zur Berechnung der Pseudoinversen. Ist A = U Σ V {\displaystyle A=U\Sigma V^{*}} die Singulärwertzerlegung von A {\displaystyle A} , dann gilt
A + = V Σ + U . {\displaystyle A^{+}=V\Sigma ^{+}U^{*}.}
Bei einer Diagonalmatrix wie Σ C m × n {\displaystyle \Sigma \in \mathbb {C} ^{m\times n}} entsteht die Pseudoinverse, indem man die Matrix transponiert und die von null verschiedenen Elemente invertiert, also Σ + C n × m {\displaystyle \Sigma ^{+}\in \mathbb {C} ^{n\times m}} bildet mit
( Σ ) i j + = { 1 σ i , falls  i = j σ i 0 , 0 , sonst . {\displaystyle (\Sigma )_{ij}^{+}={\begin{cases}{\frac {1}{\sigma _{i}}},&{\text{falls }}i=j\wedge \sigma _{i}\neq 0,\\0,&{\text{sonst}}.\end{cases}}}
  • Mit Hilfe der Ränderung von Matrizen kann die Pseudoinverse implizit dargestellt oder auch berechnet werden.[6]
  • Der Algorithmus von Greville ist eine endliche iterative Methode zur spaltenweisen Berechnung der Moore-Penrose-Inversen.[7]

Das Verfahren, bei dem man die Matrix A A {\displaystyle A^{*}A} benötigt, wird zwar bei der numerischen Berechnung der Lösung überbestimmter Gleichungssysteme der Bequemlichkeit halber öfter benutzt, ist jedoch numerisch instabil, da die Kondition der Matrix quadriert wird. Als stabile und effiziente numerische Methode gilt die Verwendung der QR-Zerlegung. Das auf der Singulärwertzerlegung aufbauende Verfahren ist das aufwendigste, aber auch das numerisch gutartigste. Das auf der Ränderung beruhende Verfahren bietet einen Kompromiss zwischen Aufwand und numerischer Stabilität.

Einen Überblick über numerischen Aufwand und Stabilität der Verfahren gibt auch [6].

Anwendungen

Ist das Gleichungssystem A x = b {\displaystyle Ax=b} nicht lösbar, so lässt sich mit der Pseudoinversen die Lösung nach der Methode der kleinsten Quadrate, also die mit kleinster euklidischer Norm A x b 2 {\displaystyle \|Ax-b\|_{2}} als x ¯ = A + b {\displaystyle {\bar {x}}=A^{+}b} berechnen.

Gibt es für das Gleichungssystem A x = b {\displaystyle Ax=b} unendlich viele Lösungen, so kann man diese über

x = A + b + ( E n A + A ) y , y C n {\displaystyle x=A^{+}b+(E_{n}-A^{+}A)y,\quad y\in \mathbb {C} ^{n}}

bestimmen. Dabei ist x {\displaystyle x} diejenige Lösung des Gleichungssystems, die von y {\displaystyle y} den kleinsten Abstand bezüglich der euklidischen Norm hat.

Ausgewählte weitere Versionen von verallgemeinerten Inversen

Drazin-Inverse

Sei A R n × n {\displaystyle A\in \mathbb {R} ^{n\times n}} eine Matrix mit Index i N {\displaystyle i\in \mathbb {N} } (der Index von A {\displaystyle A} ist die minimale ganze Zahl i 1 {\displaystyle i\geq 1} für die A i {\displaystyle A^{i}} und A ( i + 1 ) {\displaystyle A^{(i+1)}} den gleichen Kern haben). Dann ist die Drazin-Inverse diejenige eindeutig definierte n × n {\displaystyle n\times n} -Matrix A D {\displaystyle A^{D}} , die den Bedingungen

  • A i A D A = A i {\displaystyle A^{i}A^{D}\,A=A^{i}}
  • A D A A D = A D {\displaystyle A^{D}\,AA^{D}=A^{D}}
  • A A D = A D A {\displaystyle A\,A^{D}=A^{D}\,A}

genügt. Sie wurde von Michael Drazin eingeführt.

Berechnung

Zur Berechnung kann man die Zerlegung

T ( Λ 0 0 N ) T 1 = A {\displaystyle T{\begin{pmatrix}\Lambda &0\\0&N\end{pmatrix}}T^{-1}=A}

der Matrix A R n × n {\displaystyle A\in \mathbb {R} ^{n\times n}} in Jordan-Normalform nutzen, wobei Λ {\displaystyle \Lambda } der reguläre Teil der Jordan-Form sei und N {\displaystyle N} nilpotent. Die Drazin-Inverse ergibt sich dann zu

A D := T ( Λ 1 0 0 0 ) T 1 {\displaystyle A^{D}:=T{\begin{pmatrix}\Lambda ^{-1}&0\\0&0\end{pmatrix}}T^{-1}} .

Die Drazin-Inverse einer Matrix mit Index i = 1 {\displaystyle i=1} (für die also N {\displaystyle N} gleich der Nullmatrix ist) bezeichnet man auch als Gruppen-Inverse. Die Gruppen-Inverse ist eine Pseudoinverse nach der Definition von Koecher.[5]

Anwendungen

1. Eine wichtige Anwendung für die Drazin-Inverse ist die analytische Darstellung der Lösung zeitinvarianter linearer Deskriptorsysteme. Als Beispiel diene die Differenzengleichung

A x k + 1 = x k {\displaystyle Ax_{k+1}=x_{k}}

eines zeitdiskreten Deskriptorsystems mit der reellen n × n {\displaystyle n\times n} -Matrix A {\displaystyle A} . Die Lösung { x k } k = 0 , 1 , 2 , R n {\displaystyle \{x_{k}\}_{k=0,1,2,\ldots }\subset \mathbb {R} ^{n}} der Differenzengleichung erfüllt die Gleichungen A k x k = x 0 {\displaystyle A^{k}x_{k}=x_{0}} mit k = 0 , 1 , 2 , {\displaystyle k=0,1,2,\ldots } . Anfangswerte x 0 {\displaystyle x_{0}} sind also nur dann konsistent, wenn sie in allen Bildern der Matrizen A k {\displaystyle A^{k}} liegen (sonst bricht die Lösung nach endlich vielen Schritten ab). Die Lösung der Differenzengleichung ist dann x k = ( A D ) k x 0 {\displaystyle x_{k}=(A^{D})^{k}x_{0}} .

2. Für reelle oder komplexe n × n {\displaystyle n\times n} -Matrizen  A {\displaystyle A} mit Index  i {\displaystyle i} gilt die Gleichung

0 t exp ( A t ¯ ) d t ¯ = k = 0 i 1 t k + 1 ( k + 1 ) ! A k ( E n A A D ) + ( exp ( A t ) E n ) A D . {\displaystyle \int _{0}^{t}\exp(A{\bar {t}})\,d{\bar {t}}=\sum _{k=0}^{i-1}{\frac {t^{k+1}}{(k+1)!}}A^{k}(E_{n}-AA^{D})+(\exp(At)-E_{n})A^{D}.}

Damit lässt sich die Sprungantwort eines linearen zeitinvarianten dynamischen Systems

x ˙ = A x + b u y = c T x {\displaystyle {\begin{aligned}{\dot {x}}&=Ax+bu\\y&=c^{T}x\end{aligned}}}

mit Eingangssignal

u ( t ) = { 0 bei  t < 0 1 bei  t 0 , {\displaystyle u(t)={\begin{cases}0&{\text{bei }}t<0\\1&{\text{bei }}t\geq 0,\end{cases}}}

Zustandsvektor x : R R n × 1 {\displaystyle x:\mathbb {R} \rightarrow \mathbb {R} ^{n\times 1}} ( x ( 0 ) {\displaystyle x(0)} Nullvektor), Systemmatrix A R n × n {\displaystyle A\in \mathbb {R} ^{n\times n}} und Ein- beziehungsweise Ausgabevektoren b , c R n × 1 {\displaystyle b,c\in \mathbb {R} ^{n\times 1}} in der Form

y ( t ) = { 0 bei  t < 0 c T ( k = 0 i 1 t k + 1 ( k + 1 ) ! A k ( E n A A D ) + ( exp ( A t ) E n ) A D ) b bei  t 0 {\displaystyle y(t)={\begin{cases}0&{\text{bei }}t<0\\c^{T}\left(\sum _{k=0}^{i-1}{\frac {t^{k+1}}{(k+1)!}}A^{k}(E_{n}-AA^{D})+(\exp(At)-E_{n})A^{D}\right)b&{\text{bei }}t\geq 0\end{cases}}}

darstellen.

Restringierte verallgemeinerte Inversen – die Bott-Duffin-Inverse

Bei manchen praktischen Aufgabenstellungen ist die Lösung x {\displaystyle x} eines linearen Gleichungssystems

A x = b ( mit vorgegebenem  A R m × n  und  b R m ) {\displaystyle Ax=b\qquad ({\text{mit vorgegebenem }}A\in \mathbb {R} ^{m\times n}{\text{ und }}b\in \mathbb {R} ^{m})}

nur dann zulässig, wenn sie innerhalb eines gewissen linearen Teilraumes L {\displaystyle L} von R m {\displaystyle \mathbb {R} ^{m}} liegt. Man sagt auch, dass das Problem durch ein restringiertes lineares Gleichungssystem beschrieben wird (englisch constrained linear equation).

Im Folgenden werde der orthogonale Projektor auf L {\displaystyle L} mit P L {\displaystyle P_{L}} bezeichnet. Das restringierte lineare Gleichungssystem

A x = b x L {\displaystyle Ax=b\qquad x\in L}

ist genau dann lösbar, wenn das für das unrestringierte Gleichungssystem

( A P L ) x = b x R m {\displaystyle (AP_{L})x=b\qquad x\in \mathbb {R} ^{m}}

zutrifft. Ist der Unterraum L {\displaystyle L} ein echter Teilraum von R m {\displaystyle \mathbb {R} ^{m}} , so ist die Systemmatrix des unrestringierten Problems ( A P L ) {\displaystyle (AP_{L})} auch dann singulär, wenn sich die Systemmatrix A {\displaystyle A} des restringierten Problems invertieren lässt (in diesem Fall gilt m = n {\displaystyle m=n} ). Das erklärt, dass für die Lösung restringierter Probleme auch Pseudoinverse herangezogen werden. Man bezeichnet eine Pseudoinverse von ( A P L ) {\displaystyle (AP_{L})} dann auch als L {\displaystyle L} -restringierte Pseudoinverse von A {\displaystyle A} . Diese Definition scheint zunächst der Forderung 1 aus Abschnitt Allgemeine Pseudoinversen zu widersprechen. Dieser Widerspruch relativiert sich jedoch wieder, wenn man bedenkt, dass die L {\displaystyle L} -restringierte Pseudoinverse für bijektives  A {\displaystyle A} auf dem interessierenden Raum  L {\displaystyle L} injektiv ist und dass der Bildraum die gleiche Dimension wie  L {\displaystyle L} hat.

Ein Beispiel für eine Pseudoinverse mit der sich die Lösung eines restringierten Problems ermitteln lässt, ist die Bott-Duffin-Inverse von A {\displaystyle A} bzgl. L {\displaystyle L} , die durch die Gleichung

A L ( 1 ) := P L ( A P L + P L ) 1 {\displaystyle A_{L}^{(-1)}:=P_{L}(AP_{L}+P_{L^{\perp }})^{-1}}

definiert ist, falls die auf der rechten Seite auftretende gewöhnliche Inverse existiert.

Anwendungen

Die Bott-Duffin-Inverse kann zur Lösung der Gleichungen eines affin-linearen elektrischen Netzwerkes benutzt werden, wenn sich die Relation zwischen Zweigspannungsbelegungen   v R n {\displaystyle v\in \mathbb {R} ^{n}} und Zweigstrombelegungen   i R n {\displaystyle i\in \mathbb {R} ^{n}} in der Form

A i + u = u 0 i L , u L {\displaystyle Ai+u=u^{0}\qquad i\in L,\;u\in L^{\perp }}

darstellen lassen, wobei L {\displaystyle L} der Raum aller die kirchhoffschen Knotengleichungen erfüllenden Strombelegungen i {\displaystyle i} ist und u 0 {\displaystyle u^{0}} die Spaltenmatrix der in die Zweige eingespeisten unabhängigen Quellspannungen sein soll. An dieser Stelle fließt der graphentheoretische Satz von Tellegen ein, der besagt, dass die Räume der Zweigspannungsbelegungen und Zweigstrombelegungen, die die kirchhoffschen Maschen- beziehungsweise Knotengleichungen erfüllen, orthogonal komplementär zueinander sind.

Eine Eigenschaft der Bott-Duffin-Inversen ist, dass mit ihrer Hilfe die zu einer vorgegebenen Quellspannungsbelegung u 0 {\displaystyle u^{0}} zugehörigen Zweigströme

i = A L ( 1 ) u 0 {\displaystyle i=A_{L}^{(-1)}u^{0}}

und Zweigspannungen

u = ( E m A A L ( 1 ) ) u 0 {\displaystyle u=(E_{m}-AA_{L}^{(-1)})u^{0}}

berechnet werden können ( E m {\displaystyle E_{m}} steht für die Einheitsmatrix im R m × m {\displaystyle \mathbb {R} ^{m\times m}} ).

Literatur

  • W. Mackens, H. Voß: Mathematik I für Studierende der Ingenieurwissenschaften
  • A. Kielbasinski, H. Schwetlick: Numerische lineare Algebra, Deutscher Verlag der Wissenschaften, 1988

Einzelnachweise

  1. E. H. Moore: On the reciprocal of the general algebraic matrix. In: Bulletin of the American Mathematical Society 26, S. 394–395, 1920
  2. Roger Penrose: A generalized inverse for matrices. In: Proceedings of the Cambridge Philosophical Society 51, S. 406–413, 1955, doi:10.1017/S0305004100030401
  3. J. Stoer: Numerische Mathematik 1. Springer Verlag, 2002, ISBN 3-540-66154-9
  4. a b c Adi Ben-Israel, Thomas N.E. Greville: Generalized Inverses. Springer-Verlag, 2003, ISBN 0-387-00293-6
  5. a b Max Koecher: Lineare Algebra und analytische Geometrie, Springer-Verlag Berlin, 1997
  6. a b Nobuo Shinozaki, Masaaki Sibuya, and Kunio Tanabe: Numerical algorithms for the Moore-Penrose inverse of a matrix: Direct methods. Annals of the Institute of Statistical Mathematics, Springer Netherlands, Vol. 24, No. 1, Dec. 1972, pp. 193–203, doi:10.1007/BF02479751.
  7. T. N. E. Greville: Some applications of the pseudo inverse of a matrix. SIAM Rev., No. 2, 1960, pp. 15–22, doi:10.1137/1002004, JSTOR:2028054.