Eigenvalue Algorithms – Power Method, QR Iteration, and Applications
1. Motivation: Why Compute Eigenvalues Numerically?
---------------------------------------------------
• Many problems reduce to eigenvalue computations:
- Stability of dynamical systems (eigenvalues of Jacobian matrices)
- Vibrations and normal modes (eigenvalues of stiffness/mass matrices)
- Principal Component Analysis (PCA) in data science
- Markov chains and Google’s PageRank
• Exact formulas (like the quadratic formula) exist only for low degrees; in practice, we
use iterative numerical algorithms.
2. Power Method
----------------
• Goal: Approximate the dominant eigenvalue λ_1 (largest in magnitude) and corresponding
eigenvector v_1 of a matrix A.
• Algorithm (basic version):
1) Choose an initial vector x^{(0)} ≠ 0.
2) For k = 0, 1, 2, …:
y^{(k+1)} = A x^{(k)}
x^{(k+1)} = y^{(k+1)} / ||y^{(k+1)}||
3) Approximate eigenvalue:
λ^{(k)} ≈ (x^{(k)})^T A x^{(k)} / (x^{(k)})^T x^{(k)} (Rayleigh quotient).
• Convergence intuition:
- If A has eigenvalues |λ_1| > |λ_2| ≥ … and x^{(0)} has a component in the v_1
direction, then A^k x^{(0)} ≈ λ_1^k v_1.
- After normalization, x^{(k)} tends to v_1.
• Limitations:
- Only finds one eigenpair (the dominant one).
- Convergence slows when |λ_1| ≈ |λ_2|.
- Sensitive to defects and non-normality.
3. Inverse and Shifted Power Methods
------------------------------------
• To find the smallest eigenvalue or an eigenvalue near a shift σ, we apply the power
method to (A − σI)^{-1}.
• Inverse iteration:
x^{(k+1)} = (A − σI)^{-1} x^{(k)} / ||(A − σI)^{-1} x^{(k)}||
• Shift-and-invert:
- Eigenvalues of (A − σI)^{-1} are 1 / (λ_j − σ).
- If σ is close to some λ_j, then |1 / (λ_j − σ)| is large, so power method on (A −
σI)^{-1} finds the eigenvalue near σ.
• Practical note:
- Each iteration requires solving a linear system (A − σI) y = x, often via LU
decomposition.
4. Rayleigh Quotient Iteration
-------------------------------
• More advanced scheme using adaptive shifts:
1) Start with x^{(0)} normalized and σ_0 = Rayleigh quotient R(x^{(0)}).
2) For k = 0, 1, 2, …:
Solve (A − σ_k I) y^{(k+1)} = x^{(k)}
Normalize x^{(k+1)} = y^{(k+1)} / ||y^{(k+1)}||
Update σ_{k+1} = R(x^{(k+1)})
• When it converges to an eigenpair, it does so very fast (often cubically) for Hermitian
matrices.
5. QR Iteration – Basic Idea
----------------------------
• QR iteration is a general algorithm to compute all eigenvalues of a matrix A.
• One step of QR iteration:
1) Compute the QR factorization of A_k:
A_k = Q_k R_k
where Q_k is orthogonal (or unitary) and R_k is upper triangular.
2) Form
A_{k+1} = R_k Q_k
• Fundamental property:
- A_{k+1} = Q_k^T A_k Q_k (in the real case), so A_{k+1} is similar to A_k.
- Similar matrices have the same eigenvalues.
- Under mild conditions, A_k tends towards an upper triangular matrix whose diagonal
entries are the eigenvalues of A.
6. Practical QR Algorithm Enhancements
--------------------------------------
• Direct QR on A is expensive; instead, we:
- First reduce A to Hessenberg form (for general matrices) or tridiagonal form (for
symmetric ones) using orthogonal transformations.
- Perform QR iteration on the reduced form, which is cheaper per step.
• Shifts:
- Instead of factoring A_k, we factor A_k − µ_k I, where µ_k is a shift approximating
an eigenvalue.
- Then set
A_{k+1} = R_k Q_k + µ_k I.
- Proper choice of shifts dramatically accelerates convergence.
• Implicit QR with Wilkinson shift is the cornerstone of practical eigenvalue software
(like LAPACK routines).
7. Symmetric Eigenproblems
--------------------------
• For real symmetric (or complex Hermitian) matrices:
- All eigenvalues are real.
- There exists an orthogonal (unitary) matrix Q such that
A = Q Λ Q^T.
- Specialized algorithms exploit symmetry:
* Tridiagonal reduction
* Divide-and-conquer methods
* MRRR (Multiple Relatively Robust Representations)
• These methods are more stable and efficient than general eigenvalue algorithms.
8. Applications: PCA and SVD
----------------------------
• Principal Component Analysis (PCA):
- Given data matrix X, the covariance matrix C = (1/n) X^T X is symmetric positive
semi-definite.
- Eigenvectors of C give principal directions; eigenvalues measure variance along these
directions.
- Numerically, PCA is often implemented via Singular Value Decomposition (SVD) instead
of naive eigen-decomposition.
• SVD:
- Any m×n matrix A can be written as
A = U Σ V^T
where U and V are orthogonal (or unitary) and Σ is diagonal with non-negative entries
(the singular values).
- Singular values are square roots of eigenvalues of A^T A.
• Many eigenvalue algorithms can be adapted to compute the SVD or are related to it
(e.g., Golub–Kahan bidiagonalization).
9. Numerical Stability Considerations
-------------------------------------
• Orthogonal (unitary) transformations are norm-preserving:
- They do not magnify errors, which is crucial in finite precision arithmetic.
- QR algorithm is built as a sequence of such transformations, making it numerically
stable.
• Ill-conditioned eigenproblems:
- When eigenvalues are very close together or the matrix is nearly defective, small
perturbations can cause large changes in eigenvectors.
- Error analysis uses condition numbers of eigenvalues and spectral projectors.
10. Summary
-----------
• Power method is simple and useful to approximate a dominant eigenpair.
• Shifted and inverse variants, along with Rayleigh quotient iteration, refine this idea
for interior eigenvalues.
• QR iteration (with shifts and reductions) is the core method for finding all
eigenvalues in practice.
• Eigenvalue algorithms are essential for applications in data science, physics,
engineering, and beyond.