| sparseQR-class | R Documentation |
sparseQR is the class of sparse, row- and column-pivoted
QR factorizations of m \times n (m \ge n)
real matrices, having the general form
P_1 A P_2 = Q R = \begin{bmatrix} Q_1 & Q_2 \end{bmatrix} \begin{bmatrix} R_1 \\ 0 \end{bmatrix} = Q_1 R_1
or (equivalently)
A = P_1' Q R P_2' = P_1' \begin{bmatrix} Q_1 & Q_2 \end{bmatrix} \begin{bmatrix} R_1 \\ 0 \end{bmatrix} P_2' = P_1' Q_1 R_1 P_2'
where
P_1 and P_2 are permutation matrices,
Q = \prod_{j = 1}^{n} H_j
is an m \times m orthogonal matrix
(Q_1 contains the first n column vectors)
equal to the product of n Householder matrices H_j, and
R is an m \times n upper trapezoidal matrix
(R_1 contains the first n row vectors and is
upper triangular).
qrR(qr, complete = FALSE, backPermute = TRUE, row.names = TRUE)
qr |
an object of class |
complete |
a logical indicating if |
backPermute |
a logical indicating if |
row.names |
a logical indicating if |
The method for qr.Q does not return Q but rather the
(also orthogonal) product P_1' Q. This behaviour
is algebraically consistent with the base implementation
(see qr), which can be seen by noting that
qr.default in base does not pivot rows, constraining
P_1 to be an identity matrix. It follows that
qr.Q(qr.default(x)) also returns P_1' Q.
Similarly, the methods for qr.qy and qr.qty multiply
on the left by P_1' Q and Q' P_1
rather than Q and Q'.
It is wrong to expect the values of qr.Q (or qr.R,
qr.qy, qr.qty) computed from “equivalent”
sparse and dense factorizations
(say, qr(x) and qr(as(x, "matrix")) for x
of class dgCMatrix) to compare equal.
The underlying factorization algorithms are quite different,
notably as they employ different pivoting strategies,
and in general the factorization is not unique even for fixed
P_1 and P_2.
On the other hand, the values of qr.X, qr.coef,
qr.fitted, and qr.resid are well-defined, and
in those cases the sparse and dense computations should
compare equal (within some tolerance).
The method for qr.R is a simple wrapper around qrR,
but not back-permuting by default and never giving row names.
It did not support backPermute = TRUE until Matrix
1.6-0, hence code needing the back-permuted result should
call qrR if Matrix >= 1.6-0 is not known.
Dim, Dimnamesinherited from virtual class
MatrixFactorization.
betaa numeric vector of length Dim[2],
used to construct Householder matrices; see V below.
Van object of class dgCMatrix
with Dim[2] columns. The number of rows nrow(V)
is at least Dim[1] and at most Dim[1]+Dim[2].
V is lower trapezoidal, and its column vectors generate the
Householder matrices H_j that compose the orthogonal
Q factor. Specifically, H_j is constructed as
diag(Dim[1]) - beta[j] * tcrossprod(V[, j]).
Ran object of class dgCMatrix
with nrow(V) rows and Dim[2] columns.
R is the upper trapezoidal R factor.
p, q0-based integer vectors of length
nrow(V) and Dim[2], respectively,
specifying the permutations applied to the rows and columns of
the factorized matrix. q of length 0 is valid and
equivalent to the identity permutation, implying no column pivoting.
Using R syntax, the matrix P_1 A P_2
is precisely A[p+1, q+1]
(A[p+1, ] when q has length 0).
Class QR, directly.
Class MatrixFactorization, by class
QR, distance 2.
Objects can be generated directly by calls of the form
new("sparseQR", ...), but they are more typically obtained
as the value of qr(x) for x inheriting from
sparseMatrix (often dgCMatrix).
determinantsignature(from = "sparseQR", logarithm = "logical"):
computes the determinant of the factorized matrix A
or its logarithm.
expand1signature(x = "sparseQR"):
see expand1-methods.
expand2signature(x = "sparseQR"):
see expand2-methods.
qr.Qsignature(qr = "sparseQR"):
returns as a dgeMatrix either
P_1' Q or P_1' Q_1,
depending on optional argument complete. The default
is FALSE, indicating P_1' Q_1.
qr.Rsignature(qr = "sparseQR"):
qrR returns R, R_1,
R P2', or R_1 P2',
depending on optional arguments complete and
backPermute. The default in both cases is FALSE,
indicating R_1, for compatibility with base.
The class of the result in that case is
dtCMatrix. In the other three cases,
it is dgCMatrix.
qr.Xsignature(qr = "sparseQR"):
returns A as a dgeMatrix,
by default. If m > n and optional argument
ncol is greater than n, then the result
is augmented with P_1' Q J, where
J is composed of columns (n+1) through
ncol of the m \times m identity matrix.
qr.coefsignature(qr = "sparseQR", y = .):
returns as a dgeMatrix or vector
the result of multiplying y on the left by
P_2 R_1^{-1} Q_1' P_1.
qr.fittedsignature(qr = "sparseQR", y = .):
returns as a dgeMatrix or vector
the result of multiplying y on the left by
P_1' Q_1 Q_1' P_1.
qr.residsignature(qr = "sparseQR", y = .):
returns as a dgeMatrix or vector
the result of multiplying y on the left by
P_1' Q_2 Q_2' P_1.
qr.qtysignature(qr = "sparseQR", y = .):
returns as a dgeMatrix or vector
the result of multiplying y on the left by
Q' P_1.
qr.qysignature(qr = "sparseQR", y = .):
returns as a dgeMatrix or vector
the result of multiplying y on the left by
P_1' Q.
solvesignature(a = "sparseQR", b = .):
see solve-methods.
Davis, T. A. (2006). Direct methods for sparse linear systems. Society for Industrial and Applied Mathematics. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1137/1.9780898718881")}
Golub, G. H., & Van Loan, C. F. (2013). Matrix computations (4th ed.). Johns Hopkins University Press. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.56021/9781421407944")}
Class dgCMatrix.
Generic function qr from base,
whose default method qr.default “defines”
the S3 class qr of dense QR factorizations.
qr-methods for methods defined in Matrix.
Generic functions expand1 and expand2.
The many auxiliary functions for QR factorizations:
qr.Q, qr.R, qr.X,
qr.coef, qr.fitted, qr.resid,
qr.qty, qr.qy, and qr.solve.
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))
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.