alphabet_tiling

Annotate an arbitrary FASTA with a user-supplied alphabet of monomer consensi.

Given an alphabet (a FASTA of monomer consensus sequences — e.g. satellite / higher-order-repeat "letters") and an input FASTA, alphabet_tiling greedily tiles each input sequence into monomers, labels each with the nearest alphabet letter, and reports the per-monomer alignment (CIGAR), identity and strand. Sequence that does not match the alphabet is reported separately as non-alpha — the tool does not force spurious annotations (no false positives; see below).

It is alphabet-agnostic: pass any set of monomer consensi (α-satellite, HSat, a custom HOR alphabet, …).

Install

git clone https://github.com/aglabx/alphabet_tiling
cd alphabet_tiling
cargo build --release
# binary: ./target/release/alphabet_tiling

Usage

alphabet_tiling --alphabet <ALPHABET.fa> --input <INPUT.fa> [--out PREFIX] [OPTIONS]

  -a, --alphabet <FILE>     FASTA of monomer consensi (>letter\nSEQUENCE). Any alphabet.
  -i, --input <FILE>        FASTA to annotate.
  -o, --out <PREFIX>        write <PREFIX>.monomers.tsv + <PREFIX>.nonalpha.bed
                            (default: monomers TSV to stdout)
      --min-identity <F>    minimum identity (0..1) to call a monomer; below => non-alpha [default: 0.62]
      --top-k <N>           candidate letters per position from the k-mer prefilter [default: 8]
      --threads <N>         threads (0 = all cores) [default: 0]
      --min-len <N>         stop tiling within this many bp of a record's end [default: 100]
      --names <FILE>        TSV (letter<TAB>alias) -> adds a `name` column: the nearest *named*
                            class for each monomer (e.g. a standard nomenclature crosswalk)

Output

<PREFIX>.monomers.tsv — one row per monomer:

column	meaning
seq_id	input record id
start, end	0-based half-open monomer span in the record
strand	+ / - (which strand of the alphabet letter the monomer matches)
letter	best-matching alphabet letter
identity	alignment identity to that consensus (0..1)
cigar	per-monomer alignment to the consensus (= match, X mismatch, I insertion, D deletion)

<PREFIX>.nonalpha.bed — stretches with no alphabet match (seq_id start end nonalpha).

Example

alphabet_tiling -a example/alphabet.fa -i example/demo.fa -o out

example/demo.fa is a synthetic array built as 8× the higher-order pattern A B A C of the three monomers in example/alphabet.fa, followed by 600 bp of random (non-alpha) sequence. The tool recovers it:

seq_id      start  end  strand  letter  identity  cigar
demo_array  0      120  +       monA    0.958     24=2X44=1X43=1X1=1X3=
demo_array  120    238  +       monB    0.992     93=1X24=
demo_array  238    358  +       monA    0.950     38=1X22=1X15=1X7=1X2=1X11=1X19=
demo_array  358    480  +       monC    0.984     15=1X56=1X49=
...

→ 32 monomers (monA ×16, monB ×8, monC ×8 — the A-B-A-C HOR), and the random tail goes to out.nonalpha.bed (demo_array 3840 4440 nonalpha). Expected output: example/expected_output.monomers.tsv.

No false positives

The identity threshold cleanly separates real monomers from everything else. On 868-letter α-satellite testing the tool produced 0 monomers (100 % non-alpha) on:

• random DNA,
• euchromatin (non-centromeric genomic sequence),
• an unrelated tandem satellite (a real satellite that is not in the alphabet),

while recovering a real array at mean identity ≈ 0.99. A non-matching satellite is not mislabelled as a member of the alphabet.

Bundled alphabets

• alphabets/alpha_satellite_868.fa — an 868-letter pan-primate α-satellite monomer alphabet (built from human HPRC + great-ape T2T + gibbon assemblies). Use it to annotate human/primate α-satellite directly:

  alphabet_tiling -a alphabets/alpha_satellite_868.fa -i your_centromere.fa -o out

  On non-α sequence (random DNA, euchromatin, an unrelated satellite) this alphabet yields 0 monomers (100 % non-alpha) — see No false positives above.

• alphabets/alpha_satellite_868.alexandrov.tsv — crosswalk from each of the 868 letters to the nearest Alexandrov / HumAS SF1 class (581 base monomers W1..W5/R1/R2/J1/J2/D1/D2/M1, and SF1-HOR monomers S1C…H…). Pass it with --names to add the standard nomenclature to every monomer:

  alphabet_tiling -a alphabets/alpha_satellite_868.fa \
                  --names alphabets/alpha_satellite_868.alexandrov.tsv \
                  -i your_centromere.fa -o out
  # out.monomers.tsv gains a `name` column, e.g. S1C1/5/19H1L.5

Method

Greedy left-to-right tiling. At each position a 13-mer prefilter (indexed over both strands of every consensus) selects the top-K candidate letters; each candidate consensus is banded-aligned (both orientations) anchored at the position; the highest-identity hit is emitted as a monomer of the consumed length and we advance. Identity below --min-identity ⇒ non-alpha (advance and log). Because matching uses the whole consensus (not just an anchor motif), it is robust on divergent repeats where anchor-based cutters fail.

License

MIT.