7dccb248b397869dd8223dcd8f76ed2d5d4ad31a
[loomr.git] / R / loom.R
1 #' @import hdf5r
2 #' @importFrom R6 R6Class
3 NULL
4
5 #' A class for loom files
6 #'
7 #' @docType class
8 #' @name loom-class
9 #' @rdname loom-class
10 #' @return Object of \code{\link{R6::R6Class}} to generate \code{loom} objects
11 #' @format An \code{\link{R6::R6Class}} object
12 #' @seealso \code{\link{hdf5r::H5File}}
13 #'
14 #' @field version Version of loomR object was created under
15 #' @field shape Shape of \code{/matrix} in columns (cells) by rows (genes)
16 #' @field chunksize Chunks set for this dataset in columns (cells) by rows (genes)
17 #' @field matrix The main data matrix, stored as columns (cells) by rows (genes)
18 #' @field layers Additional data matricies, the same shape as \code{/matrix}
19 #' @field col.attrs Extra information about cells
20 #' @field row.attrs Extra information about genes
21 #'
22 #' @section Methods:
23 #' \describe{
24 #'   \item{\code{add.layer(layer)}}{Add a data layer to this loom file, must be in column (cells) by row (genes) orientation}
25 #'   \item{\code{add.attribute(attribute, MARGIN)}}{
26 #'     Add extra information to this loom file where
27 #'     \code{attribute} is a named list where each element is a vector that is as long as one dimension of \code{/matrix} and
28 #'     \code{MARGIN} is either 1 for cells or 2 for genes
29 #'   }
30 #'   \item{\code{add.row.attribute(attribute)}}{A wrapper for \code{add.attribute(attribute, MARGIN = 2)}}
31 #'   \item{\code{add.col.attribute(attribute)}}{A wrapper for \code{add.attribute(attribute, MARGIN = 1)}}
32 #'   \item{\code{add.meta.data(meta.data)}}{A wrapper for \code{add.attribute(attribute, MARGIN = 1)}}
33 #'   \item{\code{batch.scan(chunk.size, MARGIN, index.use, dataset.use)}, \code{batch.next()}}{
34 #'     Scan a dataset in the loom file from \code{index.use[1]} to \code{index.use[2]}, iterating by \code{chunk.size}.
35 #'     \code{dataset.use} can be the name, not \code{group/name}, unless the name is present in multiple groups.
36 #'     If the dataset is in col_attrs, pass \code{MARGIN = 1}; if in row_attrs, pass \code{MARGIN = 2}.
37 #'     Otherwise, pass \code{MARGIN = 1} to iterate on cells or \code{MARGIN = 2} to iterate on genes.
38 #'     \code{chunk.size} defaults to \code{self$chunksize}, \code{MARGIN} defaults to 1,
39 #'     \code{index.use} defaults to \code{1:self$shape[MARGIN]}, \code{dataset.use} defaults to 'matrix'
40 #'   }
41 #' }
42 #'
43 #' @importFrom iterators nextElem
44 #' @importFrom utils packageVersion
45 #' @importFrom itertools ihasNext ichunk hasNext
46 #'
47 #' @export
48 #'
49 loom <- R6Class(
50   # The loom class
51   # Based on the H5File class from hdf5r
52   # Not clonable (no loom$clone method), doesn't make sense since all data is on disk, not in memory
53   # Yes to portability, other packages can subclass the loom class
54   # Class is locked, other fields and methods cannot be added
55   classname = 'loom',
56   inherit = hdf5r::H5File,
57   cloneable = FALSE,
58   portable = TRUE,
59   lock_class = TRUE,
60   # Public fields and methods
61   # See above for documentation
62   public = list(
63     # Fields
64     version = NULL,
65     shape = NULL,
66     chunksize = NULL,
67     matrix = NULL,
68     layers = NULL,
69     col.attrs = NULL,
70     row.attrs = NULL,
71     # Methods
72     initialize = function(filename = NULL, mode = c('a', 'r', 'r+', 'w', 'w-'), ...) {
73       # If the file exists, run validation steps
74       do.validate <- file.exists(filename) && !(mode %in% c('w', 'w+'))
75       super$initialize(filename = filename, mode = mode, ...)
76       if (do.validate) {
77         # Run the validation steps
78         validateLoom(object = self)
79         # Store /matrix and the shape of /matrix
80         self$matrix <- self[['matrix']]
81         self$shape <- self[['matrix']]$dims
82         # Store the chunk size
83         chunks <- h5attr(x = self, which = 'chunks')
84         chunks <- gsub(pattern = '(', replacement = '', x = chunks, fixed = TRUE)
85         chunks <- gsub(pattern = ')', replacement = '', x = chunks, fixed = TRUE)
86         chunks <- unlist(x = strsplit(x = chunks, split = ','))
87         self$chunksize <- as.integer(x = chunks)
88         # Store version information
89         self$version <- as.character(x = tryCatch(
90           # Try getting a version
91           # If it doesn't exist, can we write to the file?
92           # If so, store the version as this version of loomR
93           # Otherwise, store the version as NA_character_
94           expr = h5attr(x = self, which = 'version'),
95           error = function(e) {
96             if (mode != 'r') {
97               version <- packageVersion(pkg = 'loomR')
98               h5attr(x = self, which = 'version') <- version
99             } else {
100               version <- NA_character_
101             }
102             return(version)
103           }
104         ))
105         # Load layers
106         private$load_layers()
107         # Load attributes
108         private$load_attributes(MARGIN = 1) # Cells (col_attrs)
109         private$load_attributes(MARGIN = 2) # Genes (row_attrs)
110       } else {
111         # Assume new HDF5 file
112         self$version <- as.character(x = packageVersion(pkg = 'loomR'))
113       }
114     },
115     add.layer = function(layers) {
116       # Value checking
117       if (!is.list(x = layers) || is.null(x = names(x = layers))) {
118         stop("'layers' must be a named list")
119       }
120       if (is.null(x = self$shape)) {
121         stop(private$err_msg)
122       }
123       # Add layers
124       for (i in 1:length(x = layers)) {
125         if (!is.matrix(x = layers[[i]])) {
126           layers[[i]] <- as.matrix(x = layers[[i]])
127         }
128         do.transpose <- FALSE
129         if (any(dim(x = layers[[i]]) != self$shape)) {
130           if (all(rev(x = dim(x = layers[[i]])) == self$shape)) {
131             do.transpose <- TRUE
132           } else {
133             stop(paste(
134               "All layers must have",
135               self$shape[1],
136               "rows for cells and",
137               self$shape[2],
138               "columns for genes"
139             ))
140           }
141         }
142         if (do.transpose) {
143           self[['layers', names(x = layers)[i]]] <- t(x = layers[[i]])
144         } else {
145           self[['layers', names(x = layers)[i]]] <- layers[[i]]
146         }
147       }
148       self$flush()
149       private$load_layers()
150       invisible(x = self)
151     },
152     add.attribute = function(attribute, MARGIN) {
153       # Value checking
154       if (!is.list(x = attribute) || is.null(x = names(x = attribute))) {
155         stop("'attribute' must be a named list")
156       }
157       if (length(x = attribute) > 1) {
158         for (i in attribute) {
159           if (!is.vector(x = attribute)) {
160             stop("All attributes must be one-dimensional vectors")
161           }
162         }
163       }
164       if (length(x = which(x = names(x = attribute) != '')) != length(x = attribute)) {
165         stop("Not all attributes had names provided")
166       }
167       if (!MARGIN %in% c(1, 2)) {
168         stop("'MARGIN' must be 1 or 2")
169       }
170       # Add the attributes as datasets for our MARGIN's group
171       if (is.null(x = self$shape)) {
172         stop(private$err_msg)
173       }
174       grp.name <- c('col_attrs', 'row_attrs')[MARGIN]
175       grp <- self[[grp.name]]
176       for (i in 1:length(x = attribute)) {
177         if (length(attribute[[i]]) != self$shape[MARGIN])
178           stop(paste(
179             "All",
180             switch(EXPR = MARGIN, '1' = 'cell', '2' = 'gene'),
181             "attributes must be of length",
182             self$shape[MARGIN]
183           ))
184         grp[[names(x = attribute)[i]]] <- attribute[[i]]
185       }
186       self$flush()
187       gc(verbose = FALSE)
188       # Load the attributes for this margin
189       private$load_attributes(MARGIN = MARGIN)
190       invisible(x = self)
191     },
192     # Add attributes for genes
193     add.row.attribute = function(attribute) {
194       self$add.attribute(attribute = attribute, MARGIN = 2)
195       invisible(x = self)
196     },
197     # Add attributes for cells
198     add.col.attribute = function(attribute) {
199       self$add.attribute(attribute = attribute, MARGIN = 1)
200       invisible(x = self)
201     },
202     # Add metadata, follows cells
203     add.meta.data = function(meta.data) {
204       self$add.col.attribute(attribute = meta.data)
205       invisible(x = self)
206     },
207     # Batch scan
208     batch.scan = function(
209       chunk.size = NULL,
210       MARGIN = 1,
211       index.use = NULL,
212       dataset.use = 'matrix'
213     ) {
214       if (is.null(x = private$it) || !grepl(pattern = dataset.use, x = private$iter.dataset)) {
215         # Check the existence of the dataset
216         private$iter.dataset <- grep(
217           pattern = dataset.use,
218           x = list.datasets(object = self),
219           value = TRUE
220         )
221         if (length(x = private$iter.dataset) != 1) {
222           stop(paste0("Cannot find dataset '", dataset.use, "' in the loom file"))
223         }
224         # Check the margin
225         if (!(MARGIN %in% c(1, 2))) {
226           stop("MARGIN must be 1 (cells) or 2 (genes)")
227         } else {
228           private$iter.margin <- MARGIN
229         }
230         if (is.null(x = chunk.size)) {
231           chunk.size <- self$chunksize[MARGIN]
232         }
233         # Set the indices to use
234         if (is.null(x = index.use)) {
235           index.use <- c(1, self$shape[MARGIN])
236         } else if (length(x = index.use) == 1) {
237           index.use <- index.use <- 1:index.use
238         } else {
239           index.use <- c(index.use[1], index.use[2])
240         }
241         index.use[1] <- max(1, index.use[1])
242         index.use[2] <- min(index.use[2], self$shape[MARGIN])
243         if (index.use[1] > index.use[2]) {
244           stop(paste0(
245             "Starting index (",
246             index.use[1],
247             ") must be lower than the ending index (",
248             index.use[2],
249             ")"
250           ))
251         }
252         # Setup our iterator
253         private$it <- ihasNext(iterable = ichunk(
254           iterable = index.use[1]:index.use[2],
255           chunkSize = chunk.size
256         ))
257         private$iter.index <- c(index.use[1], ceiling(x = index.use[2] / chunk.size))
258       }
259       # Return the times we iterate, this isn't important, we only need the length of this vector
260       return(private$iter.index[1]:private$iter.index[2])
261     },
262     batch.next = function(return.data = TRUE) {
263       # Ensure that we have a proper iterator
264       if (!'hasNext.ihasNext' %in% suppressWarnings(expr = methods(class = class(x = private$it)))) {
265         stop("Please setup the iterator with self$batch.scan")
266       }
267       # Do the iterating
268       if (hasNext(obj = private$it)) {
269         # Get the indices for this chunk
270         chunk.indices <- unlist(x = nextElem(obj = private$it))
271         if (return.data) {
272           # If we're working with a matrix dataset, ensure chunking on proper axis
273           if (private$iter.dataset == 'matrix' || grepl(pattern = 'layers', x = private$iter.dataset)) {
274             to.return <- switch(
275               EXPR = private$iter.margin,
276               '1' = self[[private$iter.dataset]][chunk.indices, ],
277               '2' = self[[private$iter.dataset]][, chunk.indices]
278             )
279           } else {
280             # Otherwise, iterating over an attribute (1 dimensional)
281             to.return <- self[[private$iter.dataset]][chunk.indices]
282           }
283         } else {
284           to.return <- chunk.indices
285         }
286         # Determine if we reset the iterator
287         if (!hasNext(obj = private$it)) {
288           private$reset_batch()
289         }
290         return(to.return)
291       } else {
292         # Just in case
293         private$reset_batch()
294         return(NULL)
295       }
296     },
297     apply = function() {}
298   ),
299   # Private fields and methods
300   # @field err_msg A simple error message if this object hasn't been created with loomR::create or loomR::connect
301   # @field it Iterator object for batch.scan and batch.next
302   # @field iter.dataset Dataset for iterating on
303   # @field iter.margin Margin to iterate over
304   # @field iter.index # Index values for iteration
305   # \describe{
306   #   \item{\code{load_attributes(MARGIN)}}{Load attributes of a given MARGIN into \code{self$col.attrs} or \code{self$row.attrs}}
307   #   \item{\code{load_layers()}}{Load layers into \code{self$layers}}
308   #   \item{\code{reset_batch()}}{Reset the batch iterator fields}
309   # }
310   private = list(
311     # Fields
312     err_msg = "This loom object has not been created with either loomR::create or loomR::connect, please use these function to create or connect to a loom file",
313     it = NULL,
314     iter.dataset = NULL,
315     iter.margin = NULL,
316     iter.index = NULL,
317     # Methods
318     load_attributes = function(MARGIN) {
319       attribute <- switch(
320         EXPR = MARGIN,
321         '1' = 'col_attrs',
322         '2' = 'row_attrs',
323         stop('Invalid attribute dimension')
324       )
325       group <- self[[attribute]]
326       attributes <- unlist(x = lapply(
327         X = names(x = group),
328         FUN = function(x) {
329           d <- list(group[[x]])
330           names(x = d) <- x
331           return(d)
332         }
333       ))
334       if (MARGIN == 1) {
335         self$col.attrs <- attributes
336       } else if (MARGIN == 2) {
337         self$row.attrs <- attributes
338       }
339     },
340     load_layers = function() {
341       self$layers <- unlist(x = lapply(
342         X = names(x = self[['layers']]),
343         FUN = function(n) {
344           d <- c(self[['layers', n]])
345           names(x = d) <- n
346           return(d)
347         }
348       ))
349     },
350     reset_batch = function() {
351       private$it <- NULL
352       private$iter.dataset <- NULL
353       private$iter.margin <- NULL
354       private$iter.index <- NULL
355     }
356   )
357 )
358
359 #' Create a loom object
360 #'
361 #' @param filename The name of the new loom file
362 #' @param data The data for \code{/matrix}, should be cells as rows and genes as columns
363 #' @param gene.attrs A named list of vectors with extra data for genes, each vector must be as long as the number of genes in \code{data}
364 #' @param cell.attrs A named list of vectors with extra data for cells, each vector must be as long as the number of cells in \code{data}
365 #' @param chunk.dims A one- or two-length integer vector of chunksizes for \code{/matrix}, defaults to 'auto' to automatically determine chunksize
366 #'
367 #' @return A connection to a loom file
368 #'
369 #' @importFrom utils packageVersion
370 #'
371 #' @seealso \code{\link{loom-class}}
372 #'
373 #' @export
374 #'
375 create <- function(
376   filename,
377   data,
378   gene.attrs = NULL,
379   cell.attrs = NULL,
380   layers = NULL,
381   chunk.dims = 'auto'
382 ) {
383   if (file.exists(filename)) {
384     stop(paste('File', file, 'already exists!'))
385   }
386   if (!is.matrix(x = data)) {
387     data <- as.matrix(x = data)
388   }
389   if (length(x = chunk.dims) > 2 || length(x = chunk.dims) < 1) {
390     stop("'chunk.dims' must be a one- or two-length integer vector or 'auto'")
391   } else if (length(x = chunk.dims == 1)) {
392     if (!grepl(pattern = '^auto$', x = chunk.dims, perl = TRUE)) {
393       chunk.dims <- rep.int(x = as.integer(x = chunk.dims), times = 2)
394     }
395   } else {
396     chunk.dims <- as.integer(x = chunk.dims)
397   }
398   new.loom <- loom$new(filename = filename, mode = 'w-')
399   # Create the matrix
400   new.loom$create_dataset(
401     name = 'matrix',
402     robj = data,
403     chunk_dims = chunk.dims
404   )
405   new.loom$matrix <- new.loom[['matrix']]
406   new.loom$shape <- new.loom[['matrix']]$dims
407   # Groups
408   new.loom$create_group(name = 'layers')
409   new.loom$create_group(name = 'row_attrs')
410   new.loom$create_group(name = 'col_attrs')
411   # Check for the existance of gene or cell names
412   if (!is.null(x = colnames(x = data))) {
413     new.loom$add.row.attribute(attribute = list('gene_names' = colnames(x = data)))
414   }
415   if (!is.null(x = rownames(x = data))) {
416     new.loom$add.col.attribute(attribute = list('cell_names' = rownames(x = data)))
417   }
418   # Store some constants as HDF5 attributes
419   h5attr(x = new.loom, which = 'version') <- new.loom$version
420   h5attr(x = new.loom, which = 'chunks') <- paste0(
421     '(',
422     paste(new.loom[['matrix']]$chunk_dims, collapse = ', '),
423     ')'
424   )
425   # Add layers
426   if (!is.null(x = layers)) {
427     new.loom$add.layer(layer = layers)
428   }
429   if (!is.null(x = gene.attrs)) {
430     new.loom$add.row.attribute(attribute = gene.attrs)
431   }
432   if (!is.null(x = cell.attrs)) {
433     new.loom$add.col.attribute(attribute = cell.attrs)
434   }
435   # Set last bit of information
436   chunks <- new.loom[['matrix']]$chunk_dims
437   chunks <- gsub(pattern = '(', replacement = '', x = chunks, fixed = TRUE)
438   chunks <- gsub(pattern = ')', replacement = '', x = chunks, fixed = TRUE)
439   chunks <- unlist(x = strsplit(x = chunks, split = ','))
440   new.loom$chunksize <- as.integer(x = chunks)
441   # Return the connection
442   return(new.loom)
443 }
444
445 #' Validate a loom object
446 #'
447 #' @param object A loom object
448 #'
449 #' @return None, errors out if object is an invalid loom connection
450 #'
451 #' @seealso \code{\link{loom-class}}
452 #'
453 #' @export
454 #'
455 validateLoom <- function(object) {
456   # A loom file is a specific HDF5
457   # We need a dataset in /matrix that's a two-dimensional dense matrix
458   root.datasets <- list.datasets(object = object, path = '/', recursive = FALSE)
459   if (length(x = root.datasets) != 1) {
460     stop("There can only be one dataset at the root of the loom file")
461   }
462   if (root.datasets != 'matrix') {
463     stop("The root dataset must be called '/matrix'")
464   }
465   # There must be groups called '/col_attrs', '/row_attrs', and '/layers'
466   required.groups <- c('row_attrs', 'col_attrs', 'layers')
467   dim.matrix <- object[[root.datasets]]$dims # Columns x Rows
468   names(x = dim.matrix) <- required.groups[c(2, 1)]
469   root.groups <- list.groups(object = object, path = '/', recursive = FALSE)
470   group.msg <- paste0(
471     "There can only be three groups in the loom file: '",
472     paste(required.groups, collapse = "', '"),
473     "'"
474   )
475   if (length(x = root.groups) != 3) {
476     stop(group.msg)
477   }
478   if (!all(required.groups %in% root.groups)) {
479     stop(group.msg)
480   }
481   unlist(x = sapply(
482     X = required.groups[1:2],
483     FUN = function(group) {
484       if (length(x = list.groups(object = object[[group]], recursive = FALSE)) > 0) {
485         stop(paste("Group", group, "cannot have subgroups"))
486       }
487       if (length(x = list.attributes(object = object[[group]])) > 0) {
488         stop(paste("Group", group, "cannot have subattributes"))
489       }
490       for (dataset in list.datasets(object = object[[group]])) {
491         if (object[[paste(group, dataset, sep = '/')]]$dims != dim.matrix[group]) {
492           stop(paste("All datasets in group", group, "must be of length", required.groups[group]))
493         }
494       }
495     }
496   ))
497   for (dataset in list.datasets(object = object[['/layers']])) {
498     if (any(object[[paste('layers', dataset, sep = '/')]]$dims != dim.matrix)) {
499       stop(paste("All datasets in '/layers' must be", dim.matrix[1], 'by', dim.matrix[2]))
500     }
501   }
502 }
503
504 #' Connect to a loom file
505 #'
506 #' @param filename The loom file to connect to
507 #' @param mode How do we connect to it? Pass 'r' for read-only or 'r+' for read/write
508 #'
509 #' @return A loom file connection
510 #'
511 #' @seealso \code{\link{loom-class}}
512 #'
513 #' @export
514 #'
515 connect <- function(filename, mode = "r") {
516   if (!(mode %in% c('r', 'r+'))) {
517     stop("'mode' must be one of 'r' or 'r+'")
518   }
519   new.loom <- loom$new(filename = filename, mode = mode)
520   return(new.loom)
521 }
522
523 CreateLoomFromSeurat <- function(object, filename) {
524   object.data=t(object@raw.data[rownames(object@data),object@cell.names])
525   object.meta.data=object@meta.data
526   row_attrs=list(); col_attrs=list()
527   gene.names=colnames(object.data)
528   object.meta.data$ident = object@ident
529   object.meta.data$CellID = object@cell.names
530   for(i in 1:ncol(object.meta.data)) {
531     col_attrs[[colnames(object.meta.data)[i]]]=object.meta.data[,i]
532   }
533   row_attrs[["Gene"]]=gene.names
534   create(filename,object.data,gene.attrs = row_attrs, cell.attrs = col_attrs)
535 }
536
537 #' Map a function or a series of functions over a loom file
538 #'
539 #' @param X A loom object
540 #' @param MARGIN The dimmension to map over, pass 1 for cells or 2 for genes
541 #' @param FUN A function to map to the loom file
542 #' @param chunk.size Chunk size to use, defaults to \code{loomfile$chunksize[MARGIN]}
543 #' @param index.use Indices of the dataset to use, defaults to \code{1:loomfile$shape[MARGIN]}
544 #' @param dataset.use Dataset to use, defauts to 'matrix'
545 #' @param display.progress Display a progress bar
546 #' @param expected ...
547 #' @param ... Extra parameters for FUN
548 #'
549 #' @return The results of the map
550 #'
551 #' @importFrom utils txtProgressBar setTxtProgressBar
552 #'
553 #' @export
554 #'
555 map <- function(
556   X,
557   MARGIN = 1,
558   FUN,
559   chunk.size = NULL,
560   index.use = NULL,
561   dataset.use = 'matrix',
562   display.progress = TRUE,
563   expected = NULL,
564   ...
565 ) {
566   if (!inherits(x = X, what = 'loom')) {
567     stop("map only works on loom objects")
568   }
569   if (!inherits(x = FUN, what = 'function')) {
570     stop("FUN must be a function")
571   }
572   # Check for existance of dataset
573   if (!any(grepl(pattern = dataset.use, x = list.datasets(object = X)))) {
574     stop(paste("Cannot find dataset", dataset.use, "in the loom file"))
575   }
576   # Figure out if we're returning a vector or matrix
577   full.dataset <- grep(
578     pattern = dataset.use,
579     x = list.datasets(object = X),
580     value = TRUE
581   )
582   results.matrix <- TRUE
583   if (grepl(pattern = 'col_attrs', x = full.dataset)) {
584     MARGIN <- 1
585     results.matrix <- FALSE
586   } else if (grepl(pattern = 'row_attrs', x = full.dataset)) {
587     MARGIN <- 2
588     results.matrix <- FALSE
589   }
590   if (!is.null(x = expected)) {
591     results.matrix <- switch(
592       EXPR = expected,
593       'vector' = FALSE,
594       'matrix' = TRUE,
595       stop("'expected' must be one of 'matrix', 'vector', or NULL")
596     )
597   }
598   # Determine the shape of our results
599   if (!(MARGIN %in% c(1, 2))) {
600     stop("MARGIN must be either 1 (cells) or 2 (genes)")
601   }
602   if (is.null(x = index.use)) {
603     index.use <- 1:X$shape[MARGIN]
604   } else if (length(x = index.use) == 1) {
605     index.use <- c(1, index.use)
606   }
607   index.use[1] <- max(1, index.use[1])
608   index.use[2] <- min(index.use[2], X$shape[MARGIN])
609   batch <- X$batch.scan(
610     chunk.size = chunk.size,
611     MARGIN = MARGIN,
612     index.use = index.use,
613     dataset.use = dataset.use
614   )
615   # Create our results holding object
616   if (results.matrix) {
617     switch(
618       EXPR = MARGIN,
619       '1' = results <- matrix(ncol = X$shape[2], nrow = length(x = batch)),
620       '2' = results <- matrix(nrow = X$shape[1], ncol = length(x = batch))
621     )
622   } else {
623     results <- vector(length = length(x = batch))
624   }
625   if (display.progress) {
626     pb <- txtProgressBar(char = '=', style = 3)
627   }
628   for (i in 1:length(x = batch)) {
629     if (results.matrix) {
630       if (MARGIN == 1) {
631         results[i, ] <- FUN(X$batch.next(), ...)
632       } else if (MARGIN == 2) {
633         results[, i] <- FUN(X$batch.next(), ...)
634       }
635     } else {
636       results[]
637     }
638     if (display.progress) {
639       setTxtProgressBar(pb = pb, value = i / length(x = batch))
640     }
641   }
642   invisible(x = NULL)
643 }
644
645 #need to comment
646 #need to add progress bar
647 #but otherwise, pretty cool
648 #for paul to try :
649 # f <- connect("~/Downloads/10X43_1.loom")
650 # mean_var = map(f,f_list = c(mean,var),chunksize = 5000)
651 # nGene <- map(f, f_list = function(x) length(which(x>0)), MARGIN = 2)
652 map <- function(self, f_list = list(mean, var), MARGIN=1, chunksize=1000, selection) {
653   n_func = length(f_list)
654   if (n_func == 1) {
655     f_list <- list(f_list)
656   }
657   if (MARGIN == 1) {
658     results <- list()
659     for (j in 1:n_func) {
660       results[[j]] <- numeric(0)
661     }
662     rows_per_chunk <- chunksize
663     ix <- 1
664     while (ix <= self@shape[1]) {
665       rows_per_chunk <- min(rows_per_chunk, self@shape[1] - ix + 1)
666       chunk <- self["matrix"][ix:(ix + rows_per_chunk - 1), ]
667       for (j in 1:n_func) {
668         new_results <- apply(chunk, 1, FUN = f_list[[j]])
669         results[[j]] <- c(results[[j]], new_results)
670       }
671       ix <- ix + chunksize
672     }
673   }
674   if (MARGIN == 2) {
675     results <- list()
676     for (j in 1:n_func) {
677       results[[j]] <- numeric(0)
678     }
679     cols_per_chunk <- chunksize
680     ix <- 1
681     while (ix <= self@shape[2]) {
682       cols_per_chunk <- min(cols_per_chunk, self@shape[2] - ix + 1)
683       chunk <- self["matrix"][, ix:(ix + cols_per_chunk - 1)]
684       for (j in 1:n_func) {
685         new_results <- apply(chunk, 2, FUN = f_list[[j]])
686         results[[j]] <- c(results[[j]], new_results)
687       }
688       ix <- ix + chunksize
689     }
690   }
691   if (n_func == 1) {
692     results <- results[[1]]
693   }
694   return(results)
695 }