The SCArray package provides large-scale single-cell RNA-seq data manipulation using Genomic Data Structure (GDS) files. It combines dense/sparse matrices stored in GDS files and the Bioconductor infrastructure framework (SingleCellExperiment and DelayedArray) to provide out-of-memory data storage and manipulation using the R programming language. As shown in the figure, SCArray provides a SingleCellExperiment
object for downstream data analyses. GDS is an alternative to HDF5. Unlike HDF5, GDS supports the direct storage of a sparse matrix without converting it to multiple vectors.
Figure 1: Workflow of SCArray
Requires R (>= v3.5.0), gdsfmt (>= v1.24.0)
Bioconductor repository
To install this package, start R and enter:
if (!requireNamespace("BiocManager", quietly=TRUE))
install.packages("BiocManager")
::install("SCArray") BiocManager
The SCArray package can convert a single-cell experiment object (SingleCellExperiment) to a GDS file using the function scConvGDS()
. For example,
suppressPackageStartupMessages(library(SCArray))
suppressPackageStartupMessages(library(SingleCellExperiment))
# load a SingleCellExperiment object
system.file("extdata", "example.rds", package="SCArray")
fn <- readRDS(fn)
sce <-
# convert to a GDS file
scConvGDS(sce, "test.gds")
## Output: test.gds
## Compression: LZMA_RA
## Dimension: 250 x 500
## Assay List:
## counts |+ counts { SparseReal32 250x500 LZMA_ra(13.2%), 77.2K }
## rowData:
## colData:
## Cell_ID
## Cell_type
## Timepoint
## Done.
# list data structure in the GDS file
scOpen("test.gds")); scClose(f) (f <-
## Object of class "SCArrayFileClass"
## File: test.gds (82.6K)
## + [ ] *
## |--+ feature.id { Str8 250 LZMA_ra(71.2%), 1.1K }
## |--+ sample.id { Str8 500 LZMA_ra(15.1%), 1.2K }
## |--+ counts { SparseReal32 250x500 LZMA_ra(13.2%), 77.2K }
## |--+ feature.data [ ]
## |--+ sample.data [ ]
## | |--+ Cell_ID { Str8 500 LZMA_ra(15.1%), 1.2K }
## | |--+ Cell_type { Str8 500 LZMA_ra(4.49%), 141B }
## | \--+ Timepoint { Str8 500 LZMA_ra(5.31%), 193B }
## \--+ meta.data [ ]
The input of scConvGDS()
can be a dense or sparse matrix for count data:
library(Matrix)
matrix(0, nrow=4, ncol=8)
cnt <-set.seed(100); cnt[sample.int(length(cnt), 8)] <- rpois(8, 4)
as(cnt, "dgCMatrix")) (cnt <-
## 4 x 8 sparse Matrix of class "dgCMatrix"
##
## [1,] . . . . . . . 6
## [2,] 3 1 . . . 4 . .
## [3,] . . . . . 3 . 4
## [4,] 4 . 3 . . . . .
# convert to a GDS file
scConvGDS(cnt, "test.gds")
## Output: test.gds
## Compression: LZMA_RA
## Dimension: 4 x 8
## Assay List:
## counts |+ counts { SparseReal32 4x8 LZMA_ra(159.4%), 109B }
## Done.
When a single-cell GDS file is available, users can use scExperiment()
to load a SingleCellExperiment object from the GDS file. The assay data in the SingleCellExperiment object are DelayedMatrix objects to avoid the memory limit.
# a GDS file in the SCArray package
system.file("extdata", "example.gds", package="SCArray")) (fn <-
## [1] "/tmp/RtmpOcrPA6/Rinst2fc8a450606ae0/SCArray/extdata/example.gds"
# load a SingleCellExperiment object from the file
scExperiment(fn)
sce <- sce
## class: SingleCellExperiment
## dim: 1000 850
## metadata(0):
## assays(1): counts
## rownames(1000): MRPL20 GNB1 ... RPS4Y1 CD24
## rowData names(0):
## colnames(850): 1772122_301_C02 1772122_180_E05 ... 1772122_180_B06
## 1772122_180_D09
## colData names(3): Cell_ID Cell_type Timepoint
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
# it is a DelayedMatrix (the whole matrix is not loaded)
assays(sce)$counts
## <1000 x 850> sparse matrix of class SC_GDSMatrix and type "double":
## 1772122_301_C02 1772122_180_E05 ... 1772122_180_B06
## MRPL20 3 2 . 0
## GNB1 11 6 . 0
## RPL22 3 5 . 6
## PARK7 1 7 . 2
## ENO1 8 19 . 7
## ... . . . .
## SSR4 0 6 . 5
## RPL10 11 4 . 1
## SLC25A6_loc1 4 5 . 3
## RPS4Y1 0 5 . 2
## CD24 18 3 . 0
## 1772122_180_D09
## MRPL20 2
## GNB1 0
## RPL22 6
## PARK7 2
## ENO1 4
## ... .
## SSR4 1
## RPL10 3
## SLC25A6_loc1 1
## RPS4Y1 4
## CD24 2
# column data
colData(sce)
## DataFrame with 850 rows and 3 columns
## Cell_ID Cell_type Timepoint
## <character> <character> <character>
## 1772122_301_C02 1772122_301_C02 eNb1 day_35
## 1772122_180_E05 1772122_180_E05 eNb1 day_35
## 1772122_300_H02 1772122_300_H02 eNb1 day_35
## 1772122_180_B09 1772122_180_B09 eNb1 day_35
## 1772122_180_G04 1772122_180_G04 eNb1 day_35
## ... ... ... ...
## 1772122_181_F11 1772122_181_F11 eRgld day_35
## 1772122_181_E02 1772122_181_E02 eRgld day_35
## 1772122_180_C03 1772122_180_C03 eRgld day_35
## 1772122_180_B06 1772122_180_B06 eRgld day_35
## 1772122_180_D09 1772122_180_D09 eRgld day_35
# row data
rowData(sce)
## DataFrame with 1000 rows and 0 columns
SCArray provides a SingleCellExperiment
object for downstream data analyses. At first, we create a log count matrix logcnt
from the count matrix. Note that logcnt
is also a DelayedMatrix without actually generating the whole matrix.
assays(sce)$counts
cnt <- log2(cnt + 1)
logcnt <-assays(sce)$logcounts <- logcnt
logcnt
## <1000 x 850> sparse matrix of class SC_GDSMatrix and type "double":
## 1772122_301_C02 1772122_180_E05 ... 1772122_180_B06
## MRPL20 2.000000 1.584963 . 0.000000
## GNB1 3.584963 2.807355 . 0.000000
## RPL22 2.000000 2.584963 . 2.807355
## PARK7 1.000000 3.000000 . 1.584963
## ENO1 3.169925 4.321928 . 3.000000
## ... . . . .
## SSR4 0.000000 2.807355 . 2.584963
## RPL10 3.584963 2.321928 . 1.000000
## SLC25A6_loc1 2.321928 2.584963 . 2.000000
## RPS4Y1 0.000000 2.584963 . 1.584963
## CD24 4.247928 2.000000 . 0.000000
## 1772122_180_D09
## MRPL20 1.584963
## GNB1 0.000000
## RPL22 2.807355
## PARK7 1.584963
## ENO1 2.321928
## ... .
## SSR4 1.000000
## RPL10 2.000000
## SLC25A6_loc1 1.000000
## RPS4Y1 2.321928
## CD24 1.584963
The DelayedMatrixStats package provides functions operating on rows and columns of DelayedMatrix objects. For example, we can calculate the mean for each column or row of the log count matrix.
suppressPackageStartupMessages(library(DelayedMatrixStats))
DelayedMatrixStats::colMeans2(logcnt)
col_mean <-str(col_mean)
## num [1:850] 1.51 1.95 2.25 1.95 1.75 ...
DelayedMatrixStats::rowMeans2(logcnt)
row_mean <-str(row_mean)
## num [1:1000] 1.27 1.51 2.62 1.98 3.75 ...
The scater package can perform the uniform manifold approximation and projection (UMAP) for the cell data, based on the data in a SingleCellExperiment object.
suppressPackageStartupMessages(library(scater))
# run umap analysis
runUMAP(sce) sce <-
plotReducedDim()
plots cell-level reduced dimension results (UMAP) stored in the SingleCellExperiment object:
plotReducedDim(sce, dimred="UMAP")
# print version information about R, the OS and attached or loaded packages
sessionInfo()
## R Under development (unstable) (2023-01-10 r83596)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 22.04.1 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.17-bioc/R/lib/libRblas.so
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## time zone: America/New_York
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] scater_1.27.2 ggplot2_3.4.0
## [3] scuttle_1.9.4 DelayedMatrixStats_1.21.0
## [5] SingleCellExperiment_1.21.0 SummarizedExperiment_1.29.1
## [7] Biobase_2.59.0 GenomicRanges_1.51.4
## [9] GenomeInfoDb_1.35.13 SCArray_1.7.5
## [11] DelayedArray_0.25.0 IRanges_2.33.0
## [13] S4Vectors_0.37.3 MatrixGenerics_1.11.0
## [15] matrixStats_0.63.0 BiocGenerics_0.45.0
## [17] Matrix_1.5-3 gdsfmt_1.35.5
##
## loaded via a namespace (and not attached):
## [1] beeswarm_0.4.0 gtable_0.3.1 xfun_0.36
## [4] bslib_0.4.2 ggrepel_0.9.2 lattice_0.20-45
## [7] vctrs_0.5.2 tools_4.3.0 bitops_1.0-7
## [10] generics_0.1.3 parallel_4.3.0 tibble_3.1.8
## [13] fansi_1.0.4 highr_0.10 BiocNeighbors_1.17.1
## [16] pkgconfig_2.0.3 sparseMatrixStats_1.11.1 assertthat_0.2.1
## [19] lifecycle_1.0.3 GenomeInfoDbData_1.2.9 farver_2.1.1
## [22] FNN_1.1.3.1 compiler_4.3.0 munsell_0.5.0
## [25] codetools_0.2-18 vipor_0.4.5 htmltools_0.5.4
## [28] sass_0.4.5 RCurl_1.98-1.9 yaml_2.3.7
## [31] pillar_1.8.1 crayon_1.5.2 jquerylib_0.1.4
## [34] uwot_0.1.14 BiocParallel_1.33.9 cachem_1.0.6
## [37] viridis_0.6.2 tidyselect_1.2.0 rsvd_1.0.5
## [40] digest_0.6.31 BiocSingular_1.15.0 dplyr_1.0.10
## [43] labeling_0.4.2 cowplot_1.1.1 fastmap_1.1.0
## [46] grid_4.3.0 colorspace_2.1-0 cli_3.6.0
## [49] magrittr_2.0.3 utf8_1.2.2 withr_2.5.0
## [52] scales_1.2.1 ggbeeswarm_0.7.1 rmarkdown_2.20
## [55] XVector_0.39.0 gridExtra_2.3 ScaledMatrix_1.7.0
## [58] beachmat_2.15.0 evaluate_0.20 knitr_1.42
## [61] viridisLite_0.4.1 irlba_2.3.5.1 rlang_1.0.6
## [64] Rcpp_1.0.10 DBI_1.1.3 glue_1.6.2
## [67] jsonlite_1.8.4 R6_2.5.1 zlibbioc_1.45.0