[수치해석 20] Givens rotation, QR algorithm and Hessenberg matrix

저번 포스트에서 QR factorization을 하는 방법 중 Projection과 Reflection 관점에서 QR로 분해하는 방법을 정리했다.

[수치해석] - [수치해석 19] Gram-schmidts, Householder transformation

이번 포스트에서는 Projection, Reflection에 이어 Rotation 관점에서 분해하는 방법인 Givens rotation과 QR을 활용하는 iterative mdethod를 소개할 예정이다.

 

먼저 소개할 Givens rotation은 이전 소개한 Gram-schmidts 방법이나 Housholder transformation과 마찬가지로, orthogonal matrix $Q'$을 곱해나가 $A$의 각 column의 원소들을 지워나가 upper triangular 형태의 행렬로 만드는 방법을 사용한다. 그렇다면 orthogonal matrix를 어떻게 구성해야 rotation이 되고, upper triangular 형태가 되도록 원소들을 지워나갈 수 있을까?

---

Givens rotation (Rotation)

 

아이디어는 그림에서 보이는 것과 같이 단순하다. $xy$라는 좌표계와 $(a, b)$라는 점이 주어져 있을 때, $(a, b)$를 $X$축 위에 위치시키도록 $XY$ 좌표계를 회전시켜 정의한다.

 

$xy$ 좌표계에서 $(a, b)$라는 점이 $XY$에서 $(a', b')$이라고 한다면 다음과 같은 관계식이 성립한다.

$$\begin {align} a &= a' \cos {\theta} - b' \sin {\theta} \\ b &= a'\sin {\theta} + b'\cos {\theta} \end {align}$$

 

이를 행렬로 표현하면

$$\begin {bmatrix} a \\ b \end {bmatrix} = \begin {bmatrix} \cos\theta & -\sin\theta \\ \sin\theta & \cos\theta \end {bmatrix} \begin {bmatrix} a' \\ b' \end {bmatrix} $$

이때 우리가 원하는 대로 $b' = 0$으로 만들고 역함수를 양변에 곱해주면

$$\begin {bmatrix}\sqrt {a^2 + b^2} \\ 0 \end {bmatrix}= \begin {bmatrix} \cos\theta & \sin\theta \\ -\sin\theta & \cos\theta \end {bmatrix} \begin {bmatrix} a \\ b\end {bmatrix}, \;\; \cos\theta = \frac {a}{\sqrt {a^2 + b^2}}, \sin\theta = \frac {b}{\sqrt {a^2+b^2}} $$

라는 식을 얻을 수 있다. 즉 이걸로 준비는 끝났다. Rotation 시켜서 벡터의 하나의 값으로 다른 하나의 값을 지울 수 있다.

 

Example)

$x = \begin {bmatrix} 1 \\ 2 \\ 3 \\ 4 \end {bmatrix}$라는 벡터가 주어졌을 때 2번째 원소를 통해 4번째 원소를 지워보자. 먼저 두 값이 2,4이기 때문에 회전시킬 각 $\theta$에 대해서 $\cos\theta = \frac {2}{\sqrt {2^2 + 4^2}} = \frac {1}{\sqrt {5}}, \;\; \sin\theta = \frac {2}{\sqrt {5}}$가 된다.

 

그럼 2번째 원소로 4번째 원소를 제거하는 방법은 

$$\begin {bmatrix}1&0&0&0\\0&c&0&s\\0&0&1&0\\0&-s&0&c\end {bmatrix} \begin {bmatrix}1\\2\\3\\4\end {bmatrix} =  \begin {bmatrix}1\\2\sqrt {5}\\3\\0\end {bmatrix} $$

 

이와 같은 방법으로 1번째 원소로 3번째 원소를 제거하면 

$$\begin {bmatrix} c&0&s&0\\0&1&0&0\\-s&0&c&0\\0&0&0&1\end {bmatrix} \begin {bmatrix}1\\2\\3\\4\end {bmatrix}= \begin {bmatrix}\sqrt {10}\\2\\0\\4\end {bmatrix}$$

 

위 결과 행렬을 보면 특이한 점이 하나 있다. 2번째 원소로 4번째 원소를 제거할 때 1,3번째 원소는 변함이 없다. 1번째 원소로 3번째 원소를 제거할 때 2,4번째 원소에 대한 변화도 없다. 동시에 진행해도 상관없는 연산이다. 즉 병렬적으로 가능한 연산이 되어 속도 측면에서 가속할 수 있는 특징을 가진다.

---

Householder tranformation과의 비교

