[Linear Analysis] 특이값(Singular Value)

특이값분해(Singular Value Decomposition)

내용

• 비정방행렬의 고유값분해
• 특이값

비정방행렬의 고유값분해

고유 분해, 스펙트럴 분해 등은 가역적인 정방 행렬을 대상으로 합니다. 이에 반해 특이값을 이용하는 특이값 분해(Singular Value Decomposition)는 비정방 행렬을 정방행렬로 변형하여 분해할 수 있는 방법으로 선형 대수의 계산에서 가장 유용하게 사용되고 있습니다.

비정방행렬은 행렬에 그 행렬의 전치행렬과의 행렬곱으로 정방행렬로 만들 수 있습니다.

import numpy as np
import numpy.linalg as la
from sympy import *

np.random.seed(1)
X=np.random.randint(0, 11, (3, 4))
print(X)

[[5 8 9 5]
 [0 0 1 7]
 [6 9 2 4]]

XTX=np.dot(X.T, X)
print(XTX)

[[ 61  94  57  49]
 [ 94 145  90  76]
 [ 57  90  86  60]
 [ 49  76  60  90]]

위 결과와 같이 전치행렬과의 행렬곱으로 생성한 정방행렬은 대칭행렬이 됩니다. 이 대칭행렬의 고유행렬은 다음과 같이 전치행렬과 역행렬이 같으므로 정규직교(Orthonormal)행렬이 됩니다. 이 결과는 정방행렬이면서 대칭행렬이므로 행렬의 대각화에 기반을 둔 분해(decomposition)가 가능합니다.

d, P=la.eig(XTX)
np.allclose(P.T, la.inv(P))

True

정방행렬 A의 고유값과 고유벡터의 관계와 고유벡터는 단위벡터인 점을 적용하면 식 1이 성립합니다. 즉, 행렬 A와 고유벡터(v)와의 내적의 크기는 고유벡터에 대응하는 고유값의 크기와 같아집니다.

$$\begin{array}{rl}A{v}_1& ={\lambda }_1{v}_1,\Vert v_1\Vert =1\\ \text{(식 1)}& \Vert A{v}_1\Vert & =\Vert {\lambda }_1{v}_1\Vert & =|{\lambda }_1|\Vert v_1\Vert \\ \\ & =|{\lambda }_1|\end{array}$$

위와 같은 정방 행렬의 특성은 비정방 행렬을 정방행렬로 변환 후 적용할 수 있습니다.

예 1)

ℝ³ 벡터 x가 행렬 A에 의해 ℝ²로 선형변환이 이루어진다면 ||Ax||가 최대가 되는 단위벡터 x를 계산해 봅니다.

$$A=\left[\begin{array}{ccc}4& 11& 14\\ 8& 7& -2\end{array}\right]$$

행렬 Ax의 최대크기를 결정하기 위해 고유값과 고유벡터를 사용할 수 있습니다. 그러나 A는 정방행렬이 아니므로 전치행렬과의 내적곱으로 정방행렬(A1)로 전환해야 합니다. 이 전환된 행렬 A1은 식 2와 같이 대칭행렬이 됩니다.

$$\begin{array}{rl}A1& =A^{T}A\\ \text{(식 2)}& (A^{T}A{)}^T& =A^T(A^T{)}^T& =A^{T}A\\ \therefore A{1}^T& =A1\end{array}$$

비정방행렬인‖Ax‖의 최대값은 식 3과 같이 결정할 수 있습니다.

$$\begin{array}{rl}\Vert Ax\Vert & =(Ax{)}^T(Ax)\\ \text{(식 3)}& & =x^T{A}^{T}Ax& =x^T(A^{T}A)x\end{array}$$

위 식의 x^T(A^TA)x는 이차 형식으로서 ||x||² = 1(x^Tx = 1)의 제한 조건에서의 A^TA의 최대값은 이 대칭 행렬의 가장 큰 고유값이 됩니다(제한된 최대 최소 참조). 그 고유값에 대응하는 단위 고유 벡터(v₁)와의 곱인 ||Av₁||은 최대가 됩니다.

A=np.array([[4,11,14], [8, 7, -2]])
A1=A.T@A
print(A1)

[[ 80 100  40]
 [100 170 140]
 [ 40 140 200]]

d, P=la.eig(A1)
print(d.round(3))

[360.   0.  90.]

대칭행렬의 고유벡터의 L₂ norm은 1입니다. 즉, 단위벡터입니다.

for i in range(3):
    print("고유벡터 v%d의 norm: %.1f" %(i, la.norm(P[:,i])))

