rsvd {rsvd}R Documentation

Randomized Singular Value Decomposition (rsvd).

Description

The randomized SVD computes the near-optimal low-rank approximation of a rectangular matrix using a fast probablistic algorithm.

Usage

rsvd(A, k = NULL, nu = NULL, nv = NULL, p = 10, q = 2, sdist = "normal")

Arguments

A

array_like;
a real/complex (m,n)(m, n) input matrix (or data frame) to be decomposed.

k

integer;
specifies the target rank of the low-rank decomposition. kk should satisfy k<<min(m,n)k << min(m,n).

nu

integer, optional;
number of left singular vectors to be returned. nunu must be between 00 and kk.

nv

integer, optional;
number of right singular vectors to be returned. nvnv must be between 00 and kk.

p

integer, optional;
oversampling parameter (by default p=10p=10).

q

integer, optional;
number of additional power iterations (by default q=2q=2).

sdist

string c(unif,normal,rademacher)c( 'unif', 'normal', 'rademacher'), optional;
specifies the sampling distribution of the random test matrix:
unif'unif' : Uniform '[-1,1]'.
normal'normal' (default) : Normal '~N(0,1)'.
rademacher'rademacher' : Rademacher random variates.

Details

The singular value decomposition (SVD) plays an important role in data analysis, and scientific computing. Given a rectangular (m,n)(m,n) matrix AA, and a target rank k<<min(m,n)k << min(m,n), the SVD factors the input matrix AA as

A=Ukdiag(dk)Vk A = U_{k} diag(d_{k}) V_{k}^\top

The kk left singular vectors are the columns of the real or complex unitary matrix UU. The kk right singular vectors are the columns of the real or complex unitary matrix VV. The kk dominant singular values are the entries of dd, and non-negative and real numbers.

pp is an oversampling parameter to improve the approximation. A value of at least 10 is recommended, and p=10p=10 is set by default.

The parameter qq specifies the number of power (subspace) iterations to reduce the approximation error. The power scheme is recommended, if the singular values decay slowly. In practice, 2 or 3 iterations achieve good results, however, computing power iterations increases the computational costs. The power scheme is set to q=2q=2 by default.

If k>(min(n,m)/4)k > (min(n,m)/4), a deterministic partial or truncated svd algorithm might be faster.

Value

rsvd returns a list containing the following three components:

d

array_like;
singular values; vector of length (k)(k).

u

array_like;
left singular vectors; (m,k)(m, k) or (m,nu)(m, nu) dimensional array.

v

array_like;
right singular vectors; (n,k)(n, k) or (n,nv)(n, nv) dimensional array.

Note

The singular vectors are not unique and only defined up to sign (a constant of modulus one in the complex case). If a left singular vector has its sign changed, changing the sign of the corresponding right vector gives an equivalent decomposition.

Author(s)

N. Benjamin Erichson, erichson@berkeley.edu

References

See Also

svd, rpca

Examples

library('rsvd')

# Create a n x n Hilbert matrix of order n,
# with entries H[i,j] = 1 / (i + j + 1).
hilbert <- function(n) { i <- 1:n; 1 / outer(i - 1, i, "+") }
H <- hilbert(n=50)

# Low-rank (k=10) matrix approximation using rsvd
k=10
s <- rsvd(H, k=k)
Hre <- s$u %*% diag(s$d) %*% t(s$v) # matrix approximation
print(100 * norm( H - Hre, 'F') / norm( H,'F')) # percentage error
# Compare to truncated base svd
s <- svd(H)
Hre <- s$u[,1:k] %*% diag(s$d[1:k]) %*% t(s$v[,1:k]) # matrix approximation
print(100 * norm( H - Hre, 'F') / norm( H,'F')) # percentage error


[Package rsvd version 1.0.5 Index]