R sparseQR-class 稀疏 QR 分解

說明

sparseQR 是 m \times n ( m \ge n ) 實數矩陣的稀疏、行和 column-pivoted QR 分解類，具有一般形式

或(等效地)

其中 P_1 和 P_2 是置換矩陣，Q = \prod_{j = 1}^{n} H_j 是 m \times m 正交矩陣(Q_1 包含第一個 n 列向量)，等於 n Householder 矩陣 H_j 的乘積, 和R 是具有非負對角線元素的 m \times n 上梯形矩陣(R_1 包含第一個 n 行向量，並且是上三角矩陣)。

用法

qrR(qr, complete = FALSE, backPermute = TRUE, row.names = TRUE)

參數

細節

qr.Q 的方法不會返回 Q ，而是返回(也是正交的)乘積 P_1' Q 。此行為在代數上與 base 實現一致(請參閱 qr )，可以通過注意到 base 中的 qr.default 不旋轉行來看出這一點，從而將 P_1 限製為單位矩陣。因此 qr.Q(qr.default(x)) 也返回 P_1' Q 。

同樣， qr.qy 和 qr.qty 的方法在左側乘以 P_1' Q 和 Q' P_1 而不是 Q 和 Q' 。

預期 qr.Q (或 qr.R 、 qr.qy 、 qr.qty )的值是根據 “equivalent” 稀疏和密集分解(例如，qr(x) 和 qr(as(x, "matrix")) 對於類的 x 計算得出的)是錯誤的dgCMatrix ) 比較相等。底層分解算法有很大不同，特別是它們采用不同的旋轉策略，並且一般來說，即使對於固定的 P_1 和 P_2 ，分解也不是唯一的。

另一方麵， qr.X 、 qr.coef 、 qr.fitted 和 qr.resid 的值是明確定義的，在這些情況下，稀疏和密集計算應該比較相等(在一定的容差範圍內)。

qr.R 的方法是 qrR 的簡單包裝，但默認情況下不是 back-permuting 並且從不給出行名稱。它不支持 backPermute = TRUE 直到 Matrix 1.6-0 ，因此如果 Matrix >= 1.6-0 未知，則需要 back-permuted 結果的代碼應調用 qrR。

插槽

Dim , Dimnames

從虛擬類MatrixFactorization繼承。

beta

長度為 Dim[2] 的數值向量，用於構造 Householder 矩陣；請參閱下麵的V。

V

具有 Dim[2] 列的類 dgCMatrix 的對象。行數 nrow(V) 至少為 Dim[1] 且最多為 Dim[1]+Dim[2] 。 V 是下梯形，其列向量生成組成正交 Q 因子的 Householder 矩陣 H_j。具體來說， H_j 被構造為 diag(Dim[1]) - beta[j] * tcrossprod(V[, j]) 。

R

類dgCMatrix 的對象，具有nrow(V) 行和Dim[2] 列。 R 是具有非負對角線元素的上梯形 R 因子。

p , q

長度從 0 開始的整數向量nrow(V)和Dim[2]，分別指定應用於因式分解矩陣的行和列的排列。q長度為 0 的值是有效的並且等價於恒等排列，意味著沒有列旋轉。使用R語法，矩陣P_1 A P_2正是A[p+1, q+1](A[p+1, ]當q長度為 0)。

擴展

直接類 QR 。類 MatrixFactorization ，按類 QR ，距離 2。

實例化

對象可以通過 new("sparseQR", ...) 形式的調用直接生成，但它們更通常作為繼承自 sparseMatrix (通常是 dgCMatrix )的 x 的 qr(x) 值獲得。

方法

determinant

signature(from = "sparseQR", logarithm = "logical") ：計算因式分解矩陣A 的行列式或其對數。

expand1

signature(x = "sparseQR") ：參見expand1-methods 。

expand2

signature(x = "sparseQR") ：參見expand2-methods 。

qr.Q

signature(qr = "sparseQR") ：返回 dgeMatrix P_1' Q 或 P_1' Q_1 ，具體取決於可選參數 complete 。默認為 FALSE ，表示 P_1' Q_1 。

qr.R

signature(qr = "sparseQR") : qrR 返回 R 、 R_1 、 R P2' 或 R_1 P2' ，具體取決於可選參數 complete 和 backPermute 。兩種情況下的默認值都是 FALSE ，表示 R_1 ，以與 base 兼容。這種情況下結果的類是 dtCMatrix 。在其他三種情況下，它是 dgCMatrix 。

qr.X

signature(qr = "sparseQR") ：默認情況下，將 A 作為 dgeMatrix 返回。如果 m > n 和可選參數 ncol 大於 n ，則結果將用 P_1' Q J 進行擴充，其中 J 由 m \times m 單位矩陣的 (n+1) 到 ncol 列組成。

qr.coef

signature(qr = "sparseQR", y = .) ：以 dgeMatrix 或向量形式返回左側 y 乘以 P_2 R_1^{-1} Q_1' P_1 的結果。

qr.fitted

signature(qr = "sparseQR", y = .) ：以 dgeMatrix 或向量形式返回左側 y 乘以 P_1' Q_1 Q_1' P_1 的結果。

qr.resid

signature(qr = "sparseQR", y = .) ：以 dgeMatrix 或向量形式返回左側 y 乘以 P_1' Q_2 Q_2' P_1 的結果。