고유벡터 v0의 norm: 1.0
고유벡터 v1의 norm: 1.0
고유벡터 v2의 norm: 1.0

A^TA의 고유값 중에 최대값에 대응하는 벡터는 첫번째 고유 벡터 입니다. 즉, 다음의 고유행렬 P의 0열에 해당합니다.

Av1=A@P[:,0]
print(Av1)

[-18.  -6.]

la.norm(Av1).round(3)

18.974

특이값

예 1로부터 m×n 형태의 행렬 A에 대한 A^TA는 대칭 행렬이므로 정규직교적으로 대각화가 가능합니다. 즉, 그 대칭행렬의 고유벡터(v_i)는 정규직교 벡터이므로 ℝⁿ 차원의 정규 직교 기저가 됩니다. λ₁, λ₂, …, λ_n을 A^TA의 고유값이라고 하면 식 4가 성립합니다. (1 ≤ i ≤ n)

$$\begin{array}{rl}{A}^{T}A{v}_i& ={\lambda }_i{v}_i\\ \Vert A{v}_i\Vert & =(A{v}_i{)}^T(A{v}_i)\\ \text{(식 4)}& & =v_i^T{A}^{T}A{v}_i& =v_I^T({\lambda }_i{v}_i)\\ & ={\lambda }_i{v}_i^T{v}_i\\ & ={\lambda }_i\\ v_i:& \text{정규직교벡터}\to v_i^T=v_i^{-1}\end{array}$$

식 4로부터 A^TA의 고유값은 음수가 아님을 나타냅니다 (식 5).

$$\tag{식 5} λ_1 \gt λ_2 \gt … λ_n \gt 0

A^TA의 고유값(λ_i)의 제곱근을 행렬 A의 특이값(singular value)이라고 하며 식 6과 같이 σ₁, σ₂, …, σ_n으로 나타냅니다.

$$\begin{array}{rl}{\sigma }_i& =\sqrt{{\lambda }_i}\\ \text{(식 6)}& & =\Vert A{v}_i\Vert 1\le & i<n\end{array}$$

예 1에서 A^TA의 최대 고유값은 360이며 이 고유값의 제곱근이 행렬 A의 최대 특이값이 됩니다.

np.sqrt(d[0]).round(3)

18.974

이 값은 ||Ax||의 최대값과 같습니다. 벡터 x는 A^TA의 최대 고유값에 대응하는 고유벡터입니다.

예 2)

다음 행렬 A의 특이값을 결정합니다.

$$A=\left[\begin{array}{ccc}3& 2& 2\\ 2& 3& -2\end{array}\right]$$

A=np.array([[3,2,2],[2,3,-2]])
A1=A.T@A
print(A1)

[[13 12  2]
 [12 13 -2]
 [ 2 -2  8]]

A^TA는 대칭행렬이지만 비가역행렬입니다. 그러나 다음 결과와 같이 정규직교행렬입니다.

la.det(A1).round(3)

0

비가역 행렬임은 피벗열이 아닌 벡터의 존재를 나타냅니다. 이 피벗열의 수는 기저벡터의 수와 같으며 행렬의 급수(rank)로 확인할 수 있습니다. 다음 결과와 같이 기저벡터의 수는 2개입니다.

la.matrix_rank(A1)

2

d,P=la.eig(A1)
np.allclose(P.T, la.inv(P))

True

print(d.round(3))

[25.  0.  9.]

고유값 0일 경우 식 7이 성립되지 않습니다. 즉, 0은 고유값이 아니므로 특이값은 2개입니다. 이는 행렬 A1의 급수와 같습니다.

$$\begin{array}{rl}\text{(식 7)}& A^{T}Av& =\lambda v\lambda ,v:& \text{고유값, 고유벡터}\end{array}$$

sigma=np.sqrt(d)
print(sigma.round(3))

[5. 0. 3.]

그러므로 0을 제외한 고유값 25, 9의 제곱근 5와 3이 행렬 A의 특이값이 됩니다.

이 특이값들의 수는 원시행렬의 열공간(기저공간)의 수와 같습니다. 열공간은 sympy객체.columnspace() 함수로 확인할 수 있습니다.

Matrix(A).columnspace()

[Matrix([
 [3],
 [2]]),
 Matrix([
 [2],
 [3]])]

특이값

비정방행렬 A의 A^TA는 대칭행렬입니다. 그러므로 0이 아닌 고유값들에 대응하는 고유벡터들은 정규직교기저입니다.

0이 아닌 고유값의 수가 행렬 A의 특이값의 수가 되며 이는 A의 열공간(Col A) 수 즉, 급수와 같습니다(식 8).

Rank A = 특이값의 수

(식 8)