Introduction to BiocDuckDB

Introduction

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, and
  • readParquet() 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.

Installation

if (!require("BiocManager"))
    install.packages("BiocManager")
BiocManager::install("BiocDuckDB")
library(BiocDuckDB)
library(SummarizedExperiment)

The round-trip

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     8

The storage layout

writeParquet() 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.

Targeting the storage contract

The datapackage.json layout is a stable Frictionless profile — BiocDuckDB’s extensions are additive and base-conformant — so any producer can target it and any consumer can attach it. These are the write and read halves of an ingest/attach contract, exposed as public functions so a store that already holds COO Parquet does not need to reconstruct an in-memory object to interoperate with BiocDuckDB.

To assemble a manifest from resources written incrementally (a dataset streamed in parts, or promoted from another store), collect the descriptors returned by the primitive writeParquet() methods and call 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__"))

The documented model values and the dimension x layout dispatch that maps each resource to its accessor are described on ?writeParquet and ?readParquet.

Operations run as SQL

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
#>            <character> <character> <character>
#> 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    1043616

Because 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 = object.size(airway),
  duckdb    = object.size(airway_ddb)) |> format(units = "MB")
#>  in_memory     duckdb 
#> "84440888" "14918848"

Genomic coordinates are preserved

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              14

Single-cell data

Single-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  249

The filter, realize, analyze pattern

For 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.

Multi-omics

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.

When to use BiocDuckDB

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.

Session information

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.1              SingleCellExperiment_1.35.2
#>  [3] SummarizedExperiment_1.43.0 Biobase_2.73.1             
#>  [5] BiocDuckDB_0.99.5           DuckDBGRanges_0.99.1       
#>  [7] GenomicRanges_1.65.1        Seqinfo_1.3.0              
#>  [9] DuckDBArray_0.99.2          DelayedArray_0.39.3        
#> [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-5               
#> [17] DuckDBDataFrame_0.99.6      IRanges_2.47.2             
#> [19] S4Vectors_0.51.5            BiocGenerics_0.59.10       
#> [21] generics_0.1.4              bit64_4.8.2                
#> [23] BiocStyle_2.41.0           
#> 
#> loaded via a namespace (and not attached):
#>  [1] tidyselect_1.2.1                   blob_1.3.0                        
#>  [3] dplyr_1.2.1                        arrow_25.0.0                      
#>  [5] fastmap_1.2.0                      bluster_1.23.0                    
#>  [7] duckdb_1.5.4.3                     digest_0.6.39                     
#>  [9] rsvd_1.0.5                         lifecycle_1.0.5                   
#> [11] cluster_2.1.8.2                    sf_1.1-1                          
#> [13] statmod_1.5.2                      magrittr_2.0.5                    
#> [15] compiler_4.6.1                     rlang_1.3.0                       
#> [17] sass_0.4.10                        tools_4.6.1                       
#> [19] igraph_2.3.3                       yaml_2.3.12                       
#> [21] knitr_1.51                         dqrng_0.4.1                       
#> [23] bit_4.6.0                          classInt_0.4-11                   
#> [25] BiocParallel_1.47.0                KernSmooth_2.23-26                
#> [27] withr_3.0.3                        purrr_1.2.2                       
#> [29] sys_3.4.3                          grid_4.6.1                        
#> [31] beachmat_2.29.0                    e1071_1.7-17                      
#> [33] edgeR_4.11.4                       MultiAssayExperiment_1.39.0       
#> [35] cli_3.6.6                          rmarkdown_2.31                    
#> [37] otel_0.2.0                         metapod_1.21.0                    
#> [39] rjson_0.2.23                       DBI_1.3.0                         
#> [41] cachem_1.1.0                       proxy_0.4-29                      
#> [43] assertthat_0.2.1                   parallel_4.6.1                    
#> [45] BiocManager_1.30.27                XVector_0.53.0                    
#> [47] vctrs_0.7.3                        jsonlite_2.0.0                    
#> [49] BiocSingular_1.29.0                BiocNeighbors_2.7.2               
#> [51] irlba_2.3.7                        maketools_1.3.2                   
#> [53] magick_2.9.1                       locfit_1.5-9.12                   
#> [55] limma_3.69.2                       jquerylib_0.1.4                   
#> [57] units_1.0-1                        glue_1.8.1                        
#> [59] codetools_0.2-20                   ScaledMatrix_1.21.0               
#> [61] tibble_3.3.1                       pillar_1.11.1                     
#> [63] htmltools_0.5.9                    MultiAssaySpatialExperiment_0.99.2
#> [65] R6_2.6.1                           dbplyr_2.6.0                      
#> [67] evaluate_1.0.5                     lattice_0.22-9                    
#> [69] SpatialExperiment_1.23.0           scran_1.41.1                      
#> [71] bslib_0.11.0                       class_7.3-23                      
#> [73] Rcpp_1.1.2                         xfun_0.60                         
#> [75] buildtools_1.0.0                   pkgconfig_2.0.3