기본 콘텐츠로 건너뛰기

[ML] 결정트리(Decision Tree) 모델

[Linear Algebra] 대칭행렬의 대각화

대칭행렬의 대각화

대칭행렬(Symmetric matrix)은 주 대각원소를 기준으로 위 부분과 아래부분의 원소들이 대칭되는 구조이며 대표적인 특성으로 식 1에 나타낸 것과 같이 원시행렬과 전치행렬은 같습니다.

$$\tag{식 1}A = A^T$$

식 1은 모든 대칭행렬에 대한 특성이지만 고유값 분해의 조건을 만족하는 가역적 대칭행렬의 경우는 다음의 특성을 가집니다.

  • 서로 다른 고유값에 대응하는 고유벡터들은 직교(orthogonal) 관계에 있습니다
  • 정규직교(orthonormal) 관계로 변형할 수 있습니다.
  • 식 2와 같이 가역적 대칭행렬을 구성하는 각 열벡터($x_1,\; x_2$)의 내적은 0가 됩니다.

$$\tag{식 2} x_1^Tx_2=0$$

식 2의 증명은 식 3과 같이 λ1x1Tx2를 정리하는 것으로 시작합니다.

\begin{align}\lambda_1 x_1^Tx_2&= (\lambda_1 x_1)^Tx_2\\ & = (Ax_1)^Tx_2\\ \tag{식 3} & = x_1^TA^Tx_2 \\& = x_1^TAx_2\\ & = x_1^T\lambda_2x_2\\& = \lambda_2 x_1^Tx_2\\ \because&\; A=A^T\end{align}

식 3을 정리하면 식 4와 같이 나타낼 수 있습니다.

\begin{align}\lambda_1 x_1^Tx_2& = \lambda_2 x_1^Tx_2\\\tag{식 4} \Rightarrow & \lambda_1 x_1^Tx_2 - \lambda_2 x_1^Tx_2=0\\ \Rightarrow& (\lambda_1 -\lambda_2)x_1^Tx_2=0 \end{align}

식 5.3.4에서 고유값 λ1과 λ2이 서로 다르기 때문에 고유벡터들 x1과 x2의 내적이 0 이어야 합니다. 이것은 두 벡터의 내적이 0인 것으로 그 고유벡터들은 서로 직교 관계에 있음을 의미합니다. 이 고유벡터들을 모두 정규화시키면 즉, 단위벡터화하면 정규 직교 관계가 됩니다.

가역적 대칭행렬인 경우 위 식 4를 성립합니다. 즉, 고유행렬은 정규직교행렬이 됩니다. 그러나 비가역 대칭행렬의 경우 고유벡터의 정규직교성을 나타내기도 합니다. 그러므로 행렬의 정규직교성은 식 5의 확인으로 결정합니다( 정규직교 참조).

$$\tag{식 5}U^TU=I \Leftrightarrow U^T=U^{-1}$$

예 1)

벡터 A의 고유벡터들의 직교성을 결정하여 봅니다.

$$A=\begin{bmatrix}6& -2 & -1\\ -2 & 6 & 1\\ -1& -1 & 5 \end{bmatrix}$$

A=np.array([[6,-2,-1], [-2,6,-1], [-1,-1,5]])
d, P=la.eig(A)
print(P.round(2))
[[ 0.58  0.71 -0.41]
 [ 0.58 -0.71 -0.41]
 [ 0.58  0.    0.82]]

고유행렬의 각 열벡터들사이의 내적을 계산합니다.

for i, j in itertools.combinations(range(3), 2):
    print(round(P[:, i-1]@P[:, j-1],3))
-0.0
-0.0
-0.0

위 결과는 고유벡터들은 서로 직교함을 나타냅니다. 이것은 행렬 A가 대칭행렬로 가역적임을 나타냅니다.

선형독립인 대칭행렬 A는식 6과 같이 대각화가 가능하므로 고유값분해가 성립됩니다.

\begin{align}\tag{식 6}A&=PDP^{-1}\\P:&\;\text{고유행렬}\\D:&\;\text{고유값으로 구성된 대각행렬} \end{align}

다시 말하면 A가 가역적이고 대칭행렬일 경우 고유벡터 행렬은 직교행렬이 되며 직교행렬을 정규직교행렬로 전환할 수 있으므로 이 행렬은 전치행렬과 역행렬이 같습니다. 그러므로 식 7이 성립합니다.

