Vectorized element-wise division on Sparse Matrices in R

Question

A/B in R performs an element-wise division on the matrix.

However, if I generate a sparse matrix from the Matrix package, and try to divide A/B, I get this error:

> class(N)
[1] "dgCMatrix"
attr(,"package")
[1] "Matrix"
> N/N
Error in asMethod(object) : 
  Cholmod error 'problem too large' at file ../Core/cholmod_dense.c, line 105
> 

Interesting. When the sparse matrix is small in total size, I get this behavior:

> m <- sparseMatrix(i=c(1,2,1,3), j=c(1,1,3,3), x=c(1,2,1,4))
> m/m
3 x 3 Matrix of class "dgeMatrix"
     [,1] [,2] [,3]
[1,]    1  NaN    1
[2,]    1  NaN  NaN
[3,]  NaN  NaN    1
> 

But when it's moderately sized (~ 20000 elements), I get the Cholmod error.

Is there a workaround or a more proper way to do element-wise division on sparse matrices in R?

Answer 1

The problem with element-wise division is that if your matrices are both sparse, then you'll have a lot of Inf and NaN in the result, and these make it dense. That's why you get the out-of-memory errors.

If you want to replace Inf and NaN with zeros in the result, then the solution is relatively easy, you just get the summary() of both matrices and work with the indices and values directly.

You'll need to restrict the A and B index vectors to their intersection and perform the division on that. To get the intersection of index pairs, one can use merge().

Here is a quick and dirty implementation:

# Some example data
A <- sparseMatrix(i=c(1,1,2,3), j=c(1,3,1,3), x=c(1,1,2,3))
B <- sparseMatrix(i=c(3,2,1), j=c(3,2,1), x=c(3,2,1))

sdiv <- function(X, Y, names=dimnames(X)) {
  sX <- summary(X)
  sY <- summary(Y)
  sRes <- merge(sX, sY, by=c("i", "j"))
  sparseMatrix(i=sRes[,1], j=sRes[,2], x=sRes[,3]/sRes[,4],
               dimnames=names)
}

sdiv(A, B)
# 3 x 3 sparse Matrix of class "dgCMatrix"
#           
# [1,] 1 . .
# [2,] . . .
# [3,] . . 1

Thanks to flodel for the suggestion about using summary and merge.