DuckDBGRanges 0.99.4
The Introduction to DuckDBGRanges vignette shows how a DuckDBGRanges behaves
like an ordinary GRanges while keeping its data on disk. This vignette asks the
follow-up question: what does that cost, and what does it buy? We compare the
DuckDB backend against an in-memory GenomicRanges object on two
range-analysis scenarios at scale:
The headline results below were produced offline (see Benchmark setup) and are rendered here from a bundled results file, so this vignette builds quickly. The A small, live comparison section runs a miniature version at build time so you can see the mechanics end to end.
The short version: DuckDBGRanges wins decisively on memory and on the filtering
that dominates the early stages of a workflow, while in-memory GRanges keeps its
edge on the interval-tree operations, findOverlaps(), subsetByOverlaps(),
nearest(). The two are complementary, which is why the idiomatic pattern is to
filter with DuckDB, then materialize to GRanges for overlap work.
DuckDBGRanges turns coordinate and metadata filters into SQL that
DuckDB pushes into the Parquet reader, so a region restriction reads only the
matching row groups instead of scanning every range. Existing Bioconductor code
picks this up automatically: calling restrict() on a DuckDBGRanges is a
WHERE clause, and chained filters compose into a single lazy query.
| Category | Operations | Backend that wins |
|---|---|---|
| Memory footprint | the R object itself | DuckDBGRanges (constant ~175 KB) |
| Coordinate / metadata filtering | restrict, logical [, pipelines |
DuckDBGRanges (predicate pushdown) |
| Nearest-feature distance | distanceToNearest |
DuckDBGRanges (SQL aggregation) |
| Inter-range transforms | reduce, disjoin |
comparable; DuckDBGRanges pulls ahead at scale |
| Overlap enumeration | findOverlaps, subsetByOverlaps |
GRanges (interval trees) |
The design vignette (Design and extension of DuckDBGRanges) shows the SQL each operation translates to.
To show the mechanics without a large dataset, we build a small range set and run
one operation on both an in-memory GRanges and a DuckDBGRanges, confirming the
results agree.
library(DuckDBGRanges)
library(GenomicRanges)
library(arrow)
set.seed(1L)
n <- 2000L
df <- data.frame(
seqnames = sample(paste0("chr", 1:5), n, replace = TRUE),
start = sample(1:1e6, n),
strand = sample(c("+", "-", "*"), n, replace = TRUE),
score = runif(n))
df$end <- df$start + sample(100:500, n, replace = TRUE)
gr <- makeGRangesFromDataFrame(df, keep.extra.columns = TRUE)
path <- tempfile(fileext = ".parquet"); write_parquet(df, path)
ddb <- DuckDBGRanges(path, seqnames = "seqnames", start = "start",
end = "end", strand = "strand", mcols = "score")
## same answer, one in memory and one queried from disk
system.time(r_gr <- restrict(gr, start = 1L, end = 500000L))
#> user system elapsed
#> 0.027 0.001 0.028
system.time(r_ddb <- restrict(ddb, start = 1L, end = 500000L))
#> user system elapsed
#> 0.153 0.002 0.155
length(r_gr) == length(r_ddb)
#> [1] TRUEAt this size the in-memory object is faster, there is nothing to gain from going to disk, and the query has fixed overhead. The advantage appears at scale, which is what the offline benchmark measures.
The full benchmark generates synthetic ranges at scale across the human autosomes: 1,000,000 peaks for the scATAC-seq scenario and 10,000,000 SNVs for the variant scenario, plus a 50,000-feature annotation set used as the subject of the overlap operations. For each scenario it times six operations on both backends and records the in-memory footprint of each object.
The backends parallelize differently, so we give each its best effort on the same
core budget: GRanges runs single-threaded C code (the in-memory baseline),
while DuckDBGRanges autotunes DuckDB’s internal threads up to the budget with
no configuration. The exact parameters are recorded with the results and shown
beneath the tables.
Rendered from the bundled offline results (inst/scripts/benchmark_results.rds);
regenerate with inst/scripts/run_vignette_benchmarks.R (see that script’s header).
Memory footprint
| Scenario | GRanges (MB) | DuckDBGRanges (MB) | Smaller by (x) |
|---|---|---|---|
| scATAC | 27.573 | 0.182 | 152 |
| variant | 275.644 | 0.182 | 1518 |
scATAC-seq scenario (1M peaks)
| Operation | GRanges (s) | DuckDBGRanges (s) | Speedup (x) | Winner |
|---|---|---|---|---|
| restrict | 0.327 | 0.138 | 2.4 | DuckDBGRanges |
| reduce | 0.417 | 0.626 | 0.7 | GRanges |
| disjoin | 0.983 | 0.820 | 1.2 | DuckDBGRanges |
| findOverlaps | 0.249 | 23.307 | 0.0 | GRanges |
| subsetByOverlaps | 0.251 | 23.865 | 0.0 | GRanges |
| distanceToNearest | 5.131 | 0.495 | 10.4 | DuckDBGRanges |
Variant scenario (10M variants)
| Operation | GRanges (s) | DuckDBGRanges (s) | Speedup (x) | Winner |
|---|---|---|---|---|
| restrict | 3.811 | 0.297 | 12.8 | DuckDBGRanges |
| reduce | 4.707 | 3.634 | 1.3 | DuckDBGRanges |
| disjoin | 9.500 | 4.322 | 2.2 | DuckDBGRanges |
| findOverlaps | 2.248 | 227.381 | 0.0 | GRanges |
| subsetByOverlaps | 2.271 | 234.302 | 0.0 | GRanges |
| distanceToNearest | 34.188 | 0.488 | 70.1 | DuckDBGRanges |
Configuration: scATAC-seq scenario = 1e+06 peaks, variant scenario = 1e+07 variants, 16-core budget. GRanges: in-memory, single-threaded. DuckDBGRanges: on-disk Parquet; DuckDB autotunes threads up to the core budget.
The two backends have complementary strengths. A DuckDBGRanges is a small,
constant-size handle (about 175 KB) regardless of how many ranges it fronts,
whereas a GRanges grows with the data (roughly 150 times smaller for a million
peaks, and over 1,500 times smaller for ten million variants). That is what makes
it practical to reference a very large interval catalogue and pay only for the
queries you actually run.
Operations that reduce data favor the DuckDB backend. restrict() and logical
subsetting compile to WHERE clauses that read only matching rows, and the
advantage grows with size (a few times faster at a million peaks, about 13 times
at ten million variants). distanceToNearest() runs as an SQL aggregation rather
than an interval-tree walk, and shows the largest single speedup measured (about
10 times at a million peaks, about 70 times at ten million variants). The
inter-range transforms reduce() and disjoin() are closer: in-memory GRanges
edges out reduce() at a million peaks, but DuckDB pulls ahead on both by ten
million ranges.
Overlap enumeration favors GRanges. findOverlaps() and subsetByOverlaps()
rely on in-memory interval trees that a set-based SQL join does not match
(GRanges is roughly 100 times faster), so those are best run on a materialized
GRanges. In practice the two compose well: use DuckDB to filter a large set down
cheaply, then hand the small result to GRanges for the interval-tree work.
| Feature | GRanges | DuckDBGRanges |
|---|---|---|
| Storage | Memory | Parquet file |
| Object footprint | Proportional to N | Constant (~175 KB) |
| Filtering | Scans all ranges | Predicate pushdown (faster) |
| Nearest-feature distance | Interval trees | SQL aggregation (faster) |
| Overlap enumeration | Interval trees (faster) | SQL join |
reduce / disjoin |
Optimized C | SQL (comparable; faster at scale) |
| Persistence | Session only | Files persist |
| Read by other languages | R only | R, Python, Julia, … |
GRanges, data fits comfortably in RAM; the fastest choice for
overlap-heavy work and small-to-mid range sets.DuckDBGRanges, larger-than-memory interval sets, filter-heavy pipelines,
keeping a huge annotation referenceable at constant memory cost, or when the
Parquet files are shared with non-R tooling. Combine with GRanges via filter,
then materialize.inst/scripts/run_vignette_benchmarks.R reproduces these numbers on your hardware
and scales to other range counts (BENCH_SCATAC_PEAKS, BENCH_VARIANTS,
BENCH_CORES; set small counts to smoke-test first). It writes
benchmark_results.rds, which this vignette renders via
inst/scripts/make_timings_table.R.
For higher-level genomic workflows built on this backend, see the BiocDuckDB package.
sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.4 LTS
#>
#> Matrix products: default
#> BLAS: /home/biocbuild/bbs-3.24-bioc/R/lib/libRblas.so
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.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] arrow_25.0.0 DuckDBGRanges_0.99.4 GenomicRanges_1.65.1
#> [4] Seqinfo_1.3.0 DuckDBDataFrame_0.99.20 IRanges_2.47.2
#> [7] S4Vectors_0.51.6 BiocGenerics_0.59.11 generics_0.1.4
#> [10] bit64_4.8.2 BiocStyle_2.41.0
#>
#> loaded via a namespace (and not attached):
#> [1] sass_0.4.10 SparseArray_1.13.2 lattice_0.22-9
#> [4] digest_0.6.39 magrittr_2.0.5 evaluate_1.0.5
#> [7] grid_4.6.1 bookdown_0.47 blob_1.3.0
#> [10] fastmap_1.2.0 jsonlite_2.0.0 Matrix_1.7-6
#> [13] DBI_1.3.0 BiocManager_1.30.27 purrr_1.2.2
#> [16] jquerylib_0.1.4 abind_1.4-8 duckdb_1.5.5
#> [19] cli_3.6.6 rlang_1.3.0 dbplyr_2.6.0
#> [22] XVector_0.53.0 withr_3.0.3 cachem_1.1.0
#> [25] DelayedArray_0.39.4 yaml_2.3.12 otel_0.2.0
#> [28] S4Arrays_1.13.0 tools_4.6.1 dplyr_1.2.1
#> [31] assertthat_0.2.1 vctrs_0.7.3 R6_2.6.1
#> [34] matrixStats_1.5.0 lifecycle_1.0.5 bit_4.6.0
#> [37] pkgconfig_2.0.3 bslib_0.12.0 pillar_1.11.1
#> [40] glue_1.8.1 xfun_0.60 tibble_3.3.1
#> [43] tidyselect_1.2.1 MatrixGenerics_1.25.0 knitr_1.51
#> [46] htmltools_0.5.9 rmarkdown_2.31 compiler_4.6.1