대칭행렬 만들고 고윳값 분해하기

 

 

 

이번 글은

2022.06.29 - [수학의 재미/행렬 이론] - 대칭행렬을 분해합시다 - 고유값 분해(eigen decomposition) #2

대칭행렬을 분해합시다 - 고유값 분해(eigen decomposition) #2

에서 이어집니다. 

정방행렬 $\mathbf{A}$ 에 대해서 eigenvalue로 이루어진 대각 행렬 $\mathbf{D}$와 eigenvector로 이루어진 행렬 $\mathbb{Q}$ 가 존재해서

$$ \mathbf{A}\mathbf{Q} = \mathbf{Q}\mathbf{D} $$

의 형태로 쓸 수 있다고 했습니다. 

만일 여기에 $\mathbf{A}$가 실수 원소이자,  대칭 행렬이라는 조건까지 가미되면

 

• $\mathbf{D}$의 원소는 모두 실수
• eigen vector의 크기가 1이 되게끔 조절하면 $\mathbf{Q}$가 직교 행렬이 되어 역행렬이 존재함. 즉 $\mathbf{Q}^{-1} = \mathbf{Q}^T$ 

$$ \mathbf{A} = \mathbf{Q} \mathbf{D} \mathbf{Q}^T $$

로 아주 멋있게 분해가 된다고 했습니다. 

 

테트리스 조각분해랑은 좀 다름~

 

그럼 대칭행렬을 하나 만들어서 분석을 해 볼까요? 

대칭 행렬 만드는 방법은 여러 가지가 있으나 대표적인 두 가지 방법을 소개해 보도록 하겠습니다.

 

1. 임의의  정방 행렬 $\mathbf{A}$에 대해 $\mathbf{A}+\mathbf{A}^T$ 은 대칭 행렬이다.

 

행렬 $\mathbf{A}+\mathbf{A}^T$ 는 대칭 행렬이 됩니다. 왜냐하면

$$ (\mathbf{A}+\mathbf{A}^T)^T = \mathbf{A}^T+\mathbf{A}$$

이기 때문이죠. 

 

 

2. 임의의 정방행렬   $\mathbf{A}$에 대해 $\mathbf{A}\cdot \mathbf{A}^T$ 은 대칭 행렬이다.

$$ ( \mathbf{A}\cdot \mathbf{A}^T)^T =( \mathbf{A}^T)^T \cdot \mathbf{A}^T = \mathbf{A}\cdot \mathbf{A}^T $$

이기 때문입니다.

여기에는 다음의 Fact가 쓰입니다.

Fact.
곱하기가 가능한 size의 행렬 $\mathbf{A}, \mathbf{B}$ 에 대해,
$$ (\mathbf{AB})^T = \mathbf{B}^T \cdot \mathbf{A}^T $$
가 성립한다.

 

---

 

이제 python을 통하여 고윳값 분해를 어떻게 하는지 알아볼까요? 

 

import numpy as np
# import scipy.stats as ss
import numpy.linalg as LA


def make_symm_matrix_decompTest():
    np.set_printoptions(suppress=True, precision=3)

    tmpA = np.random.randn(5, 5)
    A = tmpA + tmpA.T

    print('symmetric matrix A is')
    print(A)
    print('-' * 30)

    print('Check A is symmetric')
    print('A =A^t is {}'.format(np.array_equal(A, A.T)))
    print('-' * 30)

    D, Q = LA.eig(A)

    print('\n')
    print('Result: Eigen Decomposition of A : AQ = QD')
    print('1. Eigen Valus are')
    print(D)
    print('-' * 30)

    print('2. Q is')
    print(Q)
    print('-' * 30)

    print('3. D is')
    print(np.diag(D))
    print('-' * 30)

    res1 = A @ Q
    res2 = Q @ np.diag(D)

    print('\n')
    print('Check decomposition')
    print('1. A*Q is')
    print('*' * 30)
    print(res1)
    print('2. Q*D is')
    print(res2)
    print('*' * 30)

    print('Eigen Decomposition is True/False within tolerance: {}'.format(np.allclose(res1, res2)))

    print('Eigen Decomposition is exactly same? : {}'.format(np.array_equal(res1, res2)))

    print('\n')
    print('AQ - QD is')
    np.set_printoptions(suppress=False)

    print(res1 - res2)

