rsomics-seq-grep

Filter FASTA/FASTQ records by ID, full header, or sequence content — a Rust port of seqkit grep.

Install

cargo install rsomics-seq-grep

Usage

# Keep records by exact ID (default target is the ID — the header up to the
# first space/tab, not the whole header line)
rsomics-seq-grep -p seq1,seq2 in.fa

# Match against the full header line instead of just the ID
rsomics-seq-grep -n -p "seq1 description text" in.fa

# 10k-ID allow-list from a file (a pattern file always overrides -p entirely)
rsomics-seq-grep -f ids.txt in.fq.gz

# Regex against ID/name (unanchored — partial match, unlike the default
# exact-match mode)
rsomics-seq-grep -r -p 'sample_[0-9]+$' in.fa

# Search the sequence itself; both strands are searched by default
rsomics-seq-grep -s -p GAATTC in.fa

# IUPAC degenerate motif (implies -s)
rsomics-seq-grep -s -d -p RYSW in.fa

# Restrict sequence search to a 1-based region (implies -s)
rsomics-seq-grep -s -R 1:12 -p ATG in.fa

# Invert: drop matching records, keep the rest
rsomics-seq-grep -f ids.txt --invert-match in.fa

# Count matches without writing records
rsomics-seq-grep -C -f ids.txt in.fa

# stdin
cat in.fa | rsomics-seq-grep -p seq1 -

Semantics (verified against seqkit v2.13.0, both by reading

shenwei356/seqkit/shenwei356/bio source and by black-box comparison)

• Default target is the ID, defined the same way as seqkit's DefaultIDRegexp (^(\S+)\s?): the header up to the first space or tab, whichever comes first. -n/--by-name switches the target to the entire header line (ID + description) instead.
• Default matching is whole-string equality, not substring — -p seq1 does not match a record named seq1_dup. This is seqkit's own documented deviation from POSIX/GNU grep. -r/--use-regexp switches to an unanchored (partial-match) regex — the regex source is used exactly as given, not escaped.
• -s/--by-seq switches the target to the sequence itself and switches matching to substring search (not equality). Both the plus and minus (reverse-complement) strand are searched by default — -P/--only-positive-strand restricts to the plus strand only. Strand search short-circuits: the first strand that hits wins, matching seqkit's for strand := range strands { if hit { break } }.
• -d/--degenerate expands an IUPAC nucleotide motif (ACGTURYMKSWHBVDN, case-preserved) into the exact regex character classes from seqkit's bio/seq/seq.go DegenerateBaseMapNucl table, then matches as an unanchored regex. Implies -s. Protein degenerate codes (DegenerateBaseMapProt) are not implemented — -d here is nucleotide-only, the overwhelmingly common use (restriction sites, primer motifs, degenerate probes).
• -R/--region <start:end> restricts -s search to a 1-based window, applied independently per strand (so -R 1:12 on the minus strand means the first 12 bases of the reverse complement). An out-of-range region skips that strand rather than erroring. Negative indices count from the end (-12:-1 = last 12 bases); mixed-sign spans (start<0, end>0) are rejected as ambiguous, matching seqkit. Implies -s.
• An empty pattern (-p "") is legal and specific: in ID/name mode it matches only records with an empty ID/name; in -s mode (bytes.Contains semantics upstream) it matches only records with an empty sequence rather than "matches everywhere" — seqkit special-cases this rather than relying on the substring-search library's usual empty-needle-always-matches behaviour, and so do we. In regex/degenerate mode an empty pattern compiles to ^$ (exact-empty match) rather than the empty regex (which would match every position).
• -i/--ignore-case lowercases both sides for equality/substring search and prefixes (?i) for regex/degenerate search — same as upstream.
• A pattern file (-f) overrides -p entirely when both are given — this is seqkit's actual behaviour (verified empirically, not just from reading the help text): the two sources are never merged. Pattern-file lines are read verbatim except for a trailing \r/\n strip; blank lines are skipped.
• FASTA output line width defaults to 60 (-w/--line-width, 0 = no wrap), matching seqkit's default. FASTQ output always ignores -w — seqkit grep forces line width to 0 for FASTQ regardless of the flag (fastx.ForcelyOutputFastq / the per-file LineWidth = 0 override in grep.go), and so do we.
• Header lines are reproduced verbatim — grep never reconstructs or reflows the header, only the sequence/quality body is wrapped.
• Input: plain or gzip FASTA/FASTQ, content-sniffed (not by file extension); - reads stdin. Output order always matches input order.

Flags implemented vs. deferred

Implemented to real feature-completeness: -p/--pattern, -f/--pattern-file, -n/--by-name, -r/--use-regexp, -i/--ignore-case, -s/--by-seq, -d/--degenerate, --invert-match, -R/--region, -P/--only-positive-strand, -w/--line-width, -C/--count, -o/--output.

Deferred — not defined on the CLI at all, so clap itself rejects them with "unexpected argument" (exit 2) rather than silently ignoring them:

One further documented divergence: upstream auto-detects the input alphabet and forces -P (positive-strand-only) automatically for protein/unclassified sequences, since reverse-complementing a protein sequence is meaningless. This crate always honours -s's "both strands by default" rule regardless of alphabet — it only matters if -s is used without -P on protein input, an unusual combination for a nucleotide-first tool.

-v is not used for invert-match here (unlike upstream) because -v is already CommonFlags' --verbose short flag across every rsomics-* tool in this campaign; --invert-match is long-only by the established conflict-resolution convention.

Origin

Independent Rust reimplementation of seqkit grep, informed by:

• The seqkit paper: Shen, W. et al. SeqKit2: A Swiss Army Knife for Sequence and Alignment Processing. iMeta (2024) [doi:10.1002/imt2.191].
• Reading the upstream source (seqkit is MIT-licensed, so reading and citing is allowed and standard practice in this project): seqkit/cmd/grep.go (flag wiring, matching precedence, strand loop, region handling), shenwei356/bio seqio/fastx/reader.go (parseHeadIDAndDesc, DefaultIDRegexp), seqio/fastx/records.go (Format/FormatToWriter — the exact FASTA wrap / FASTQ-always-unwrapped output rule), seq/seq.go (DegenerateBaseMapNucl, RevCom/Complement), and seq/alphabet.go (the DNAredundant complement pairing table).
• Black-box behaviour comparison against seqkit v2.13.0 across ID/name exact match, regex, ignore-case, invert-match, by-seq substring, both strands, only-positive-strand, IUPAC degenerate motifs, region restriction, line-width variants, gzip input, stdin, pattern files, and count-only mode (tests/compat.rs, frozen goldens in tests/golden/).

License: MIT OR Apache-2.0. Upstream credit: seqkit (MIT).

External-dep quadrant classification

• needletail — Quadrant ① (pure Rust FASTA/FASTQ parser; content-sniffed gzip via its flate2/zlib-rs backend).
• regex — Quadrant ① (pure Rust, SIMD-accelerated literal prefilter internally).
• memchr — Quadrant ① (pure Rust, SIMD substring search for -s literal matching).
• rayon — Quadrant ① (pure Rust parallelism; drives the by-seq path's per-record matching across threads, honouring -t/--threads via rsomics-common's default-on rayon-pool wiring).
• rsomics-common, rsomics-help, clap, serde — Quadrant ④.

JSON output schema (--json)

{
  "schema_version": "1.0",
  "tool": "rsomics-seq-grep",
  "tool_version": "0.1.0",
  "status": "ok",
  "result": {
    "input": "in.fa",
    "output": "-",
    "pattern_count": 3,
    "total_records": 200000,
    "matched": 51994
  }
}