Principal component analysis (PCA) is a popular method for dimensional reduction, and has been invented independently many times in different fields, resulting in various definitions. This note attempts to unify 2 of the most popular definitions of principal components and illustrate how PCA can be done correspondingly.
Alternative definitions PCA from the SVD of the centered matrix The singular value decomposition (SVD) of an \(n\times d\) matrix \(X\) has the form
Defining utility functions burd = colorRampPalette(colors = c("blue", "white", "red"))(n = 499) blues = colorRampPalette(colors = c('#deebf7', '#08306b'))(n = 256) plot.matrix = function(m, col = burd, asp=1) { m %>% apply(MARGIN = 2, rev) %>% t() %>% image(useRaster = TRUE, axes = FALSE, col = col, asp = asp) } parse_timing_output = function(output_raw) { sapply(output_raw, function(x) { str = stringr::str_split(x,":\\s+")[[1]] return(as.numeric(str[2])) }) } An arbitrary matrix sin2d = function(a, b) { sin((a/ 500 - b / 15) * pi) } start = proc.
Singular value decomposition is an expensive operation. For rectangular matrices with significant different dimensions, i.e. very “fat” or “thin” matrices, there is a trick to make the computation cheaper. This trick is implemented in fast.svd() of the R package corpcor.
Calculate SVD The singular value decomposition of a matrix \(M\) of size \(m \times n\).
\[ M = UDV^T \]
\[ \begin{align} MM^T &= (UDV^T)(UDV^T)^T \\ &= (UDV^T)V(UD)^T \\ &= UD (V^TV) (UD)^T \\ &= UD(UD)^T \quad (V\text{ is orthogonal}) \\ &= UDD^TU^T \\ \end{align} \] Thus the decomposition of \(MM^T\) gives \(U\) and \(D^2\).
\[ \newcommand{\matrix}[1]{\mathbf{#1}} \] Let \(\matrix{A}\) be a data set of \(m\) points in \(\mathbb{R}^d\). One application of SVD is to create a compressed representation of \(\matrix{A}\). Rank-\(k\) approximation of \(A\) is created by calculating the singular value decomposition of \(\matrix{A}\) \[ \matrix{A} = \matrix{U}\matrix{\Sigma}{\matrix{V}} \] and reconstruct it with \(k \leq d\) first singular values. \[ \matrix{A_k} = \matrix{U_k}\matrix{\Sigma_k}\matrix{V_k^T} \]
SVD in R Each implementation of SVD has some varieties in the output representation.