bitnuc

A library for efficient nucleotide sequence manipulation using 2-bit encoding.

Summary

This a SIMD-accelerated two-bit encoding library for nucleotide sequences. It is meant to have fast encode/decode routines for small and large sequences and to provide a fairly unstructured interface for working with nucleotide sequences in memory.

It provides:

• 2-bit nucleotide encoding (A=00, C=01, G=10, T=11) with SIMD dispatched at runtime via fearless_simd
• A little-endian-pinned u64 kmer boundary (as_2bit / from_2bit) for hashing and fixed-width integer storage of sequences up to 32 bases
• Ambiguous base detection (ambiguous_bases) for tracking non-ACGTacgt bases
• Extraction routines for extracting aligned and unaligned ranges from packed sequences

Encoded Format

Sequences encode to plain bytes: byte k holds bases 4k..4k+4, base j at bits 2*(j % 4) of its byte. A sequence of n bases occupies exactly n.div_ceil(4) bytes and trailing pad bits are always zero.

Input Sequence: [A][C][G][T]
Output Buffer: [11 10 01 00]

[]: byte boundary
Bit pairs are MSB-first: base j occupies bits 2j..2j+1,
so the first base is the rightmost pair (A=00).

Note: Encoding is lossy for bytes outside ACGTacgt: invalid bases map to an unspecified code rather than an error. If you need to preserve ambiguous bases, detect and track them separately (see ambiguous_bases).

Encoding and Decoding

The core functions operate on byte slices, with _resize variants that manage the buffer length for you:

use bitnuc::{encode_resize, decode_resize};

fn main() -> Result<(), bitnuc::BitnucError> {
    let seq = b"ACGTACGTAC"; // 10 bases -> 3 encoded bytes

    let mut ebuf = Vec::new();
    encode_resize(seq, &mut ebuf);
    assert_eq!(ebuf.len(), 3);

    let mut dbuf = Vec::new();
    decode_resize(&ebuf, seq.len(), &mut dbuf)?;
    assert_eq!(&dbuf, seq);
    Ok(())
}

The slice-based variants write into caller-provided buffers, which lets consumers control allocation and padding (e.g. file formats that pad encoded sequences to 8-byte words):

use bitnuc::{encode, decode};

fn main() -> Result<(), bitnuc::BitnucError> {
    let seq = b"ACGTACGTAC";

    // Pad the encoded buffer to an 8-byte multiple: the layout is identical
    // to the legacy u64 packing serialized little-endian
    let mut ebuf = vec![0u8; seq.len().div_ceil(4).next_multiple_of(8)];
    encode(seq, &mut ebuf)?;

    let mut dbuf = vec![0u8; seq.len()];
    decode(&ebuf, seq.len(), &mut dbuf)?;
    assert_eq!(&dbuf, seq);
    Ok(())
}

u64 Kmer Packing

For hashing and fixed-width storage of short sequences (barcodes, UMIs, k-mers up to 32 bases), as_2bit and from_2bit pack to and from a u64.

Note: These are pinned to little-endian internally.

use bitnuc::{as_2bit, from_2bit};
use std::collections::HashMap;

fn main() -> Result<(), bitnuc::BitnucError> {
    let packed = as_2bit(b"ACGT")?;
    assert_eq!(packed, 0b11100100);

    // Efficient k-mer counting
    let mut kmer_counts = HashMap::new();
    for window in b"ACGTACGT".windows(4) {
        *kmer_counts.entry(as_2bit(window)?).or_insert(0) += 1;
    }
    assert_eq!(kmer_counts.get(&packed), Some(&2));

    // Unpacking returns a stack array of all 32 bases; slice to your length
    let unpacked = from_2bit(packed);
    assert_eq!(&unpacked[..4], b"ACGT");
    Ok(())
}

Packed kmers can be compared with hdist_scalar:

use bitnuc::{as_2bit, hdist_scalar};

fn main() -> Result<(), bitnuc::BitnucError> {
    let u = as_2bit(b"ACGT")?;
    let v = as_2bit(b"ACGA")?;
    assert_eq!(hdist_scalar(u, v, 4)?, 1);
    Ok(())
}

Identifying ambiguous bases

Ambiguous bases (non-ACGT) bases are unable to be represented with this two-bit encoding scheme. The encoding algorithm also remaps lowercase to upper case (so acgt -> ACGT internally). Use ambiguous_bases to track all positions of unrepresentable nucleotides. The position buffer is generic over its element type (usize, u64, or u32):

use bitnuc::ambiguous_bases;

fn main() {
    let seq = b"ACgTNACYAaTH"; // has unrepresentable bases (N/Y/H)

    let mut pos: Vec<usize> = Vec::default();
    ambiguous_bases(seq, &mut pos);

    assert_eq!(
        pos,
        vec![4, 7, 11],
    );
}

Note: ambiguous_bases only tracks non-ACGTacgt bases. It does not identify lowercase letters which are also not representable but which are remapped to their uppercase variants through encoding/decoding.

Extraction

It is oftentimes useful to extract subsequences from packed sequences, e.g. for hashing or extracting ranges. Use extract to pull out ranges of packed bytes into an aligned buffer.

The range is specified in the basepair range (e.g. 0..10 for the first 10 bases) which is converted to a byte range internally. The bits relevant to the range are extracted and packed into the output buffer so that the output buffer is correctly aligned.

use bitnuc::{encode_resize, extract_resize, decode_resize, BitnucError};

fn main() -> Result<(), BitnucError> {
    let seq = b"ACCAAGGTTACATGAAGTTAACCAAGAGAC";

    // encode the sequence
    let mut packed = Vec::new();
    encode_resize(seq, &mut packed);

    // extract some subsequence
    let range = 4..11; 
    let mut extracted = Vec::new();
    extract_resize(&packed, range.clone(), &mut extracted)?;

    // decode the extracted sequence
    let mut decoded = Vec::new();
    decode_resize(&extracted, range.len(), &mut decoded)?;

    // validate expected output
    assert_eq!(decoded, seq[range]);

    Ok(())
}

Memory Usage

The 2-bit encoding provides significant memory savings:

Standard encoding: 1 byte per base
ACGT = 4 bytes = 32 bits

2-bit encoding: 2 bits per base
ACGT = 1 byte = 8 bits

Performance

Throughput by sequence length, measured on an Apple M3 Pro with target-cpu=native (criterion mean, 1 byte per base):

To regenerate the table on your machine:

RUSTFLAGS="-C target-cpu=native" cargo bench --bench simd_comparison -- coding_2bit
uv run scripts/perf_table.py

SIMD Acceleration

The 2-bit encode and decode are SIMD accelerated via fearless_simd, with the instruction set (NEON, SSE, AVX2, AVX-512) selected at runtime.

Related Work

I highly recommend checking out packed-seq. They are currently the highest performance 2-bit encoding library in Rust as far as I can tell. They follow a different bit-packing scheme than this library and they can shave off a few instructions in their SIMD routines.

If you're interested in 2-bit encoding in general make sure to check out cute-nucleotides which has an excellent overview of different algorithms and their performance characteristics.