Decomposição de matrizes

Decomposição por Valores Singulares (SVD)

Foi mostrado e continuaremos mostrando a importância dos Autovalores e Autovetores, entretanto esse conceito só se aplica para matrizes quadradas, o que é um impedimento já que a grande maioria das matrizes reais não são quadradas.

Nem toda matriz possui vetores entrada que continuam na mesma direção na saída, nem toda matriz possui um conjunto ortonormal de entrada que produz um conjunto ortonormal na saída, mas será que toda matriz possui um conjunto ortonormal de entrada que produz um conjunto ortogonal na saída?

Tomando essa hipótese como verdadeira, então temos:

\[ AV = U\Sigma \]

Onde V é o conjunto ortonormal de entrada, U é conjunto ortonormal da saída e \(\Sigma\) é a matriz diagonal que escala cada coluna com uma escalar \(\sigma_n\).

Assim, toda matriz poderia ser escrita com esses termos:

\[A = U\Sigma V^T\]

Porém como garantir que U e V existem? Vamos tentar dividir o problema em duas equações com apenas U ou apenas V, uma ideia seria:

Relembre

Transposta

\[A^TA = V \Sigma U^T U \Sigma V^T\ = V \Sigma^2 V^T\]

\[AA^T = U \Sigma V^T V \Sigma U^T = U \Sigma^2 U^T\]

Note como o resultado parece a Decomposição Espectral que é garantido por essas matrizes, como \(A^TA\) e \(AA^T\) são matrizes de Gram garantimos que essa decomposição é possível.

Relembre

matriz de Gram

Então, garantimos que para toda matriz existe um conjunto entrada que é definido pelos Autovetores de \(A^TA\) que chamamos de Vetores Singulares a Direita, além de que garantimos que a saída para esse conjunto será ortogonal e será os Autovetores de \(AA^T\) que chamamos de Vetores Singulares a Esquerda.

As matrizes \(A^TA\) e \(AA^T\) também são Similares, então possuem os mesmos Autovalores, logo chamamos a raiz quadrada desses Autovetores de Valores Singulares. Cada Valor Singular está associado a 2 Vetores Singulares.

Geometricamente toda matriz pode ser vista como uma composição de rotação, multiplicação por escalar e rotação, como mestrado no vídeo:

Dica

Compute o SVD

Norma matricial

Como mostrado no vídeo, toda matriz pode ser vista como uma transformação de uma esfera para um elipsoide. Se definirmos essa esfera como todos os vetores unitários, então o ponto(vetor) mais distante do elipsoide até a origem será o vetor que for escalado pelo maior Valor Singular, assim como o ponto menos distante com o menor Valor Singular.

Essa informação é usada como uma forma de medir o "tamanho" da matriz, definimos a Norma Euclidiana da matriz como:

\[||A|| = max_{||\mathbf{x}||=1}\{||A\mathbf{x}||\}\]

Que pode ser lido como: a norma da matriz é a maior norma de um vetor saída, dado que o vetor entrada tinha norma 1.

E como vimos:

\[||A|| = \sigma_{max}\]

E como a inversa tem reverter esse elipsoide para uma esfera, então o vetor entrada que for colinear com o vetor de menor distância do elipsoide terá o maior escalamento:

\[||A^{-1}|| = \frac{1}{\sigma_{min}}\]

Forma truncada

Relembre

Base para os Subespaços com Autovetores na matriz de Gram

Note que o autovalor associado aos autovetores do Espaço nulo e Espaço Nulo Esquerdo é 0, ao colocarmos \(\sqrt{0}\) na diagonal de \(\Sigma\) acabamos por anular os Vetores Singulares associados.

Então, na prática precisamos só colocar os Valores Singulares diferentes de 0 e os Vetores Singulares associados a eles.

Podemos visualizar essa anulação ao tratar essa multiplicação de matrizes como uma soma de Produtos Externos e com a diagonal apenas como uma escalar:

