SCM

SCM Repository

[matrix] Diff of /pkg/R/Matrix.R
ViewVC logotype

Diff of /pkg/R/Matrix.R

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 2112, Mon Feb 18 08:24:46 2008 UTC revision 2115, Sat Feb 23 09:23:17 2008 UTC
# Line 45  Line 45 
45            function(x, ...) as.logical(as.vector(x)))            function(x, ...) as.logical(as.vector(x)))
46    
47  setMethod("cov2cor", signature(V = "Matrix"),  setMethod("cov2cor", signature(V = "Matrix"),
48            function(V) as(cov2cor(as(V, "matrix")), "dpoMatrix"))            function(V) { ## was as(cov2cor(as(V, "matrix")), "dpoMatrix"))
49                  r <- V
50                  p <- (d <- dim(V))[1]
51                  if(p != d[2]) stop("'V' is not a square matrix")
52                  Is <- sqrt(1/diag(V)) # diag( 1/sigma_i )
53                  if(any(!is.finite(Is)))
54                      warning("diag(.) had 0 or NA entries; non-finite result is doubtful")
55                  Is <- Diagonal(x = Is)
56                  r <- Is %*% V %*% Is
57                  r[cbind(1L:p,1L:p)] <- 1 # exact in diagonal
58                  as(forceSymmetric(r), "dpoMatrix")
59              })
60    
61  ## "base" has an isSymmetric() S3-generic since R 2.3.0  ## "base" has an isSymmetric() S3-generic since R 2.3.0
62  setMethod("isSymmetric", signature(object = "symmetricMatrix"),  setMethod("isSymmetric", signature(object = "symmetricMatrix"),
# Line 109  Line 120 
120  ##        "!" is in ./not.R  ##        "!" is in ./not.R
121    
122    
123  Matrix <-  Matrix <- function (data = NA, nrow = 1, ncol = 1, byrow = FALSE,
124      function (data = NA, nrow = 1, ncol = 1, byrow = FALSE, dimnames = NULL,                      dimnames = NULL, sparse = NULL, forceCheck = FALSE)
               sparse = NULL, forceCheck = FALSE)  
125  {  {
126      sparseDefault <- function(m) prod(dim(m)) > 2*sum(isN0(as(m, "matrix")))      sparseDefault <- function(m) prod(dim(m)) > 2*sum(isN0(as(m, "matrix")))
127    
# Line 120  Line 130 
130          class(data) <- "matrix" # "matrix" first for S4 dispatch          class(data) <- "matrix" # "matrix" first for S4 dispatch
131      if(is.null(sparse1 <- sparse) && (i.M || is(data, "matrix")))      if(is.null(sparse1 <- sparse) && (i.M || is(data, "matrix")))
132          sparse <- sparseDefault(data)          sparse <- sparseDefault(data)
133        sM <- FALSE
134      doDN <- TRUE      doDN <- TRUE
135      if (i.M) {      if (i.M) {
136          if(!missing(nrow) || !missing(ncol)|| !missing(byrow))          if(!missing(nrow) || !missing(ncol)|| !missing(byrow))
# Line 139  Line 149 
149              ## Matrix(0, ...) : always sparse unless "sparse = FALSE":              ## Matrix(0, ...) : always sparse unless "sparse = FALSE":
150              if(is.null(sparse)) sparse1 <- sparse <- TRUE              if(is.null(sparse)) sparse1 <- sparse <- TRUE
151              i.M <- sM <- TRUE              i.M <- sM <- TRUE
152                isSym <- nrow == ncol
153              ## will be sparse: do NOT construct full matrix!              ## will be sparse: do NOT construct full matrix!
154              data <- new(if(is.numeric(data)) "dgTMatrix" else              data <- new(paste(if(is.numeric(data)) "d" else
155                          if(is.logical(data)) "lgTMatrix" else                                if(is.logical(data)) "l" else
156                          stop("invalid 'data'"),                          stop("invalid 'data'"),
157                                  if(isSym) "s" else "g", "CMatrix", sep=''),
158                            p = rep.int(0L, ncol+1L),
159                          Dim = as.integer(c(nrow,ncol)),                          Dim = as.integer(c(nrow,ncol)),
160                          Dimnames = if(is.null(dimnames)) list(NULL,NULL)                          Dimnames = if(is.null(dimnames)) list(NULL,NULL)
161                          else dimnames)                          else dimnames)
# Line 170  Line 183 
183          isTri <- isTriangular(data)          isTri <- isTriangular(data)
184      isDiag <- isSym # cannot be diagonal if it isn't symmetric      isDiag <- isSym # cannot be diagonal if it isn't symmetric
185      if(isDiag)      if(isDiag)
186          isDiag <- isDiagonal(data)          isDiag <- !isTRUE(sparse1) && isDiagonal(data)
187    
188      ## Find proper matrix class 'cl'      ## try to coerce ``via'' virtual classes
189      cl <-      if(isDiag) { ## diagonal is preferred to sparse !
190          if(isDiag && !isTRUE(sparse1))          data <- as(data, "diagonalMatrix")
191              "diagonalMatrix" # -> will automatically check for type          isSym <- FALSE
192          else {      } else if(sparse && !sM)
             ## consider it's type  
             ctype <-  
                 if(is(data,"Matrix")) class(data)  
                 else {  
                     if("complex" == (ctype <- typeof(data)))  
                         "z" else ctype  
                 }  
             ctype <- substr(ctype, 1,1) # "d", "l", "i" or "z"  
             if(ctype == "z")  
                 stop("complex matrices not yet implemented in Matrix package")  
             if(ctype == "i") {  
                 warning("integer matrices not yet implemented in 'Matrix'; ",  
                         "using 'double' ones'")  
                 ctype <- "d"  
             }  
             paste(ctype,  
                   if(sparse) {  
                       if(isSym) "sCMatrix" else  
                       if(isTri) "tCMatrix" else "gCMatrix"  
                   } else { ## dense  
                       if(isSym) "syMatrix" else  
                       if(isTri) "trMatrix" else "geMatrix"  
                   }, sep="")  
         }  
   
     ## Can we coerce and be done?  
     if(!canCoerce(data,cl)) { ## try to coerce ``via'' virtual classes  
         if(sparse && !sM)  
193              data <- as(data, "sparseMatrix")              data <- as(data, "sparseMatrix")
194          else if(!sparse && !is(data, "denseMatrix"))      else if(!sparse) {
195            if(i.M) { ## data is 'Matrix'
196                if(!is(data, "denseMatrix"))
197              data <- as(data, "denseMatrix")              data <- as(data, "denseMatrix")
198          if(isTri && !is(data, "triangularMatrix"))          } else { ## data is "matrix" (and result "dense" -> go via "general"
199              data <- as(data, "triangularMatrix")              ctype <- typeof(data)
200          else if(isSym && !is(data, "symmetricMatrix"))              if (ctype == "complex")
201              data <- as(data, "symmetricMatrix")                  stop("complex matrices not yet implemented in Matrix package")
202                if (ctype == "integer") ## integer Matrices not yet implemented
203                    storage.mode(data) <- "double"
204                data <- new(paste(.M.kind(data), "geMatrix", sep=''),
205                            Dim = dim(data),
206                            Dimnames = .M.DN(data),
207                            x = c(data))
208            }
209      }      }
210      ## now coerce in any case .. maybe producing sensible error message:  
211      as(data, cl)      if(isTri && !is(data, "triangularMatrix")) {
212            data <- if(attr(isTri,"kind") == "L") tril(data) else triu(data)
213                                            #was as(data, "triangularMatrix")
214        } else if(isSym && !is(data, "symmetricMatrix"))
215            data <- forceSymmetric(data) #was as(data, "symmetricMatrix")
216    
217        data
218  }  }
219    
220  ## Methods for operations where one argument is numeric  ## Methods for operations where one argument is numeric
# Line 480  Line 480 
480              ans              ans
481    
482          } else { ## non-sparse : dense          } else { ## non-sparse : dense
483                ##---- NEVER happens:  'denseMatrix' has its own setMethod(.) !
484              message("m[ <ij-matrix> ]: inefficiently indexing single elements")              message("m[ <ij-matrix> ]: inefficiently indexing single elements")
485              i1 <- ij[,1]              i1 <- ij[,1]
486              i2 <- ij[,2]              i2 <- ij[,2]
487              ## very inefficient for large m -- FIXME! --              ## very inefficient for large m
488              unlist(lapply(seq_len(m), function(j) x[i1[j], i2[j]]))              unlist(lapply(seq_len(m), function(j) x[i1[j], i2[j]]))
489          }          }
490      } else { # 1 <= m <= 3      } else { # 1 <= m <= 3

Legend:
Removed from v.2112  
changed lines
  Added in v.2115

root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge