DuckDBGRanges 0.99.4
This vignette is for developers: it documents how DuckDBGRanges represents genomic ranges on top of DuckDBDataFrame, how range operations become SQL, and how the classes are structured so you can maintain or extend the package. For day-to-day use see Introduction to DuckDBGRanges; for performance see Benchmarking DuckDBGRanges.
library(DuckDBGRanges)
library(GenomicRanges)DuckDBGRanges layers two classes on the DuckDBDataFrame foundation, each extending its in-memory GenomicRanges counterpart:
| Class | Extends | Role |
|---|---|---|
DuckDBGRanges |
GRanges |
genomic intervals with a DuckDB backend |
DuckDBGRangesList |
GRangesList |
grouped features (transcripts, per-gene peaks) |
Each holds its coordinates in a DuckDBDataFrame (from
DuckDBDataFrame) describing a Parquet file, rather than in memory.
Because a DuckDBGRanges is a GRanges, the full GenomicRanges
API applies unchanged; the package overrides the operations that benefit from
running as SQL, and falls back to materialization for the rest.
showClass("DuckDBGRanges")
#> Class "DuckDBGRanges" [package "DuckDBGRanges"]
#>
#> Slots:
#>
#> Name: frame seqinfo elementMetadata elementType
#> Class: DuckDBDataFrame Seqinfo DataFrame character
#>
#> Name: metadata
#> Class: list
#>
#> Extends:
#> Class "GenomicRanges", directly
#> Class "Ranges", by class "GenomicRanges", distance 2
#> Class "GenomicRanges_OR_missing", by class "GenomicRanges", distance 2
#> Class "GenomicRanges_OR_GenomicRangesList", by class "GenomicRanges", distance 2
#> Class "GenomicRanges_OR_GRangesList", by class "GenomicRanges", distance 2
#> Class "List", by class "GenomicRanges", distance 3
#> Class "Vector", by class "GenomicRanges", distance 4
#> Class "list_OR_List", by class "GenomicRanges", distance 4
#> Class "Annotated", by class "GenomicRanges", distance 5
#> Class "vector_OR_Vector", by class "GenomicRanges", distance 5A DuckDBGRanges carries three slots:
frame: a DuckDBDataFrame with the coordinate columns (see The five
required columns). Extra columns become metadata columns.seqinfo: a Seqinfo with sequence lengths and genome build, kept
materialized (it is small and needed for validation).elementMetadata: an optional DuckDBDataFrame of metadata columns
(mcols) sharing the same backend as frame.DuckDBGRangesList has the same slots, but its coordinate columns are DuckDB
LIST[] columns (one row per group) and its elementMetadata is per group
rather than per range (see Grouped ranges as LIST columns).
Every DuckDBGRanges resolves to five coordinate columns, mirroring the
components of a GRanges:
| Column | Meaning |
|---|---|
seqnames |
sequence/chromosome name ("chr1", …) |
start |
start position (1-based, inclusive) |
end |
end position (1-based, inclusive) |
width |
interval width (end - start + 1) |
strand |
"+", "-", or "*" |
The constructor names the columns backing these components; anything else named in
mcols= becomes a metadata column. Only the columns an operation actually needs
are ever read from the Parquet file.
gr <- GRanges(c("chr1:100-150:+", "chr1:200-275:-", "chr2:300-360:+"),
gene_id = paste0("GENE", 1:3))
gr_df <- as.data.frame(gr)
gr_df$range_id <- paste0("range", seq_len(nrow(gr_df)))
gr_path <- tempfile(fileext = ".parquet"); arrow::write_parquet(gr_df, gr_path)
gr_ddb <- DuckDBGRanges(gr_path,
seqnames = "seqnames", start = "start",
end = "end", strand = "strand",
keycol = list(range_id = gr_df$range_id),
mcols = "gene_id")
gr_ddb
#> DuckDBGRanges object with 3 ranges and 1 metadata column:
#> seqnames start end width strand | gene_id
#> <character> <integer> <integer> <integer> <character> | <character>
#> range1 chr1 100 150 51 + | GENE1
#> range2 chr1 200 275 76 - | GENE2
#> range3 chr2 300 360 61 + | GENE3
#> -------
#> seqinfo: 2 sequences from an unspecified genome; no seqlengthsThe performance story rests on one idea: rather than materialize the ranges and
run the in-memory algorithm, DuckDBGRanges translates the
operations that map cleanly onto set-based queries into SQL on the underlying
DuckDBDataFrame, so DuckDB does the scan and the arithmetic and only the result
crosses back into R.
Coordinate filters (restrict(), logical subsetting) become WHERE
clauses, which DuckDB pushes into the Parquet reader, so a region filter reads
only the matching row groups instead of scanning everything:
-- restrict(x, start = s, end = e) becomes
SELECT * FROM ranges WHERE start >= s AND end <= eOverlap joins (findOverlaps() and friends) become a self- or cross-join on
the interval-overlap predicate:
SELECT a.id, b.id
FROM ranges_a a JOIN ranges_b b
ON a.seqnames = b.seqnames
AND a.start <= b.end AND a.end >= b.start
AND (a.strand = b.strand OR a.strand = '*' OR b.strand = '*')Distance becomes a CASE expression over the gap between two intervals:
SELECT CASE
WHEN a.end < b.start THEN b.start - a.end
WHEN b.end < a.start THEN a.start - b.end
ELSE 0
END AS distance
FROM ranges_a a, ranges_b bFilters compose lazily, a chain of restrict() and logical subsets becomes a
single query, evaluated only when results are pulled. Operations without a natural
set-based form fall back to materializing the (already filtered) ranges to an
in-memory GRanges, which is correct for any input and cheap once the data has
been cut down. Benchmarking DuckDBGRanges shows which operations win in SQL and
which are better left to the in-memory interval-tree implementations.
DuckDBGRangesList stores grouped features using DuckDB LIST[] columns: each
row is one group, and the coordinate columns hold lists of per-range values.
transcript_id | seqnames | start | width | strand
ENST001 | [chr1,chr1,chr1]| [100,200,300]| [50,75,60] | [+,-,+]
ENST002 | [chr2,chr2] | [150,250] | [80,45] | [-,+]This is compact (one row per group instead of one per range), matches the
biological grouping directly, and lets DuckDB’s list functions operate on groups
without unnesting. elementNROWS() reads the list lengths; extracting a group
unnests it back to a GRanges.
grl_data <- data.frame(
transcript_id = c("ENST001", "ENST002"),
seqnames = I(list(rep("chr1", 3), rep("chr2", 2))),
start = I(list(c(100L, 200L, 300L), c(150L, 250L))),
width = I(list(c(50L, 75L, 60L), c(80L, 45L))),
strand = I(list(c("+", "-", "+"), c("-", "+"))))
grl_path <- tempfile(fileext = ".parquet"); arrow::write_parquet(grl_data, grl_path)
grl_ddb <- DuckDBGRangesList(grl_path,
seqnames = "seqnames", start = "start",
width = "width", strand = "strand",
keycol = list(transcript_id = grl_data$transcript_id))
elementNROWS(grl_ddb)
#> ENST001 ENST002
#> 3 2A DuckDBGRanges coerces to an ordinary GRanges with as(x, "GRanges"), which
runs the accumulated query and reads the result into memory, the materialize
half of the filter, then materialize pattern. Coercion is how overlap-heavy
work is handed off to the in-memory interval-tree algorithms.
as(gr_ddb[seqnames(gr_ddb) == "chr1"], "GRanges")
#> GRanges object with 2 ranges and 1 metadata column:
#> seqnames ranges strand | gene_id
#> <Rle> <IRanges> <Rle> | <character>
#> range2 chr1 200-275 - | GENE2
#> range1 chr1 100-150 + | GENE1
#> -------
#> seqinfo: 2 sequences from an unspecified genome; no seqlengthsBecause the on-disk representation is plain Parquet, the same files are readable outside R (Python, Julia, DuckDB itself, cloud query engines), the property the BiocDuckDB package relies on for its language-neutral on-disk experiment format. A new range-bearing backend mostly needs to define how its data maps onto the five coordinate columns and which operations to push into SQL; everything else follows from the GenomicRanges contract.
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