\[U\Sigma V^T = \sigma_1 \mathbf{u_1} \mathbf{v_1}^T + ... + \sigma_p \mathbf{u_p} \mathbf{v_p}^T + \underbrace{ 0 \mathbf{u_{p+1}} \mathbf{v_{p+1}}^T + ...}_{anulados}\]

Onde p corresponde ao Posto da matriz A.

Relembre

Multiplicação como Soma de Produtos Externos

---

A = QR

Condição: matriz deve conter colunas Linearmente Independentes.

Decompomos A em uma matriz Q Ortogonal e em uma matriz R que é escalonada.

Relembre

Matriz Ortogonal

Escalonamento

Existem dois processos para isso: Gram-Schimdt e Reflexões de Householder.

O compromisso entre esses dois processos é que Gram-Schimdt é mais rápido, porém Reflexões de Householder é mais numericamente estável.

Dica

Essa discussão de estabilidade é vista em Álgebra Linear Computacional

Gram-Schimdt

Nesse método achamos Q e o subproduto é R.

Considerando A como um conjunto de vetores Linearmente Independentes, podemos dizer que esse conjunto é uma Base para algum Subespaço, porém esse conjunto não é necessariamente ortonormal.

O processo consiste em normalizar as colunas e tornar cada coluna independente uma da outra, assim Q será a Base ortonormal derivada da Base atual.

Relembre

Independência

Primeiramente, definimos o primeiro vetor da nossa base como a primeira coluna normalizada:

\[ \mathbf{b_1} = \frac{\mathbf{c_1}}{||\mathbf{c_1}||} \]

Agora ao definirmos \(\mathbf{b_2}\), devemos tirar a informação que pertence a \(\mathbf{b_1}\), assim decompomos o vetor e retiramos sua parte na direção de \(\mathbf{b_1}\), depois normalizamos.

Relembre

Decomposição de vetor

\[ \mathbf{b_2} = (\mathbf{c_2} - (\mathbf{c_2} \cdot \mathbf{b_1})\mathbf{b_1}) \frac{1}{||\mathbf{c_2} - (\mathbf{c_2} \cdot \mathbf{b_1})\mathbf{b_1}||} \]

E para cada coluna nova repetimos esse processo com as bases já criadas.

\[ \mathbf{b_3} = (\mathbf{c_3} - (\mathbf{c_3} \cdot \mathbf{b_1})\mathbf{b_1} - (\mathbf{c_3} \cdot \mathbf{b_2})\mathbf{b_2}) \frac{1}{||\mathbf{c_3} - (\mathbf{c_3} \cdot \mathbf{b_1})\mathbf{b_1} - (\mathbf{c_3} \cdot \mathbf{b_2})\mathbf{b_2}||} \]

Quando terminarmos Q, então achamos R:

\[ R = Q^TA \]

Reflexões de Householder

Nesse método achamos R e o subproduto é Q.

Considerando a matriz como um conjunto ordenado vetores colunas, então sempre podemos fazer uma rotação em todos esses vetores para que o primeiro vetor caia no primeiro eixo.

Assim, no primeiro vetor da matriz resultante, o primeiro componente será um número e todos os componentes abaixo serão 0.

A reflexão é simplesmente a forma mais fácil de achar uma rotação que satisfaça esse requerimento (existe uma fórmula pronta), essa reflexão é definida como uma matriz H:

Onde \(\mathbf{\vec{u}}\) é a normal do eixo de reflexão, então o vetor \(\mathbf{\vec{u}}\) é basicamente o vetor normalizado que conecta \(\mathbf{v}\) com sua reflexão.

Relembre

Conectamos vetores com substração

Note que a reflexão \(H\mathbf{v} = ||\mathbf{v}|| \mathbf{\vec{e_1}}\), onde \(\mathbf{e_1}\) é o primeiro eixo \(\begin{pmatrix} 1 \\ 0 \\0 \\ ...\end{pmatrix}\), pois uma rotação não altera a norma dos vetores.

