Symmetric Matrix classes
On Nov 26, 2009, at 11:00 PM, Gad Abraham wrote:
On Fri, Nov 27, 2009 at 2:22 PM, David Winsemius <dwinsemius at comcast.net
wrote: On Nov 26, 2009, at 9:45 PM, Gad Abraham wrote:
Hi,
I'd like to store large covariance matrices using Matrix classes.
dsyMatrix seems like the right one, but I want to specify just the
upper/lower triangle and diagonal and not have to instantiate a huge
n^2 vector just for the sake of having half of it ignored:
Dumb example:
M <- new("dsyMatrix", uplo="U", x=rnorm(1e4),
Dim=as.integer(c(100, 100)))
diag(M) <- 1
This doesn't work:
M <- new("dsyMatrix", uplo="U", x=0, Dim=as.integer(c(100, 100)))
Per help
Slots
uplo:A character object indicating if the upper triangle ("U") or
the lower
triangle ("L") is stored. At present only the lower triangle form is
allowed.
This "works" ( at least to the extent of not producing an error
message)
after addressing a couple of other error messages that were helpful.
M <- new("dsyMatrix", uplo="L", x=as.double(1:(100*100)),
Dim=as.integer(c(100, 100)))
I'm not sure that I would have expected the result, though it makes a
certain amount of sense after considering the contradictory
requirements:
M <- new("dsyMatrix", uplo="L", x=as.double(1:(10*10)),
Dim=as.integer(c(10, 10)))
M
10 x 10 Matrix of class "dsyMatrix"
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,] 1 2 3 4 5 6 7 8 9 10
[2,] 2 12 13 14 15 16 17 18 19 20
[3,] 3 13 23 24 25 26 27 28 29 30
[4,] 4 14 24 34 35 36 37 38 39 40
[5,] 5 15 25 35 45 46 47 48 49 50
[6,] 6 16 26 36 46 56 57 58 59 60
[7,] 7 17 27 37 47 57 67 68 69 70
[8,] 8 18 28 38 48 58 68 78 79 80
[9,] 9 19 29 39 49 59 69 79 89 90
[10,] 10 20 30 40 50 60 70 80 90 100
But I don't want to specify x=foo where foo is of length n^2, because that uses a lot of memory (I have 20000 by 20000 matrices). I just want to be able to fill in one triangle and diagonal, like you can with regular non-Matrix matrices: M <- matrix(0, m, n) gdata::lowerTriangle(M) <- x diag(M) <- y
Then maybe you should pick a Matrix class that suits your desires. You are the one who said he wanted a symmetric class. There is a triangular sparse class. I doubt you will get the lowerTriangle function to work on these instances of Matrix. When I use the base function lower.tri to index that 10 x 10 sparse matrix I get an error message that makes me think you will need to populate these Matrices element by element. M[lower.tri] <- 0 Error in .local(x, i, j, ..., value) : not-yet-implemented 'Matrix[<-' method diag(M) <- 1 # did work.
-- Gad Abraham PhD Student, Dept. CSSE and NICTA The University of Melbourne Parkville 3010, Victoria, Australia email: gabraham at csse.unimelb.edu.au web: http://www.csse.unimelb.edu.au/~gabraham
David Winsemius, MD Heritage Laboratories West Hartford, CT