SCM

SCM Repository

[matrix] Annotation of /pkg/R/sparseMatrix.R
ViewVC logotype

Annotation of /pkg/R/sparseMatrix.R

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1911 - (view) (download)

1 : bates 684 ### Define Methods that can be inherited for all subclasses
2 :    
3 : maechler 925 ### Idea: Coercion between *VIRTUAL* classes -- as() chooses "closest" classes
4 :     ### ---- should also work e.g. for dense-triangular --> sparse-triangular !
5 : maechler 868
6 : maechler 1472 ##-> see als ./dMatrix.R, ./ddenseMatrix.R and ./lMatrix.R
7 : maechler 868
8 : maechler 1472 setAs("ANY", "sparseMatrix", function(from) as(from, "CsparseMatrix"))
9 : maechler 868
10 : maechler 1845 setAs("sparseMatrix", "generalMatrix", as_gSparse)
11 : maechler 1472
12 : maechler 1845 setAs("sparseMatrix", "symmetricMatrix", as_sSparse)
13 :    
14 :     setAs("sparseMatrix", "triangularMatrix", as_tSparse)
15 :    
16 : maechler 1852 spMatrix <- function(nrow, ncol, i,j,x) {
17 :     ## Author: Martin Maechler, Date: 8 Jan 2007, 18:46
18 :     dim <- c(as.integer(nrow), as.integer(ncol))
19 :     ## The conformability of (i,j,x) with itself and with 'dim'
20 :     ## is checked automatically by internal "validObject()" inside new(.):
21 :     kind <- .M.kind(x)
22 :     new(paste(kind, "gTMatrix", sep=''), Dim = dim,
23 :     x = if(kind == "d") as.double(x) else x,
24 :     ## our "Tsparse" Matrices use 0-based indices :
25 :     i = as.integer(i - 1L),
26 :     j = as.integer(j - 1L))
27 :     }
28 :    
29 :    
30 : maechler 871 ## "graph" coercions -- this needs the graph package which is currently
31 :     ## ----- *not* required on purpose
32 :     ## Note: 'undirected' graph <==> 'symmetric' matrix
33 :    
34 : maechler 1271 ## Add some utils that may no longer be needed in future versions of the 'graph' package
35 :     graph.has.weights <- function(g) "weight" %in% names(edgeDataDefaults(g))
36 :    
37 :     graph.wgtMatrix <- function(g)
38 :     {
39 :     ## Purpose: work around "graph" package's as(g, "matrix") bug
40 :     ## ----------------------------------------------------------------------
41 :     ## Arguments: g: an object inheriting from (S4) class "graph"
42 :     ## ----------------------------------------------------------------------
43 :     ## Author: Martin Maechler, based on Seth Falcon's code; Date: 12 May 2006
44 :    
45 :     ## MM: another buglet for the case of "no edges":
46 :     if(numEdges(g) == 0) {
47 :     p <- length(nd <- nodes(g))
48 :     return( matrix(0, p,p, dimnames = list(nd, nd)) )
49 :     }
50 :     ## Usual case, when there are edges:
51 :     has.w <- "weight" %in% names(edgeDataDefaults(g))
52 :     if(has.w) {
53 :     w <- unlist(edgeData(g, attr = "weight"))
54 :     has.w <- any(w != 1)
55 :     } ## now 'has.w' is TRUE iff there are weights != 1
56 :     m <- as(g, "matrix")
57 :     ## now is a 0/1 - matrix (instead of 0/wgts) with the 'graph' bug
58 :     if(has.w) { ## fix it if needed
59 :     tm <- t(m)
60 :     tm[tm != 0] <- w
61 :     t(tm)
62 :     }
63 :     else m
64 :     }
65 :    
66 :    
67 :     setAs("graphAM", "sparseMatrix",
68 : bates 862 function(from) {
69 : maechler 1271 symm <- edgemode(from) == "undirected" && isSymmetric(from@adjMat)
70 :     ## This is only ok if there are no weights...
71 :     if(graph.has.weights(from)) {
72 :     as(graph.wgtMatrix(from),
73 :     if(symm) "dsTMatrix" else "dgTMatrix")
74 :     }
75 :     else { ## no weights: 0/1 matrix -> logical
76 :     as(as(from, "matrix"),
77 : maechler 1548 if(symm) "nsTMatrix" else "ngTMatrix")
78 : maechler 1271 }
79 : bates 862 })
80 : maechler 1271
81 : bates 1476 setAs("graph", "CsparseMatrix",
82 : bates 1479 function(from) as(as(from, "graphNEL"), "CsparseMatrix"))
83 : maechler 687
84 : bates 1474 setAs("graphNEL", "CsparseMatrix",
85 : maechler 1565 function(from) as(as(from, "TsparseMatrix"), "CsparseMatrix"))
86 :    
87 :     setAs("graphNEL", "TsparseMatrix",
88 : maechler 1271 function(from) {
89 : maechler 1565 nd <- nodes(from)
90 : bates 1474 dm <- rep.int(length(nd), 2)
91 : maechler 1271 symm <- edgemode(from) == "undirected"
92 : bates 1474
93 : maechler 1565 if(graph.has.weights(from)) {
94 :     eWts <- edgeWeights(from)
95 :     lens <- unlist(lapply(eWts, length))
96 :     i <- rep.int(0:(dm[1]-1), lens) # column indices (0-based)
97 :     To <- unlist(lapply(eWts, names))
98 : maechler 1832 j <- as.integer(match(To,nd) - 1L) # row indices (0-based)
99 : maechler 1565 ## symm <- symm && <weights must also be symmetric>: improbable
100 :     ## if(symm) new("dsTMatrix", .....) else
101 :     new("dgTMatrix", i = i, j = j, x = unlist(eWts),
102 :     Dim = dm, Dimnames = list(nd, nd))
103 :     }
104 :     else { ## no weights: 0/1 matrix -> logical
105 :     edges <- lapply(from@edgeL[nd], "[[", "edges")
106 :     lens <- unlist(lapply(edges, length))
107 :     ## nnz <- sum(unlist(lens)) # number of non-zeros
108 :     i <- rep.int(0:(dm[1]-1), lens) # column indices (0-based)
109 :     j <- as.integer(unlist(edges) - 1) # row indices (0-based)
110 :     if(symm) { # symmetric: ensure upper triangle
111 :     tmp <- i
112 :     flip <- i > j
113 :     i[flip] <- j[flip]
114 :     j[flip] <- tmp[flip]
115 :     new("nsTMatrix", i = i, j = j, Dim = dm,
116 :     Dimnames = list(nd, nd), uplo = "U")
117 :     } else {
118 :     new("ngTMatrix", i = i, j = j, Dim = dm,
119 :     Dimnames = list(nd, nd))
120 :     }
121 : bates 1474 }
122 : maechler 1271 })
123 : maechler 687
124 : maechler 871 setAs("sparseMatrix", "graph", function(from) as(from, "graphNEL"))
125 :     setAs("sparseMatrix", "graphNEL",
126 : maechler 1271 function(from) as(as(from, "TsparseMatrix"), "graphNEL"))
127 : maechler 908
128 : maechler 1348 Tsp2grNEL <- function(from) {
129 :     d <- dim(from)
130 :     if(d[1] != d[2])
131 : maechler 1513 stop("only square matrices can be used as incidence matrices for graphs")
132 : maechler 1348 n <- d[1]
133 :     if(n == 0) return(new("graphNEL"))
134 :     if(is.null(rn <- dimnames(from)[[1]]))
135 :     rn <- as.character(1:n)
136 :     from <- uniq(from) ## Need to 'uniquify' the triplets!
137 : maechler 908
138 : maechler 1348 if(isSymmetric(from)) { # either "symmetricMatrix" or otherwise
139 :     ##-> undirected graph: every edge only once!
140 :     if(!is(from, "symmetricMatrix")) {
141 :     ## a general matrix which happens to be symmetric
142 :     ## ==> remove the double indices
143 :     from <- tril(from)
144 :     }
145 : maechler 1565 eMode <- "undirected"
146 :     } else {
147 :     eMode <- "directed"
148 : maechler 1348 }
149 : maechler 1565 ## every edge is there only once, either upper or lower triangle
150 : maechler 1832 ft1 <- cbind(rn[from@i + 1L], rn[from@j + 1L])
151 : maechler 1565 ## not yet: graph::ftM2graphNEL(.........)
152 :     ftM2graphNEL(ft1, W = from@x, V= rn, edgemode= eMode)
153 : maechler 871
154 : maechler 1348 }
155 :     setAs("TsparseMatrix", "graphNEL", Tsp2grNEL)
156 : maechler 871
157 : maechler 1348
158 : maechler 868 ### Subsetting -- basic things (drop = "missing") are done in ./Matrix.R
159 : maechler 687
160 : maechler 925 ### FIXME : we defer to the "*gT" -- conveniently, but not efficient for gC !
161 : maechler 687
162 : maechler 925 ## [dl]sparse -> [dl]gT -- treat both in one via superclass
163 :     ## -- more useful when have "z" (complex) and even more
164 : maechler 687
165 : maechler 925 setMethod("[", signature(x = "sparseMatrix", i = "index", j = "missing",
166 : maechler 868 drop = "logical"),
167 : maechler 925 function (x, i, j, drop) {
168 : maechler 1725 cld <- getClassDef(class(x))
169 : maechler 1751 if(!extends(cld, "generalMatrix")) x <- as(x, "generalMatrix")
170 :     viaCl <- paste(.M.kind(x, cld), "gTMatrix", sep='')
171 :     x <- callGeneric(x = as(x, viaCl), i=i, drop=drop)
172 :     ## try_as(x, c(cl, sub("T","C", viaCl)))
173 :     if(is(x, "Matrix") && extends(cld, "CsparseMatrix"))
174 :     as(x, "CsparseMatrix") else x
175 :     })
176 : maechler 687
177 : maechler 925 setMethod("[", signature(x = "sparseMatrix", i = "missing", j = "index",
178 : maechler 868 drop = "logical"),
179 : maechler 925 function (x, i, j, drop) {
180 : maechler 1725 cld <- getClassDef(class(x))
181 : maechler 1751 if(!extends(cld, "generalMatrix")) x <- as(x, "generalMatrix")
182 :     viaCl <- paste(.M.kind(x, cld), "gTMatrix", sep='')
183 :     x <- callGeneric(x = as(x, viaCl), j=j, drop=drop)
184 :     if(is(x, "Matrix") && extends(cld, "CsparseMatrix"))
185 :     as(x, "CsparseMatrix") else x
186 :     })
187 : maechler 868
188 : maechler 925 setMethod("[", signature(x = "sparseMatrix",
189 : maechler 886 i = "index", j = "index", drop = "logical"),
190 : maechler 925 function (x, i, j, drop) {
191 : maechler 1725 cld <- getClassDef(class(x))
192 : maechler 1665 ## be smart to keep symmetric indexing of <symm.Mat.> symmetric:
193 : maechler 1725 doSym <- (extends(cld, "symmetricMatrix") &&
194 : maechler 1665 length(i) == length(j) && all(i == j))
195 : maechler 1751 if(!doSym && !extends(cld, "generalMatrix"))
196 :     x <- as(x, "generalMatrix")
197 : maechler 1725 viaCl <- paste(.M.kind(x, cld),
198 : maechler 1665 if(doSym) "sTMatrix" else "gTMatrix", sep='')
199 :     x <- callGeneric(x = as(x, viaCl), i=i, j=j, drop=drop)
200 : maechler 1725 if(is(x, "Matrix") && extends(cld, "CsparseMatrix"))
201 : maechler 1751 as(x, "CsparseMatrix") else x
202 : maechler 1665 })
203 : maechler 868
204 :    
205 : maechler 1673 ## setReplaceMethod("[", .........)
206 :     ## -> ./Tsparse.R
207 :     ## & ./Csparse.R
208 :     ## FIXME: also for RsparseMatrix
209 : maechler 1226
210 :    
211 : maechler 1714 ## Group Methods
212 : maechler 1226
213 : maechler 1737 setMethod("Math",
214 :     signature(x = "sparseMatrix"),
215 :     function(x) callGeneric(as(x, "CsparseMatrix")))
216 : maechler 868
217 : maechler 1737 ## further group methods -> see ./Ops.R
218 : maechler 1472
219 :    
220 : maechler 1737
221 : maechler 687 ### --- show() method ---
222 :    
223 : maechler 1389 ## FIXME(?) -- ``merge this'' (at least ``synchronize'') with
224 :     ## - - - prMatrix() from ./Auxiliaries.R
225 : maechler 1903 prSpMatrix <- function(x, digits = getOption("digits"),
226 : maechler 1705 maxp = getOption("max.print"), zero.print = ".",
227 : maechler 1903 col.names = FALSE, note.dropping.colnames = TRUE,
228 : maechler 1870 col.trailer = '', align = c("fancy", "right"))
229 : maechler 1705 ## FIXME: prTriang() in ./Auxiliaries.R should also get align = "fancy"
230 : maechler 687 {
231 : maechler 1903 cl <- getClassDef(class(x))
232 : maechler 1737 stopifnot(extends(cl, "sparseMatrix"))
233 : maechler 1903 d <- dim(x)
234 : maechler 1389 if(prod(d) > maxp) { # "Large" => will be "cut"
235 :     ## only coerce to dense that part which won't be cut :
236 :     nr <- maxp %/% d[2]
237 : maechler 1903 m <- as(x[1:max(1, nr), ,drop=FALSE], "Matrix")
238 : maechler 1389 } else {
239 : maechler 1903 m <- as(x, "matrix")
240 : maechler 1389 }
241 : maechler 1737 logi <- extends(cl,"lsparseMatrix") || extends(cl,"nsparseMatrix")
242 : maechler 1315 if(logi)
243 : maechler 1903 cx <- array("N", # or as.character(NA),
244 : maechler 1572 dim(m), dimnames=dimnames(m))
245 : maechler 1673 else { ## numeric (or --not yet-- complex):
246 : maechler 1903 cx <- apply(m, 2, format)
247 :     if(is.null(dim(cx))) {# e.g. in 1 x 1 case
248 :     dim(cx) <- dim(m)
249 :     dimnames(cx) <- dimnames(m)
250 : maechler 1315 }
251 : maechler 687 }
252 : maechler 1903 if(!col.names)
253 :     cx <- emptyColnames(cx, msg.if.not.empty = note.dropping.colnames)
254 : maechler 687 if(is.logical(zero.print))
255 :     zero.print <- if(zero.print) "0" else " "
256 : maechler 1315 if(logi) {
257 : maechler 1903 cx[!m] <- zero.print
258 :     cx[m] <- "|"
259 : maechler 1315 } else { # non logical
260 :     ## show only "structural" zeros as 'zero.print', not all of them..
261 :     ## -> cannot use 'm'
262 : maechler 1903 d <- dim(cx)
263 :     ne <- length(iN0 <- 1L + encodeInd(non0ind(x, cl), nr = d[1]))
264 : maechler 1705 if(0 < ne && ne < prod(d)) {
265 :     align <- match.arg(align)
266 :     if(align == "fancy") {
267 :     fi <- apply(m, 2, format.info) ## fi[3,] == 0 <==> not expo.
268 :     ## now 'format' the zero.print by padding it with ' ' on the right:
269 :     ## case 1: non-exponent: fi[2,] + as.logical(fi[2,] > 0)
270 :     ## the column numbers of all 'zero' entries -- (*large*)
271 : maechler 1832 cols <- 1L + (0:(prod(d)-1L))[-iN0] %/% d[1]
272 : maechler 1705 pad <-
273 :     ifelse(fi[3,] == 0,
274 :     fi[2,] + as.logical(fi[2,] > 0),
275 :     ## exponential:
276 :     fi[2,] + fi[3,] + 4)
277 : maechler 1737 ## now be efficient ; sprintf() is relatively slow
278 :     ## and pad is much smaller than 'cols'; instead of "simply"
279 :     ## zero.print <- sprintf("%-*s", pad[cols] + 1, zero.print)
280 :     if(any(doP <- pad > 0)) {#
281 :     ## only pad those that need padding - *before* expanding
282 :     z.p.pad <- rep.int(zero.print, length(pad))
283 :     z.p.pad[doP] <- sprintf("%-*s", pad[doP] + 1, zero.print)
284 :     zero.print <- z.p.pad[cols]
285 :     }
286 :     else
287 :     zero.print <- rep.int(zero.print, length(cols))
288 : maechler 1705 } ## else "right" : nothing to do
289 :    
290 : maechler 1903 cx[-iN0] <- zero.print
291 : maechler 1705 } else if (ne == 0)# all zeroes
292 : maechler 1903 cx[] <- zero.print
293 : maechler 1315 }
294 : maechler 1870 if(col.trailer != '')
295 : maechler 1903 cx <- cbind(cx, col.trailer, deparse.level = 0)
296 : maechler 1705 ## right = TRUE : cheap attempt to get better "." alignment
297 : maechler 1903 print(cx, quote = FALSE, right = TRUE, max = maxp)
298 :     invisible(x)
299 : maechler 687 }
300 :    
301 : maechler 1903 setMethod("print", signature(x = "sparseMatrix"), prSpMatrix)
302 :    
303 : maechler 687 setMethod("show", signature(object = "sparseMatrix"),
304 :     function(object) {
305 :     d <- dim(object)
306 :     cl <- class(object)
307 :     cat(sprintf('%d x %d sparse Matrix of class "%s"\n', d[1], d[2], cl))
308 :     maxp <- getOption("max.print")
309 :     if(prod(d) <= maxp)
310 : maechler 1389 prSpMatrix(object, maxp = maxp)
311 : maechler 687 else { ## d[1] > maxp / d[2] >= nr : -- this needs [,] working:
312 : maechler 1737
313 : maechler 687 nR <- d[1] # nrow
314 : maechler 1737 useW <- getOption("width") - (format.info(nR)[1] + 3+1)
315 :     ## space for "[<last>,] "
316 : maechler 1845
317 :     ## --> suppress rows and/or columns in printing ...
318 :    
319 : maechler 1737 suppCols <- (d[2] * 2 > useW)
320 :     nc <- if(suppCols) (useW - (1 + 6)) %/% 2 else d[2]
321 : maechler 1870 ## sp+ col.trailer
322 :     col.trailer <- if(suppCols) "......" else ""
323 : maechler 1737 nr <- maxp %/% nc
324 :     suppRows <- (nr < nR)
325 :     if(suppRows) {
326 :     if(suppCols)
327 :     object <- object[ , 1:nc, drop = FALSE]
328 :     n2 <- ceiling(nr / 2)
329 :     prSpMatrix(object[seq_len(min(nR, max(1, n2))), , drop=FALSE],
330 : maechler 1870 col.trailer = col.trailer)
331 : maechler 1737 cat("\n ..............................",
332 :     "\n ..........suppressing rows in show(); maybe adjust 'options(max.print= *)'",
333 :     "\n ..............................\n\n", sep='')
334 :     ## tail() automagically uses "[..,]" rownames:
335 :     prSpMatrix(tail(object, max(1, nr-n2)),
336 : maechler 1870 col.trailer = col.trailer)
337 : maechler 1737 }
338 :     else if(suppCols) {
339 :     prSpMatrix(object[ , 1:nc , drop = FALSE],
340 : maechler 1870 col.trailer = col.trailer)
341 : maechler 1737
342 :     cat("\n .....suppressing columns in show(); maybe adjust 'options(max.print= *)'",
343 :     "\n ..............................\n", sep='')
344 :     }
345 :     else stop("logic programming error in prSpMatrix(), please report")
346 :    
347 : maechler 687 invisible(object)
348 :     }
349 :     })
350 : maechler 886
351 :    
352 : maechler 1852 ## For very large and very sparse matrices, the above show()
353 :     ## is not really helpful; Use summary() as an alternative:
354 :    
355 :     setMethod("summary", signature(object = "sparseMatrix"),
356 :     function(object, ...) {
357 :     d <- dim(object)
358 :     T <- as(object, "TsparseMatrix")
359 :     ## return a data frame (int, int, {double|logical|...}) :
360 :     r <- data.frame(i = T@i + 1L, j = T@j + 1L, x = T@x)
361 :     attr(r, "header") <-
362 : maechler 1855 sprintf('%d x %d sparse Matrix of class "%s", with %d entries',
363 :     d[1], d[2], class(object), nnzero(object))
364 : maechler 1852 ## use ole' S3 technology for such a simple case
365 :     class(r) <- c("sparseSummary", class(r))
366 :     r
367 :     })
368 :    
369 :     print.sparseSummary <- function (x, ...) {
370 :     cat(attr(x, "header"),"\n")
371 :     print.data.frame(x, ...)
372 :     invisible(x)
373 :     }
374 :    
375 : maechler 1108 setMethod("isSymmetric", signature(object = "sparseMatrix"),
376 : maechler 973 function(object, tol = 100*.Machine$double.eps) {
377 : maechler 886 ## pretest: is it square?
378 :     d <- dim(object)
379 :     if(d[1] != d[2]) return(FALSE)
380 :     ## else slower test
381 : maechler 973 if (is(object, "dMatrix"))
382 : maechler 886 ## use gC; "T" (triplet) is *not* unique!
383 : maechler 1659 isTRUE(all.equal(.as.dgC.0.factors( object),
384 :     .as.dgC.0.factors(t(object)), tol = tol))
385 : maechler 973 else if (is(object, "lMatrix"))
386 : maechler 886 ## test for exact equality; FIXME(?): identical() too strict?
387 :     identical(as(object, "lgCMatrix"),
388 :     as(t(object), "lgCMatrix"))
389 : maechler 1548 else if (is(object, "nMatrix"))
390 :     ## test for exact equality; FIXME(?): identical() too strict?
391 :     identical(as(object, "ngCMatrix"),
392 :     as(t(object), "ngCMatrix"))
393 : maechler 886 else stop("not yet implemented")
394 :     })
395 : maechler 1174
396 : maechler 1238
397 : maechler 1659 ## These two are not (yet?) exported:
398 : maechler 1174 setMethod("isTriangular", signature(object = "sparseMatrix"),
399 : maechler 1238 function(object, upper = NA)
400 :     isTriC(as(object, "CsparseMatrix"), upper))
401 : maechler 1174
402 :     setMethod("isDiagonal", signature(object = "sparseMatrix"),
403 :     function(object) {
404 : maechler 1799 d <- dim(object)
405 :     if(d[1] != d[2]) return(FALSE)
406 :     ## else
407 : maechler 1174 gT <- as(object, "TsparseMatrix")
408 :     all(gT@i == gT@j)
409 :     })
410 :    
411 : maechler 1290
412 : maechler 1472 setMethod("diag", signature(x = "sparseMatrix"),
413 :     function(x, nrow, ncol = n) diag(as(x, "CsparseMatrix")))
414 :    
415 : maechler 1845 setMethod("dim<-", signature(x = "sparseMatrix", value = "ANY"),
416 :     function(x, value) {
417 :     if(!is.numeric(value) || length(value) != 2)
418 :     stop("dim(.) value must be numeric of length 2")
419 :     if(prod(dim(x)) != prod(value <- as.integer(value)))
420 :     stop("dimensions don't match the number of cells")
421 :     ## be careful to keep things sparse
422 :     as(spV2M(as(x, "sparseVector"), nrow=value[1], ncol=value[2]),
423 :     class(x))
424 :     })
425 :    
426 :    
427 : bates 1895 lm.fit.sparse <-
428 :     function(x, y, offset = NULL, method = c("qr", "cholesky"),
429 :     tol = 1e-7, singular.ok = TRUE, transpose = FALSE, ...)
430 :     ### Fit a linear model using a sparse QR or a sparse Cholesky factorization
431 :     {
432 :     stopifnot(is(x, "dsparseMatrix"))
433 :     yy <- as.numeric(y)
434 :     if (!is.null(offset)) {
435 :     stopifnot(length(offset) == length(y))
436 :     yy <- yy - as.numeric(offset)
437 :     }
438 :     ans <- switch(as.character(method)[1],
439 :     cholesky = .Call(dgCMatrix_cholsol,
440 :     as(if (transpose) x else t(x), "dgCMatrix"), yy),
441 :     qr = .Call(dgCMatrix_qrsol,
442 :     as(if (transpose) t(x) else x, "dgCMatrix"), yy),
443 :     stop(paste("unknown method", dQuote(method)))
444 :     )
445 :     ans
446 :     }
447 :    
448 : maechler 1911 fac2sparse <- function(from, to = c("d","i","l","n","z"))
449 :     {
450 :     ## factor(-like) --> sparseMatrix {also works for integer, character}
451 :     levs <- levels(fact <- factor(from)) # drop unused levels
452 :     n <- length(fact)
453 :     to <- match.arg(to)
454 :     res <- new(paste(to, "gCMatrix", sep=''))
455 :     res@i <- as.integer(fact) - 1L # 0-based
456 :     res@p <- 0:n
457 :     res@Dim <- c(length(levs), n)
458 :     res@Dimnames <- list(levs, NULL)
459 :     if(to != "n")
460 :     res@x <- rep.int(switch(to,
461 :     "d" = 1., "i" = 1L, "l" = TRUE, "z" = 1+0i),
462 :     n)
463 :     res
464 :     }
465 : bates 1895
466 : maechler 1911 setAs("factor", "sparseMatrix", function(from) fac2sparse(from, to = "d"))

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