R rankMatrix 矩陣的秩

說明

計算‘the’矩陣秩，理論上是一個定義明確的函數(*)，但在實踐中有些模糊。我們提供了幾種方法，默認對應Matlab的定義。

(*) n \times m 矩陣 A 、 rk(A) 的秩是線性獨立列(或行)的最大數量；因此rk(A) \le 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	比如說，維度為 n \times m 的數字矩陣。
tol	非負數，指定用於測試 “practically zero” 的(相對，“scalefree”)容差，具體含義取決於 method ；默認情況下， max(dim(x)) * .Machine$double.eps 是根據 Matlab 的默認值(因為它唯一的方法是我們的 method="tolNorm2" )。
method	指定排名計算方法的字符串，可以縮寫為： "tolNorm2"： 奇異值的數量>= tol * max(sval)； "qrLINPACK"： 對於稠密矩陣，這是qr(x, tol, LAPACK=FALSE)(這是qr(...)$rank); 這個("qr*"，密集)版本曾經是 過去一段時間推薦的計算矩陣秩的方法。 對於稀疏 x ，這相當於 "qr.R" 。 "qr.R"： 這是三角矩陣 R 的秩，其中 qr() 使用 LAPACK 或 "sparseQR" 方法(參見 qr-methods )來計算分解 QR 。 R 的排名然後定義為 R 的 “non-zero” 對角線條目 d_i 的數量，並且 “non-zero” 滿足 |d_i| \ge \mathrm{tol}\cdot\max(|d_i|) 。 "qr"： 用於向後兼容；對於密集 x ，它對應於 "qrLINPACK" ，而對於稀疏 x ，它使用 "qr.R" 。 對於所有 "qr*" 方法，不使用奇異值 sval，這對於大型稀疏矩陣 x 可能至關重要，在這種情況下，當未指定 sval 時，默認計算 svd() 目前將 x 強製轉換為密集矩陣。 "useGrad"： 考慮到(遞減的)奇異值的“gradient”，最小間隙的索引。 "maybeGrad"： 僅當看起來合理時才選擇方法"useGrad"；否則使用"tolNorm2"。
sval	x 的非遞增奇異值的數值向量；通常未指定，並在需要時根據 x 計算，即除非 method = "qr" 。
warn.t	邏輯指示 rankMatrix() 在需要 t(x) 而不是 x 時是否應發出警告。目前，僅對於method = "qr"，默認情況下會發出警告，因為調用者通常可以更有效地直接傳遞t(x)。
warn.qr	在 QR 情況下(即，如果 method 以 "qr" 開頭)， rankMatrix() 調用 qr2rankMarix(.., do.warn = warn.qr) ，請參見下文。
qr	一個R對象產生於qr(x,..)，即通常繼承自class "qr"或者"sparseQR".
isBqr	logical 指示 qr 是否由 base qr() 產生。 (否則，它通常來自 Matrix 包稀疏 qr 。)
do.warn	邏輯性；如果為 true，則警告 QR 分解的 R 矩陣中的非有限(或在 sparseQR 情況下為負)對角線條目。不要輕易改變！

細節

qr2rankMatrix() 通常從 rankMatrix() 調用 "qr" * method ，但可以直接使用 - 如果 qr 分解可用，則效率更高。

值

如果x是所有0(或零維度)的矩陣，則秩為零；否則，通常是 1:min(dim(x)) 中的正整數，其屬性詳細說明所使用的方法。

在極少數情況下，稀疏QR分解“fails”在R、d_i(見上文)的對角線條目中以非有限(通常是NaN條目)結束。然後，發出警告信號(除非 warn.qr /do.warn 不為 true)並返回 NA (具體來說， NA_integer_ )。

注意

對於大型稀疏矩陣 x ，除非您可以自己指定 sval ，否則當前 method = "qr" 可能是唯一可行的矩陣，因為其他矩陣需要 sval 並調用 svd() ，當前將 x 強製為 denseMatrix 這可能非常慢或不可能，具體取決於矩陣維度。

請注意，在稀疏 x 、 method = "qr" 的情況下，計算所有非嚴格零對角線條目 d_i ，直至包括 Matrix 版本 1.1-0，即該方法隱式使用 tol = 0 ，請參閱還有下麵的set.seed(42) 示例。

例子

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)

作者

Martin Maechler; for the "*Grad" methods building on suggestions by Ravi Varadhan.

也可以看看

qr、svd。