--- title: "Design and extension of DuckDBGRanges" author: - name: Patrick Aboyoun email: aboyounp@gene.com affiliation: Genentech, Inc. date: "Compiled: `r BiocStyle::doc_date()`; Modified: 6 July 2026" package: DuckDBGRanges vignette: > %\VignetteIndexEntry{3. Design and extension of DuckDBGRanges} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} output: BiocStyle::html_document: number_sections: true toc: true toc_depth: 3 --- ```{r setup, include=FALSE} library(BiocStyle) knitr::opts_chunk$set(collapse = TRUE, comment = "#>", error = FALSE, warning = FALSE, message = FALSE) ``` # Scope This vignette is for developers: it documents how `r Biocpkg("DuckDBGRanges")` represents genomic ranges on top of `r Biocpkg("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*. ```{r load} library(DuckDBGRanges) library(GenomicRanges) ``` # Architecture `r Biocpkg("DuckDBGRanges")` layers two classes on the `r Biocpkg("DuckDBDataFrame")` foundation, each extending its in-memory `r Biocpkg("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 `r Biocpkg("DuckDBDataFrame")`) describing a Parquet file, rather than in memory. Because a `DuckDBGRanges` *is* a `GRanges`, the full `r Biocpkg("GenomicRanges")` API applies unchanged; the package overrides the operations that benefit from running as SQL, and falls back to materialization for the rest. ```{r classes} showClass("DuckDBGRanges") ``` ## Slots A `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]). # The five required 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. ```{r construct} 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 ``` # From range operations to SQL The performance story rests on one idea: rather than materialize the ranges and run the in-memory algorithm, `r Biocpkg("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: ```sql -- restrict(x, start = s, end = e) becomes SELECT * FROM ranges WHERE start >= s AND end <= e ``` **Overlap joins** (`findOverlaps()` and friends) become a self- or cross-join on the interval-overlap predicate: ```sql 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: ```sql 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 b ``` Filters 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. # Grouped ranges as LIST columns `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`. ```{r grl} 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) ``` # Materialization and interoperability A `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. ```{r materialize} as(gr_ddb[seqnames(gr_ddb) == "chr1"], "GRanges") ``` Because the on-disk representation is plain Parquet, the same files are readable outside R (Python, Julia, DuckDB itself, cloud query engines) --- the property the `r Biocpkg("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 `r Biocpkg("GenomicRanges")` contract. # Session information ```{r sessioninfo} sessionInfo() ```