30.38 고유값 분해의 계산 복잡도와 대규모 행렬에 대한 근사 기법

30.38 고유값 분해의 계산 복잡도와 대규모 행렬에 대한 근사 기법

1. 정확한 고유값 분해의 계산 복잡도

1.1 밀집 행렬(Dense Matrix)의 경우

n \times n 밀집 행렬의 모든 고유값을 계산하는 표준 알고리즘은 다음의 복잡도를 갖는다.

단계일반 행렬대칭 행렬
하이젠베르크/삼대각 환원O(\frac{10}{3}n^3)O(\frac{4}{3}n^3)
QR 반복 (시프트 포함)O(n^3)O(n^2) (삼대각)
고유벡터 계산 (선택적)O(n^3)O(n^3)
총 비용O(n^3)O(n^3) (고유벡터 포함)

대칭 삼대각 행렬에 대한 QR 반복은 O(n^2)이므로, 고유값만 필요한 경우 대칭 행렬의 전체 비용은 환원 단계의 O(n^3)이 지배적이다.

1.2 메모리 요구량

n \times n 밀집 행렬의 저장에 O(n^2)의 메모리가 필요하다. n = 10^5이면 약 80 GB(배정밀도)가 소요되어 현실적으로 밀집 행렬 알고리즘의 적용이 불가능하다.

2. 대규모 행렬에서의 문제 설정

딥러닝과 데이터 과학에서 발생하는 행렬은 흔히 n \gg 10^4이다. 이 규모에서는 모든 고유값을 계산하는 것이 불필요하거나 불가능하므로, 다음의 문제 설정이 적절하다:

  • 소수의 최대 (또는 최소) 고유값과 고유벡터만 계산: PCA에서 상위 k개 주성분, 스펙트럼 반지름 추정 등.
  • 행렬이 희소(sparse)하거나 암묵적(implicit): 행렬-벡터 곱 Ax만 계산 가능한 경우.
  • 근사 해가 허용: 정확한 고유값이 아닌 근사 고유값으로 충분한 경우.

3. 크릴로프 부분 공간 방법(Krylov Subspace Methods)

3.1 크릴로프 부분 공간의 정의

초기 벡터 b \in \mathbb{R}^n에 대하여 **크릴로프 부분 공간(Krylov subspace)**은

\mathcal{K}_k(A, b) = \operatorname{span}\{b, Ab, A^2 b, \ldots, A^{k-1} b\}

으로 정의된다. 크릴로프 부분 공간은 행렬-벡터 곱만으로 구성되며, A의 원소에 직접 접근할 필요가 없다.

3.2 란초스 알고리즘(Lanczos Algorithm)

란초스 알고리즘은 대칭 행렬 A에 대하여 크릴로프 부분 공간 위에서 A를 삼대각 행렬로 사영하는 방법이다.

알고리즘:

초기 벡터 q_1 (\lVert q_1 \rVert = 1)을 선택한다. q_0 = 0, \beta_0 = 0으로 설정한다. j = 1, 2, \ldots, k에 대하여:

w_j = Aq_j - \beta_{j-1} q_{j-1}

\alpha_j = q_j^T w_j

w_j \leftarrow w_j - \alpha_j q_j

\beta_j = \lVert w_j \rVert

q_{j+1} = w_j / \beta_j

k단계 후, 삼대각 행렬 T_k의 대각 성분은 \alpha_1, \ldots, \alpha_k이고 부대각 성분은 \beta_1, \ldots, \beta_{k-1}이다. T_k의 고유값은 A의 극값 고유값(최대, 최소)의 좋은 근사이며, k가 증가함에 따라 근사가 개선된다.

계산 비용: 반복당 행렬-벡터 곱 1회(O(n^2) 밀집, O(\text{nnz}) 희소)와 O(n)의 벡터 연산. k번 반복의 총 비용은 O(kn^2) (밀집) 또는 O(k \cdot \text{nnz}) (희소)이다. k \ll n이면 전체 고유값 분해 O(n^3)보다 현저히 효율적이다.

3.3 아놀디 알고리즘(Arnoldi Algorithm)

비대칭 행렬에 대한 크릴로프 방법은 **아놀디 알고리즘(Arnoldi algorithm)**이다. 란초스 알고리즘의 일반화로서, A를 상부 하이젠베르크 행렬로 사영한다. 아놀디 알고리즘은 비대칭 행렬의 극값 고유값을 근사하는 데 사용된다.

4. 암묵적 재시작 란초스/아놀디(IRLM/IRAM)

기본 란초스/아놀디 알고리즘은 수치적 불안정성(직교성 상실)과 수렴 속도의 문제를 갖는다. 암묵적 재시작(implicit restart) 기법은 원하는 고유값에 대한 수렴을 가속하면서 메모리를 제한하는 방법이다.

ARPACK 라이브러리와 SciPy의 scipy.sparse.linalg.eigsh (대칭), scipy.sparse.linalg.eigs (비대칭)는 암묵적 재시작 란초스/아놀디 알고리즘을 구현한 것이다.