qr.qty

signature(qr = "sparseQR", y = .) ：以 dgeMatrix 或向量形式返回左側 y 乘以 Q' P_1 的結果。

qr.qy

signature(qr = "sparseQR", y = .) ：以 dgeMatrix 或向量形式返回左側 y 乘以 P_1' Q 的結果。

solve

signature(a = "sparseQR", b = .) ：參見solve-methods 。

例子

showClass("sparseQR")
set.seed(2)

m <- 300L
n <- 60L
A <- rsparsematrix(m, n, 0.05)

## With dimnames, to see that they are propagated :
dimnames(A) <- dn <- list(paste0("r", seq_len(m)),
                          paste0("c", seq_len(n)))

(qr.A <- qr(A))
str(e.qr.A <- expand2(qr.A, complete = FALSE), max.level = 2L)
str(E.qr.A <- expand2(qr.A, complete =  TRUE), max.level = 2L)

t(sapply(e.qr.A, dim))
t(sapply(E.qr.A, dim))

## Horribly inefficient, but instructive :
slowQ <- function(V, beta) {
    d <- dim(V)
    Q <- diag(d[1L])
    if(d[2L] > 0L) {
        for(j in d[2L]:1L) {
            cat(j, "\n", sep = "")
            Q <- Q - (beta[j] * tcrossprod(V[, j])) %*% Q
        }
    }
    Q
}

ae1 <- function(a, b, ...) all.equal(as(a, "matrix"), as(b, "matrix"), ...)
ae2 <- function(a, b, ...) ae1(unname(a), unname(b), ...)

## A ~ P1' Q R P2' ~ P1' Q1 R1 P2' in floating point
stopifnot(exprs = {
    identical(names(e.qr.A), c("P1.", "Q1", "R1", "P2."))
    identical(names(E.qr.A), c("P1.", "Q" , "R" , "P2."))
    identical(e.qr.A[["P1."]],
              new("pMatrix", Dim = c(m, m), Dimnames = c(dn[1L], list(NULL)),
                  margin = 1L, perm = invertPerm(qr.A@p, 0L, 1L)))
    identical(e.qr.A[["P2."]],
              new("pMatrix", Dim = c(n, n), Dimnames = c(list(NULL), dn[2L]),
                  margin = 2L, perm = invertPerm(qr.A@q, 0L, 1L)))
    identical(e.qr.A[["R1"]], triu(E.qr.A[["R"]][seq_len(n), ]))
    identical(e.qr.A[["Q1"]],      E.qr.A[["Q"]][, seq_len(n)] )
    identical(E.qr.A[["R"]], qr.A@R)
 ## ae1(E.qr.A[["Q"]], slowQ(qr.A@V, qr.A@beta))
    ae1(crossprod(E.qr.A[["Q"]]), diag(m))
    ae1(A, with(e.qr.A, P1. %*% Q1 %*% R1 %*% P2.))
    ae1(A, with(E.qr.A, P1. %*% Q  %*% R  %*% P2.))
    ae2(A.perm <- A[qr.A@p + 1L, qr.A@q + 1L], with(e.qr.A, Q1 %*% R1))
    ae2(A.perm                               , with(E.qr.A, Q  %*% R ))
})

## More identities
b <- rnorm(m)
stopifnot(exprs = {
    ae1(qrX <- qr.X     (qr.A   ), A)
    ae2(qrQ <- qr.Q     (qr.A   ), with(e.qr.A, P1. %*% Q1))
    ae2(       qr.R     (qr.A   ), with(e.qr.A, R1))
    ae2(qrc <- qr.coef  (qr.A, b), with(e.qr.A, solve(R1 %*% P2., t(qrQ)) %*% b))
    ae2(qrf <- qr.fitted(qr.A, b), with(e.qr.A, tcrossprod(qrQ) %*% b))
    ae2(qrr <- qr.resid (qr.A, b), b - qrf)
    ae2(qrq <- qr.qy    (qr.A, b), with(E.qr.A, P1. %*% Q %*% b))
    ae2(qr.qty(qr.A, qrq), b)
})

## Sparse and dense computations should agree here
qr.Am <- qr(as(A, "matrix")) # <=> qr.default(A)
stopifnot(exprs = {
    ae2(qrX, qr.X     (qr.Am   ))
    ae2(qrc, qr.coef  (qr.Am, b))
    ae2(qrf, qr.fitted(qr.Am, b))
    ae2(qrr, qr.resid (qr.Am, b))
})

參考

Davis, T. A. (2006). Direct methods for sparse linear systems. Society for Industrial and Applied Mathematics. doi:10.1137/1.9780898718881

Golub, G. H., & Van Loan, C. F. (2013). Matrix computations (4th ed.). Johns Hopkins University Press. doi:10.56021/9781421407944

也可以看看

類dgCMatrix。

來自 base 的通用函數 qr ，其默認方法 qr.default “defines” 是密集 QR 分解的 S3 類 qr。

qr-methods 用於 Matrix 中定義的方法。

通用函數 expand1 和 expand2 。

QR 分解的許多輔助函數： qr.Q , qr.R , qr.X , qr.coef , qr.fitted , qr.resid , qr.qty , qr.qy 和 qr.solve 。