How do I obtain the -true- QR decomposition of a sparse matrix in R?

I have a simple question: say I have a sparse matrix A in R. If I code

QR<-qr(A)
Q<-qr.Q(QR)
R<-qr.R(QR)

I do not obtain the actual QR decomposition of A, but rather the matrices are such that PAP_=QR where P and P_ are some permutation matrices (that are given). In most instances this is probably not an issue, but for my purposes I believe that I absolutely need the 'real' matrices Q and R such that A=QR, with R being upper triangular. As far as I'm aware, this is not available for sparse matrices? At the moment, I am considering converting my matrices to dense matrices, and then applying the qr() function (because for dense matrices, there appear to be no such permutations); however, that seems a bit stupid. I suppose hand-coding a QR factorization is far too complex, no? Can anyone help with this issue?

P.S. I'm new to stack overflow, I'm usually in the math forums; can I not write in Latex here?

Edit: I should mention that the above applies to the qr decomposition in the matrix package. Apparently, I can use the base qr function, but numerically this seems to not be very nice. (It may be my best option though?)

Solution

The qr method for class dgCMatrix factorizes an m-by-n matrix A (with m >= n) as:

P1 * A * P2 = Q * R    <==>    A = P1' * Q * R * P2'

where P1 and P2 are permutation matrices. At C level, there is an option to disable column pivoting, so that P2 is an identity matrix. Then you have A = Q~ * R, where Q~ = P1' * Q is orthogonal (because P1 and Q are orthogonal), and you would obtain Q~ with qr.Q and R with qr.R.

Unfortunately, the option is not exposed at R level in the current Matrix version 1.5-4. To disable column pivoting from R, you must call the C code directly (which outside of this example I strongly discourage):

library(Matrix)
set.seed(1)
x <- rsparsematrix(10L, 10L, 0.5)
stopifnot(is(x, "dgCMatrix"), validObject(x),
          packageVersion("Matrix") == "1.5-4")

qr.x <- qr(x)
all.equal(as(qr.Q(qr.x) %*% qr.R(qr.x), "matrix"), as(x, "matrix"))
## [1] "Mean relative difference: 1.059398"

qr.x <- .Call(Matrix:::dgCMatrix_QR, x, 0L, FALSE)
all.equal(as(qr.Q(qr.x) %*% qr.R(qr.x), "matrix"), as(x, "matrix"))
## [1] TRUE

It so happens that I've spent quite a bit of time this past month refactoring the factorization code in Matrix. In Matrix 1.6-0 (not yet released), you will have the option to disable column pivoting, with qr(x, order = 0L).

A side effect is that the documentation for all matrix factorization classes and methods will be very much improved. So, people can look forward to that ...