Prismalign

N-color nucleotide-conversion alignment engine with pluggable backends.

Prismalign maps sequencing reads from any nucleotide-conversion chemistry (bisulfite-seq C→T, SLAM-seq T→C, m6A / A-to-I A→G, MK/KM dual-base, or a custom 3rd channel) using a HISAT-3N-style strategy:

1. build a converted reference index (scheme.ref_from → ref_to)
2. transform each read per color channel and align it to the converted index via a pluggable backend (bwamem by default; a built-in native C k-mer backend; WFA2-lib; minimap2/mappy optional)
3. re-score every hit against the original reference so that real conversions are rewarded (not counted as mismatches), emitting a color-correct MD plus per-channel Y/Z counts in BAM tags.

All per-read kernels are C (bwamem / WFA2-lib / the built-in k-mer aligner); the Python layer is a thin, friendly wrapper.

Install

pip install -e .              # bwamem + built-in C k-mer backends
pip install -e "./[mappy]"    # + minimap2 backend

Usage — Python (clean wrapper)

import prismalign as ps

# one-shot mapping -> BAM (builds indexes, maps, cleans up)
ps.map_reads("reads.fq", "ref.fa", "out.bam",
             scheme="MK", backend="bwamem", threads=4)

# object API / reuse
with ps.NColorMapper(scheme=ps.BS, backend="python") as mapper:
    mapper.map_file("reads.fq", ref_files=["ref.fa"], output_files=["bs.bam"])

Usage — CLI

# classic two-color (MK: A->G + C->T) on bwamem
prismalign map -s MK --backend bwamem -r ref.fa -o out.bam reads.fq

# bisulfite-seq (3-nt single channel C->T)
prismalign map -s BS -r genome.fa -o bs.bam --index-dir idx reads.fq

# parallel (2 copies of the reads, byte-identical output to -t 1)
prismalign map -s MK -r ref.fa -o out.bam -t 4 reads.fq

# list built-in schemes
prismalign schemes

Schemes

Python API

from prismalign import NColorMapper, BS

mapper = NColorMapper(scheme=BS, backend="bwamem", index_dir="idx")
mapper.map_file(r1_file="reads.fq", ref_files=["genome.fa"],
                output_files=["out.bam"])

Backends

Prismalign's engine only needs align() -> [RawHit] from a backend (re-scoring against the original reference is engine-side), so adding one is easy:

Direct vs CLI backends. bwamem, minibwa, mappy, wfa2 and python are direct/in-process (native bindings / compiled C). The only CLI (subprocess) backends are sam (generic) and strobealign (ultra-fast short reads — no Rust→Python binding). wfa2 reuses PythonBackend's k-mer seeding to anchor a diagonal and runs WFA2's exact gap-affine alignment for true I/D CIGARs — the same "one core algorithm" as wfmash/gem3, minus the CLI layer. On exact / simple-mismatch reads every backend's output is identical; on gapped reads WFA2 may pick a different-but-equally-valid split of the M-runs around an indel than BWA (same position and I/D set), so byte-identity applies to the mapping, not to the exact CIGAR representation.

Full inventory — including where each Python wrapper lives — is in docs/backends.md.

Speed & IO

• Parallel mapping: -t/--threads N maps reads in an ordered fork+COW process pool (any backend); batches are drained in read order so the BAM is byte-identical to threads=1. --batch-size tunes reads per worker.
• Reduced repeated IO: references are copy+converted once even when reused across layers (cache keyed by path+scheme); per-hit reference fetch is cached in memory for small contigs (RNA/transcript references), so only one indexed read per contig.

Limitations (v0.0.1)

• paired-end needs the bwamem backend (python/mappy backends are SE-only)