Merge branch 'develop'

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: Signac
 Title: Analysis of Single-Cell Chromatin Data
-Version: 1.2.0
-Date: 2021-04-28
+Version: 1.2.1
+Date: 2021-05-11
 Authors@R: c(
   person(given = 'Tim', family = 'Stuart', email = 'tstuart@nygenome.org', role = c('aut', 'cre'), comment = c(ORCID = '0000-0002-3044-0897')),
   person(given = 'Avi', family = 'Srivastava', email = 'asrivastava@nygenome.org', role = 'aut', comment = c(ORCID = '0000-0001-9798-2079')),

NAMESPACE

@@ -73,6 +73,7 @@ S3method(SetMotifData,Seurat)
 S3method(as.ChromatinAssay,Assay)
 S3method(dim,Motif)
 S3method(dimnames,Motif)
+S3method(head,Fragment)
 S3method(merge,ChromatinAssay)
 S3method(subset,ChromatinAssay)
 S3method(subset,Motif)

NEWS.md

@@ -1,3 +1,17 @@
+# Signac 1.2.1
+
+New functionality:
+
+* Added `head()` method for `Fragment`-class objects.
+
+Bug fixes:
+
+* Fixed bug in `ChromatinAssay` merging ([#596](https://github.com/timoast/signac/pull/596))
+
+Other changes:
+
+* Added support for fragment files containing headers (cellranger-atac v2; [#609](https://github.com/timoast/signac/issues/609))
+
 # Signac 1.2.0
 
 New functionality:

R/fragments.R

@@ -2,6 +2,32 @@
 #'
 NULL
 
+#' Return the first rows of a fragment file
+#'
+#' Returns the first \code{n} rows of a fragment file. This allows the content
+#' of a fragment file to be inspected.
+#'
+#' @param x a \code{Fragment} object
+#' @param n an integer specifying the number of rows to return from the fragment
+#' file
+#' @param ... additional arguments passed to \code{\link[utils]{read.table}}
+#'
+#' @return The first \code{n} rows of a fragment file as a \code{data.frame}
+#' with the following columns: chrom, start, end, barcode, readCount.
+#'
+#' @export
+#' @method head Fragment
+#' @concept fragments
+#'
+head.Fragment <- function(x, n = 6L, ...) {
+  fpath <- GetFragmentData(object = x, slot = "path")
+  df <- read.table(file = fpath, nrows = n, ...)
+  if (ncol(x = df) == 5) {
+    colnames(x = df) <- c("chrom", "start", "end", "barcode", "readCount")
+  }
+  return(df)
+}
+
 #' Count fragments
 #'
 #' Count total fragments per cell barcode present in a fragment file.

R/objects.R

@@ -1147,7 +1147,7 @@ merge.ChromatinAssay <- function(
       if (all(count_nonzero > 0)) {
         merged.counts <- RowMergeSparseMatrices(
           mat1 = all.counts[[1]],
-          mat2 = all.counts[[2:length(x = all.counts)]]
+          mat2 = all.counts[2:length(x = all.counts)]
         )
         reduced.ranges <- StringToGRanges(regions = rownames(x = merged.counts))
       } else {
@@ -1157,7 +1157,7 @@ merge.ChromatinAssay <- function(
       if (all(data_nonzero > 0)) {
         merged.data <- RowMergeSparseMatrices(
           mat1 = all.data[[1]],
-          mat2 = all.data[[2:length(x = all.data)]]
+          mat2 = all.data[2:length(x = all.data)]
         )
         reduced.ranges <- SetIfNull(
           x = reduced.ranges,

R/preprocessing.R

@@ -156,6 +156,13 @@ CreateMotifMatrix <- function(
     }
   }
   rownames(motif.matrix) <- GRangesToString(grange = features, sep = sep)
+  if (is.null(x = names(x = pwm))) {
+    warning("No 'names' attribute found in PFMatrixList. ",
+            "Extracting names from individual entries.")
+    colnames(x = motif.matrix) <- vapply(
+      X = pwm, FUN = slot, FUN.VALUE = "character", "name"
+    )
+  }
   return(motif.matrix)
 }
 
@@ -752,7 +759,7 @@ TSSEnrichment <- function(
 
   # exclude chrM
   sn <- seqnames(x = tss.positions)
-  tss.positions <- tss.positions[!as.character(sn) %in% c("chrM", "Mt")]
+  tss.positions <- tss.positions[!as.character(sn) %in% c("chrM", "Mt", "MT")]
 
   if (fast) {
     # just compute the TSS enrichment score without storing the full matrix

src/filter.cpp

@@ -47,13 +47,27 @@ int filterCells(
   cb_seq.reserve(32);
   line_seq.reserve(buffer_length);
 
-  // looping over the fragments file
-  while(gzgets(ifileHandler, buffer, buffer_length) !=0 ){
+  bool eof_check;
+  while ((eof_check = gzgets(ifileHandler, buffer, buffer_length)) !=0) {
     line_seq.clear();
     line_seq.append(buffer);
 
-    cb_char = strtok ( buffer, "\t" );
+    if (line_seq.at(0) != '#') {
+      break;
+    }
+  }
+
+  if (!eof_check) {
+    Rcpp::Rcerr << "Error: fragment file contains header only\n" << std::flush;
+    gzclose(ifileHandler);
+    return 1;
+  }
 
+  // looping over the fragments file
+  do {
+    line_seq.append(buffer);
+
+    cb_char = strtok ( buffer, "\t" );
     for (auto i=1; i<=3; i++) {
       cb_char = strtok (NULL, "\t");
 
@@ -82,7 +96,9 @@ int filterCells(
     if (is_ten_mil) {
       Rcpp::checkUserInterrupt();
     }
-  }
+
+    line_seq.clear();
+  } while(gzgets(ifileHandler, buffer, buffer_length) !=0 );
 
   //Cleanup
   gzclose(ifileHandler);

src/group.cpp

@@ -50,8 +50,9 @@ SEXP groupCommand(
   }
 
   // char * to string extraction
-  std::string cb_seq;
+  std::string cb_seq, line_seq;
   cb_seq.reserve(32);
+  line_seq.reserve(buffer_length);
 
   // metadata from fragments file
   size_t start, end, reads;
@@ -64,8 +65,25 @@ SEXP groupCommand(
     read_count.assign(num_whitelist_cells, 0);
   }
 
+  // skip header if present
+  bool eof_check;
+  while ((eof_check = gzgets(fileHandler, buffer, buffer_length)) !=0) {
+    line_seq.clear();
+    line_seq.append(buffer);
+
+    if (line_seq.at(0) != '#') {
+      break;
+    }
+  }
+
+  if (!eof_check) {
+    Rcpp::Rcerr << "Error: fragment file contains header only" << std::flush;
+    gzclose(fileHandler);
+    return (Rcpp::DataFrame::create());
+  }
+
   // looping over the fragments file
-  while(gzgets(fileHandler, buffer, buffer_length) !=0 ){
+  do {
     cb_char = strtok ( buffer, "\t" );
 
     for (auto i=1; i<=4; i++) {
@@ -136,7 +154,8 @@ SEXP groupCommand(
     if (line_counter % 2000000 == 0) {
       Rcpp::checkUserInterrupt();
     }
-  }
+    line_seq.clear();
+  } while(gzgets(fileHandler, buffer, buffer_length) !=0 );
 
   //Cleanup
   gzclose(fileHandler);

src/split.cpp

@@ -44,7 +44,7 @@ int splitFragments(
     } else {
       out->open( fileName.c_str() );
     }
- 
+
     streams.push_back(out);
   }
 
@@ -80,8 +80,25 @@ int splitFragments(
   cb_seq.reserve(32);
   line_seq.reserve(buffer_length);
 
+  // skip header if present
+  bool eof_check;
+  while ((eof_check = gzgets(ifileHandler, buffer, buffer_length)) !=0) {
+    line_seq.clear();
+    line_seq.append(buffer);
+
+    if (line_seq.at(0) != '#') {
+      break;
+    }
+  }
+
+  if (!eof_check) {
+    Rcpp::Rcerr << "Error: fragment file contains header only\n" << std::flush;
+    gzclose(ifileHandler);
+    return 1;
+  }
+
   // looping over the fragments file
-  while(gzgets(ifileHandler, buffer, buffer_length) !=0 ){
+  do {
     line_seq.clear();
     line_seq.append(buffer);
 
@@ -124,7 +141,7 @@ int splitFragments(
     if (is_ten_mil) {
       Rcpp::checkUserInterrupt();
     }
-  }
+  } while(gzgets(ifileHandler, buffer, buffer_length) !=0 );
 
   // Cleanup
   gzclose(ifileHandler);

src/validate.cpp

@@ -46,11 +46,29 @@ bool validateCells(
   }
 
   // char * to string extraction
-  std::string cb_seq;
+  std::string cb_seq, line_seq;
   cb_seq.reserve(32);
+  line_seq.reserve(buffer_length);
+
+  // skip header if present
+  bool eof_check;
+  while ((eof_check = gzgets(fileHandler, buffer, buffer_length)) !=0) {
+    line_seq.clear();
+    line_seq.append(buffer);
+
+    if (line_seq.at(0) != '#') {
+      break;
+    }
+  }
+
+  if (!eof_check) {
+    Rcpp::Rcerr << "Error: fragment file contains header only\n" << std::flush;
+    gzclose(fileHandler);
+    return (false);
+  }
 
   // looping over the fragments file
-  while(gzgets(fileHandler, buffer, buffer_length) !=0 ){
+  do {
     cb_char = strtok ( buffer, "\t" );
 
     for (auto i=1; i<=3; i++) {
@@ -83,7 +101,7 @@ bool validateCells(
     if (line_counter % 2000000 == 0) {
       Rcpp::checkUserInterrupt();
     }
-  }
+  } while(gzgets(fileHandler, buffer, buffer_length) !=0 );
 
   //Cleanup
   gzclose(fileHandler);