(original) (raw)

--- title: "Eigenvalues: Spectral Decomposition" author: "Michael Friendly" date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Eigenvalues: Spectral Decomposition} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, echo = FALSE} knitr::opts_chunk$set( warning = FALSE, message = FALSE ) options(digits=4) ``` ```{r} library(matlib) # use the package ``` ## Setup This vignette uses an example of a 3times33 \times 33times3 matrix to illustrate some properties of eigenvalues and eigenvectors. We could consider this to be the variance-covariance matrix of three variables, but the main thing is that the matrix is **square** and **symmetric**, which guarantees that the eigenvalues, lambdai\lambda_ilambdai are real numbers, and non-negative, lambdaige0\lambda_i \ge 0lambdaige0. ```{r} A <- matrix(c(13, -4, 2, -4, 11, -2, 2, -2, 8), 3, 3, byrow=TRUE) A ``` Get the eigenvalues and eigenvectors using `eigen()`; this returns a named list, with eigenvalues named `values` and eigenvectors named `vectors`. We call these `L` and `V` here, but in formulas they correspond to a diagonal matrix, mathbfLambda=diag(lambda1,lambda2,lambda3)\mathbf{\Lambda} = diag(\lambda_1, \lambda_2, \lambda_3)mathbfLambda=diag(lambda1,lambda2,lambda3), and a (orthogonal) matrix mathbfV\mathbf{V}mathbfV. ```{r} ev <- eigen(A) # extract components (L <- ev$values) (V <- ev$vectors) ``` ## Matrix factorization 1. Factorization of A: A = V diag(L) V'. That is, the matrix mathbfA\mathbf{A}mathbfA can be represented as the product mathbfA=mathbfVmathbfLambdamathbfV′\mathbf{A}= \mathbf{V} \mathbf{\Lambda} \mathbf{V}'mathbfA=mathbfVmathbfLambdamathbfV′. ```{r} V %*% diag(L) %*% t(V) ``` 2. V diagonalizes A: L = V' A V. That is, the matrix mathbfV\mathbf{V}mathbfV transforms mathbfA\mathbf{A}mathbfA into the diagonal matrix mathbfLambda\mathbf{\Lambda}mathbfLambda, corresponding to orthogonal (uncorrelated) variables. ```{r} diag(L) zapsmall(t(V) %*% A %*% V) ``` ## Spectral decomposition The basic idea here is that each eigenvalue--eigenvector pair generates a rank 1 matrix, lambdaimathbfvimathbfvi′\lambda_i \mathbf{v}_i \mathbf{v}_i 'lambdaimathbfvimathbfvi′, and these sum to the original matrix, mathbfA=sumilambdaimathbfvimathbfvi′\mathbf{A} = \sum_i \lambda_i \mathbf{v}_i \mathbf{v}_i 'mathbfA=sumilambdaimathbfvimathbfvi′. ```{r} A1 = L[1] * V[,1] %*% t(V[,1]) A1 A2 = L[2] * V[,2] %*% t(V[,2]) A2 A3 = L[3] * V[,3] %*% t(V[,3]) A3 ``` Then, summing them gives `A`, so they do decompose `A`: ```{r} A1 + A2 + A3 all.equal(A, A1+A2+A3) ``` ### Further properties 1. Sum of squares of A = sum of sum of squares of A1, A2, A3 ```{r} sum(A^2) c( sum(A1^2), sum(A2^2), sum(A3^2) ) sum( sum(A1^2), sum(A2^2), sum(A3^2) ) #' same as tr(A' A) tr(crossprod(A)) ``` 2. Each squared eigenvalue gives the sum of squares accounted for by the latent vector ```{r} L^2 cumsum(L^2) # cumulative ``` 3. The first iii eigenvalues and vectors give a rank iii approximation to `A` ```{r} R(A1) R(A1 + A2) R(A1 + A2 + A3) # two dimensions sum((A1+A2)^2) sum((A1+A2)^2) / sum(A^2) # proportion ```