SVD Demonstration |
|
Not sure what SVD’s purpose is in a world of better compression algorithms - stuffing this in here after pulling it out of the PCA write up where it used to serve as a conceptual stepping stone - i don’t think it works as a seperate document and think it needs to be recombined
Singular value decompositions uses the raw data matrix. When performing SVD we use the unaltered data matrix with no centering or scaling. Using the modern convention the SVD takes the form \(A=U\:\Sigma\:V^T\) where:
\[AA^T=U\Sigma V^T(U\Sigma V^T)^T=U\Sigma^2U^T\]
if \(U\) is orthogonal eigenvector matrix this is the diagonalisation of the positive definite \(AA^T\) matrix: \[U\Sigma^2U^T\sim Q\Lambda Q^{T} \sim S\Lambda S^{-1}\\\:\\ \therefore\quad AA^T\:S=S\Lambda\;\sim\;AA^T\:U=U\Sigma^2\]
So, the left singular vectors \(\vec{\mathbb{u}_i}\in\mathbb{R}^n\) span observation space and correspond to both the eigenvectors of \(AA^T\) and the square roots of the eigenvalues \(\lambda_i\) known as the singular values: \(\sigma_i=\sqrt{\lambda_i}\). \[ AA^T\:\vec{\mathbb{u}_i}=\sigma_i^2\vec{\mathbb{u}_i} \]
\[A^TA=(U\Sigma V^T)^TU\Sigma V^T=V\Sigma^2V^T\]
if \(V\) is orthogonal eigenvector matrix this is the diagonalisation of the positive definite \(A^TA\) matrix: \[V\Sigma^2V^T\sim Q\Lambda Q^{T} \sim S\Lambda S^{-1}\\\:\\ \therefore\quad A^TA\:S=S\Lambda\;\sim\;A^TA\:V=V\Sigma^2\]
So the right singular vectors \(\vec{\mathbb{v}}_i\in\mathbb{R}^m\) span the feature space and they correspond to both the eigenvectors of \(A^TA\) and the singular values \(\sigma_i\) which are the square root of the eigenvalues. \(V\) was orthogonal - lucky that! \[ A^TA\:\vec{\mathbb{v}_i}=\sigma_i^2\vec{\mathbb{v}_i} \]
If using a colour image either select a channel or take the mean of each colour channel to get your data matrix.
img <- read.bmp("./seg001.bmp")
dim(img) # Check the structure of the image
## [1] 206 206 3
green_channel <- img[, , 2]
A <- apply(img, c(1, 2), mean)
Make the call to svd
and collect the components if you
like… (Note how the diagonal singlar values matrix has to be built from
the diagonal values vector)
mySVDObject<-svd(A)
U <- mySVDObject$u
Sigma <- diag(mySVDObject$d)
V <- mySVDObject$v
Check to see when singular values level off - looks about dimension #75 to me.
plot(mySVDObject$d, type = "b", pch = 19, col = "blue",
xlab = "Index of Singular Value", ylab = "Singular Values")
abline(h = 0, col = "red", lty = 2)
Perform reconstructions
A_approx <- U[,1:75] %*% Sigma[1:75,1:75] %*% t(V[,1:75])
Checkout the reconstruction and compare with original…
par(mfrow = c(1, 2))
image(A,main="Actual: Rank 206")
image(A_approx,main=paste("Approx: Rank", 75))