# SCM Repository

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

# Diff of /pkg/R/Matrix.R

revision 871, Fri Aug 26 17:26:49 2005 UTC revision 1174, Mon Jan 16 20:03:48 2006 UTC
# Line 1  Line 1
1  #### Toplevel ``virtual'' class "Matrix"  #### Toplevel ``virtual'' class "Matrix"
2
3  ## probably not needed eventually:  ## ## probably not needed eventually:
4  setAs(from = "ddenseMatrix", to = "matrix",  ## setAs(from = "ddenseMatrix", to = "matrix",
5        function(from) {  ##       function(from) {
6            if(length(d <- dim(from)) != 2) stop("dim(.) has not length 2")  ##        if(length(d <- dim(from)) != 2) stop("dim(.) has not length 2")
7            array(from@x, dim = d, dimnames = dimnames(from))  ##        array(from@x, dim = d, dimnames = dimnames(from))
8        })  ##       })

## private function to be used as show() method possibly more than once
prMatrix <- function(object) {
d <- dim(object)
cl <- class(object)
cat(sprintf('%d x %d Matrix of class "%s"\n', d[1], d[2], cl))
m <- as(object, "matrix")
maxp <- getOption("max.print")
if(prod(d) <= maxp) print(m)
else { ## d[1] > maxp / d[2] >= nr :
nr <- maxp %/% d[2]
n2 <- ceiling(nr / 2)
cat("\n ..........\n\n")
print(tail(m, max(1, nr - n2)))
}
## DEBUG: cat("str(.):\n") ; str(object)
invisible(object)# as print() S3 methods do
}

setMethod("show", signature(object = "ddenseMatrix"), prMatrix)

