Randomized SVD decomposes a matrix by extracting the first k singular values/vectors using k+p random projections. This works surprisingly well for large matrices.
My question concerns the singular values that are output from the algorithm. Why aren't the values equal to the first k-singular values if you do the full SVD?
Below I have a simple implementation in R. Any suggestions on improving the performance would be appreciated.
rsvd = function(A, k=10, p=5){
n = nrow(A)
y = A %*% matrix(rnorm(n * (k+p)), nrow=n)
q = qr.Q(qr(y))
b = t(q) %*% A
svd = svd(b)
list(u=q %*% svd$u, d=svd$d, v=svd$v)
}
> set.seed(10)
> A <- matrix(rnorm(500*500),500,500)
> svd(A)$d[1:15]
[1] 44.94307 44.48235 43.78984 43.44626 43.27146 43.15066 42.79720 42.54440 42.27439 42.21873 41.79763 41.51349 41.48338 41.35024 41.18068
> rsvd.o(A,10,5)$d
[1] 34.83741 33.83411 33.09522 32.65761 32.34326 31.80868 31.38253 30.96395 30.79063 30.34387 30.04538 29.56061 29.24128 29.12612 27.61804
B <- matrix(rnorm(500*50),500,500) # rank 50
> rsvd(B,10,5)$d
[1] 86.48035 83.02114 81.03988 80.04358 77.24979 76.10945 74.47357 74.08382
[9] 72.85898 72.06897 69.59526 67.70750 66.53867 62.96446 61.50838
> svd(B)$d[1:15]
[1] 92.44779 91.47689 88.71948 88.08170 87.24533 85.13312 84.14741 83.71757
[9] 82.80832 81.43005 80.73903 79.92959 78.87421 78.33509 77.38431
As Joris pointed out, I have this posted on stackoverflow as well. you can find the revelant conversation here
https://stackoverflow.com/questions/4224031/randomized-svd-singular-values
Also see the relevant paper by Martinsson et al: A randomized algorithm for the decomposition of matrices