Householder transformation에 비해 알고리즘 자체가 매우 간단하다. 그럼 Givens만 쓰면 되는 거 아닌가? Givens rotation의 가장 큰 단점은 한 번에 하나의 원소만 지울 수 있다. 즉 행렬이 매우 크고 dense 한 행렬이라면 Householder transformation이 훨씬 빠르고 효율적이다. Givens rotation은 반대로 하나씩의 원소만 지울 수 있다는 점에서 sparse 한 행렬에서 효율적이다.

 

현실 세계의 문제에서는 sparse 한 행렬을 다루는 경우가 많기 때문에 Givens rotation을 활용하는 경우가 많다. 하지만 행렬의 크기가 1000 by 1000이고 dense 한 행렬일 때 첫 번째 column을 $\begin {bmatrix}1&0&0&\cdots&0 \end {bmatrix}^\top$로 만들기 위해서는 999회의 rotation이 필요해 버리기 때문에, 한 번에 변환할 수 있는 Householder transformation을 무시할 수 없다.

 

여러 가지 방법을 통해서 QR factorization을 하는 방법을 알아봤으니 이제 활용하는 방법에 대해서 정리해 보자

---

QR algorithm

먼저 알고리즘을 소개하고 설명하자.

 

For given $A$, take orthogonal matrix $Q_0 \in \mathbb {R}^{n \times n}$ as initial

$T_0 = Q_0^\top A Q_0$

for $k = 1,2,\cdots$

    $T_{k-1} = Q_k R_k$ (QR factorization)

    $T_k = R_k Q_k$

 

단순히 QR factorization을 하고, QR의 순서를 바꿔서 곱해주는 것의 반복이다. 이를 통해 뭘 얻을 수 있을까?

 

$T_1 = R_1Q_1 = Q_1^\top T_0 Q_1 = Q_1^\top Q_0^\top A Q_0 Q_1$

$T_k = Q_k^\top \cdots Q_0^\top A Q_0 \cdots Q_k = Q^\top A Q$

즉 반복을 통해서 얻는 $T_k = Q^\top A Q$가 되고, 이러한 orthogonal 한 행렬에 의한 $Q\top A Q$ 변환은 eigenvalue를 변화시키지 않는다. 이때 $T_k$가 upper triangular와 같이 주대각선 원소들이 정직하게 eigenvalue가 된다면, 이러한 QR factorization과 순서를 바꾼 곱만으로 Eigenvalue들을 구할 수 있게 된다.

 

그럼 $Q_0$를 upper triangular인 $T_0$로 한 번에 만들 수 있다면 정말 베스트인데, 안타깝게도 그러한 방법은 아직 없다. upper triangular에 최대한 가까운 형태로 변환하는 방법은 있을까?

---

Hessenberg matrix

$H \in \mathbb {R}^{n \times n}, h_{ij} = 0 \;\; \forall i > j+1$

Shape wise $\begin {bmatrix} *&*&*&*&*\\*&*&*&*&*\\0&*&*&*&*\\0&0&*&*&*\\0&0&0&*&* \end {bmatrix}$

 

이러한 형태의 행렬을 Hessenberg matrix라고 하고, 임의의 행렬 $A$에 대해서 한 번에 Hessenberg 형태를 만들어낼 수 있다. 그리고 Hessenberg matrix의 특징은 Givens rotation을 통해 $H = QR$로 QR factorization 후 $RQ$를 연산하면, 여전히 Hessenberg matrix이다.

 

즉 $A$라는 행렬에 대해서 Householder transformation을 통해 Hessen reduction (Hessenberg matrix 형태로 만드는 것)을 한 번에 진행하고, Givens transformation을 통해 QR 분해 및 RQ 연산으로 iteration을 진행하면, 굉장히 빠르게 upper triangular 형태의 $T_k$를 구할 수 있고, 연산량이 너무 많아 구하기 힘들었던 Eigenvalue problem을 이전에 Power method처럼 하나씩 구하지 않고 한 번에 모든 Eigenvalue를 얻을 수 있다.

---

이번 포스트에서는 Givens rotation을 통해 sparse 한 행렬에 대해서 QR factorization을 할 수 있는 방법에 대해서 정리했고, QR algorithm을 통해 Eigenvalue problem을 해결할 수 있는 방법에 대해서 정리했다.

 

다음 포스트에서 정리할 부분은 Hessenberg matrix를 한 번에 구할 수 있다고 했는데, 그러한 Hessen reduction은 어떻게 진행할지, 다음으로 Givens rotation을 통해 QR algorithm을 하면 왜 upper triangular matrix로 수렴할지에 대해서 알아볼 예정이다.