rankMatrix | R Documentation |
Compute ‘the’ matrix rank, a well-defined functional in theory(*), somewhat ambiguous in practice. We provide several methods, the default corresponding to Matlab's definition.
(*) The rank of a n x m matrix A, rk(A), is the maximal number of linearly independent columns (or rows); hence rk(A) <= min(n,m).
rankMatrix(x, tol = NULL, method = c("tolNorm2", "qr.R", "qrLINPACK", "qr", "useGrad", "maybeGrad"), sval = svd(x, 0, 0)$d, warn.t = TRUE, warn.qr = TRUE) qr2rankMatrix(qr, tol = NULL, isBqr = is.qr(qr), do.warn = TRUE)
x |
numeric matrix, of dimension n x m, say. |
tol |
nonnegative number specifying a (relative,
“scalefree”) tolerance for testing of
“practically zero” with specific meaning depending on
|
method |
a character string specifying the computational method for the rank, can be abbreviated:
|
sval |
numeric vector of non-increasing singular values of
|
warn.t |
logical indicating if |
warn.qr |
in the QR cases (i.e., if |
qr |
an R object resulting from |
isBqr |
|
do.warn |
logical; if true, warn about non-finite (or in the
|
qr2rankMatrix()
is typically called from rankMatrix()
for
the "qr"
* method
s, but can be used directly - much more
efficiently in case the qr
-decomposition is available anyway.
If x
is a matrix of all 0
(or of zero dimension), the rank
is zero; otherwise, typically a positive integer in 1:min(dim(x))
with attributes detailing the method used.
There are rare cases where the sparse QR decomposition
“fails” in so far as the diagonal entries of R, the
d_i (see above), end with non-finite, typically NaN
entries. Then, a warning is signalled (unless warn.qr
/
do.warn
is not true) and NA
(specifically,
NA_integer_
) is returned.
For large sparse matrices x
, unless you can specify
sval
yourself, currently method = "qr"
may
be the only feasible one, as the others need sval
and call
svd()
which currently coerces x
to a
denseMatrix
which may be very slow or impossible,
depending on the matrix dimensions.
Note that in the case of sparse x
, method = "qr"
, all
non-strictly zero diagonal entries d_i where counted, up to
including Matrix version 1.1-0, i.e., that method implicitly
used tol = 0
, see also the set.seed(42)
example below.
Martin Maechler; for the "*Grad" methods building on suggestions by Ravi Varadhan.
qr
, svd
.
rankMatrix(cbind(1, 0, 1:3)) # 2 (meths <- eval(formals(rankMatrix)$method)) ## a "border" case: H12 <- Hilbert(12) rankMatrix(H12, tol = 1e-20) # 12; but 11 with default method & tol. sapply(meths, function(.m.) rankMatrix(H12, method = .m.)) ## tolNorm2 qr.R qrLINPACK qr useGrad maybeGrad ## 11 11 12 12 11 11 ## The meaning of 'tol' for method="qrLINPACK" and *dense* x is not entirely "scale free" rMQL <- function(ex, M) rankMatrix(M, method="qrLINPACK",tol = 10^-ex) rMQR <- function(ex, M) rankMatrix(M, method="qr.R", tol = 10^-ex) sapply(5:15, rMQL, M = H12) # result is platform dependent ## 7 7 8 10 10 11 11 11 12 12 12 {x86_64} sapply(5:15, rMQL, M = 1000 * H12) # not identical unfortunately ## 7 7 8 10 11 11 12 12 12 12 12 sapply(5:15, rMQR, M = H12) ## 5 6 7 8 8 9 9 10 10 11 11 sapply(5:15, rMQR, M = 1000 * H12) # the *same* ## "sparse" case: M15 <- kronecker(diag(x=c(100,1,10)), Hilbert(5)) sapply(meths, function(.m.) rankMatrix(M15, method = .m.)) #--> all 15, but 'useGrad' has 14. sapply(meths, function(.m.) rankMatrix(M15, method = .m., tol = 1e-7)) # all 14 ## "large" sparse n <- 250000; p <- 33; nnz <- 10000 L <- sparseMatrix(i = sample.int(n, nnz, replace=TRUE), j = sample.int(p, nnz, replace=TRUE), x = rnorm(nnz)) (st1 <- system.time(r1 <- rankMatrix(L))) # warning+ ~1.5 sec (2013) (st2 <- system.time(r2 <- rankMatrix(L, method = "qr"))) # considerably faster! r1[[1]] == print(r2[[1]]) ## --> ( 33 TRUE ) ## another sparse-"qr" one, which ``failed'' till 2013-11-23: set.seed(42) f1 <- factor(sample(50, 1000, replace=TRUE)) f2 <- factor(sample(50, 1000, replace=TRUE)) f3 <- factor(sample(50, 1000, replace=TRUE)) D <- t(do.call(rbind, lapply(list(f1,f2,f3), as, 'sparseMatrix'))) dim(D); nnzero(D) ## 1000 x 150 // 3000 non-zeros (= 2%) stopifnot(rankMatrix(D, method='qr') == 148, rankMatrix(crossprod(D),method='qr') == 148) ## zero matrix has rank 0 : stopifnot(sapply(meths, function(.m.) rankMatrix(matrix(0, 2, 2), method = .m.)) == 0)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.