\begin{align}A&=PDP^{-1}\\\tag{식 7} & = PDP^T\\\because &\; P^TP=I \Rightarrow P^T =P^{-1} \end{align}

식 7에서 D는 서로 다른 고유값을 대각 요소로하는 대각행렬, P는 각 고유값에 대응하는 고유행렬입니다. 행렬 A에 대응하는 대각행렬은 다음과 같습니다.

D=np.diag(d)
print(D)
[[3. 0. 0.]
 [0. 8. 0.]
 [0. 0. 6.]]

다음 결과와 같이 A의 고유행렬 P의 각 벡터의 크기는 1이므로 행렬 P는 정규직교행렬입니다. 이 행렬의 경우 전치행렬과 역행렬은 같습니다.

for i in range(P.shape[1]):
    print(F"Norm of {'P'+str(i)}:{la.norm(P[:,i]).round(3)}")
Norm of P0:1.0
Norm of P1:1.0
Norm of P2:1.0
print(np.isclose(la.inv(P), P.T))
[[ True  True  True]
 [ True  True  True]
 [ True  True  True]]

예 2)

다음 행렬 B의 대각화 과정에서의 p는 직교행렬인지를 확인해 봅니다.

$$B=\begin{bmatrix}3& -2& 4\\-2& 6& 2\\4& 2& 3 \end{bmatrix}$$

B=np.array([[3,-2,4],[-2,6,2],[4,2,3]])
la.det(B).round(3)
-98.0

대칭행렬 B는 가역적입니다.

d,P=la.eig(B)
print(d)
[ 7. -2.  7.]

위 결과와 같이 두개의 고유값이 같으므로 행렬 B의 열차원과 고유값의 수에서 차이가 납니다. 이는 가역적 대칭행렬의 특성(식 2)에 위배되므로 식 7은 성립하지 않습니다.

D=np.diag(d)
np.allclose(B, P@D@la.inv(P))
True
np.allclose(B, P@D@P.T)
False

즉, 대각화가 이루어지지 않으며 고유행렬은 직교행렬이 아닙니다.

for i, j in itertools.combinations(range(3), 2):
    print("%d와 %d 열벡터의 내적: %.3f" %(i+1, j+1, P[:, i]@P[:,j]))
1와 2 열벡터의 내적: -0.000
1와 3 열벡터의 내적: -0.111
2와 3 열벡터의 내적: 0.000

위 결과와 같이 고유벡터 P1과 P3는 직교관계가 성립하지 않습니다. 그러므로 PT = P-1의 관계가 성립하지 않습니다. 이 관계를 성립시키기 위해 P[:,0]과 P[:,2] 벡터 사이의 직교벡터들로 고유행렬을 대체할 수 있습니다.

Gram-Schmit 과정으로 그 직교벡터들을 계산할 수 있지만 이 결과는 QR 분해의 Q와 같습니다. linalg.qr()함수를 적용하여 정규직교벡터들을 계산하면 다음과 같습니다.

Q,R=la.qr(P[:,[0,2]])
print(Q.round(3))
[[-0.745  0.   ]
 [ 0.298 -0.894]
 [-0.596 -0.447]]

위 결과들을 P[:,0]과 P[:,1]을 대체하여 새로운 고유행렬(P2)을 생성합니다.

P2=np.c_[Q[:,1],P[:,1], Q[:,0]]
print(P2.round(3))
[[ 0.    -0.667 -0.745]
 [-0.894 -0.333  0.298]
 [-0.447  0.667 -0.596]]

행렬 P2는 대칭행렬의 직교성을 확보합니다. 즉, P2T = P2-1관계가 성립합니다.

np.allclose(P2.T, la.inv(P2))
True

그러므로 식 5.3.5가 성립합니다.

np.allclose(B, P2@D@la.inv(P2))
True
np.allclose(B, P2@D@P2.T)
True

댓글

이 블로그의 인기 게시물

[Linear Algebra] 유사변환(Similarity transformation)

