행렬로 바라보는 피보나치 수열
다음의 책을 보니 재미있는 내용이 있습니다. ^^
프로그래머를 위한 선형대수
; http://www.yes24.com/24/goods/39446808
(평을 보시면 아시겠지만, 저 역시 추천하고 싶은 책입니다. ^^)
249페이지에 보면 "자기회귀모델(AR: AutoRegressive)"의 이산시간에 대한 예로,
오늘의 ζ(t)는 어제의 ζ(t - 1), 이틀 전의 ζ(t - 2), 사흘 전의 ζ(t - 3)과 오늘의 u(t)에 따라 다음과 같이 정해진다.
ζ(t) = -0.5ζ(t - 1) + 0.34ζ(t - 2) + 0.08ζ(t - 3) + 2u(t)
초기 조건 ζ(0) = 0.78, ζ(-1) = 0.8, ζ(-2) = 1.5
소개가 되면서 다음과 같이 행렬 표현을 합니다.
저걸 보니, 피보나치 수열이 생각났습니다.
황금비율 증명 - 피보나치 수와 연분수의 관계
; https://www.sysnet.pe.kr/2/0/1312
역시 초깃값이 주어지고 x(t)는 x(t - 1)에 의해 결정되니까요. 따라서 위와 같은 기준으로 피보나치 수열을 바라보면 다음과 같이 정리가 됩니다.
ζ(t) = 1ζ(t - 1) + 1ζ(t - 2)
초기 조건 ζ(0) = 1, ζ(-1) = 0
간단하게 t = 1 ~ 4까지 테스트하면 이렇게 되고,
t = 1일 때, ζ(1) = ζ(1 - 1) + ζ(1 - 2) = ζ(0) + ζ(-1) = 1 + 0 = 1
t = 2일 때, ζ(2) = ζ(2 - 1) + ζ(2 - 2) = ζ(1) + ζ(0) = 1 + 1 = 2
t = 3일 때, ζ(3) = ζ(3 - 1) + ζ(3 - 2) = ζ(2) + ζ(1) = 2 + 1 = 3
t = 4일 때, ζ(4) = ζ(4 - 1) + ζ(4 - 2) = ζ(3) + ζ(2) = 3 + 2 = 5
이를 행렬로 표현하면 다음과 같습니다.
마찬가지로 t = 1 ~ 4까지에 대해 행렬로 계산하면 이렇게 됩니다.
따라서 (초깃값 2개를 넘어) n 번째 피보나치 수열은,
보는 바와 같이 행렬 [1 1; 1 0]에 대해 n 승을 하고 그 값을 [1 0] 행렬에 곱하면 n 번째 피보나치 수열이 구해지는 것입니다. 실제로 octave 같은 도구를 이용해 다음과 같이 행렬 계산을 바로 해볼 수 있습니다.
function fib_1()
a = [1 1; 1 0]
b = [1;0]
a ^ 1 * b
a ^ 2 * b
a ^ 3 * b
a ^ 4 * b
a ^ 5 * b
endfunction
위의 함수를 실행하면 2*1 행렬이 5개가 출력되는 데 그것의 첫 번째 원소들을 보면 1, 2, 3, 5, 8로 피보나치 수열이 나옵니다.
행렬로 표현된 피보나치 계산에서 고윳값/고유벡터를 이용해 풀어보면 재미있는 결과가 나옵니다.
[선형대수학 #3] 고유값과 고유벡터 (eigenvalue & eigenvector)
; http://darkpgmr.tistory.com/105
행렬 [1 1; 1 0]에 대한 고윳값, 고유벡터를 계산해 보면,
위의 행렬식을 구하면,
= (1 - λ)(0 - λ) - 1
= λ2 - λ -1
위와 같이 구한 특성 다항식을 특성 방정식에 따라 0 값이 나오는 해를 구하면,
det(A - λ E) = 0
λ2 - λ -1 = 0
근의 공식에 따라,
와 같이 계산됩니다. 고윳값을 구했으니 고유벡터까지 구해볼까요? ^^
연립 방정식으로 풀으면,
(1 - λ)vx + vy = 0
vx - λvy = 0
vx = λvy
따라서, vx가 vy의 λ배로 이뤄진 무수히 많은 벡터 = [λt, t]
그럼 고유 벡터를 아무거나 다음과 같이 선정할 수 있습니다.
따라서 고윳값 λ의 2가지 값에 대해,
이 중에서 고유 벡터를 [(1 + sqrt(5)) / 2, 1]인 쌍으로 골라 보겠습니다. 이를 다시 Gram-Schmidt 정규 직교로 바꾸면,
Matlab/Octave로 Gram-Schmidt 정규 직교 집합 구하는 방법
; https://www.sysnet.pe.kr/2/0/11235
(0.52573, -0.85065), (-0.85065, -0.52573)로 구할 수 있습니다. 즉, 이 2개의 벡터 각각에 대응하는 λ배의 모든 벡터들이 고유 벡터들이 됩니다.
실제로 위의 과정들을 간단하게 octave로 구할 수 있습니다.
a = [1 1; 1 0]
[ev, ei] = eig(a)
ev =
0.52573 -0.85065
-0.85065 -0.52573
ei =
Diagonal Matrix
-0.61803 0
0 1.61803
또한, Av = λv인 것도 다음과 같이 쉽게 계산해볼 수 있습니다.
a * [0.52573, -0.85065]'
ans =
-0.32492
0.52573
-0.61803 * [0.52573 -0.85065]'
ans =
-0.32492
0.52573
피보나치 수열의 고윳값과 고유벡터를 구했으니 n 번째 값을 구하는 방법에 대해 행렬의 성질로 다시 살펴보겠습니다.
"
[선형대수학 #3] 고유값과 고유벡터 (eigenvalue & eigenvector)" 글에 보면 다음과 같은 공식이 나옵니다.
A = 행렬
P = 행렬 A의 고유벡터들을 열벡터로 하는 행렬
Λ = 교윳값들을 대각 원소로 하는 대각 행렬
AP = PΛ
A = PΛP-1
이를 기반으로 A의 n 승을 다음과 같이 쉽게 구할 수 있는 방법을 포함하고 있습니다.
Ak = (PΛP-1)k
= (PΛP-1)(PΛP-1)......(PΛP-1)
= PΛkP-1
= Pdiag(λk1,......,λkn)P-1
따라서, 가령 5번째 피보나치 수를 구하고 싶다면 고유 벡터와 그것의 역행렬만 구한 후 고윳값 2개를 대각 행렬로 갖는 것만 5 승을 해주면 되는 것입니다. 이것을 octave로 다음과 같이 테스트할 수 있습니다.
ev * ei ^ 5 * inverse(ev)
ans =
8.0000 5.0000
5.0000 3.0000
즉, 고윳값을 알기 전에는 다음과 같은 행렬 계산이었지만,
고윳값을 알게 된 이상, 그것은 대각행렬의 n 승으로 바뀌었기 때문에 단순히 스칼라 값인 고윳값 2개만 n 승을 해주면 되는 문제로 바뀐 것입니다.
[이 글에 대해서 여러분들과 의견을 공유하고 싶습니다. 틀리거나 미흡한 부분 또는 의문 사항이 있으시면 언제든 댓글 남겨주십시오.]