if __name__ == '__main__':
    make_symm_matrix_decompTest()

 

code를 간략히 살펴봅시다.

 

우선 python에서 eigen decomposition을 제공해주는 툴은 

numpy.linalg 모듈 하의 eig라는 함수

입니다. 이미 numpy는 많이 썼던 라이브러리라, 따로 인스톨하지 않고 numpy.lialg를 import 하여 씁니다.

 

 

    np.set_printoptions(suppress=True, precision=3)

○ 행렬을 print 하는데, e+ 형태의 숫자가 안 나오게끔 suppress를 True로, 또 소수 셋째 자리로 보이게끔 precision은 3으로 설정합니다.

 

 

 

    tmpA = np.random.randn(5, 5)
    A = tmpA + tmpA.T

○ 대칭 행렬을 만들기 위해 우선, 정규분포 난수를 원소로 가지는 5 by 5 행렬을 만듭니다.

○ 위의 방법 1을 써서 대칭행렬 A를 만듭니다.

 

 

 

    print('Check A is symmetric')
    print('A =A^t is {}'.format(np.array_equal(A, A.T)))

○  A가 대칭 행렬임을 보이려면 $A =A^T$ 를 체크하면 됩니다.

○  A.T 는 A의 transpose matrix를 뜻합니다. 

○ numpy.array_equal 은 두 array가 동일한 놈인지를 check 하는 함수입니다. 같은 배열이면 True, 다르면 False를 반환하죠.

 

 

    D, Q = LA.eig(A)

○ LA는 import 부분에서 numpy.linalg를 LA라는 별칭으로 부르기로 정의해 놨습니다. 즉, LA.eig는 numpy.linalg밑의 고윳값 분해를 도와주는 함수입니다.

○ eigen value를 배열인(아직 대각 행렬 아님!) D와 eigen vector로 이루어진 행렬 Q를 반환합니다.

 

 

 

    print('3. D is')
    print(np.diag(D))
    print('-' * 30)

○  D를  대각 행렬로 만들기 위해 np.diag라는 함수를 씁니다. np.diag를 하면, D에 있는 원소를 대각선으로 배열시켜 대각 행렬을 만들어 줍니다.

 

 

 

    res1 = A @ Q
    res2 = Q @ np.diag(D)

    print('\n')
    print('Check decomposition')
    print('1. A*Q is')
    print('*' * 30)
    print(res1)
    print('2. Q*D is')
    print(res2)
    print('*' * 30)

    print('Eigen Decomposition is True/False within tolerance: {}'.format(np.allclose(res1, res2)))

    print('Eigen Decomposition is exactly same? : {}'.format(np.array_equal(res1, res2)))

○  AQ = QD 가 만족하는 지를 보기 위해, AQ를 계산하여 res1,  QD를 계산하여 res2에 저장합니다.

○ 행렬의 곱을 수행하는 python operator는 @입니다.

○ np.allclose는 두 배열(행렬)이 허용오차 안에서 같은지, 다른지를 판단하는 함수입니다.

○ np.array_equal은 두 배열(행렬)이 정확하게 같은지 다른지를 판단하는 함수입니다.

 

 

 

    print('AQ - QD is')
    np.set_printoptions(suppress=False)

    print(res1 - res2)

○ res1과 res2의 차이를 출력합니다.

○ printoption에서 suppress를 false로 조정해서 과학적 표기(scientific notation)를 허용합니다. 여기를 보면 다음과 같이 설명하고 있네요.