유사변환(Similarity transformation) n×n 차원의 정방 행렬 A, B 그리고 가역 행렬 P 사이에 식 1의 관계가 성립하면 행렬 A와 B는 유사행렬(similarity matrix)이 되며 행렬 A를 가역행렬 P와 B로 분해하는 것을 유사 변환(similarity transformation) 이라고 합니다. $$\tag{1} A = PBP^{-1} \Leftrightarrow P^{-1}AP = B $$ 식 2는 식 1의 양변에 B의 고유값을 고려한 것입니다. \begin{align}\tag{식 2} B - \lambda I &= P^{-1}AP – \lambda P^{-1}P\\ &= P^{-1}(AP – \lambda P)\\ &= P^{-1}(A - \lambda I)P \end{align} 식 2의 행렬식은 식 3과 같이 정리됩니다. \begin{align} &\begin{aligned}\textsf{det}(B - \lambda I ) & = \textsf{det}(P^{-1}(AP – \lambda P))\\ &= \textsf{det}(P^{-1}) \textsf{det}((A – \lambda I)) \textsf{det}(P)\\ &= \textsf{det}(P^{-1}) \textsf{det}(P) \textsf{det}((A – \lambda I))\\ &= \textsf{det}(A – \lambda I)\end{aligned}\\ &\begin{aligned}\because \; \textsf{det}(P^{-1}) \textsf{det}(P) &= \textsf{det}(P^{-1}P)\\ &= \textsf{det}(I)\end{aligned}\end{align} 유사행렬의 특성 유사행렬인 두 정방행렬 A와 B는 'A ~ B' 와 같

[matplotlib] 히스토그램(Histogram)

히스토그램(Histogram) 히스토그램은 확률분포의 그래픽적인 표현이며 막대그래프의 종류입니다. 이 그래프가 확률분포와 관계가 있으므로 통계적 요소를 나타내기 위해 많이 사용됩니다. plt.hist(X, bins=10)함수를 사용합니다. x=np.random.randn(1000) plt.hist(x, 10) plt.show() 위 그래프의 y축은 각 구간에 해당하는 갯수이다. 빈도수 대신 확률밀도를 나타내기 위해서는 위 함수의 매개변수 normed=True로 조정하여 나타낼 수 있다. 또한 매개변수 bins의 인수를 숫자로 전달할 수 있지만 리스트 객체로 지정할 수 있다. 막대그래프의 경우와 마찬가지로 각 막대의 폭은 매개변수 width에 의해 조정된다. y=np.linspace(min(x)-1, max(x)+1, 10) y array([-4.48810153, -3.54351935, -2.59893717, -1.65435499, -0.70977282, 0.23480936, 1.17939154, 2.12397372, 3.0685559 , 4.01313807]) plt.hist(x, y, normed=True) plt.show()

R 미분과 적분

내용 expression 미분 2차 미분 mosaic를 사용한 미분 적분 미분과 적분 R에서의 미분과 적분 함수는 expression()함수에 의해 생성된 표현식을 대상으로 합니다. expression expression(문자, 또는 식) 이 표현식의 평가는 eval() 함수에 의해 실행됩니다. > ex1<-expression(1+0:9) > ex1 expression(1 + 0:9) > eval(ex1) [1] 1 2 3 4 5 6 7 8 9 10 > ex2<-expression(u, 2, u+0:9) > ex2 expression(u, 2, u + 0:9) > ex2[1] expression(u) > ex2[2] expression(2) > ex2[3] expression(u + 0:9) > u<-0.9 > eval(ex2[3]) [1] 0.9 1.9 2.9 3.9 4.9 5.9 6.9 7.9 8.9 9.9 미분 D(표현식, 미분 변수) 함수로 미분을 실행합니다. 이 함수의 표현식은 expression() 함수로 생성된 객체이며 미분 변수는 다음 식의 분모의 변수를 의미합니다. $$\frac{d}{d \text{변수}}\text{표현식}$$ 이 함수는 어떤 함수의 미분의 결과를 표현식으로 반환합니다. > D(expression(2*x^3), "x") 2 * (3 * x^2) > eq<-expression(log(x)) > eq expression(log(x)) > D(eq, "x") 1/x > eq2<-expression(a/(1+b*exp(-d*x))); eq2 expression(a/(1 + b * exp(-d * x))) > D(eq2, "x") a * (b * (exp(-d * x) * d))/(1 + b