Zarr (AnnData) sparse matrices

.Rbuildignore

@@ -0,0 +1,2 @@
+^.*\.Rproj$
+^\.Rproj\.user$

.gitignore

@@ -0,0 +1,5 @@
+.Rproj.user
+.Rhistory
+.RData
+.Ruserdata
+*Rproj

DESCRIPTION

@@ -29,14 +29,22 @@ Authors@R: c(
 		comment=c(ORCID="0000-0003-2725-0694")))
 Depends: R (>= 3.4), methods, SparseArray, DelayedArray
 Imports: stats, tools, BiocGenerics, S4Vectors, IRanges, S4Arrays,
-	Rarr (>= 1.11.33)
-Suggests: paws.storage, HDF5Array, testthat, knitr, rmarkdown, BiocStyle
+	Rarr (>= 1.11.33), jsonlite
+Suggests: paws.storage, HDF5Array, testthat, knitr, rmarkdown, BiocStyle, 
+          anndataR
 VignetteBuilder: knitr
 Collate:
 	utils.R
 	options.R
+	zarr_utils.R
+	zarr_mread.R
 	ZarrArraySeed-class.R
 	ZarrArray-class.R
+	ZarrSparseMatrixSeed-class.R
+	ZarrSparseMatrix-class.R
+	ZarrADMatrix-class.R
+  ZarrADMatrixSeed-class.R
 	writeZarrArray-auto-args.R
 	writeZarrArray.R
 	zzz.R
+Config/roxygen2/version: 8.0.0

NAMESPACE

@@ -18,9 +18,33 @@ importFrom(Rarr, read_zarr_array, create_empty_zarr_array, update_zarr_array)
 exportClasses(
     "ZarrArraySeed",
     "ZarrArray", "ZarrMatrix",
-    "ZarrRealizationSink"
+    "ZarrRealizationSink",
+    "ZarrSparseMatrixSeed", "CSC_ZarrSparseMatrixSeed", "CSR_ZarrSparseMatrixSeed",
+    "ZarrSparseMatrix",
+    "ZarrADMatrixSeed",
+    "Dense_ZarrADMatrixSeed", "CSC_ZarrADMatrixSeed", "CSR_ZarrADMatrixSeed",
+    "ZarrADMatrix"
 )
 
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Export S3 methods
+###
+
+S3method(t, CSC_ZarrSparseMatrixSeed)
+S3method(t, CSR_ZarrSparseMatrixSeed)
+S3method(t, CSC_ZarrADMatrixSeed)
+S3method(t, CSR_ZarrADMatrixSeed)
+
+### We also export them thru the export() directive so that (a) they can be
+### called directly, (b) tab-completion on the name of the generic shows them,
+### and (c) methods() doesn't asterisk them.
+
+export(
+    t.CSC_ZarrSparseMatrixSeed,
+    t.CSR_ZarrSparseMatrixSeed,
+    t.CSC_ZarrADMatrixSeed,
+    t.CSR_ZarrADMatrixSeed
+)
 
 ### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 ### Export S4 methods for generics not defined in ZarrArray
