Power Iteration

Power iteration is an extremely simple iterative method to approximate the eigenvector associated with the largest eigenvalue (the dominant eigenvector) of a matrix. The key idea: if you repeatedly multiply a vector by a matrix $A$, the component along the eigenvector with the largest eigenvalue $\lambda_{\max}$ gets amplified the most, so the direction of the vector converges to that eigenvector (after normalization each step). This works well for large, sparse matrices where full eigen-decomposition (e.g., via SVD or QR) is too expensive.

In applications like PageRank and PCA, we care primarily about this dominant eigenvector: PageRank’s ranking vector is the dominant eigenvector of a stochastic matrix; in PCA, the first principal component is the dominant eigenvector of the covariance matrix. Power iteration is attractive here because it uses only matrix–vector products and normalizations, which are easy to implement and scale to large dimensions.

---

1. Background Knowledge (concepts & theory)

• Eigenvalues/eigenvectors For a square matrix $A \in \mathbb{R}^{n \times n}$, a nonzero vector $v$ is an eigenvector if:

$$A v = \lambda v$$

for some scalar $\lambda$, the eigenvalue. Geometrically, applying $A$ to $v$ only scales it (by $\lambda$), not change its direction.

• Dominant eigenvalue/eigenvector Let the eigenvalues be $\lambda_1, \dots,$\lambda$_n$ ordered so that $|$\lambda$_1| > |$\lambda$_2| \ge \dots \ge |$\lambda$_n|$.

• $\lambda_1$ is the dominant eigenvalue.

• Its eigenvector $v_1$ is the dominant eigenvector. Power iteration converges to $v_1$ under mild conditions (e.g., $|$\lambda$_1| > |$\lambda$_2|$ and the starting vector has a nonzero component in direction $v_1$).

• Why normalization? If you repeatedly compute $x_{k+1} = A x_k$, the norm $\|x_k\|$ tends to blow up or decay (depending on $|$\lambda$_1|$). Normalizing:

$$v_{k+1} = \frac{A v_k}{\|A v_k\|}$$

keeps vectors on the unit sphere so the sequence converges in direction to the dominant eigenvector.

---

2. Algorithm / Approach

The power iteration pattern is:

1. Pick an initial nonzero vector $v_0$ (often random).
2. Repeat:

• Multiply by matrix: $w = A v_k$
• Normalize: $v_{k+1} = \dfrac{w}{\|w\|}$

3. Stop when convergence: e.g., $\|v_{k+1} - v_k\|$ is below a tolerance or max iterations reached.
4. (Optionally) approximate the dominant eigenvalue with the Rayleigh quotient:

$$\hat{\lambda} = v_k^\top A v_k$$

This is a simple fixed-point iteration: you are iterating a function $f(v) = \dfrac{A v}{\|A v\|}$ until it stops changing much.

---

3. Step-by-Step Strategy (for implementation)

Assuming you have:

• Matrix $A$ as an $n \times n$ array
• Parameters: max_iters, tol (tolerance)

A typical strategy:

1. Initialize vector

• Create a random vector:

v = np.random.randn(n)
v = v / np.linalg.norm(v)

• Ensure it is not the zero vector.

2. Main iteration loop

for k in range(max_iters):
# 1) Matrix-vector product
w = A @ v # or A.dot(v)

# 2) Compute norm (e.g., L2)
norm_w = np.linalg.norm(w)

# 3) Normalize
v_new = w / norm_w

# 4) Check convergence
if np.linalg.norm(v_new - v) < tol:
v = v_new
break

v = v_new

3. Return the dominant eigenvector

• v is your approximate dominant eigenvector.

4. Optional: estimate eigenvalue

lambda_est = float(v.T @ (A @ v))

When coding for a platform problem, you’ll wrap this logic in a function and use the platform’s provided matrix/vector types.

---

4. Common Pitfalls

• Not normalizing each step Without normalization, the vector may overflow or underflow numerically; convergence in direction is then unstable.

• Bad stopping criteria

• Using only k == max_iters without checking convergence can waste iterations or return a poor approximation.

• A robust choice is ||v_{k+1} - v_k|| < tol or ||A v_k - λ v_k|| < tol.

• Starting vector issues

• If the initial vector is orthogonal to the dominant eigenvector, convergence fails. With a random vector, this almost never happens in practice, but avoid starting with the zero vector or a very special structure.

• Non-dominant or close eigenvalues

• If $|$\lambda$_1|$ is close to $|$\lambda$_2|$, convergence can be slow.

• If there are multiple eigenvalues with the same largest magnitude, the method may converge to some vector in the subspace they span, not a unique eigenvector.

• Non-symmetric matrices

• Power iteration still works for general matrices, but:

• The dominant eigenvalue might be complex.

• Implementation with real arithmetic usually assumes $A$ has a real dominant eigenvalue/eigenvector.

• Numerical stability in norm

• Use a stable norm function (e.g., np.linalg.norm); don’t manually compute sqrt(sum(x[i]**2)) if you can avoid it.

---

5. Time & Space Complexity

Let $n$ be the dimension of the matrix (i.e., $A$ is $n \times n$) and let T be the number of iterations until convergence.

• Time complexity

• Each iteration does a matrix–vector product $A v$.

• Dense $A$: $O(n^2)$ per iteration.

• Sparse $A$: $O($\text{nnz}$(A))$ per iteration, where $\text{nnz}$ is number of nonzeros.

• Total: $O(T \cdot n^2)$ for dense, or $O(T \cdot$\text{nnz}$(A))$ for sparse.

• Space complexity

• You store:

• The matrix $A$: $O(n^2)$ for dense, or $O($\text{nnz}$(A))$ for sparse.

• A few vectors of size $n$: $O(n)$.

• Additional overhead is constant, so asymptotically it is dominated by the storage of $A$.

This complexity profile explains why power iteration is popular in large-scale ML/vision tasks: it uses simple, repeated matrix–vector products that can leverage sparsity and GPU acceleration.