\[ \mathbf{\vec{u}} = (\mathbf{v} - ||\mathbf{v}|| \mathbf{\vec{e_1}}) \frac{1}{||\mathbf{v} - ||\mathbf{v}|| \mathbf{\vec{e_1}}||} \]

Essa fórmula da matriz H pode ser vista nesse vídeo:

Exemplo

\[ A = \begin{pmatrix} 2 & 6 \\ 2 & -9 \\ 1 & 12 \end{pmatrix} \]

Após aplicarmos o primeiro H, manteremos a primeira coluna \(\mathbf{a_1}\) com a mesma norma, porém na direção do eixo:

\[||\mathbf{a_1}|| = \sqrt{2^2 + 2^2 + 1^2} = 3\]

\[ \mathbf{\vec{u}} = \left( \begin{pmatrix} 2\\ 2\\ 1 \end{pmatrix} - \begin{pmatrix} 3\\ 0\\ 0 \end{pmatrix} \right) \frac{1}{\sqrt{(-1)^2 + 2^2 + 1^2}} = \frac{1}{\sqrt{6}} \begin{pmatrix} -1\\ 2\\ 1 \end{pmatrix} \]

\[ H_1 = I - 2 \cdot \frac{1}{6} \begin{pmatrix} -1\\ 2\\ 1 \end{pmatrix} \begin{pmatrix} -1 & 2 & 1 \end{pmatrix} \]

\[ H_1 = \frac{1}{3} \begin{pmatrix} 2 & 2 & 1\\ 2 & -1 & -2\\ 1 & -2 & 2 \end{pmatrix} \]

\[ H_1A = \begin{pmatrix} 3 & 10 \\ 0 & 3\\ 0 & 4 \end{pmatrix} \]

Note que iremos alterar todas as colunas, mas chegamos no resultado esperado na primeira coluna

Agora devemos zerar o terceiro elemento da segunda coluna, o que podemos fazer é ignorar a primeira coluna e linha, assim temos uma matriz independente.

\[ B = \begin{pmatrix} 3\\ 4 \end{pmatrix} \]

\[||\mathbf{b_1}|| = \sqrt{ 3^2 + 4^2} = 5\]

\[\mathbf{\vec{u}} = \left( \begin{pmatrix} 3\\ 4 \end{pmatrix} - \begin{pmatrix} 5\\ 0 \end{pmatrix} \right) \frac{1}{\sqrt{(-2)^2 + 4^2}} = \frac{1}{\sqrt{20}} \begin{pmatrix} -2\\ 4 \end{pmatrix} \]

\[\begin{aligned} H' &= I - 2 \cdot 20 \cdot \begin{pmatrix} -2\\ 4 \end{pmatrix} \begin{pmatrix} -2 & 4 \end{pmatrix}\\ H' &= \begin{pmatrix} -159 & 320\\ 320 & -639 \end{pmatrix} \end{aligned}\]

Para encadearmos as reflexões devemos pegar apenas a submatriz do resultado da direita, então:

\[ H_2 = \begin{pmatrix} 1 & 0 & 0 \\ 0 & -159 & 320 \\ 0 & 320 & -639 \end{pmatrix} \]

Assim, achamos R:

\[ R = H_2 H_1 A \]

\[ R = \begin{pmatrix} 3 & 10\\ 0 & 5\\ 0 & 0\\ \end{pmatrix} \]

E usamos a propriedade que os H's são matrizes Ortogonais (não normalizados) e simétricas para achar Q:

\[ Q = H_1 H_2 A \]

---

A = LU

Decompomos A em uma matriz Triangular Inferior (Lower) e em uma matriz Triangular Superior (Upper).

Essa decomposição é um resultado direto da Eliminação de Gauss, onde U é a forma Escalonada de A e L reverte o escalonamento.

Relembre

