Approximate string matching inside DuckDB.
DuckHTS exposes FASTA and FASTQ records as DuckDB relations. Sassy supplies the approximate matching kernel. ducksassy connects them: file searches return matches with edit cost, coordinates, strand and CIGAR alongside the input record’s columns. Scalar functions search sequence values from any relation. They also search ordinary ASCII text columns, so fuzzy grep can be a SQL query.
This follows DuckHTS’s composability model: the file reader and native matcher handle their respective mechanics; DuckDB handles the join, filtering and aggregation. For example, give each guide its own PAM length and search a FASTA relation:
WITH guides(guide_id, guide, pam_length) AS (
VALUES ('g1', 'ACGTNGG', 3), ('g2', 'ACGTAGG', 2)
)
SELECT r.name AS reference_name, g.guide_id,
hit.text_start, hit.text_end, hit.strand
FROM guides AS g
CROSS JOIN read_fasta('test/data/references.fasta', scan_mode := 'sequential') AS r
CROSS JOIN LATERAL unnest(sassy_crispr_matches(g.guide, r.sequence, 0,
pam_length := g.pam_length)) AS matches(hit)
ORDER BY reference_name, guide_id, hit.text_start, hit.text_end, hit.strand;| reference_name | guide_id | text_start | text_end | strand |
|---|---|---|---|---|
| forward | g1 | 2 | 9 | + |
| forward | g2 | 2 | 9 | + |
| masked | g1 | 2 | 9 | + |
| masked | g2 | 2 | 9 | + |
| reverse | g1 | 2 | 9 | - |
| reverse | g2 | 2 | 9 | - |
pam_length comes from each guide row. The lateral join expands its matches;
the output can be joined to other annotations or grouped by guide without
another file format. Use LEFT JOIN LATERAL ... ON true to retain rows with no
hit.
- Search with up to
ksubstitutions, insertions or deletions. IUPAC ambiguity codes work in patterns and targets; reverse-complement search is on by default. - Search one pattern, a panel, or a guide with a PAM. CRISPR search applies Sassy’s exact PAM endpoint filter and excludes targets above the configured N fraction. The CRISPR semantics below spell out the filter.
- Search files through DuckHTS, or call the scalar functions on any sequence column. Results stay in DuckDB for joins, counts and downstream queries.
The library and adapter comparison measures direct upstream Sassy calls, the Rust C ABI and SQL over identical resident inputs, checking the complete hit multisets. The buffer reuse measurements report the current SQL overhead and the effect of reusing native result storage.
The measured workloads include eight 23-base
guides against a 4.6 Mb E. coli reference and 262,144 sequence rows. The timed
comparison includes Ducksassy’s scalar and AVX2 backends and an upstream
Sassy CLI baseline for the same CRISPR workload. The report records startup,
threading, output work and cache conditions: persistent SQL aggregates hits,
while the upstream CLI starts a process and writes TSV.
The 4.6 Mb input is one bacterial record; these results do not establish
whole-human-genome throughput or peak memory use.
DuckDB v1.5.5 can load the stable v1.2.0 C API artifact without DuckHTS or an external SQL bootstrap:
make sdk-v1 release-v1
.deps-v1/cli/duckdb -unsigned -no-initLOAD 'build-v1/ducksassy.duckdb_extension';
SELECT sassy_count('ACGT', 'TTACGT', 0, 'dna', false);
SELECT * FROM sassy_grep('error', 'error: disk full', 0);
SELECT * FROM sassy_backend_info();Download the CLI separately as described in
host architecture, setup and measurements. Every v1 function,
with signatures and examples, is listed in the function reference.
V1 LOAD writes no catalog entries and works on a read-only database.
Scalars use positional trailing options: alphabet, rc, all_endpoints, cigar_format for matches, and pam_length, allow_pam_edits, max_n_frac, rc, cigar_format for CRISPR. Defaults and VARCHAR/BLOB support match v2. V1 has no file or
relation-name helpers; compose native calls with a table, CTE or explicitly
loaded DuckHTS reader:
LOAD '/absolute/path/duckhts.duckdb_extension';
SELECT r.*, hit
FROM read_fasta('reference.fa', scan_mode := 'sequential') r
CROSS JOIN LATERAL unnest(
sassy_matches('ACGTAGG', r.sequence, 1, 'iupac', false, false, 'both')
) AS matches(hit);Use read_fastq(...) for reads, _many with a pattern list for panels, or
sassy_crispr_matches(guide, r.sequence, k) for CRISPR. The v1 build uses
DUCKSASSY_HOST=v1; standalone CMake defaults to v2. Keep their build directories separate.
Distribution CI uses make configure release test_release, which produces the
v1 artifact under build/release/; see the CI contract.
The distribution matrix covers Linux and macOS x86-64/ARM64, Windows x86-64
MinGW/Rtools, and the default DuckDB-Wasm MVP/EH bundles. MSVC is unsupported.
Wasm uses scalar search and a fail-stop Rust panic ABI: an internal panic
invalidates that worker, which the caller terminates before starting a fresh
worker that can load and query the extension. The optional
COI/threads bundle requires an atomics-enabled Rust standard library. See the
portable build evidence.
The v2 artifact targets the DuckDB C API v2 preview.
make setup fetches the matching CLI and checks the bundled SDK. On Linux x86-64 you need C/C++
compilers, CMake ≥ 3.20, Python ≥ 3.11, Git, R ≥ 4.1 and rustup.
make setup JOBS=4
make release-v2
.deps/duckdb-build/duckdb -unsigned -no-initThen, in the DuckDB CLI:
SET GLOBAL lambda_syntax='ENABLE_SINGLE_ARROW';
INSTALL '.deps/duckhts.duckdb_extension';
INSTALL 'build/ducksassy.duckdb_extension';
.read sql/ducksassy.sqlThe lambda setting is required by DuckHTS 1.5.2 and applies to the whole
database instance, so use a dedicated one. V2 keeps its connection-scoped
macros, := scalar options and file/relation helpers. Positional trailing
options are v1-only: a same-named macro shadows a native function in v2
(binding probe).
Supported host versions are pinned in
ducksassy-package.json. When upgrading, move an
outdated .deps/sassy-source checkout aside before rerunning make setup.
Rducksassy builds the stable C API v1 extension
from bundled sources and provides rducksassy_connect() and
rducksassy_load() for the CRAN duckdb package.
install.packages("Rducksassy")
con <- Rducksassy::rducksassy_connect()
DBI::dbGetQuery(con, "SELECT * FROM sassy_grep('timeout', 'request timedout', 1)")
DBI::dbDisconnect(con, shutdown = TRUE)Install and load the signed DuckHTS community extension when FASTA, FASTQ or BAM readers are needed; the Ducksassy scalar functions accept its sequence columns.
make vendor-rust refreshes the shared Rust source archive using the committed
lockfile. Native and R package builds use it offline. make r-package stages
the R package’s sources and builds its source tarball.
Every example below runs against the bundled test fixtures and is rendered live from this README.
Search a text column with one allowed edit. ASCII mode requires rc := false;
matches return byte offsets within each value.
WITH lines(line_id, message) AS (
VALUES (1, 'error: disk full'),
(2, 'errot: disk full'),
(3, 'ready')
)
SELECT line_id, message, hit.text_start, hit.text_end, hit.cost
FROM lines
CROSS JOIN LATERAL unnest(sassy_matches('error', message, 1,
alphabet := 'ascii', rc := false)) AS matches(hit)
ORDER BY line_id, hit.text_start;| line_id | message | text_start | text_end | cost |
|---|---|---|---|---|
| 1 | error: disk full | 0 | 5 | 0 |
| 2 | errot: disk full | 0 | 5 | 1 |
The input can also be a CSV, Parquet file or any other DuckDB relation with a text column. ASCII matching is case-insensitive; it is not Unicode text search.
For one long ASCII value, sassy_grep returns rows from a stateful scan. It
searches overlapping 1 KiB regions and gives each match endpoint to one region.
DuckDB can stop requesting regions once an outer LIMIT is satisfied:
SELECT text_start, text_end, cost, cigar
FROM sassy_grep('error', 'error: ' || repeat('ready ', 100000), 1)
LIMIT 1;This table function accepts constant expressions, as DuckDB requires for table
function arguments. It returns all qualifying endpoints, including nearby
alternative alignments; sassy_matches defaults to rightmost local minima.
The scalar function above still accepts a different text value on every row.
DuckDB constructs the constant text argument before the scan starts, so
LIMIT saves search work but does not avoid holding that argument in memory.
Find a pattern in each FASTA record, with coordinates, strand and CIGAR:
SELECT name AS reference_name, hit.text_start, hit.text_end, hit.strand, hit.cigar
FROM sassy_search_fasta('test/data/references.fasta', 'ACGTAGG', 0,
alphabet := 'iupac', rc := false)
ORDER BY name, hit.text_start, hit.text_end;| reference_name | text_start | text_end | strand | cigar |
|---|---|---|---|---|
| ambiguous | 2 | 9 | + | 7= |
| forward | 2 | 9 | + | 7= |
| masked | 2 | 9 | + | 7= |
With alphabet := 'iupac', ambiguous bases in the reference (including N) match
any base they could stand for. sassy_panel_search_fasta() takes a list of
patterns instead. Both stream the file sequentially and need no index.
Allow one edit over the full guide and apply the exact PAM endpoint filter. This synthetic sequence illustrates the API; it is not an experimental guide.
SELECT name AS reference_name, hit.text_start, hit.text_end, hit.cost, hit.strand, hit.cigar
FROM sassy_crispr_search_fasta('test/data/references.fasta', 'ACGTNGG', 1)
ORDER BY name, hit.text_start, hit.text_end, hit.strand;| reference_name | text_start | text_end | cost | strand | cigar |
|---|---|---|---|---|---|
| deletion | 2 | 8 | 1 | + | 1=1I5= |
| forward | 2 | 9 | 0 | + | 7= |
| insertion | 2 | 10 | 1 | + | 2=1D5= |
| masked | 2 | 9 | 0 | + | 7= |
| reverse | 2 | 9 | 0 | - | 7= |
Use sassy_crispr_panel_search_fasta() for a guide panel, or
sassy_crispr_matches() / sassy_crispr_matches_many() on sequence values, as in
the join at the top.
SELECT name AS read_name, hit.text_start, hit.text_end, hit.cost, hit.strand, hit.cigar
FROM sassy_search_fastq('test/data/reads.fastq', 'ACGA', 0, rc := false)
ORDER BY name, hit.text_start, hit.text_end;| read_name | text_start | text_end | cost | strand | cigar |
|---|---|---|---|---|---|
| r1 | 2 | 6 | 0 | + | 4= |
| r3 | 0 | 4 | 0 | + | 4= |
Count hits for a whole panel at once. Hits carry the zero-based pattern_idx,
so duplicate patterns remain distinct:
SELECT hit.pattern_idx, count(*) AS observations
FROM sassy_panel_search_fastq('test/data/reads.fastq', ['ACGA', 'TTGC'], 0, rc := false)
GROUP BY hit.pattern_idx
ORDER BY hit.pattern_idx;| pattern_idx | observations |
|---|---|
| 0 | 2 |
Anything with a sequence column works, including views:
CREATE TEMP VIEW reads AS
SELECT * FROM read_fastq('test/data/reads.fastq', scan_mode := 'sequential');SELECT name AS read_name, hit.text_start, hit.text_end, hit.cigar
FROM sassy_search_table('reads', 'ACGA', 0, rc := false)
ORDER BY name, hit.text_start, hit.text_end;| read_name | text_start | text_end | cigar |
|---|---|---|---|
| r1 | 2 | 6 | 4= |
| r3 | 0 | 4 | 4= |
For other column names or computed sequences, use the scalar functions directly:
SELECT r.name AS read_name, hit.text_start, hit.text_end, hit.cigar
FROM reads AS r
CROSS JOIN LATERAL unnest(sassy_matches('ACGA', r.sequence, 0,
rc := false)) AS matches(hit)
ORDER BY read_name, hit.text_start, hit.text_end;| read_name | text_start | text_end | cigar |
|---|---|---|---|
| r1 | 2 | 6 | 4= |
| r3 | 0 | 4 | 4= |
SELECT sassy_contains('ACGA', 'TTACGATT', 0, rc := false) AS contains;| contains |
|---|
| true |
SELECT unnest(sassy_matches('ACGA', 'TTACGATT', 0, rc := false), recursive := true);| pattern_idx | text_start | text_end | pattern_start | pattern_end | cost | strand | cigar | cigar_ops |
|---|---|---|---|---|---|---|---|---|
| 0 | 2 | 6 | 0 | 4 | 0 | + | 4= | NULL |
| Input | One pattern | A panel of patterns |
|---|---|---|
| FASTA file | sassy_search_fasta |
sassy_panel_search_fasta |
| FASTQ file | sassy_search_fastq |
sassy_panel_search_fastq |
| Table / view | sassy_search_table |
— |
| Values | sassy_matches, sassy_count, sassy_contains |
sassy_matches_many, sassy_count_many, sassy_contains_many |
| CRISPR, FASTA | sassy_crispr_search_fasta |
sassy_crispr_panel_search_fasta |
| CRISPR, table | sassy_crispr_search_table |
sassy_crispr_panel_search_table |
| CRISPR, values | sassy_crispr_matches |
sassy_crispr_matches_many |
Value functions accept VARCHAR or BLOB (panels: VARCHAR[] or BLOB[]);
pattern and text must use the same type. BLOBs hold plain sequence bytes, so
decode DuckHTS nt16/nt4 packed symbols first. sassy_count and sassy_contains
skip building CIGAR strings.
Set cigar_format := 'packed' on the general or CRISPR match functions for
typed cigar_ops UINTEGER[] without formatting text CIGAR (cigar is NULL).
Set cigar_format := 'both' for both representations. The default
'text' returns text CIGAR and NULL cigar_ops. All formats share a fixed
result struct with both nullable fields. Packed ops include soft clips for
any partial pattern coverage returned by the search.
SELECT hit.text_start, hit.cigar_ops
FROM UNNEST(sassy_matches('ACGTTGCA', 'GGTGCAAACGTCC', 1,
alphabet := 'dna', cigar_format := 'packed')) AS matches(hit)
WHERE hit.strand = '-';When a DuckHTS build provides cigar_aligned_blocks(cigar, pos), the typed
result can be passed directly (this snippet is not an executed example):
SELECT cigar_aligned_blocks(hit.cigar_ops, hit.text_start)
FROM UNNEST(sassy_matches('ACGTTGCA', 'GGTGCAAACGTCC', 1,
alphabet := 'dna', cigar_format := 'packed')) AS matches(hit)
WHERE hit.strand = '-';Matching semantics: coordinates, strands, alphabets, defaults
Matching uses Sassy 0.2.6. Coordinates are zero-based, half-open
in the input text. + and - identify the pattern’s strand; reverse-strand
CIGAR follows the pattern direction, not SAM direction. Packed cigar_ops
follow SAM’s forward-reference order: query = pattern, reference = text,
POS = text_start (zero-based). A reverse hit reverses the structured Sassy
operation list, retaining I (pattern-only) and D (text-only). The codes are
BAM’s I=1, D=2, S=4, =7, X=8; Sassy emits explicit =/X, not M.
Unaligned query ends are soft-clipped, with clips swapped on reverse hits;
pattern_start/pattern_end describe the aligned segment in the original
pattern. Excluding soft clips, the packed query span is
pattern_end - pattern_start; its reference span is text_end - text_start.
Each run is limited to the BAM 28-bit length field; larger runs fail with a
limit error. Text CIGAR retains Sassy’s pattern-oriented notation without
clips, including on reverse hits.
General searches default to alphabet := 'iupac', rc := true and
all_endpoints := false (Sassy’s rightmost local minima). Set
all_endpoints := true to report every qualifying endpoint; this still does not
enumerate every alignment. For literal ASCII searches, use
alphabet := 'ascii', rc := false.
DNA accepts A/C/G/T; IUPAC also accepts R/Y/S/W/K/M/B/D/H/V/N. Both match case-insensitively. FASTA wrappers default to DNA and preserve the reader’s sequence representation: DuckHTS 1.5.2 uppercases FASTA sequences and may return NULL descriptions.
Patterns contain 1–4096 bytes, panels at most 4096 patterns, and k must be
smaller than every pattern length. NULL arguments produce NULL; empty text or an
empty panel produces no hits. Empty patterns and NULL panel entries are errors.
Use ORDER BY when result order matters. Reference-region offsets must be added
by the caller using the sequence’s coordinate map.
CRISPR semantics: what a guide hit means
- IUPAC matching on both strands, with
pam_length := 3by default. - Unit edit distance
kover the full guide, including PAM. Insertions and deletions count as edits and appear in the CIGAR. allow_pam_edits := false: Sassy’s exact IUPAC PAM endpoint filter, not a separately constrained PAM alignment or Cas-specific score.max_n_frac := 0.2: N/n content over the whole target match, including PAM, compared using Sassy’s float32 rule.- All qualifying endpoints are reported. Guides in one panel must share identical PAM suffix bytes; separate guide rows can use different PAMs and PAM lengths.
Strand identifies a guide occurrence, not an assigned biological binding strand. Independent DNA/RNA-bulge limits and positional penalties are not provided.
Large inputs and SQL LIMIT: current execution model
DuckHTS streams records, but Sassy receives each whole sequence or text value
and can allocate RAM for search state and matches. DuckDB’s memory_limit does
not fully count that RAM or move it to disk. A long input is still one large
search, and parallel searches can use more RAM at once. Input length and hit
count have no user-set caps; peak memory is unbounded.
An outer SQL LIMIT limits rows returned by a query. The scalar functions first
compute all matches for one text value as a LIST, then DuckDB expands that LIST
into rows. LIMIT may avoid work on later input rows, but it cannot stop the
native search or bound its memory within one value. sassy_grep is a separate
ASCII table scan that searches one bounded region at a time; its output batches
can stop after LIMIT. Table function arguments are constant expressions, so
it cannot take a column from read_fasta or a general relation.
The scalar adapter borrows DuckDB’s input string bytes for each synchronous
search call. DuckHTS’s read_fasta currently copies each parsed sequence into a
DuckDB VARCHAR vector before this adapter sees it. Avoiding that materialization
would require a reader/search scan that passes a record buffer directly to the
matcher while keeping it alive through the search.
Each worker’s searcher reuses one result allocation, retaining hit and CIGAR
buffer capacity between calls until the searcher is destroyed. Successful calls
do not allocate an error string. C callers can return consumed results with
sassy_c_result_recycle() or release them with sassy_c_result_free().
The DuckDB adapter caps LIST output at 1,048,576 hits across a
chunk and raises an error if it cannot materialize them. sassy_count and
sassy_contains avoid LIST output and CIGAR strings, but currently still run the
full native match. Splitting a long reference into regions can reduce memory per
search, but the caller must restore reference coordinates and reconcile hits in
overlapping regions. sassy_grep does this for its ASCII all-endpoints scan.
SIMD backends: inspecting and forcing a kernel
The native library selects the best available backend automatically:
SELECT * FROM sassy_backend_info() ORDER BY name;| name | compiled | supported | selected |
|---|---|---|---|
| avx2 | true | true | true |
| avx512 | true | false | false |
| neon | false | false | false |
| scalar | true | true | false |
| wasm128 | false | false | false |
Set SASSY_C_BACKEND to auto, scalar, avx2, avx512, neon or wasm128
before the first search to force one. The choice is fixed for the loaded
library; an unavailable request returns an error, and changing it requires a
fresh process.
scalar follows Sassy’s baseline and includes SSE2 on x86-64; SSE4.1 has no
separate tier. Linux x86-64 baseline and AVX2 are tested; AVX-512 is compiled
but unexecuted, and aarch64/NEON is unverified here. Windows MinGW selects AVX2
under Wine with DuckDB R 1.5.5. Emscripten MVP/EH builds use only scalar;
wasm128 is reserved for the optional threaded build. Browser tests gate the
scalar catalog, examples, CRISPR packed CIGAR, invalid-input recovery and
worker restart after an internal Rust panic; see the platform checks.
Using the C library directly
The build also produces build/libsassy_c.a.
include/sassy_c.h defines its ABI, ownership rules and
limits; test/c/test_abi.c is a working example.
make test-v2 sql-test oracle-test readmeoracle-test compares complete CRISPR hit multisets against the Sassy
0.2.6 CLI across 36 profiles and 864 guide/record pairs.
make r-test runs DBI integration against the R preview host pinned in
ducksassy-package.json, installed separately under .deps/Rlib/.
make benchmarks re-renders the benchmark report.
ducksassy is licensed under the GNU GPL, version 2 or later. Sassy is by Rick Beeloo and Ragnar Groot Koerkamp (MIT). File reading comes from DuckHTS.