O que é a decomposição de Cholesky?

Question

Explique e desenvolva com todos os detalhes. Dê um exemplo em Python como gerar uma matriz de covariância para retornos normais com uma matriz de covariância pre-especificada.

Answer 1

Una matriz simétrica definida positiva pode ser descomposta como o produto de uma matriz triangular inferior e a transposta da matriz triangular inferior. Uma condição necessária e suficiente para que uma matriz \(A\) qualquer admita feitorização de Cholesky é que seja simétrica e definida positiva. Então, \(A\) pode ser decomposta como
\[A=LL^T,\]
onde \(L\) é una matriz triangular inferior com entradas diagonais estritamente positivas e \(L^T\) representa la conjugada transposta de \(L\).

Algoritmo para a decomposição de Cholesky

Existem vários métodos para gerar a decomposição de Cholesky, uns deles é o seguinte.

Particionar a matriz \(A=LL^T\) da seguinte forma
\[\left[\begin{array}{cc} a_{11} & A_{21}^T \\ A_{21} & A_{22} \\ \end{array}\right] = \left[\begin{array}{cc} l_{11} & 0 \\ L_{21}^T & L_{22}^T \\ \end{array}\right] \left[\begin{array}{cc} l_{11} & L_{21}^T \\ 0 & L_{22}^T \\ \end{array}\right] = \left[\begin{array}{cc} l_{11}^2 & l_{11}L_{21}^T \\ l_{11}L_{21} & L_{21} L_{21}^T + L_{22} L_{22}^T\\ \end{array}\right] \]

1) Determine \(l_{11}\) e \(L_{21}\):
\[l_{11} = \sqrt{a_{11}}, \ L_{21}= \frac{A_{21}}{l_{11}}. \]

2) A partir da igualdade acima,
\[ A_{22} - L_{21} L_{21}^T = L_{22} L_{22}^T, \]
determine \(L_{22} \).

Note que esta é a decomposição de Cholesky de ordem \(n-1\); o tempo de execução do algoritmo é \(\mathcal{O}(n^3)\).

Vamos demonstrar que este algoritmo funciona com matrizes \(A\) positiva definida de ordem \(n\).
Se a é positiva definida, logo

1) \( \ a_{11}>0\)

2) \( \ A_{22} - L_{21} L_{21}^T = A_{22} - \frac{1}{a_{11}} A_{21} A_{2_1}^T \) é positiva definida pelo complemento de Schur.

Se o algoritmo funciona para \(n=m\), funciona para \(n=m-1\); claramente, funciona para \(n=1\).

O complemento de Schur estabelece que se \(A\) é positiva definida, então o complemento de Schur \(S = A_{22} - \frac{1}{a_{11}} A_{21} A_{21}^T\) é positiva definida, onde
\[ A = \left[\begin{array}{cc} a_{11}^2 & A_{21}^T \\ A_{21} & A_{22} ^T\\ \end{array}\right] \]

Exemplo

Vamos utilizar a decomposição de Cholesky para criar uma distribuição normal multivariada definido por um vetor de médias \(\mu\) e uma matriz de covariância \(\Sigma\). Obter \( L\) a partir de \(\Sigma\) usando algum algoritmo de decomposição de Cholesky. Para gerar amostras de esta distribuição \(K\)-variada, geramos uma amostra de k observações independentes a partir de uma distribuição normal padrão. Suponha que armazenamos a amostra num vetor \(z\), logo o resultado desejado será
\[ x = \mu + Cz.\]

Por que o vetor \(x\) tem correlação \(\Sigma\)?
Observe que,
\[ E[(x-u)(x-u)^T] = E[ Cz (Cz)^T] = CE[zz^T]C^T = CIC^T = \Sigma,\]
a penúltima igualdade ocorre por que as observações no vetor \(z\) são independentes; por tanto, a covariância é zero, e a variância (autocovariância) é um.

Algoritmo

O exemplo vai ser implementado em Python.

import numpy as np
from scipy.linalg import cholesky
from scipy.stats import norm

tamano_amostra = 10000

# matriz de covariancia desejada
Sigma = np.array([
        [ 1.00,  0.75,  0.25],
        [ 0.75,  1.00,  0.50],
        [ 0.25,  0.50,  1.00]
    ])

#vetor de medias desejado
u =np.array([[1.],[2.],[3.] ])



#gerar um amostragem a partir de uma variavel normalmente distribuid
#padronizada(media = 0, desvio padrao = 1)
z = norm.rvs(size=(3, tamano_amostra))


# calcular a decomposicao de Cholesky
c = cholesky(Sigma, lower=True)

#convertir a amostra numa distribuicao correlacionada segundo a definicao de sigma
#e media segundo o vetor u
x = np.dot(c, z)
x += u

print('covariancia')
print(np.cov(x))
print('media')
print(np.mean(x, axis=1))

Saida

covariancia
[[ 1.02262798  0.7712547   0.25287991]
 [ 0.7712547   1.02521253  0.51647226]
 [ 0.25287991  0.51647226  1.01653229]]
media
[ 0.99318368  1.99785094  3.01186588]

Comments

O método de decomposição de Cholesky, foi muito bem especificado, tanto teoricamente quanto computacionalmente pelo exemplo apresentado nesse post. Apenas a título de complementação, acredito ser enriquecedor salientar alguns pontos históricos e aplicações uteis dessa decomposição:
1)    Seu nome se deve a uma homenagem a André-Louis Cholesky que estabeleceu que uma matriz simétrica e positiva definida pode ser decomposta em uma matriz triangular inferior e sua transposta;
2)    A solução numérica de sistemas lineares de grande porte, em geral envolve a decomposição de Cholesky, pois esse método minimiza problemas de instabilidade numérica causada por álgebras matriciais utilizadas para encontrar matrizes inversas;
3)    Em contextos de “big data”, um conjunto cada vez maior de variáveis em projeções lineares é utilizado com frequência, a compreensão desse novo conjunto de dados envolve sempre o tratamento de matrizes de covariância, o que por sua vez utilizará a decomposição de Cholesky;
4)    No contexto de modelos de regressão linear, uma hipótese fundamental para finalidade de inferência, estabelece que a matriz de covariância seja esférica, ou seja, homocedástica e não autocorrelacionada. Quando essa hipótese é violada, a decomposição de Cholesky permite a utilização do Método de mínimos quadrados generalizados;

---

André, muito interessantes as observações. Obrigado!