Textbook: Elementary Linear Algebra with Supplemental Applications, Eleventh Edition, International student version, Howard Anton, Chris Rorres
7.2 Orthogonal Diagonalization
$n \times n$인 정사각행렬 $A$가 $A^{-1} = A^T$이면 직교행렬이라하고, 따라서 $A^TA = AA^T = I$이다.
- $n \times n$ 인 직교행렬 $A$에 대하여 아래 문장들은 뜻이다.
- $A$가 직교행렬이다.
- $A$의 행 벡터는 $R^n$ 유클리드 공간의 정규직교기저를 이룬다. orthonormal set이다.
- $A$의 열 벡터는$R^n$ 유클리드 공간의 정규직교기저를 이룬다. orthonormal set이다.
- $\| A\mathbf{x} \| = \| \mathbf{x} \|$
- $A \mathbf{x} \cdot A \mathbf{y} = \mathbf{x} \cdot \mathbf{y}$
- 직교행렬의 역행렬도 직교행렬이다.
- 직교행렬의 곱셈도 직교행렬이다.
- $A$가 직교행렬이면, $\mathrm{det}(A) = \pm 1$ 이다.
두 정사각행렬 $A$, $B$에 대하여 $P^TAP = B$를 만족하는 행렬 $P$가 존재하면, $A$와 $B$는 직교닯음(orthogonally similar)이라 한다.
이때, $A$를 직교대각화 가능하다(orthogonally diagonalizable)고 하며, $P$가 $A$를 직교대각화한다라고 한다.
- $n \times n$ 인 정사각행렬 $A$에 대하여 아래 문장들은 같은 뜻이다
- $A$는 직교대각화가능하다.
- $A$는 $n$개의 고윳값을 갖는 정규직교기저를 갖는다.
- $A$는 대칭행렬이다.
- 대칭행렬이면, 모든 교윳값은 실수이다.
- 대칭행렬이면, 서로 다른 고유벡터들은 서로 직교이다.
$n \times n$ 대칭행렬의 직교대각화
- $A$의 eigenspace의 basis를 구한다.
- Gram-Schumidt를 이용하여 각 eigenspace의 orthonormal basis를 구한다.
- $P$의 열벡터를 (2)에서 구한 벡터로 구성한다. 이 행렬은 $A$를 직교대각화한 것으로, $D = P^TAP$는 주대각성분이 $A$의 고윳값인 행렬과 같다.
Spectral Decomposition (스펙트럼 분해)
$P = \begin{bmatrix} \mathbf{u}_1 & \dots \mathbf{u}_n \end{bmatrix}$가 대칭행렬 $A$를 직교대각화할 때, ($\mathbf{u}_i$는 $A$의 고윳값 $\lambda_i$에 대응하는 단위고유벡터), $D = P^TAP$는 고윳값을 갖는다는 것을 안다.
이를 전개하면
$A = P^TDP = \begin{bmatrix} \mathbf{u}_1 & \dots \mathbf{u}_n \end{bmatrix} \begin{bmatrix} \lambda_1 & 0 & \cdots & 0 \\ 0 & \lambda_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & \lambda_n \end{bmatrix} \begin{bmatrix} \mathbf{u}_1^T \\ \mathbf{u}_2^T \\ \vdots \\ \mathbf{u}_n^T \end{bmatrix} = \displaystyle\sum_{i=1}^{n}\lambda_i \mathbf{u}_1 \mathbf{u}_1^T$
9.5 Singular Value Decomposition
$n \times n$ symmetric matrix의 대각화 이론을 $m \times n$ 행렬로 확장해보자.
그 결과는 디지털 정보의 압축, 저장, 전송(compression, storage, transmission)에 적용 가능하고, 많은 computational algorithm에 적용된다.
모든 대칭 행렬 $A$에 대하여
$$A = PDP^T$$
가 성립한다. 이때 $P$는 $A$의 고유벡터를 모은 $n \times n$ 의 직교행렬이고, $D$는 $P$의 열벡터에 각각 대응되는 고윳값들을 대각성분으로 갖는 대각행렬이다. 이를 고윳값분해(eigenvalue decomposition)이라 부르고 줄여서 EVD라 부른다.
$A$가 $n \times n$이지만 대칭행렬이 아니라면, EVD는 존재하지 않는다. 대신에 Hessenberg decomposition이 존재하여
$$A = PHP^T$$
를 만족하는 직교행렬 $P$, 상헤센베르크 행렬(upper Hessenberg form) $H$가 존재한다.
$A$가 실수 고윳값을 갖는 경우, 슈어 분해(Schur decomposition)가 존재하여
$$A = PSP^T$$
를 만족하는 직교행렬 $P$와 상삼각행렬(upper triangular) S가 존재한다. (슈어분해는 복소수 정사각행렬에서도 분해가 가능하고 이때는 켤레전치를 이용한다)
eigenvalue, Hessenberg, Schur decomposition은 $D, H, S$가 $A$보다 단순한 형태인 것도 중요하지만, 반올림오차(roundoff error) 없이 직교행렬을 나타낼 수 있기 때문이다.
일반적인 정사각행렬의 분해는 $A = PJP^{-1}$의 형태인 것을 알 수 있다. ($P$는 invertible하지만 반드시 직교행렬일 필요는 없다.) 혹은 $A = U \Sigma V^{T}$의 형태라고 생각할 수 있다. ($U$, $V$는 직교행렬이지만 반드시 같은 형태일 필요는 없다). 첫번째 형태의 경우 Jordan canonical form이라 부른다. 특정 분야에서 중요하지만 수치적으로 중요도가 떨어지는데 그 이유는 roundoff difficulty가 존재하기 때문이다. 따라서 이 책에서는 두번째 형태에 집중하기로 한다.
임의의 $m \times n$ 행렬 $A$에 대하여 다음이 성립한다.
- $A$와 $A^TA$는 같은 null space를 갖는다.
- $A$와 $A^TA$는 같은 row space를 갖는다.
- $A^T$와 $A^TA$는 같은 column space를 갖는다.
- $A^T$와 $A^TA$는 같은 rank를 갖는다.
- $A^TA$는 직교대각화가 가능하다. (Orthogonally diagonalizable)
- $A^TA$의 고윳값들은 음수가 아니다.
$m \times n$ 행렬 $A$에 대하여, $A^TA$의 고윳값들을 $\lambda_i$ ($i = 1, 2, ..., n$)라 할 때
$$\sigma_i = \sqrt{\lambda_i} (i=1, 2, ..., n)$$
을 $A$의 특이값(singular value)이라 한다.
- 일반적으로 고윳값을 표현할 때 $\lambda_1 \ge \lambda_2 \ge \dots \lambda_n \ge 0$이므로 특이값도 $\sigma_1 \ge \sigma_2 \ge \dots \sigma_n \ge 0$이다.
Example 1. Singular Values
$A = \begin{bmatrix}1 & 1 \\ 0 & 1 \\ 1 & 0 \end{bmatrix}$의 특이값을 구해보자.
$A^TA = \begin{bmatrix} 2 & 1 \\ 1 & 2 \end{bmatrix}$이고, 특성방정식을 구하면 $(\lambda - 3)(\lambda -1) = 0$이므로 고윳값은 $\lambda_1 = 3, \lambda_2 = 1$이다. 따라서 특이값은 $\sigma_1 = \sqrt{3}, \sigma_2 = 1$이다.
특이값 분해(Singular Value Decomposition, SVD)
$A$가 $m \times n$ 행렬일 때,
$$A = U \Sigma V^T$$
의 형태로 나타낼 수 있다. 이 때 $U, V$는 직교행렬이고, $\Sigma$는 주대각선의 원소들이 $A$의 특이값을 갖고 나머지는 0인 $m \times n$인 대각행렬이다.
$k$가 고윳값의 개수일 때, 아래처럼 분해할 수 있다.
$$
A = U \Sigma V^T =
\begin{matrix}
\underbrace{\left[\begin{matrix}\mathbf{u_1} & \mathbf{u_2} & \dots & \mathbf{u_k}
\end{matrix}
\right.}& \underbrace{\left.\begin{matrix}\mathbf{u_{k+1}} & \dots & \mathbf{u_m}
\end{matrix}\right]}\\ \text{col} (A) & \text{null} (A^T)
\end{matrix}
\begin{bmatrix}
\sigma_1 & 0 & \dots & 0 & 0 & \dots & 0 \\
0 & \sigma_2 & \dots & 0 & 0 & \dots & 0 \\
\dots& & & & & \\
0 & 0 & \dots & \sigma_k & 0 & \dots & 0 \\
0 & 0 & \dots & 0 & 0 & \dots & 0 \\
\dots& & & & & \\
0 & 0 & \dots & 0 & 0 & \dots & 0
\end{bmatrix}
\begin{bmatrix} \mathbf{v_1^T} \\ \mathbf{v_2^T} \\ \dots \\ \mathbf{v_k^T} \\ \mathbf{v_{k+1}^T} \\ \dots \\ \mathbf{v_n^T} \end{bmatrix}
\begin{matrix} \left.\vphantom{\begin{bmatrix} \vect v_1^T \\ \vect v_2^T \\ \dots \\ \vect v_k^T \end{bmatrix}}\right\}\text{row} (A) \\ \left.\vphantom{\begin{bmatrix} \vect v_{k+1}^T \\ \dots \\ \vect v_n^T \end{bmatrix}}\right\}\text{null} (A)
\end{matrix}
$$
$V = [\mathbf{v}_1 \mathbf{v}_2 \dots \mathbf{v}_n]$은 $A^TA$의 직교대각화이다.
$\mathbf{u}_i = \cfrac{A\mathbf{v}_i}{\| A\mathbf{v}_i\|} = \cfrac{1}{\sigma_i}A\mathbf{v}_i$
$\{ \mathbf{u}_1, ..., \mathbf{u}_k \}$는 $\text{col}(A)$의 정규직교기저(orthonormal basis)이다.
Example 2. SVD if $A$ is not square
Example 1의 $A = \begin{bmatrix} 1 & 1 \\ 0 & 1 \\ 1 & 0 \end{bmatrix}$의 SVD를 구해보자.
SVD를 이용하면 데이터를 압축하여 적은 공간으로 저장하거나 빠르게 전송할 수 있다.
고윳값과 특이값은 $\lambda_1 = 3, \lambda_2 = 1, \sigma_1 = \sqrt{3}, \sigma_2 = 1$임을 예제1에서 구했다.
eigenvector를 구하면
$$\mathbf{v}_1 = \begin{bmatrix} \frac{\sqrt{2}}{2} \\ \frac{\sqrt{2}}{2} \end{bmatrix}, \mathbf{v}_2 = \begin{bmatrix} \frac{\sqrt{2}}{2} \\ -\frac{\sqrt{2}}{2} \end{bmatrix}$$
$U$의 열벡터를 구해보면
$$\mathbf{u}_1 = \cfrac{1}{\sigma_1}A\mathbf{v}_1 = \begin{bmatrix} \frac{\sqrt{6}}{3} \\ \frac{\sqrt{6}}{6} \\ \frac{\sqrt{6}}{6} \end{bmatrix}$$
$$\mathbf{u}_2 = \cfrac{1}{\sigma_2}A\mathbf{v}_2 = \begin{bmatrix} 0 \\ -\frac{\sqrt{2}}{2} \\ \frac{\sqrt{2}}{2} \end{bmatrix}$$
$\{ \mathbf{u}_1, \mathbf{u}_2 \}$을 $R^3$로 확장하자. 계산의 편의를 위해 적당한 상수배를 하여 $\sqrt{6}\mathbf{u}_1 = \begin{bmatrix} 2 & 1 & 1 \end{bmatrix}^T, \sqrt{2}\mathbf{u}_2 = \begin{bmatrix} 0 & -1 & 1 \end{bmatrix}^T$ 을 이용하여 $\mathbf{u}_3$을 구해보면
$$\begin{bmatrix} 2 & 1 & 1 \\ 0 & -1 & 1 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \end{bmatrix}$$
$$\begin{bmatrix} x_1 \\ x_2 \\ x_3 \end{bmatrix} = t \begin{bmatrix} -1 \\ 1 \\ 1 \end{bmatrix}$$
따라서 $\mathbf{u}_3 = \begin{bmatrix} -\frac{1}{\sqrt{3}} & \frac{1}{\sqrt{3}} & \frac{1}{\sqrt{3}} \end{bmatrix} ^T$ 이다.
따라서 $A$의 SVD는
$$A = U \Sigma V^T$$
이때,
$U = \begin{bmatrix} \frac{\sqrt{2}}{\sqrt{3}} & 0 & -\frac{1}{\sqrt{3}} \\ \frac{1}{\sqrt{6}} & -\frac{1}{\sqrt{2}} & \frac{1} {\sqrt{3}} \\ \frac{1}{\sqrt{6}} & \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{3}} \end{bmatrix}$, $\Sigma = \begin{bmatrix} \sqrt{3} & 0 \\ 0 & 1 \\ 0 & 0 \end{bmatrix}$, $V^T = \begin{bmatrix} \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} \\ \frac{1}{\sqrt{2}} & -\frac{1}{\sqrt{2}} \end{bmatrix}$
이다.
9.6 Data Compression Using Singular Value Decomposition
수식적으로, 0 행과 0열을 불필요하다. zero blocks를 제거한 형태를 reduced singular value decomposition이라 부른다. 이 책에서는 $A = U_1 \Sigma_1 V_1^T$로 표기하고
$$A = \begin{matrix} \begin{bmatrix} \mathbf{u}_1 & \mathbf{u}_2 & \dots \mathbf{u}_k \end{bmatrix} \begin{bmatrix} \sigma_1 & 0 & \dots & 0 \\ 0 & \sigma_2 & \dots & 0 \\ 0 & 0 & \dots & \sigma_k \end{bmatrix} \begin{bmatrix} \mathbf{v}_1^T \\ \mathbf{v}_2^T \\ \vdots \\ \mathbf{v}_k^T \end{bmatrix} \end{matrix} $$
이다.
예제 이미지는 아래 사이트에서 무료로 구할 수 있다.
Digital Image Processing
Digital Image Processing 2009/2010, Uni-Konstanz
cms.uni-konstanz.de
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
def compress_image(image_path):
image = Image.open(image_path).convert('L')
print('Image size:', image.size)
col, row = image.size
original_image_size = row * col
# PIL -> numpy.array
A = np.array(image, 'uint8')
# SVD
U, S, Vt = np.linalg.svd(A, full_matrices=False)
print('The shape of U, S, Vt:', U.shape, S.shape, Vt.shape)
ranks = [4, 10, 20, 50, 128]
plt.figure(figsize=(18, 10))
for i, rank in enumerate(ranks, 1):
# The @ operator can be used as a shorthand for np.matmul on ndarrays.
new_image = U[:, :rank] @ np.diag(S[:rank]) @ Vt[:rank, :]
compressed_image_size = rank * (row + col + 1)
ratio = compressed_image_size / original_image_size * 100
plt.subplot(1, len(ranks), i)
plt.imshow(new_image, cmap='gray')
plt.title(f'Rank: {rank}, Ratio: {ratio:.2f}')
new_image, _ = image_path.split('.')
new_image_path = f'{new_image}_compressed.png'
plt.savefig(f'{new_image_path}')
print(f'{new_image_path} saved successfully.\n')
if __name__ == '__main__':
compress_image('baboon.png')
compress_image('runner1.jpg')
RGB 이미지 압축
RGB 이미지 역시 각 색상별 2차원 행렬에 대하여 SVD를 수행하고 다시 3차원으로 합쳐 재구성하면 된다.
이때 rank가 작은 경우, 복원이 완전하지 않아 색상이 불완전하게 출력될 수 있다.
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
import sys
def load_image(image_path):
print('load image... ', end='')
img = Image.open(image_path).convert('RGB')
print('Done.')
return img
def image_to_SVD(image, path):
img = np.array(image)
img_r, img_g, img_b = img[:, :, 0], img[:, :, 1], img[:, :, 2]
print('Apply SVD... ', end='')
U_r, S_r, V_r = np.linalg.svd(img_r, full_matrices=False)
U_g, S_g, V_g = np.linalg.svd(img_g, full_matrices=False)
U_b, S_b, V_b = np.linalg.svd(img_b, full_matrices=False)
print('Done.')
print('##### Shape of USV from SVD(image) #####')
print('R:', U_r.shape, S_r.shape, V_r.shape)
print('G:', U_g.shape, S_g.shape, V_g.shape)
print('B:', U_b.shape, S_b.shape, V_b.shape)
print()
print('Plot 6 samples... ', end='')
comps = [3648, 1, 5, 10, 15, 20]
plt.figure(figsize=(12, 6))
for i in range(len(comps)):
tmp = np.zeros((img.shape[0], img.shape[1], 3))
tmp[:, :, 0] = U_r[:, :comps[i]] @ np.diag(S_r[:comps[i]]) @ V_r[:comps[i], :]
tmp[:, :, 1] = U_g[:, :comps[i]] @ np.diag(S_g[:comps[i]]) @ V_g[:comps[i], :]
tmp[:, :, 2] = U_b[:, :comps[i]] @ np.diag(S_b[:comps[i]]) @ V_b[:comps[i], :]
tmp = tmp.astype(np.uint8)
plt.subplot(2, 3, i+1)
plt.imshow(tmp)
if(i == 0):
plt.title(f'Actual Image with n_components = {comps[i]}')
else:
plt.title(f'n_components = {comps[i]}')
plt.savefig(path)
print('Done.')
if __name__ == '__main__':
if len(sys.argv) != 2:
print('please pass the image path')
exit(0)
print('#### Config #####')
print('image path: {}'.format(sys.argv[1]))
image = load_image(sys.argv[1])
image_to_SVD(image, 'result.png')
수정
2023.03.05) 3차원 RGB 이미지 압축도 소스 코드 추가. 예시는 모나리자
'스터디 > 선형대수학' 카테고리의 다른 글
| [선형대수] Vector, Matrix, Gaussian Elimination, LU-decomposition, LDU-decomposition, inverse (0) | 2023.02.20 |
|---|---|
| [선형대수학] LU-Decompositions (0) | 2023.02.07 |