##- ## FIXME: The following is only for the "dMatrix" objects that are not
##- ##        "dense" nor "sparse" -- i.e. "packed" ones :
##- ## But these could be printed better -- "." for structural zeros.
##- setMethod("show", signature(object = "dMatrix"), prMatrix)
##- ## and improve this as well:
##- setMethod("show", signature(object = "pMatrix"), prMatrix)
##- ## this should now be superfluous [keep for safety for the moment]:
setMethod("show", signature(object = "Matrix"), prMatrix)
9
10  ## should propagate to all subclasses:  ## should propagate to all subclasses:
11  setMethod("as.matrix", signature(x = "Matrix"), function(x) as(x, "matrix"))  setMethod("as.matrix", signature(x = "Matrix"), function(x) as(x, "matrix"))
12    ## for 'Matrix' objects, as.array() should be equivalent:
13    setMethod("as.array",  signature(x = "Matrix"), function(x) as(x, "matrix"))
14
15    ## slow "fall back" method {subclasses should have faster ones}:
16    setMethod("as.vector", signature(x = "Matrix", mode = "missing"),
17              function(x) as.vector(as(x, "matrix")))
18
19
20    ## Note that isSymmetric is *not* exported
21    ## but that "base" has an isSymmetric() S3-generic since R 2.3.0
22    setMethod("isSymmetric", signature(object = "symmetricMatrix"),
23              function(object,tol) TRUE)
24    setMethod("isSymmetric", signature(object = "triangularMatrix"),
25              ## TRUE iff diagonal:
26              function(object,tol) isDiagonal(object))
27    if(paste(R.version\$major, R.version\$minor, sep=".") < "2.3")
28        ## need a "matrix" method as in R 2.3 and later
29        setMethod("isSymmetric", signature(object = "matrix"),
30                  function(object, tol = 100*.Machine\$double.eps, ...)
31              {
32                  ## pretest: is it square?
33                  d <- dim(object)
34                  if(d[1] != d[2]) return(FALSE)
35                  test <-
36                      if(is.complex(object))
37                          all.equal.numeric(object, Conj(t(object)), tol = tol, ...)
38                      else              # numeric, character, ..
39                          all.equal(object, t(object), tol = tol, ...)
40                  isTRUE(test)
41              })
42
43
44    setMethod("isTriangular", signature(object = "triangularMatrix"),
45              function(object,tol) TRUE)
46    setMethod("isTriangular", signature(object = "matrix"),
47              .is.triangular)
48
49    setMethod("isDiagonal", signature(object = "matrix"), .is.diagonal)
50
51
52
53  setMethod("dim", signature(x = "Matrix"),  setMethod("dim", signature(x = "Matrix"),
54            function(x) x@Dim, valueClass = "integer")            function(x) x@Dim, valueClass = "integer")
# Line 62  Line 72
72            function(obj) { obj@Dimnames <- list(NULL,NULL); obj})            function(obj) { obj@Dimnames <- list(NULL,NULL); obj})
73
74  Matrix <-  Matrix <-
75      function (data = NA, nrow = 1, ncol = 1, byrow = FALSE, dimnames = NULL)      function (data = NA, nrow = 1, ncol = 1, byrow = FALSE, dimnames = NULL,
76                  sparse = NULL)
77  {  {
78      if (is(data, "Matrix")) return(data)      sparseDefault <- function(m)
79      if (is.matrix(data)) { val <- data }          prod(dim(m)) > 2*sum(as(m, "matrix") != 0)
80      else { ## cut & paste from "base::matrix" :
81        i.M <- is(data, "Matrix")
82        if(is.null(sparse) && (i.M || is(data, "matrix")))
83            sparse <- sparseDefault(data)
84
85        if (i.M) {
86            sM <- is(data,"sparseMatrix")
87            if((sparse && sM) || (!sparse && !sM))
88                return(data)
89            ## else : convert  dense <-> sparse -> at end
90        }
91        else if (!is.matrix(data)) { ## cut & paste from "base::matrix" :
92          if (missing(nrow))          if (missing(nrow))
93              nrow <- ceiling(length(data)/ncol)              nrow <- ceiling(length(data)/ncol)
94          else if (missing(ncol))          else if (missing(ncol))
95              ncol <- ceiling(length(data)/nrow)              ncol <- ceiling(length(data)/nrow)
96          val <- .Internal(matrix(data, nrow, ncol, byrow))          data <- .Internal(matrix(data, nrow, ncol, byrow))
97          dimnames(val) <- dimnames          if(is.null(sparse))
98                sparse <- sparseDefault(data)
99            dimnames(data) <- dimnames
100        }
101
102        ## 'data' is now a "matrix" or "Matrix"
103
104        ## check for symmetric / triangular / diagonal :
105        isSym <- isSymmetric(data)
106        if((isTri <- !isSym))
107            isTri <- isTriangular(data)
108        isDiag <- isSym # cannot be diagonal if it isn't symmetric
109        if(isDiag)
110            isDiag <- isDiagonal(data)
111
112    ### TODO: Compare with as.Matrix() and its tests in ./dgeMatrix.R
113
114        ## Find proper matrix class 'cl'
115        cl <-
116            if(isDiag)
117                "diagonalMatrix" # -> will automatically check for type
118            else {
119                ## consider it's type
120                ctype <-
121                    if(is(data,"Matrix")) class(data)
122                    else {
123                        if("complex" == (ctype <- typeof(data)))
124                            "z" else ctype
125                    }
126                ctype <- substr(ctype, 1,1) # "d", "l", "i" or "z"
127                if(ctype == "z")
128                    stop("complex matrices not yet implemented in Matrix package")
129                if(ctype == "i") {
130                    warning("integer matrices not yet implemented in 'Matrix'; ",
131                            "using 'double' ones'")
132                    ctype <- "d"
133      }      }
134      as(val, "dgeMatrix")              paste(ctype,
135                      if(sparse) {
136                          if(isSym) "sCMatrix" else
137                          if(isTri) "tCMatrix" else "gCMatrix"
138                      } else { ## dense
139                          if(isSym) "syMatrix" else
140                          if(isTri) "trMatrix" else "geMatrix"
141                      }, sep="")
142            }
143
144        ## Now coerce and return
145        as(data, cl)
146  }  }
147
148  ## Methods for operations where one argument is numeric  ## Methods for operations where one argument is numeric
# Line 97  Line 165
165  setMethod("solve", signature(a = "Matrix", b = "numeric"),  setMethod("solve", signature(a = "Matrix", b = "numeric"),
166            function(a, b, ...) callGeneric(a, as.matrix(b)))            function(a, b, ...) callGeneric(a, as.matrix(b)))
167
168    ## bail-out methods in order to get better error messages
169    setMethod("%*%", signature(x = "Matrix", y = "Matrix"),
170              function (x, y)
171              stop(gettextf('not-yet-implemented method for <%s> %%*%% <%s>',
172                            class(x), class(y))))
173
174    setMethod("crossprod", signature(x = "Matrix", y = "ANY"),
175              function (x, y = NULL) .bail.out.2(.Generic, class(x), class(y)))
176    setMethod("crossprod", signature(x = "ANY", y = "Matrix"),
177              function (x, y = NULL) .bail.out.2(.Generic, class(x), class(y)))
178
179    ## There are special sparse methods; this is a "fall back":
180    setMethod("kronecker", signature(X = "Matrix", Y = "ANY",
181                                     FUN = "ANY", make.dimnames = "ANY"),
182              function(X, Y, FUN, make.dimnames, ...) {
183                  X <- as(X, "matrix") ; Matrix(callGeneric()) })
184    setMethod("kronecker", signature(X = "ANY", Y = "Matrix",
185                                     FUN = "ANY", make.dimnames = "ANY"),
186              function(X, Y, FUN, make.dimnames, ...) {
187                  Y <- as(Y, "matrix") ; Matrix(callGeneric()) })
188
189
190    setMethod("t", signature(x = "Matrix"),
191              function(x) .bail.out.1(.Generic, class(x)))
192
193    ## Group Methods
194    setMethod("+", signature(e1 = "Matrix", e2 = "missing"), function(e1) e1)
195    ## "fallback":
196    setMethod("-", signature(e1 = "Matrix", e2 = "missing"),
197              function(e1) {
198                  warning("inefficient method used for \"- e1\"")
199                  0-e1
200              })
201
202    ## bail-outs:
203    setMethod("Compare", signature(e1 = "Matrix", e2 = "Matrix"),
204              function(e1, e2) {
205                  d <- dimCheck(e1,e2)
206                  .bail.out.2(.Generic, class(e1), class(e2))
207              })
208    setMethod("Compare", signature(e1 = "Matrix", e2 = "ANY"),
209              function(e1, e2) .bail.out.2(.Generic, class(e1), class(e2)))
210    setMethod("Compare", signature(e1 = "ANY", e2 = "Matrix"),
211              function(e1, e2) .bail.out.2(.Generic, class(e1), class(e2)))
212
213
214
215  ### --------------------------------------------------------------------------  ### --------------------------------------------------------------------------
216  ###  ###
217  ### Subsetting "["  and  ### Subsetting "["  and
218  ### SubAssign  "[<-" : The "missing" cases can be dealt with here, "at the top":  ### SubAssign  "[<-" : The "missing" cases can be dealt with here, "at the top":
219
220    ## Using "index" for indices should allow
221    ## integer (numeric), logical, or character (names!) indices :
222
223  ## "x[]":  ## "x[]":
224  setMethod("[", signature(x = "Matrix",  setMethod("[", signature(x = "Matrix",
225                           i = "missing", j = "missing", drop = "ANY"),                           i = "missing", j = "missing", drop = "ANY"),
# Line 109  Line 227
227  ## missing 'drop' --> 'drop = TRUE'  ## missing 'drop' --> 'drop = TRUE'
228  ##                     -----------  ##                     -----------
229  ## select rows  ## select rows
230  setMethod("[", signature(x = "Matrix", i = "numeric", j = "missing",  setMethod("[", signature(x = "Matrix", i = "index", j = "missing",
231                           drop = "missing"),                           drop = "missing"),
232            function(x,i,j, drop) callGeneric(x, i=i, drop= TRUE))            function(x,i,j, drop) callGeneric(x, i=i, drop= TRUE))
233  ## select columns  ## select columns
234  setMethod("[", signature(x = "Matrix", i = "missing", j = "numeric",  setMethod("[", signature(x = "Matrix", i = "missing", j = "index",
235                           drop = "missing"),                           drop = "missing"),
236            function(x,i,j, drop) callGeneric(x, j=j, drop= TRUE))            function(x,i,j, drop) callGeneric(x, j=j, drop= TRUE))
237  setMethod("[", signature(x = "Matrix", i = "numeric", j = "numeric",  setMethod("[", signature(x = "Matrix", i = "index", j = "index",
238                           drop = "missing"),                           drop = "missing"),
239            function(x,i,j, drop) callGeneric(x, i=i, j=j, drop= TRUE))            function(x,i,j, drop) callGeneric(x, i=i, j=j, drop= TRUE))
240
241    ## bail out if any of (i,j,drop) is "non-sense"
242    setMethod("[", signature(x = "Matrix", i = "ANY", j = "ANY", drop = "ANY"),
243              function(x,i,j, drop)
244              stop("invalid or not-yet-implemented 'Matrix' subsetting"))
245
246  ## "FIXME:"  ## "FIXME:"
247  ## How can we get at   A[ ij ]  where ij is (i,j) 2-column matrix?  ## How can we get at   A[ ij ]  where ij is (i,j) 2-column matrix?
248  ##  and                A[ LL ]  where LL is a logical *vector*  ##  and                A[ LL ]  where LL is a logical *vector*
249    ## -> [.data.frame uses nargs() - can we do this in the *generic* ?
250
251
252  ### "[<-" : -----------------  ### "[<-" : -----------------
253
254  ## x[] <- value :  ## x[] <- value :
255  setReplaceMethod("[", signature(x = "Matrix", i = "missing", j = "missing",  setReplaceMethod("[", signature(x = "Matrix", i = "missing", j = "missing",
256                                  value = "vector"),##  double/logical/...                                  value = "index"),##  double/logical/...
257            function (x, value) { x@x <- value ; validObject(x); x })            function (x, value) { x@x <- value ; validObject(x); x })
258
259  ## Otherwise (value is not "vector"): bail out  ## Otherwise (value is not "index"): bail out
260  setReplaceMethod("[", signature(x = "Matrix", i = "ANY", j = "ANY",  setReplaceMethod("[", signature(x = "Matrix", i = "ANY", j = "ANY",
261                                  value = "ANY"),                                  value = "ANY"),
262            function (x, i, j, value) stop("RHS 'value' must be of class \"vector\""))            function (x, i, j, value)
263                     if(!is(value,"index"))
264                     stop("RHS 'value' must be of class \"index\"")
265                     else stop("not-yet-implemented 'Matrix[<-' method"))
266
267
268
269    ## NOTE: the following only works for R 2.2.x (and later) ---
270    ## ----  *and* 'Matrix' must have been *installed* by R >= 2.2.x
271
272    if(paste(R.version\$major, R.version\$minor, sep=".") >= "2.2") {
273
274        ## The trivial methods :
275        setMethod("cbind2", signature(x = "Matrix", y = "NULL"),
276                  function(x, y) x)
277        setMethod("cbind2", signature(x = "Matrix", y = "missing"),
278                  function(x, y) x)
279        setMethod("cbind2", signature(x = "NULL", y="Matrix"),
280                  function(x, y) x)
281
282        setMethod("rbind2", signature(x = "Matrix", y = "NULL"),
283                  function(x, y) x)
284        setMethod("rbind2", signature(x = "Matrix", y = "missing"),
285                  function(x, y) x)
286        setMethod("rbind2", signature(x = "NULL", y="Matrix"),
287                  function(x, y) x)
288
289        ## Makes sure one gets x decent error message for the unimplemented cases:
290        setMethod("cbind2", signature(x = "Matrix", y = "Matrix"),
291                  function(x, y) {
292                      rowCheck(x,y)
293                      stop(gettextf("cbind2() method for (%s,%s) not-yet defined",
294                                    class(x), class(y)))
295                  })
296
297  if(FALSE) ## The following can't work as long as cbind is function(..., *)      ## Use a working fall back {particularly useful for sparse}:
298  setMethod("cbind", signature(a = "Matrix", b = "Matrix"),      ## FIXME: implement rbind2 via "cholmod" for C* and Tsparse ones
299            function(a, b, ...) {      setMethod("rbind2", signature(x = "Matrix", y = "Matrix"),
300                da <- Dim(a)                function(x, y) {
301                db <- Dim(b)                    colCheck(x,y)
302                if(da[1] != db[1])                    t(cbind2(t(x), t(y)))
stop("Matrices must have same number of rows for cbind()ing")
303            })            })
304