Merge branch 'master' of https://github.com/cole-trapnell-lab/cicero-…

…release

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: cicero
 Type: Package
-Title: Precict cis-co-accessibility from single-cell chromatin accessibility data
-Version: 1.11.1
+Title: Predict cis-co-accessibility from single-cell chromatin accessibility data
+Version: 1.13.4
 Authors@R: c(
     person("Hannah", "Pliner", email = "hpliner@uw.edu", role = c("aut", "cre")), 
     person("Cole", "Trapnell", email = "coletrap@uw.edu", role = c("aut")))
@@ -40,10 +40,11 @@ Imports:
     tidyr,
     VGAM (>= 1.0-5),
     utils
-RoxygenNote: 7.1.0
+RoxygenNote: 7.1.2
 Suggests:
     AnnotationDbi (>= 1.38.2),
     knitr,
+    markdown,
     rmarkdown,
     rtracklayer (>= 1.36.6),
     testthat,

R/aggregate.R

@@ -186,8 +186,8 @@ aggregate_by_cell_bin <- function(cds, group_col) {
                                      lowerDetectionLimit=0))
 
     compart_cds <- detectGenes(compart_cds, min_expr=0.1)
-    compart_cds <- estimateSizeFactors(compart_cds)
-    compart_cds <- estimateDispersions(compart_cds)
+    compart_cds <- estimateSizeFactorsSimp(compart_cds)
+    compart_cds <- estimateDispersionsSimp(compart_cds)
 
     fData(compart_cds)$use_for_ordering <- FALSE
 

R/runCicero.R

@@ -19,14 +19,19 @@
 #' @param size_factor_normalize Logical, should accessibility values be
 #'   normalized by size factor?
 #' @param silent Logical, should warning and info messages be printed?
+#' @param return_agg_info Logical, should a list of the assignments of cells to
+#'   aggregated bins be output? When \code{TRUE}, this function returns a list 
+#'   of two items, first, the aggregated CDS object and second, a data.frame 
+#'   with the binning information.
 #'
 #' @details Aggregation of similar cells is done using a k-nearest-neighbors
 #'   graph and a randomized "bagging" procedure. Details are available in the
 #'   publication that accompanies this package. Run \code{citation("cicero")}
 #'   for publication details. KNN is calculated using
 #'   \code{\link[FNN]{knn.index}}
 #'
-#' @return Aggregated CDS object.
+#' @return Aggregated CDS object. If return_agg_info is \code{TRUE}, a list 
+#'   of the aggregated CDS object and a data.frame of aggregation info.
 #' @export
 #'
 #' @examples
@@ -46,7 +51,8 @@ make_cicero_cds <- function(cds,
                             k=50,
                             summary_stats = NULL,
                             size_factor_normalize = TRUE,
-                            silent = FALSE) {
+                            silent = FALSE, 
+                            return_agg_info = FALSE) {
 
   assertthat::assert_that(is(cds, "CellDataSet"))
   assertthat::assert_that(is.data.frame(reduced_coordinates) |
@@ -73,6 +79,7 @@ make_cicero_cds <- function(cds,
   }
   assertthat::assert_that(is.logical(size_factor_normalize))
   assertthat::assert_that(is.logical(silent))
+  assertthat::assert_that(is.logical(return_agg_info))
 
   reduced_coordinates <- as.data.frame(reduced_coordinates)
 
@@ -128,12 +135,24 @@ make_cicero_cds <- function(cds,
                    "\nMean shared cells bin-bin: ", mean(shared),
                    "\nMedian shared cells bin-bin: ", median(shared)))
 
-    if (mean(shared)/k > .1) warning("On average, more than 10% of cells are shared between paired bins.")
+    if (mean(shared)/k > .1)
+      warning("On average, more than 10% of cells are shared between paired bins.")
   }
 
   exprs_old <- exprs(cds)
 
-  mask <- sapply(seq_len(nrow(cell_sample)), function(x) seq_len(ncol(exprs_old)) %in% cell_sample[x,,drop=FALSE])
+  mask <- sapply(seq_len(nrow(cell_sample)), 
+                 function(x) seq_len(ncol(exprs_old)) %in% cell_sample[x,,drop=FALSE])
+
+  if (return_agg_info) {
+    row.names(mask) <- colnames(exprs_old)
+    colnames(mask) <- paste0("agg_", chosen)
+    agg_map <- reshape2::melt(mask)
+    agg_map <- agg_map[agg_map$value,]
+    agg_map$value <- NULL
+    names(agg_map) <- c("cell", "agg_cell")
+  }
+
   mask <- Matrix::Matrix(mask)
   new_exprs <- exprs_old %*% mask
 
@@ -176,7 +195,7 @@ make_cicero_cds <- function(cds,
                                                  lowerDetectionLimit=0))
 
   cicero_cds <- monocle::detectGenes(cicero_cds, min_expr = .1)
-  cicero_cds <- BiocGenerics::estimateSizeFactors(cicero_cds)
+  cicero_cds <- estimateSizeFactorsSimp(cicero_cds)
   #cicero_cds <- suppressWarnings(BiocGenerics::estimateDispersions(cicero_cds))
 
   if (any(!c("chr", "bp1", "bp2") %in% names(fData(cicero_cds)))) {
@@ -192,6 +211,9 @@ make_cicero_cds <- function(cds,
       t(t(Biobase::exprs(cicero_cds))/Biobase::pData(cicero_cds)$Size_Factor)
   }
 
+  if (return_agg_info) {
+    return(list(cicero_cds, agg_map))
+  }
   cicero_cds
 }