BiocDuckDB 0.99.22
Bioconductor’s experiment containers, SummarizedExperiment,
SingleCellExperiment, MultiAssayExperiment, hold their assays, annotations,
and reduced dimensions in memory. For modern datasets that is increasingly
impractical: a million-cell single-cell experiment is tens of gigabytes, and a
multi-omics study multiplies that across assays.
BiocDuckDB is the integration layer of the BiocDuckDB suite. It provides two high-level functions,
writeParquet() serializes a Bioconductor experiment object to a structured
directory of columnar Parquet files, andreadParquet() reads that directory back as a DuckDB-backed experiment
object.So an experiment keeps its assays on disk while presenting the same
Bioconductor API. Assays become DuckDBArray matrices, rowData and
colData become DuckDBDataFrames, and rowRanges becomes a
DuckDBGRanges list, and operations on them run as SQL against the
Parquet files. The on-disk format is a self-describing
Frictionless Data Package, readable by other tools
and languages, not just R.
This vignette covers the round-trip workflow and the operations it enables. For the SQL-optimized scran/scuttle analysis methods and how they perform at scale, see Benchmarking BiocDuckDB.
if (!require("BiocManager"))
install.packages("BiocManager")
BiocManager::install("BiocDuckDB")library(BiocDuckDB)
library(SummarizedExperiment)We write a SummarizedExperiment to Parquet and read it back as a DuckDB-backed
object. The airway bulk RNA-seq dataset is a RangedSummarizedExperiment,
so it also exercises the genomic-coordinate path.
data(airway, package = "airway")
se_path <- file.path(tempdir(), "airway_se")
writeParquet(airway, se_path)
airway_ddb <- readParquet(se_path)
class(airway_ddb)
#> [1] "RangedSummarizedExperiment"
#> attr(,"package")
#> [1] "SummarizedExperiment"
dim(airway_ddb)
#> [1] 63677 8writeParquet() lays out each component as its own Parquet under a
Frictionless-style directory, described by a datapackage.json:
head(list.files(se_path, recursive = TRUE), 10)
#> [1] "assay_counts/part-0.parquet" "datapackage.json"
#> [3] "features/part-0.parquet" "samples/part-0.parquet"Assays are stored in coordinate (COO) form, one row per stored entry with
__feature__/__sample__ index columns and a value, which is compact for the
sparse matrices typical of single-cell data, and lets DuckDB scan only what a query
touches.
The rest of this section is for interoperating with other tools or writing data
incrementally; skip ahead to “Operations run as SQL” below if you only need the
writeParquet()/readParquet() round trip above.
The datapackage.json that writeParquet() writes follows the
Frictionless Data Package specification, a widely used,
language-agnostic way to describe a directory of data files. BiocDuckDB adds a few
of its own fields on top (dimension and layout, used below), but the file
otherwise reads like any other Frictionless package. In practice that means two
things: another tool, or another language, can write COO Parquet files with a
matching datapackage.json and readParquet() will read them directly, with no R
object involved in producing them; and R code that already has COO Parquet on disk
from some other process can attach to it directly (below), instead of writing it
out again through writeParquet().
To build a datapackage.json from pieces written separately, for example a
dataset streamed in one part at a time, or migrated from another store, collect the
descriptor each writeParquet() call returns and pass them all to
writeDatapackage():
resources <- c(
writeParquet(features, file.path(dir, "features"), dimension = "feature"),
writeParquet(assays, dir, dimension = "crossed"))
writeDatapackage("summarized_experiment", resources, dir)Going the other way, the DuckDBMatrix() / DuckDBArray() / DuckDBTable()
constructors attach an existing COO Parquet in place, wrapping it as a lazy,
out-of-core object without copying or materializing it:
mat <- DuckDBMatrix(coo_path, datacol = "value",
keycols = c("__feature__", "__sample__"))coo_path is the directory holding the COO Parquet; datacol names the column
holding the stored values, and keycols names the two index columns that locate
each value in the matrix. The documented model values and the dimension x
layout dispatch that maps each resource to its accessor are described on
?writeParquet and ?readParquet.
The DuckDB-backed object supports the standard API; component access, subsetting, and summaries are pushed to DuckDB:
assay(airway_ddb, "counts")
#> <63677 x 8> sparse DuckDBMatrix object of type "integer":
#> SRR1039508 SRR1039509 ... SRR1039520 SRR1039521
#> ENSG00000000003 679 448 . 770 572
#> ENSG00000000005 0 0 . 0 0
#> ENSG00000000419 467 515 . 417 508
#> ENSG00000000457 260 211 . 233 229
#> ENSG00000000460 60 55 . 76 60
#> ... . . . . .
#> ENSG00000273489 0 0 . 0 0
#> ENSG00000273490 0 0 . 0 0
#> ENSG00000273491 0 0 . 0 0
#> ENSG00000273492 0 0 . 0 0
#> ENSG00000273493 0 0 . 0 0
colData(airway_ddb)[, 1:3]
#> DuckDBDataFrame with 8 rows and 3 columns
#> SampleName cell dex
#> <factor> <factor> <factor>
#> SRR1039508 GSM1275862 N61311 untrt
#> SRR1039509 GSM1275863 N61311 trt
#> SRR1039512 GSM1275866 N052611 untrt
#> SRR1039513 GSM1275867 N052611 trt
#> SRR1039516 GSM1275870 N080611 untrt
#> SRR1039517 GSM1275871 N080611 trt
#> SRR1039520 GSM1275874 N061011 untrt
#> SRR1039521 GSM1275875 N061011 trt
## subset on disk, then summarize
sub <- airway_ddb[1:1000, 1:4]
colSums(assay(sub, "counts"))
#> __sample__
#> SRR1039508 SRR1039509 SRR1039512 SRR1039513
#> 1352811 1216995 1792134 1043616Because only stored values and the touched columns are read, the in-memory footprint of the object stays small relative to the full data:
c(in_memory = format(object.size(airway), units = "MB"),
duckdb = format(object.size(airway_ddb), units = "MB"))
#> in_memory duckdb
#> "80.5 Mb" "14.2 Mb"writeParquet() stores a RangedSummarizedExperiment’s rowRanges as Parquet
LIST[] columns, so they come back as a DuckDBGRanges list with the
genomic API intact:
rr <- rowRanges(airway_ddb)
class(rr)
#> [1] "DuckDBGRangesList"
#> attr(,"package")
#> [1] "DuckDBGRanges"
elementNROWS(rr)[1:8] # exons per gene, queried from disk
#> ENSG00000000003 ENSG00000000005 ENSG00000000419 ENSG00000000457 ENSG00000000460
#> 17 10 29 30 72
#> ENSG00000000938 ENSG00000000971 ENSG00000001036
#> 26 45 14Single-cell experiments are where the on-disk representation pays off most: the
matrices are large and sparse, so COO storage compresses well and the DuckDB-backed
object holds essentially only metadata. A small sparse SingleCellExperiment
illustrates the round-trip:
library(SingleCellExperiment)
library(Matrix)
set.seed(1L)
counts <- Matrix(rpois(2000 * 500, lambda = 0.3), nrow = 2000, ncol = 500,
sparse = TRUE)
rownames(counts) <- paste0("Gene", seq_len(nrow(counts)))
colnames(counts) <- paste0("Cell", seq_len(ncol(counts)))
sce <- SingleCellExperiment(assays = list(counts = counts))
sce_path <- file.path(tempdir(), "demo_sce")
writeParquet(sce, sce_path)
sce_ddb <- readParquet(sce_path)
## QC filtering stays on disk (SQL), then summarize
totals <- colSums(assay(sce_ddb, "counts"))
keep <- sce_ddb[, totals > median(totals)]
dim(keep)
#> [1] 2000 249For data larger than memory, the idiomatic workflow uses DuckDB for the parts it does well, filtering and summarizing on disk, and only then realizes the manageable subset into memory for the standard analysis pipeline:
library(scran)
library(scater)
sce_ddb <- readParquet("data/raw_counts") # DuckDB-backed, low RAM
## 1. filter on disk (SQL-optimized)
keep <- colSums(assay(sce_ddb, "counts")) > 1000
sce_ddb <- sce_ddb[, keep]
detected <- rowSums(assay(sce_ddb, "counts") > 0)
sce_ddb <- sce_ddb[detected >= 10, ]
## 2. realize the filtered subset into memory
sce <- as(sce_ddb, "SingleCellExperiment")
## 3. standard in-memory pipeline
sce <- logNormCounts(sce)
dec <- modelGeneVar(sce)
sce <- runPCA(sce[getTopHVGs(dec, n = 2000), ])
## 4. persist results back to Parquet
writeParquet(sce, "data/processed")Many scran/scuttle steps can run directly on the DuckDB-backed object as SQL, before any realization, that is the subject of Benchmarking BiocDuckDB.
For a MultiAssayExperiment, write each experiment to its own directory and read
them back as DuckDB-backed experiments:
library(MultiAssayExperiment)
writeParquet(rna_sce, file.path(mae_path, "rna"))
writeParquet(protein_se, file.path(mae_path, "protein"))
mae <- MultiAssayExperiment(experiments = list(
rna = readParquet(file.path(mae_path, "rna")),
protein = readParquet(file.path(mae_path, "protein"))))Each modality stays on disk until accessed. Spatial experiments
(MultiAssaySpatialExperiment) are supported too, with their geometry written as
GeoParquet via DuckDBSpatial.
A good fit for experiments larger than memory, for exploration, QC, and
filtering of large datasets, and for a self-describing, cross-language on-disk
format (the Parquet files are readable by Python, Julia, and DuckDB directly).
Parquet is also a compact archival format, typically several-fold smaller than
.rds for sparse data, and self-describing.
Realize a subset to an ordinary in-memory object for iterative algorithms with random access patterns or dense linear algebra; and small datasets that fit comfortably in memory need no backend at all.
sessionInfo()
#> R version 4.6.1 Patched (2026-06-24 r90190)
#> Platform: x86_64-apple-darwin20
#> Running under: macOS Ventura 13.7.8
#>
#> Matrix products: default
#> BLAS: /Library/Frameworks/R.framework/Versions/4.6-x86_64/Resources/lib/libRblas.0.dylib
#> LAPACK: /Library/Frameworks/R.framework/Versions/4.6-x86_64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.1
#>
#> locale:
#> [1] C/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
#>
#> time zone: America/New_York
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats4 stats graphics grDevices utils datasets methods
#> [8] base
#>
#> other attached packages:
#> [1] scuttle_1.23.2 SingleCellExperiment_1.35.2
#> [3] SummarizedExperiment_1.43.0 Biobase_2.73.2
#> [5] BiocDuckDB_0.99.22 DuckDBGRanges_0.99.8
#> [7] GenomicRanges_1.65.4 Seqinfo_1.3.2
#> [9] DuckDBArray_0.99.8 DelayedArray_0.39.6
#> [11] SparseArray_1.13.2 S4Arrays_1.13.0
#> [13] abind_1.4-8 MatrixGenerics_1.25.0
#> [15] matrixStats_1.5.0 Matrix_1.7-6
#> [17] DuckDBDataFrame_0.99.26 IRanges_2.47.5
#> [19] S4Vectors_0.51.9 BiocGenerics_0.59.12
#> [21] generics_0.1.4 bit64_4.8.6
#> [23] BiocStyle_2.41.0
#>
#> loaded via a namespace (and not attached):
#> [1] DBI_1.3.0 rlang_1.3.0
#> [3] magrittr_2.0.5 otel_0.2.0
#> [5] e1071_1.7-17 compiler_4.6.1
#> [7] vctrs_0.7.3 pkgconfig_2.0.3
#> [9] SpatialExperiment_1.23.0 fastmap_1.2.0
#> [11] dbplyr_2.6.0 magick_2.9.1
#> [13] XVector_0.53.0 rmarkdown_2.32
#> [15] purrr_1.2.2 bit_4.6.0
#> [17] xfun_0.60 MultiAssayExperiment_1.39.1
#> [19] bluster_1.23.1 cachem_1.1.0
#> [21] beachmat_2.29.2 jsonlite_2.0.0
#> [23] blob_1.3.0 BiocParallel_1.47.0
#> [25] irlba_2.3.7 parallel_4.6.1
#> [27] cluster_2.1.8.3 R6_2.6.1
#> [29] MultiAssaySpatialExperiment_0.99.12 bslib_0.12.0
#> [31] limma_3.99.0 jquerylib_0.1.4
#> [33] Rcpp_1.1.2 bookdown_0.48
#> [35] assertthat_0.2.1 knitr_1.52
#> [37] igraph_2.3.3 tidyselect_1.2.1
#> [39] yaml_2.3.12 codetools_0.2-20
#> [41] lattice_0.23-1 tibble_3.3.1
#> [43] withr_3.0.3 evaluate_1.0.5
#> [45] sf_1.1-2 units_1.0-1
#> [47] proxy_0.4-29 pillar_1.11.1
#> [49] BiocManager_1.30.27 KernSmooth_2.23-27
#> [51] class_7.3-24 glue_1.8.1
#> [53] metapod_1.21.0 tools_4.6.1
#> [55] BiocNeighbors_2.7.3 ScaledMatrix_1.21.0
#> [57] locfit_1.5-9.12 scran_1.41.1
#> [59] grid_4.6.1 edgeR_4.99.4
#> [61] duckdb_1.5.5 BiocSingular_1.29.1
#> [63] cli_3.6.6 rsvd_1.0.5
#> [65] arrow_25.0.1 dplyr_1.2.1
#> [67] airway_1.33.2 sass_0.4.10
#> [69] digest_0.6.39 classInt_0.4-11
#> [71] dqrng_0.4.1 rjson_0.2.23
#> [73] htmltools_0.5.9 lifecycle_1.0.5
#> [75] statmod_1.5.2