Verstehen, wie Numpy SVD macht

  • Ich habe verschiedene Methoden angewandt, um sowohl den Rang einer Matrix als auch die Lösung eines Matrixsystems von Gleichungen zu berechnen. Ich bin auf die Funktion linalg.svd gestoßen. Vergleicht man dies mit meinen eigenen Bemühungen, das System mit Gaussian Elimination zu lösen, scheint es sowohl schneller als auch präziser zu sein. Ich versuche zu verstehen, wie das möglich ist.

    Soweit ich weiß, verwendet die Funktion linalg.svd einen QR-Algorithmus, um die Eigenwerte meiner Matrix zu berechnen. Ich weiß, wie das mathematisch funktioniert, aber ich weiß nicht, wie Numpy es so schnell und ohne viel Präzision schafft.

    Also meine Frage: Wie funktioniert die numpy.svd Funktion funktionieren, und genauer gesagt, wie gelingt es, dies schnell und genau zu tun (im Vergleich zur Gaußschen Eliminierung)?

    05 April 2012
    Rob Allen
1 answer
  • Aufgrund Ihrer Fragestellung gehe ich davon aus, dass Ihre Matrix quadratisch ist. Die SVD-Routinen von LAPACK, wie z. B. zgesvd , gehen für quadratische Matrizen im Wesentlichen in drei Stufen vor:

    1. Implizite Berechnung der unitären Matrizen $ U_A $ und $ V_A $ als Produkte von Householder-Transformationen, so dass die allgemeine Matrix $ A $ auf eine echte, obere bidiagonale Matrix $ B reduziert wird: = U_A ^ HA V_A $. Beim Verlassen dieser Subroutine werden die Householder-Vektoren (normalisiert, so dass ihre erste Eingabe implizit eins ist) für $ U_A $ und $ V_A $ jeweils in den Abschnitten von $ B $ unterhalb und rechts von Haupt- und Superdiagonale gespeichert . Dieser Schritt erfordert $ O (n ^ 3) $ work.
    2. Eine Variation eines Algorithmus zur Berechnung der Eigenwertzerlegung einer reellen symmetrischen tridiagonalen Matrix wird zur Berechnung einer bidiagonalen SVD verwendet. Obwohl der bekannteste Algorithmus für echte symmetrische tridiagonale EVPs (MRRR) meines Wissens noch nicht stabil für die bidiagonale SVD ist, gibt es eine interessante Diskussion hier . LAPACK verwendet derzeit einen Divide and Conquer -Ansatz für die bidiagonale SVD. Die bidiagonale SVD ergibt $ \ {U_B, V_B, \ Sigma \} $, so dass $ B = U_B \ Sigma V_B ^ H $. Für diesen Schritt ist $ O (n ^ 2) $ work erforderlich, wenn MRRR für SVD stabil gemacht wird, derzeit jedoch so viel wie $ O (n ^ 3) $ work erforderlich ist.
    3. Dieser Schritt ist trivial. Da $ U_A B V_A ^ H = A $, $ A = (U_A U_B) \ Sigma (V_A V_B) ^ H $ ist, so dass die SVD nach dem Anwenden der Householder-Transformationen explizit gebildet werden kann, die $ U_A $ und $ V_A implizit definieren $ bis $ U_B $ bzw. $ V_B $. Für diesen Schritt war $ O (n ^ 3) $ erforderlich.
    10 April 2012
    Rex MFrank Schwieterman