[R] Sparse matrix no longer sparse (Matrix Package)

John Kane jrkrideau at inbox.com
Thu Jul 11 16:56:55 CEST 2013


The message got through but not the attachment. The R help list tends to strip off attachements for security reasons.  Files of types  txt, png, & pdf should get through.

In most cases the accepted method of sending data is to use the dput() function to output a file in the console and then copy and paste the results into your email.

So for file "dat1" one would just use dput(dat1) and paste the results into an email. 

John Kane
Kingston ON Canada


> -----Original Message-----
> From: hdoran at air.org
> Sent: Thu, 11 Jul 2013 09:53:40 +0000
> To: r-help at r-project.org
> Subject: [R] Sparse matrix no longer sparse (Matrix Package)
> 
> I sent this message yesterday with an attachment allowing for
> reproduction
> of the issue. But I think the attachment is preventing the message from
> coming through. If anyone is interested I will forward the attachment
> directly allowing for you to reproduce the issue I observe.
> 
> On 7/10/13 2:38 PM, "Doran, Harold" <HDoran at air.org> wrote:
> 
> >I have zero'd in on what appears to be the issue. This seems to be a bug
> >in Matrix, but I am not sure yet. I am attaching files that would allow
> >others to replicate this with my toy data.
>> 
> >Notice the elements of D1 in the attached data are all integers. It is a
> >sparse, diagonal matrix.
>> 
>>> library(Matrix)
>>> class(D1)
> >[1] "ddiMatrix"
> >attr(,"package")
> >[1] "Matrix"
>> 
> >Now, I find the inverse of the matrix A as follows:
>>> A <- Ir + ZtZ %*% D1
>>> A.Inv <- solve(A, Ir)
>> 
> >Notice now the inverse of A remains a dgCMatrix and it is relatively
> >small in size, only 33424 bytes.
>>> class(A.Inv)
> >[1] "dgCMatrix"
> >attr(,"package")
> >[1] "Matrix"
>> 
>>> object.size(A.Inv)
> >33424 bytes
>> 
> >Now, if I change an element of the matrix D1 to be non-integer, D1 still
> >has the same class as it did before
>> 
>>> D1[1] <- 1.2
>> 
>>> class(D1)
> >[1] "ddiMatrix"
> >attr(,"package")
> >[1] "Matrix"
>> 
> >Now, if I use this new version of D1 in the same calculations as above,
> >notice that A.Inv is no longer a dgCMatrix but instead becomes a
> >dgeMatrix. It then increases from an object of size 33424 bytes to an
> >object of size 2001112 bytes!
>> 
>>> A <- Ir + ZtZ %*% D1
>>> A.Inv <- solve(A, Ir)
>>> class(A.Inv)
> >[1] "dgeMatrix"
> >attr(,"package")
> >[1] "Matrix"
>>> object.size(A.Inv)
> >2001112 bytes
>> 
> >What I desire is that the object A.Inv remain sparse at all times and
> not
> >become dense. But, perhaps there is a reason this change occurs that I
> >don't fully understand.
>> 
> >I can of course coerce it back to a sparse matrix and it reduces back in
> >size.
>>>  object.size(as(A.Inv, 'sparseMatrix'))
> >33424 bytes
>> 
> >I of course recognize it requires more memory to store floating points
> >than integers, but is this large increase on the order of magnitude that
> >seems about right?
>> 
> >Is there a reason the floating point in D1 causes for A.Inv to no longer
> >remain sparse?
>> 
> >Thank you for your help,
> >Harold
>> 
>> 
>> 
>> 
>> 
> >-----Original Message-----
> >From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
> >On Behalf Of Doran, Harold
> >Sent: Wednesday, July 10, 2013 11:42 AM
> >To: r-help at r-project.org
> >Cc: dmbates at gmail.com; 'maechler at stat.math.ethz.ch'
> >Subject: [R] Sparse matrix no longer sparse (Matrix Package)
>> 
> >I have a large function computing an iterative algorithm for fitting
> >mixed linear models. Almost all code relies on functions from the Matrix
> >package. I've come across an issue that I do not believe previously
> >occurred in earlier versions of R or Matrix.
>> 
> >I have a large, sparse matrix, A as
>> 
>>> class(A);dim(A)
> >[1] "dgCMatrix"
> >attr(,"package")
> >[1] "Matrix"
> >[1] 12312 12312
>> 
> >I am in a position where I must find its inverse.  I realize this is
> less
> >than ideal, and I have two ways of doing this
>> 
> >A.Inv <- solve(A, Ir) or just solve(A)
>> 
> >Where Ir is an identity matrix with the same dimensions as A and it is
> >also sparse
>> 
>>> class(Ir)
> >[1] "ddiMatrix"
> >attr(,"package")
> >[1] "Matrix"
>> 
> >The issue, however, is that the inverse of A is converted into a dense
> >matrix and this becomes a huge memory hog, causing the rest of the
> >algorithm to fail. In prior versions this remained as a sparse matrix.
>> 
>>> A.Inv[1:5, 1:5]
> >5 x 5 Matrix of class "dgeMatrix"
>>          [,1]      [,2]      [,3]      [,4]      [,5]
> >[1,] 0.6878713 0.0000000 0.0000000 0.0000000 0.0000000 [2,] 0.0000000
> >0.6718767 0.0000000 0.0000000 0.0000000 [3,] 0.0000000 0.0000000
> >0.5076945 0.0000000 0.0000000 [4,] 0.0000000 0.0000000 0.0000000
> >0.2324122 0.0000000 [5,] 0.0000000 0.0000000 0.0000000 0.0000000
> 0.2139975
>> 
> >I could coerce this matrix to become sparse such as
>> 
>>> AA <- as(A.Inv, 'sparseMatrix')
>>> class(AA)
> >[1] "dgCMatrix"
> >attr(,"package")
> >[1] "Matrix"
>> 
>>> AA[1:5, 1:5]
> >5 x 5 sparse Matrix of class "dgCMatrix"
>> 
> >[1,] 0.6878713 .         .         .         .
> >[2,] .         0.6718767 .         .         .
> >[3,] .         .         0.5076945 .         .
> >[4,] .         .         .         0.2324122 .
> >[5,] .         .         .         .         0.2139975
>> 
> >But I don't think this is best.
>> 
> >So, my question is why is a matrix that is sparse turning into a dense
> >matrix? Can I avoid that and keep it sparse without having to coerce it
> >to be sparse after it is created?
>> 
> >Thank you very much
> >Harold
>> 
>> 
>>> sessionInfo()
> >R version 3.0.1 (2013-05-16)
> >Platform: x86_64-w64-mingw32/x64 (64-bit)
>> 
> >locale:
> >[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United
> >States.1252 [3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C [5]
> >LC_TIME=English_United States.1252
>> 
> >attached base packages:
> >[1] stats     graphics  grDevices utils     datasets  methods   base
>> 
> >other attached packages:
> >[1] lme4_0.999999-2 Matrix_1.0-12   lattice_0.20-15
>> 
> >loaded via a namespace (and not attached):
> >[1] grid_3.0.1   nlme_3.1-109 stats4_3.0.1 tools_3.0.1
>> 
> >	[[alternative HTML version deleted]]
>> 
> >______________________________________________
> >R-help at r-project.org mailing list
> >https://stat.ethz.ch/mailman/listinfo/r-help
> >PLEASE do read the posting guide
> >http://www.R-project.org/posting-guide.html
> >and provide commented, minimal, self-contained, reproducible code.
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

____________________________________________________________
FREE 3D EARTH SCREENSAVER - Watch the Earth right on your desktop!



More information about the R-help mailing list