5. 무작위화 알고리즘(Randomized Algorithms)

5.1 무작위화 SVD/고유값 분해

Halko, Martinsson, Tropp (2011)의 “Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions“에서 제안된 무작위화 알고리즘은 대규모 행렬의 저계수 근사를 효율적으로 계산한다.

알고리즘 (무작위화 범위 탐색):

  1. 가우시안 랜덤 행렬 \Omega \in M_{n \times (k+p)}(\mathbb{R})을 생성한다 (p는 과표본 매개변수, 통상 p = 5 \sim 10).
  2. 표본 행렬 Y = A \Omega를 계산한다. 비용: k+p번의 행렬-벡터 곱.
  3. Y의 QR 분해 Y = QR을 수행하여 A의 열 공간의 근사 기저 Q를 구한다.
  4. 소규모 행렬 B = Q^T A Q (k \times k)를 구성하고, B의 고유값 분해를 수행한다.

계산 비용: 행렬-벡터 곱 O(k)회, 총 O(kn^2) (밀집) 또는 O(k \cdot \text{nnz}) (희소). 소규모 행렬의 고유값 분해 O(k^3).

오차 보장: 높은 확률로

\lVert A - Q Q^T A \rVert_2 \leq (1 + \epsilon) \sigma_{k+1}(A)

여기서 \sigma_{k+1}A(k+1)번째 특이값이다. 특이값이 급속히 감소하면 근사가 정밀하다.

5.2 거듭제곱 반복에 의한 개선

특이값의 감소가 완만한 경우, 거듭제곱 반복(power iteration)을 결합하여 정밀도를 향상시킨다:

Y = (AA^T)^q A \Omega

q번의 거듭제곱은 상위 k개 특이값의 우세를 증폭시켜 근사의 질을 개선한다. 비용은 (2q+1)(k+p)번의 행렬-벡터 곱이다.

6. 분할 정복법(Divide and Conquer)

대칭 삼대각 행렬의 고유값 분해에 대한 **분할 정복법(divide and conquer algorithm)**은 Cuppen (1981)에 의하여 제안되었으며, 다음 단계로 구성된다:

  1. 삼대각 행렬을 두 개의 작은 삼대각 행렬과 계수 1의 수정(rank-one modification)으로 분할한다.
  2. 각 부분 행렬의 고유값 분해를 재귀적으로 구한다.
  3. 세속 방정식(secular equation)을 풀어 전체 고유값을 결정한다.

분할 정복법의 평균 비용은 O(n^2)이며 (고유벡터 포함 시 O(n^3)의 최악 경우가 있으나 실용적으로 드물다), 높은 병렬성을 가져 GPU 구현에 유리하다. LAPACK의 dstevd와 cuSOLVER의 GPU 고유값 분해가 이 알고리즘을 사용한다.

7. 알고리즘 선택 지침

상황권장 알고리즘비용
밀집, 모든 고유값QR 알고리즘O(n^3)
대칭 밀집, 모든 고유값분할 정복법O(n^2) ~ O(n^3)
희소, 소수의 극값 고유값란초스/아놀디 (ARPACK)O(k \cdot \text{nnz})
대규모, 저계수 근사무작위화 SVDO(k \cdot \text{nnz})
최대 고유값 1개멱승법반복당 O(\text{nnz})
특정 고유값 1개시프트 역멱승법반복당 O(n^2)

8. 딥러닝에서의 대규모 고유값 근사

대규모 공분산 행렬의 주성분. 고차원 데이터의 PCA에서 공분산 행렬을 명시적으로 구성하지 않고, 데이터 행렬과의 곱으로 무작위화 SVD를 수행하는 것이 표준적이다. Facebook의 fbpca 라이브러리와 scikit-learn의 TruncatedSVD가 이 접근을 채택한다.

헤시안의 상위 고유값 추정. 심층 신경망의 손실 함수 헤시안은 매개변수 수가 10^6 이상으로서 밀집 행렬의 저장 자체가 불가능하다. 그러나 헤시안-벡터 곱 Hv는 자동 미분(automatic differentiation)으로 O(n)에 계산할 수 있으므로, 란초스 알고리즘을 적용하여 헤시안의 상위 고유값과 고유벡터를 추정할 수 있다.

스펙트럼 정규화의 최대 특이값 추정. GAN의 스펙트럼 정규화에서 가중치 행렬의 최대 특이값은 멱승법 1~2회의 반복으로 추정한다. 이는 O(n \cdot d)의 비용으로 (n \times d 가중치 행렬), 학습 루프 내에서 실시간으로 수행 가능하다.

그래프 신경망의 라플라시안 스펙트럼. 그래프 라플라시안의 하위 고유값과 고유벡터는 그래프의 클러스터링 구조와 주파수 영역 표현을 결정한다. 수백만 노드의 그래프에서 란초스 알고리즘으로 하위 k개의 고유쌍을 효율적으로 계산한다.