15

I can not find any function or package to calculate the null space or (QR decomposition) of a bigmatrix (from library(bigmemory)) in R. For example:

library(bigmemory)

a <- big.matrix(1000000, 1000, type='double', init=0)

I tried the following but got the errors shown. How can I find the null space of a bigmemory object?

a.qr <- Matrix::qr(a)
# Error in as.vector(data) : 
#   no method for coercing this S4 class to a vector
q.null <- MASS::Null(a)
# Error in as.vector(data) : 
#   no method for coercing this S4 class to a vector
9
  • Do any of these work ?qr, or ?Matrix::qr, or ?MASS::Null
    – user20650
    Commented Sep 16, 2017 at 11:55
  • Yes. I do, but these functions don't work for a bigmatrix (S4 class) or I couldn't use them for big matrices. I can use these functions only for regular matrices, not for bigmatrices.
    – Mahin
    Commented Sep 16, 2017 at 12:12
  • okay, I wasnt sure if you had a big <space> matrix or bigmatrix ;). Currently, your question is off topic as it directly ask for a package recommendation, and in its present state it may get closed. But it is interesting. Could you edit your question with further details please. For example, could you add a small example of a bigmatrix (including any packages used), illustrate how standard tools dont work, and maybe ask for an alternative. Thanks
    – user20650
    Commented Sep 16, 2017 at 12:22
  • 1
    The bigalgebra package has started some methods, but the QR functions are incomplete, however, this fork , github.com/cdeterman/bigalgebra , has added QR functionality. It does give a warning Warning: This is not advised. - you could ask author why
    – user20650
    Commented Sep 16, 2017 at 14:52
  • 1
    I have a bigmatrix. For example: library(bigmemory) a<-big.matrix(1000000, 1000, type='double', init=5) options(bigmemory.allow.dimnames=TRUE) But I can not find its nullspace or QR decompostion. Thanks.
    – Mahin
    Commented Sep 18, 2017 at 7:53

2 Answers 2

10
+200

If you want to compute the full SVD of the matrix, you can use package bigstatsr to perform computations by block. A FBM stands for a Filebacked Big Matrix and is an object similar to a filebacked big.matrix object of package bigmemory.

library(bigstatsr)
options(bigstatsr.block.sizeGB = 0.5)

# Initialize FBM with random numbers
a <- FBM(1e6, 1e3)
big_apply(a, a.FUN = function(X, ind) {
  X[, ind] <- rnorm(nrow(X) * length(ind))
  NULL
}, a.combine = 'c')

# Compute t(a) * a
K <- big_crossprodSelf(a, big_scale(center = FALSE, scale = FALSE))

# Get v and d where a = u * d * t(v) the SVD of a
eig <- eigen(K[])
v <- eig$vectors
d <- sqrt(eig$values)

# Get u if you need it. It will be of the same size of u
# so that I store it as a FBM.
u <- FBM(nrow(a), ncol(a))
big_apply(u, a.FUN = function(X, ind, a, v, d) {
  X[ind, ] <- sweep(a[ind, ] %*% v, 2, d, "/")
  NULL
}, a.combine = 'c', block.size = 50e3, ind = rows_along(u),
a = a, v = v, d = d)

# Verification
ind <- sample(nrow(a), 1000)
all.equal(a[ind, ], tcrossprod(sweep(u[ind, ], 2, d, "*"), v))

This takes approximately 10 minutes on my computer.

6
  • Hi. Thanks for your answer, but I need to calculate the "full" SVD for an m×n matrix A (A=u.d.t(v)), where u is an m×m matrix, d is an m×n matrix and v is an n×n matrix. I really need the matrix v to get into the null space, but the matrix elements get from your proposed method are different from what I want. For example see: math.stackexchange.com/questions/1771013/null-space-from-svd
    – Mahin
    Commented Sep 25, 2017 at 13:53
  • 1
    @Mahin From my answer, you get u is m x n, d is diagonal n x n (here just the diagonal is given) and v is n x n like you want. It is full in the sense that you can perfectly reconstruct a from these 3 matrices.
    – F. Privé
    Commented Sep 25, 2017 at 15:53
  • 3
    @Mahin From what I understand, the full SVD that you talk about is just adding some (useless) columns to u. v should be the same and you can verify it by comparing svd(mat)$v and svd(mat, nu = nrow(mat), nv = ncol(mat))$v for any matrix with more rows than columns.
    – F. Privé
    Commented Sep 26, 2017 at 9:44
  • 2
    @ F. Privé You are right. I was wrong. In both cases, svd(mat)$v and svd(mat, nu = nrow(mat), nv = ncol(mat))$v, the matrix elements are the same, and I can now calculate the nullspace from v. Thank you very much.
    – Mahin
    Commented Sep 26, 2017 at 10:57
  • Is it possible to have a 1e6x1e6 bigmatrix and to calculate its SVD by these commands in R, on your computer? I do not know how I can work for this size of matrix. Thanks.
    – Mahin
    Commented Oct 2, 2017 at 16:07
2

@Mahon @user20650 @F.Privė For clarity I pinged the bigmemory team and asked

Essentially, is there an implementation of the QR function (QR Decomposition) that works with big memory matrixes?

I felt it useful to get clarity on the original question asked. @F.Privė - nice answer. Hopefully your answer, and their response will help guide people in the future. Their response below:

Thanks for the note. There is not currently an implementation of the qr decomposition. Ideally, you would implement this using Householder reflections (if the matrix is dense) or Givens rotations (if it is sparse).

The irlba package is compatible with bigmemory. It provides a truncated singular value decomposition. So, if your matrix is relatively sparse, you could truncate at the rank of the matrix. This is probably your best option. If you don't know the rank then you can use the package to update the truncation iteratively.

Please note that if your matrix is (tall and skinny or short and fat) then the SO solution is OK. However, anytime you resort to calculating the cross-product you lose some numerical stability. This can be an issue if you are planning on inverting the matrix.

4
  • Hello. Thank you very much for your guide. But I do not want to make a Linear Model Fitting. I just want to do a complete QR decomposition (or full SVD) of a big.matrix and then get into the null space of it. I run the 'RcppEigen' package but this package did not give me the QR decomposition of mat.obj and only gave the Linear Model Fitting. Thanks.
    – Mahin
    Commented Sep 23, 2017 at 8:06
  • 2
    Thanks techno, useful comments. So the answer really is that there is no QR decomposition, or null space calculation out of the box, and so, c++ code will need to be written.
    – user20650
    Commented Sep 27, 2017 at 18:37
  • @user20650 - Yes, that is the current situation - unfortunate. Michael K. was super helpful here (A.k.a BigMemory). That said I do like the F.Privė alternate approach, a package I plan to take a look at in more detail. Commented Sep 27, 2017 at 18:44
  • @ techno Thanks for your attention. I'm currently doing the F.Privė's solution, but I hope that the relevant functions will soon be written.
    – Mahin
    Commented Sep 30, 2017 at 14:59

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.