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.
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:
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:
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.
Step 3’s runPCA() is a partial exception to needing step
2 first: BiocSingular’s
runSVD() (and so scran::fixedPCA() /
scater::runPCA()) already compute a correct PCA directly on
sce_ddb, through the SQL-pushdown
%*%/crossprod methods DuckDBArray
implements, no realization required for correctness. Speed is a separate
question: the ordinary BSPARAM = BiocSingular::IrlbaParam()
calls %*% once per solver iteration, so on a
DuckDBMatrix that cost is paid repeatedly rather than once,
covered in detail (with measured numbers) in Benchmarking
BiocDuckDB. Passing BSPARAM = DuckDBIrlbaParam()
instead closes that gap: it materializes just the matrix
runPCA() touches (already filtered and HVG-selected here,
not the whole SingleCellExperiment the way step 2’s
as(sce_ddb, "SingleCellExperiment") would) and drives
irlba directly on it, falling back to the ordinary path
whenever that materialization doesn’t apply. It is a lighter-weight,
PCA-specific alternative to full realization, not just a slower
substitute for it.
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 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04 LTS
#>
#> Matrix products: default
#> BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so; LAPACK version 3.12.0
#>
#> locale:
#> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
#> [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: Etc/UTC
#> tzcode source: system (glibc)
#>
#> 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.23 DuckDBGRanges_0.99.8
#> [7] GenomicRanges_1.65.4 Seqinfo_1.3.2
#> [9] DuckDBArray_0.99.9 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 assertthat_0.2.1
#> [35] knitr_1.52 igraph_2.3.3
#> [37] tidyselect_1.2.1 yaml_2.3.12
#> [39] codetools_0.2-20 lattice_0.23-1
#> [41] tibble_3.3.1 withr_3.0.3
#> [43] evaluate_1.0.5 sf_1.1-2
#> [45] units_1.0-1 proxy_0.4-29
#> [47] pillar_1.11.1 BiocManager_1.30.27
#> [49] KernSmooth_2.23-27 class_7.3-24
#> [51] glue_1.8.1 metapod_1.21.0
#> [53] maketools_1.3.2 tools_4.6.1
#> [55] BiocNeighbors_2.7.3 sys_3.4.3
#> [57] ScaledMatrix_1.21.0 locfit_1.5-9.12
#> [59] buildtools_1.0.0 scran_1.41.1
#> [61] grid_4.6.1 edgeR_4.99.5
#> [63] duckdb_1.5.5 BiocSingular_1.29.1
#> [65] cli_3.6.6 rsvd_1.0.5
#> [67] arrow_25.0.1 dplyr_1.2.1
#> [69] airway_1.33.2 sass_0.4.10
#> [71] digest_0.6.39 classInt_0.4-11
#> [73] dqrng_0.4.1 rjson_0.2.23
#> [75] htmltools_0.5.9 lifecycle_1.0.5
#> [77] statmod_1.5.2