@@ -35,6 +59,12 @@ exportMethods(
 
     ## Methods for generics defined in the BiocGenerics package:
     path, type,
+
+    ## Methods for generics defined in the S4Arrays package:
+    extract_array, is_sparse, 
+
+    ## Methods for generics defined in the SparseArray package:
+    nzcount, extract_sparse_array,
 
     ## Methods for generics defined in the S4Arrays package:
     extract_array, write_block,
@@ -64,6 +94,31 @@ export(
     get_writeZarrArray_auto_chunkdim,
 
     ## writeZarrArray.R:
-    ZarrRealizationSink, writeZarrArray
+    ZarrRealizationSink, writeZarrArray,
+
+    ## ZarrSparseMatrixSeed-class.R
+    ZarrSparseMatrixSeed,
+
+    ## ZarrSparseMatrix-class.R
+    ZarrSparseMatrix,
+
+    ## ZarrADMatrixSeed-class.R
+    ZarrADMatrixSeed,
+
+    ## ZarrADMatrix-class.R
+    ZarrADMatrix
 )
 
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Export S4 generics defined in HDF5Array, and corresponding methods
+###
+
+export(
+    ## ZarrSparseMatrixSeed-class.R:
+    extractNonzeroDataByCol, extractNonzeroDataByRow
+)
+
+### Exactly the same list as above.
+exportMethods(
+    extractNonzeroDataByCol, extractNonzeroDataByRow
+)

R/ZarrADMatrix-class.R

@@ -0,0 +1,52 @@
+### =========================================================================
+### ZarrADMatrix objects
+### -------------------------------------------------------------------------
+###
+
+setClass("ZarrADMatrix",
+         contains="DelayedMatrix",
+         representation(seed="ZarrADMatrixSeed")
+)
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Constructor
+###
+
+setMethod("DelayedArray", "ZarrADMatrixSeed",
+          function(seed) new_DelayedArray(seed, Class="ZarrADMatrix")
+)
+
+### Works directly on an ZarrADMatrixSeed derivative, in which case it must
+### be called with a single argument.
+ZarrADMatrix <- function(filepath, layer=NULL)
+{
+  if (is(filepath, "ZarrADMatrixSeed")) {
+    if (!is.null(layer))
+      stop(wmsg("ZarrADMatrix() must be called with a single argument ",
+                "when passed an ZarrADMatrixSeed derivative"))
+    seed <- filepath
+  } else {
+    seed <- ZarrADMatrixSeed(filepath, layer=layer)
+  }
+  DelayedArray(seed)
+}
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Taking advantage of sparsity
+###
+
+### Will work only if the seed is an H5SparseMatrixSeed derivative, that is,
+### if it's a CSC_ZarrADMatrixSeed or CSR_ZarrADMatrixSeed object.
+setMethod("nzcount", "ZarrADMatrix", function(x) nzcount(x@seed))
+
+### Will work only if the seed is a CSC_ZarrADMatrixSeed object.
+setMethod("extractNonzeroDataByCol", "ZarrADMatrix",
+          function(x, j) extractNonzeroDataByCol(x@seed, j)
+)
+
+### Will work only if the seed is a CSR_ZarrADMatrixSeed object.
+setMethod("extractNonzeroDataByRow", "ZarrADMatrix",
+          function(x, i) extractNonzeroDataByCol(x@seed, i)
+)

R/ZarrADMatrixSeed-class.R

@@ -0,0 +1,134 @@
+### =========================================================================
+### ZarrADMatrixSeed objects
+### -------------------------------------------------------------------------
+
+
+setClass("ZarrADMatrixSeed",
+         contains=c("Array", "OutOfMemoryObject"),
+         representation("VIRTUAL")
+)
+
+setClass("Dense_ZarrADMatrixSeed",
+         contains=c("ZarrADMatrixSeed", "ZarrArraySeed"),
+         representation(dimnames="list"),
+         prototype(dimnames=list(NULL, NULL))
+)
+setClass("CSC_ZarrADMatrixSeed",
+         contains=c("ZarrADMatrixSeed", "CSC_ZarrSparseMatrixSeed")
+)
+setClass("CSR_ZarrADMatrixSeed",
+         contains=c("ZarrADMatrixSeed", "CSR_ZarrSparseMatrixSeed")
+)
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### dimnames() method for Dense_ZarrADMatrixSeed objects
+###
+
+### We overwrite the method for HDF5ArraySeed objects with a method that
+### accesses the slot, not the store
+setMethod("dimnames", "Dense_ZarrADMatrixSeed",
+          function(x) S4Arrays:::simplify_NULL_dimnames(x@dimnames)
+)
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Transposition
+###
+
+### S3/S4 combo for t.CSC_ZarrADMatrixSeed
+t.CSC_ZarrADMatrixSeed <- function(x)
+{
+  new2("CSR_ZarrADMatrixSeed", callNextMethod())
+}
+setMethod("t", "CSC_ZarrADMatrixSeed", t.CSC_ZarrADMatrixSeed)
+
+### S3/S4 combo for t.CSR_ZarrADMatrixSeed
+t.CSR_ZarrADMatrixSeed <- function(x)
+{
+  new2("CSC_ZarrADMatrixSeed", callNextMethod())
+}
+setMethod("t", "CSR_ZarrADMatrixSeed", t.CSR_ZarrADMatrixSeed)
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Constructor
+###
+
+.load_zarr_ad_rownames <- function(filepath, name="var")
+{
+  ok <- try(zarrisdataset(filepath, name), silent=TRUE)
+  if (isTRUE(ok)) {
+    ## Must use rhdf5::h5read() for now, until h5mread() knows how
+    ## to read COMPOUND datasets.
+    ans <- h5read(filepath, name)$index
+    if (!is.null(ans))
+      ans <- as.character(ans)
+    return(ans)
+  }
+  ok <- try(zarrisgroup(filepath, name), silent=TRUE)
+  if (!isTRUE(ok))
+    return(NULL)
+  ROWNAMES_DATASET <- paste0(name, "/_index")
+  ok <- try(zarrisdataset(filepath, ROWNAMES_DATASET), silent=TRUE)
+  if (!isTRUE(ok))
+    return(NULL)
+  zarr_mread(filepath, ROWNAMES_DATASET)
+}
+
+### Must return a list of length 2.
+.load_zarr_ad_dimnames <- function(filepath)
+{
+  ans_rownames <- .load_zarr_ad_rownames(filepath)
+  ans_colnames <- .load_zarr_ad_rownames(filepath, name="obs")
+  if (is.null(ans_rownames) && is.null(ans_colnames))
+    warning(wmsg("could not find dimnames in this anndata-zarr store"))
+  list(ans_rownames, ans_colnames)
+}
+
+### Returns an ZarrADMatrixSeed derivative (can be either a Dense_ZarrADMatrixSeed,
+### or a CSC_ZarrSparseMatrixSeed, or a CSR_ZarrSparseMatrixSeed object).
+ZarrADMatrixSeed <- function(filepath, layer=NULL)
+{
+  if (!isSingleString(filepath))
+    stop(wmsg("'filepath' must be a single string specifying the ",
+              "path to the anndata-zarr store"))
+  filepath <- file_path_as_absolute(filepath)
+  if (is.null(layer)) {
+    name <- "/X"
+  } else {
+    if (!isSingleString(layer) || layer == "")
+      stop(wmsg("'layer' must be NULL or a single non-empty string"))
+    name <- paste0("/layers/", layer)
+  }
+  if (!zarrexists(filepath, name)) {
+    msg <- c("Zarr object \"", name, "\" does not exist ",
+             "in this Zarr store")
+    if (is.null(layer))
+      msg <- c(msg, " Is this a valid anndata-zarr store?")
+    stop(wmsg(msg))
+  }
+  dimnames <- .load_zarr_ad_dimnames(filepath)
+
+  if (zarrisdataset(filepath, name)) {
+    ans0 <- HDF5ArraySeed(filepath, name)
+    if (length(dim(ans0)) != 2L)
+      stop(wmsg("Zarr dataset \"", name, "\" in store \"", filepath, "\" ",
+                "does not have exactly 2 dimensions. Please consider ",
+                "using the HDF5Array() constructor to access this ",
+                "dataset."))
+    ans <- new2("Dense_ZarrADMatrixSeed", ans0, dimnames=dimnames)
+  } else if (zarrisgroup(filepath, name)) {
+    ans0 <- ZarrSparseMatrixSeed(filepath, name)
+    if (is(ans0, "CSC_ZarrSparseMatrixSeed"))
+      ans_class <- "CSC_ZarrADMatrixSeed"
+    else
+      ans_class <- "CSR_ZarrADMatrixSeed"
+    ans <- new2(ans_class, ans0, dimnames=dimnames)
+  } else {
+    stop(wmsg("Zarr object \"", name, "\" in store \"", filepath, "\" ",
+              "is neither a dataset or a group. Is this a valid ",
+              "anndata-zarr store?"))
+  }
+  ans
+}

R/ZarrSparseMatrix-class.R

@@ -0,0 +1,49 @@
+### =========================================================================
+### ZarrSparseMatrix objects
+### -------------------------------------------------------------------------
+###
+
+
+setClass("ZarrSparseMatrix",
+    contains="DelayedMatrix",
+    representation(seed="ZarrSparseMatrixSeed")
+)
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Constructor
+###
+
+setMethod("DelayedArray", "ZarrSparseMatrixSeed",
+    function(seed) new_DelayedArray(seed, Class="ZarrSparseMatrix")
+)
+
+### Works directly on an ZarrSparseMatrixSeed derivative, in which case it must
+### be called with a single argument.
+ZarrSparseMatrix <- function(filepath, group)
+{
+    if (is(filepath, "ZarrSparseMatrixSeed")) {
+        if (!missing(group))
+            stop(wmsg("ZarrSparseMatrix() must be called with a single argument ",
+                      "when passed an ZarrSparseMatrixSeed object"))
+        seed <- filepath
+    } else {
+        seed <- ZarrSparseMatrixSeed(filepath, group)
+    }
+    DelayedArray(seed)
+}
+
+
+### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+### Taking advantage of sparsity
+###
+
+setMethod("nzcount", "ZarrSparseMatrix", function(x) nzcount(x@seed))
+
+setMethod("extractNonzeroDataByCol", "ZarrSparseMatrix",
+    function(x, j) extractNonzeroDataByCol(x@seed, j)
+)
+
+setMethod("extractNonzeroDataByRow", "ZarrSparseMatrix",
+    function(x, i) extractNonzeroDataByCol(x@seed, i)
+)