Cholesky Decomposition

조건

• symmetric matrix : 대칭행렬
• positive definite matrix : 양의 정부호 행렬 - 즉 모든 고윳값이 양수인 에르미트 행렬

방법

lower triangular matrix 란 행렬의 diagonal 축을 기준으로 아래 부분에만 행렬의 원소가 0이 아닌 값을 가질 수 있는 경우를 말합니다.
ex) $A=\begin{bmatrix} 2 & 0 & 0\\ 3 & 4 & 0\\ 6 & 7 & 8 \end{bmatrix}$

이때 어떤 행렬 $A$가 위의 두 조건인 1) symmetric matrix, 2) positive definite matrix 조건을 만족한다면 다음과 같이 표현할 수 있습니다.
$A=LL^T$

이때 $L$을 다음의 방법으로 구하는 것을 cholesky decomposition이라고 합니다.

$L_{j,j}=\sqrt{A_{j,j}-\Sigma_{k=1}^{j-1}L^2_{j,k}}$
$L_{i,j}={1\over L_{j,j}}(A_{i,j}-\Sigma_{k=1}^{j-1}L_{i,k}L_{j,k})\space for\space i>j$

장점

• LU decomposition 보다 속도가 훨씬 빠르다.
• 연산이 비교적 간단하여 다른 방법들에 비해 누적오차가 적다.
• 시간복잡도가 $A$가 nxn 행렬일 때 $O(n^3)$로 행렬 연산치고 준수하다.
• 구현이 쉽다.

단점

• 더 시간복잡도가 적은 방법들이 존재한다.

코드

from math import sqrt
def cholesky(A):
    n = len(A)

    # initialize L
    L = [[0.0] * n for i in range(n)]

    # Cholesky Decomposition
    for j in range(n):
        for i in range(j+1):
            tmp_sum = sum(L[j][k] * L[i][k] for k in range(i))
            if (i == j): # Diagonal elements
                L[j][i] = sqrt(A[j][j] - tmp_sum)
            else:
                L[j][i] = (1.0 / L[i][i] * (A[j][i] - tmp_sum))
    return L