suppress
If True, always print floating point numbers using fixed point notation, in which case numbers equal to zero in the current precision will print as zero. If False, then scientific notation is used when absolute value of the smallest number is < 1e-4 or the ratio of the maximum absolute value to the minimum is > 1e3. The default is False.

즉, 가장 작은 값의 크기가 1만 분의 1 보다 작거나, 최솟값에 대한 최대 절댓값의 비율이 1000배가 넘으면 scientific notation, 즉 e+ 로 표현합니다.

 

 

결과를 보실까요?

 

symmetric matrix A is
[[-3.573  0.977  0.927 -0.977  2.17 ]
 [ 0.977 -1.41  -2.149  1.037  0.947]
 [ 0.927 -2.149 -3.021 -0.872  1.241]
 [-0.977  1.037 -0.872 -4.193 -0.319]
 [ 2.17   0.947  1.241 -0.319  1.662]]
------------------------------
Check A is symmetric
A =A^t is True
------------------------------


Result: Eigen Decomposition of A : AQ = QD
1. Eigen Valus are
[ 3.003 -5.722 -3.489 -4.798  0.47 ]
------------------------------
2. Q is
[[ 0.359  0.636 -0.673 -0.116 -0.004]
 [ 0.154 -0.452 -0.298 -0.248 -0.789]
 [ 0.198 -0.407 -0.165 -0.683  0.549]
 [-0.09   0.473  0.513 -0.655 -0.276]
 [ 0.894 -0.04   0.41   0.174 -0.012]]
------------------------------
3. D is
[[ 3.003  0.     0.     0.     0.   ]
 [ 0.    -5.722  0.     0.     0.   ]
 [ 0.     0.    -3.489  0.     0.   ]
 [ 0.     0.     0.    -4.798  0.   ]
 [ 0.     0.     0.     0.     0.47 ]]
------------------------------


Check decomposition
1. A*Q is
******************************
[[ 1.079 -3.641  2.347  0.556 -0.002]
 [ 0.463  2.586  1.039  1.19  -0.37 ]
 [ 0.594  2.33   0.577  3.276  0.258]
 [-0.271 -2.704 -1.791  3.141 -0.13 ]
 [ 2.686  0.23  -1.43  -0.835 -0.006]]
2. Q*D is
[[ 1.079 -3.641  2.347  0.556 -0.002]
 [ 0.463  2.586  1.039  1.19  -0.37 ]
 [ 0.594  2.33   0.577  3.276  0.258]
 [-0.271 -2.704 -1.791  3.141 -0.13 ]
 [ 2.686  0.23  -1.43  -0.835 -0.006]]
******************************
Eigen Decomposition is True/False within tolerance: True
Eigen Decomposition is exactly same? : False


AQ - QD is
[[-1.998e-15 -4.441e-15 -1.332e-15 -4.441e-16 -1.698e-16]
 [ 4.996e-16  1.332e-15 -3.331e-15  2.220e-16  0.000e+00]
 [ 1.332e-15  2.665e-15 -2.554e-15 -2.665e-15  5.551e-16]
 [-2.220e-16 -1.332e-15  4.663e-15 -1.332e-15  4.441e-16]
 [ 0.000e+00  4.663e-15  2.442e-15  6.661e-16 -7.338e-16]]

Process finished with exit code 0

○  A는 대칭 행렬이 맞습니다.

○  1. Eigen Valus are [ 3.003 -5.722 -3.489 -4.798  0.47 ] 로서 모두 실수가 나왔네요.

○ AQ와 QD는 정확히 일치하진 않지만, 아주 극미의 차이만 보입니다. 이론적으로는 일치하는 숫자이지만, 파이썬 계산 과정에서 지지고 볶고 하다가 조금씩 missing 되는 숫자가 있는 것으로 보입니다. AQ-QD의 원소들을 suppress = False로 과학적 표기해 보면 차이가 e-15, 즉 소수점 15자리까지는 가야 그 차이를 보이죠.