Eliminação de Gauss

Se representarmos cada Operação Elementar da Eliminação de Gauss como operação entre linhas (multiplica pela esquerda), então chegamos a:

\[E_3E_2E_1A = U\]

Como os \(E_n\) são Triangulares Inferiores, então essa decomposição não é possível se houver troca de linha na Eliminação de Gauss.

Logo, podemos calcular L:

\[A = E_1^{-1}E_2^{-1}E_3^{-1}U\]

\[L = E_1^{-1}E_2^{-1}E_3^{-1}\]

Para achar a inversa dessas matrizes temos apenas que trocar o sinal dos elementos abaixo da diagonal, pois isso reverte as Operações Elementares.

Exemplo

\[ E_1 = \begin{pmatrix} 1 & 0 & 0\\ -2 & 1 & 0\\ 0 & 0 & 1\\ \end{pmatrix}; E_2 = \begin{pmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ -3 & 0 & 1\\ \end{pmatrix}; E_3 = \begin{pmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 0 & -2 & 1\\ \end{pmatrix} \]

\[ E_1^{-1} = \begin{pmatrix} 1 & 0 & 0\\ 2 & 1 & 0\\ 0 & 0 & 1\\ \end{pmatrix}; E_2^{-1} = \begin{pmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 3 & 0 & 1\\ \end{pmatrix}; E_3^{-1} = \begin{pmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 0 & 2 & 1\\ \end{pmatrix} \]

Como o encadeamento das matrizes \(E_n\) foi invertido, então sempre estamos no estado certo para aplicar a inversa, não precisamos acumular as operações passadas.

\[ L = \begin{pmatrix} 1 & 0 & 0\\ 2 & 1 & 0\\ 3 & 2 & 1\\ \end{pmatrix} \]

PA = LU

Essa é uma variação que torna a decomposição A = LU possível quando há troca de linha na Eliminação de Gauss.

P é uma matriz que muda a ordem das linhas de A e torna possível a Eliminação de Gauss da matriz \(PA\) sem troca de linhas.

---

Cholesky

A decomposição Cholensky é uma extensão da decomposição A=LU quando A é uma matriz Simétrica Definida Positiva.

O resultado esperado é achar \(A=GG^T\), onde G é uma matriz Triangular Inferior.

Exemplo

Podemos achar alguns G's diferentes que não são Triangular, exemplo:

Na Decomposição Espectral:

\[A = XDX^T\]

Como todos os Autovalores de A são positivos, então podemos tirar a raiz quadrada da diagonal, então:

\[A = X \sqrt{D} \sqrt{D} X^T\]

Se definirmos \(G=X \sqrt{D}\), então \(G^T = \sqrt{D} X^T\)

O resultado especial é que na decomposição A=LU em uma matriz Simétrica Definida Positiva, U pode ser decomposto em \(DL^T\) onde D é uma matriz Diagonal com os pivôs:

\[A = LDL^T\]

Como os pivôs são positivos, então podemos tirar a raiz quadrada:

\[A = L \sqrt{D} \sqrt{D} L^T\]

Se definirmos \(G=L \sqrt{D}\), então \(G^T=\sqrt{D} L^T\).

Aplicação

Uma aplicação da decomposição Cholesky é na simulação de Monte Carlo.

---

A = CR

Se definirmos C como a Base para o Espaço Coluna de A, então C será composto apenas das colunas Linearmente Independentes de A.

Relembre

Espaço Coluna

R será uma matriz de Combinação Linear das colunas de C para montar as colunas de A, note que a quantidade de linhas de R será igual a quantidade de colunas de C.

R também será o resultado da Eliminação Gauss (excluindo as linhas nulas), mostrando mais uma vez que a quantidade de linhas Linearmente Independentes é igual a quantidade de Colunas Linearmente Independentes, como dito pelo Teorema Fundamental da Álgebra Linear.

---

<Anterior | Próximo>