Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
47 changes: 47 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,53 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0

## [Unreleased]

### Added

- **`leech.crf`: CTC-CRF sequence models — encoder, training objective and
decode.** A second task alongside the chunk classifiers: where a classifier
maps a signal window to a label, a CRF maps one to a *sequence*, over
`n_base ** state_len` states whose Viterbi traceback emits one base per move.
Ported unchanged from `escapepod_models.crf`, which now imports it from here.
The formulation is ONT's, introduced in bonito; the architecture is
SeqTagger's published parameters (Genome Res 35:956). The equivalence checks
that prove it correct stay in escapepod-models, where the reference
implementation is available to run as the oracle. The port is pinned by an A/B against
the pre-move copy: encoder forward, loss forward *and* backward, the analytic
forward-backward's posteriors, and the decoded strings are all bit-identical
on CPU and on GPU (Triton kernels included).

Contents: `CrfEncoder` (3 convs -> 5 alternating LSTMs -> linear -> tanh/scale
-> per-state blank splice), `CtcCrfLoss` with an analytic forward-backward
that replaces a 300-step autograd replay with one elementwise multiply,
optional Triton lattice kernels, the two-pass Viterbi `decode_batch`, and the
packaged architecture config `leech/crf/configs/crf_ctc.toml`.

The subpackage imports **only torch and numpy, and only on demand** — the
config path is eager and the rest is lazy (PEP 562, as in `leech.models`), so
`from leech.crf import DEFAULT_CONFIG` costs no torch import. That is what
lets escapepod-models install leech `--no-deps` into a conda-forge pixi
environment, and it is enforced by `tests/test_crf_package.py`.

Not included yet: the trainer, the ONNX export, and the corpus paths.

- **`leech.crf.manifest`: the seam between a corpus's vocabulary and its signal.**
One table, one row per read — `read_id, pod5, anchor_end, target` plus optional
`label`/`group`/`batch`/`quality_score`/`quality_margin`/`split` — and nothing
about where those facts came from. `target` is the *resolved* sequence: a class
name may ride along in `label` for reporting, but nothing here looks one up.
Everything above the manifest is vocabulary (panels, codes, oligos, gates) and
belongs to whatever project defines it; everything below is signal ML.

Two rules it exists to enforce, both silent when broken. Label quality travels
as **numbers, never a `keep` boolean**, so the gate stays sweepable at training
time — gating one panel's labels moved accuracy from 0.875 to 0.97, and a
boolean decided at extraction would mean re-cutting an 8 GB corpus per
threshold; `quality_coverage()` reports partial scoring, because an unscored
read cannot pass a gate and is dropped without a word (this once cut a corpus
to a non-random 13.5%). And `check_geometry` **raises** on a window too short
to hold its target rather than warning, because a short window trains,
converges, and quietly discriminates on fewer bases than designed.

## [0.7.0] - 2026-08-25

### Fixed
Expand Down
99 changes: 99 additions & 0 deletions CLAUDE.md
Original file line number Diff line number Diff line change
Expand Up @@ -391,6 +391,92 @@ and a training loop contends with its own forward/backward for that pool, while
grid search runs N such processes. `LeechDataset`'s block fill and
`__getitems__` gather go through numpy for this reason.

### CTC-CRF: a second task, and five rules that do not announce themselves

`leech.crf` maps a signal window to a *sequence* rather than a label — a CRF
over `n_base ** state_len` states whose Viterbi traceback emits one base per
move. The formulation is ONT's, introduced in bonito; the architecture is
SeqTagger's published parameters. Ported here from `escapepod_models.crf`
(rnabioco/escapepod-models#40), which now imports it from leech.

**The equivalence checks that prove this code correct live in escapepod-models**
(`scripts/ldx/analysis/verify_crf_*.py`), because they run the reference
implementation as the oracle and leech carries no dependency on it. Keep them
there, and do not "simplify" them by dropping that comparison — it is the
evidence.

**`import leech.crf` must pull only torch and numpy, and only on demand.** Not
pysam, polars, escapepod, sklearn or click — ever. escapepod-models installs
leech `--no-deps` into a conda-forge pixi environment precisely so its solver
never has to reconcile leech's POD5/BAM stack against conda's pytorch, and the
CRF path is all it needs; one convenience import at the top of `encoder.py`
turns that into an install that breaks at first use. The config path is eager
and everything else is lazy (PEP 562, as in `leech.models`), so
`from leech.crf import DEFAULT_CONFIG` costs no torch import: escapepod-models'
`ldxlib` exposes it as a module constant and is imported by two dozen scripts
that want only edit distances and panel lookups. `tests/test_crf_package.py`
fails if either property regresses.

**The model cannot emit the first `state_len` bases of its target.** They fix
the initial state and nothing else, so a `target_len` target decodes to
`target_len - state_len` bases *at any window width* — widening the signal
window recovers exactly nothing. Match decodes against `target[state_len:]`;
matching the full-length target calls the same sequence but inflates every edit
distance and compresses the confidence margin that ranking depends on. Size
targets so the sacrificial bases come from a constant prefix.

**Blank is entry 0 of each state group, and the score width is 1280.**
`score_index = state * (n_base + 1) + label`, `label == 0` meaning stay. 1024 is
the *linear layer's* width (`n_base ** (state_len + 1)`); the blank is spliced in
per state afterwards. This is the layout `escapepod-demux`'s Rust decoder
assumes — move the blank to the end of each group and every shape still lines
up while every call is wrong. Output is also **time-major** `(T, N, n_score)`,
the opposite of the boundary CNN's batch-major `[B, 2, L]` in the same stack.

**The loss runs in fp32, outside autocast, and `_UNREACHABLE` is -1e30, not
-inf.** The lattice scan accumulates over `chunk // stride` timesteps and fp16
loses the tail of that sum; autocast the encoder, where the matmuls are, and
cast back before the loss. And `logaddexp(-inf, -inf)` is `-inf` forward but
differentiates to `nan`, which poisons every upstream gradient — the loss looks
perfect and training silently does nothing. A large finite floor underflows to
zero weight against any real path while keeping the backward pass finite.

**The decode is two passes and both are load-bearing.** Log-semiring
forward/backward for per-timestep edge posteriors, then max-semiring over
`log(post + 1e-8)` for the argmax edge. A one-pass Viterbi over the raw encoder
scores is a different and worse decode, and is the obvious thing to simplify
away. The floor goes on the probability, not the log, so it cannot be folded
into the softmax.

**The manifest is the seam, and vocabulary stays on the far side of it.**
`leech.crf.manifest` takes one table — `read_id, pod5, anchor_end, target` plus
optional `label`/`group`/`batch`/`quality_score`/`quality_margin`/`split` — and
nothing about where those facts came from. Which reads belong to which barcode,
which flowcell, how a label's trustworthiness was scored: that is the producing
project's business. `target` is the **resolved sequence**; a class name may ride
along in `label` for reporting, but nothing here looks one up. escapepod-models'
extractor took a `--panel` argument purely to turn a class name into a target
string, and that one thread is what kept it tied to a single assay.

Two rules the manifest exists to enforce, both silent when broken:

- **Label quality travels as numbers, never as a `keep` boolean.** The gate is
applied at training time so it stays sweepable — gating the ldx labels moved
accuracy from 0.875 to 0.97, and a boolean decided at extraction would mean
re-cutting an 8 GB corpus per threshold. `quality_coverage()` is there because
an *unscored* read cannot pass a gate and is therefore dropped without a word:
a partially scored table once cut a corpus from 56% to 13.5% of its reads,
non-randomly.
- **`anchor_end` and `target` are coupled.** `check_geometry` refuses a window
too short to hold its target rather than warning, because a short window
trains, converges, and quietly discriminates on fewer bases than designed.
Pass a *measured* `samples_per_base` — leech has dwell times — not a constant.

The analytic forward-backward in `_analytic.py` is the loss path; the plain
scans in `loss.py` are the readable reference the tests check it against, and
the Triton kernels check against those. Keep all three — the fallback chain is
what makes a wrong kernel visible.

### Key Classes and Functions

**`MoveTable` (features.py)**
Expand Down Expand Up @@ -499,6 +585,15 @@ src/leech/ # Main package source
│ ├── spec.py # ReleaseSpec: YAML model-release specification
│ ├── notes.py # Release-note rendering
│ └── github.py # gh CLI wrapper
├── crf/ # CTC-CRF sequence models (torch + numpy ONLY)
│ ├── encoder.py # CrfEncoder: signal -> transition scores
│ ├── loss.py # CtcCrfLoss + the readable reference scans
│ ├── _analytic.py # Analytic forward-backward (autograd Functions)
│ ├── _triton.py # Optional CUDA lattice kernels
│ ├── decode.py # Two-pass Viterbi -> sequences
│ ├── config.py # Architecture TOML reader
│ ├── _flags.py # LEECH_* / ESCAPEPOD_* switches
│ └── configs/ # Packaged crf_ctc.toml (the shipped geometry)
├── features.py # MoveTable, dwell times, signal levels, normalization
├── dataset.py # PyTorch Dataset classes, collate_fn, DataLoader sizing
├── training.py # Training loop with Trainer class
Expand Down Expand Up @@ -661,6 +756,10 @@ The codebase is feature-complete (v0.7.0):
training, validation and eval: auto on GPU (capped by the job's CPU
allocation), serial on CPU, never workers inside a daemonic pool worker;
`eval test` takes `--num-workers`
- ✓ CTC-CRF sequence models (`leech.crf`): encoder, training objective with an
analytic forward-backward, optional Triton lattice kernels, and the two-pass
Viterbi decode — ported from escapepod-models, torch + numpy only. The
trainer, ONNX export and corpus paths are not here yet.

All core functionality is implemented and ready for use.

Expand Down
176 changes: 176 additions & 0 deletions docs/api/crf.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,176 @@
# CRF Module

CTC-CRF sequence models: encoder, training objective, and decode.

## Overview

A second task alongside leech's chunk classifiers. Where a classifier maps a
signal window to a label, `leech.crf` maps one to a **sequence**: a CRF over
`n_base ** state_len` states whose Viterbi traceback emits one base per move.
It is what the barcode basecallers in
[escapepod-models](https://github.com/rnabioco/escapepod-models) are trained
with, and what `escapepod-demux`'s Rust decoder runs in production.

```python
import torch
from leech.crf import CrfEncoder, CtcCrfLoss, decode_batch, encoder_config_from_toml, load_config

cfg = encoder_config_from_toml(load_config()) # packaged default geometry
model = CrfEncoder(cfg)
criterion = CtcCrfLoss(cfg.n_base, cfg.state_len)

scores = model(signal) # (N, 1, chunk) -> (T, N, n_score)
loss = criterion(scores.float(), targets, target_lengths)
sequences = decode_batch(scores, cfg.n_base, cfg.state_len)
```

Importing the subpackage costs **nothing but torch and numpy**, and only when a
symbol that needs them is touched — never pysam or escapepod. That is
deliberate: escapepod-models installs leech into a conda-forge environment with
`--no-deps` and needs the CRF path alone. (`load_manifest` pulls polars, but
only when you actually read a manifest.)

## Four things to know before using it

### The model cannot emit the first `state_len` bases of its target

They fix the initial state and nothing else, so a `target_len`-base target
decodes to `target_len - state_len` bases — **at any window width**. Widening
the signal window does not lengthen the decode. Size targets so the sacrificial
bases come from a constant prefix, and match decodes against
`target[state_len:]`, never the full-length target. Matching against the full
target still calls the right sequence, but inflates every edit distance and
compresses the confidence margin that ranking depends on.

### Blank is entry 0 of each state's group

`score_index = state * (n_base + 1) + label`, with `label == 0` meaning stay.
The score width is therefore `n_states * (n_base + 1)` = **1280** for the
default geometry, not the linear layer's 1024. This is the layout
`escapepod-demux`'s Rust decoder assumes; moving the blank to the end of each
group keeps every shape valid and makes every call wrong.

### Output is time-major

`(T, N, n_score)`, not `(N, T, n_score)`. The boundary CNN in the same stack is
batch-major `[B, 2, L]`, so the two contracts sit next to each other and a
consumer that assumes the wrong one silently transposes rather than failing.

### The loss runs in fp32, outside autocast

The lattice scan accumulates over `chunk // stride` timesteps and fp16 loses the
tail of that sum. Autocast the encoder — that is where the matmuls and the speed
are — and cast the scores back before the loss.

## The manifest seam

`leech.crf` cuts a corpus from a **manifest**: one row per read, naming what
leech needs and nothing about where it came from.

| column | required | meaning |
|---|---|---|
| `read_id` | yes | the read, as POD5 and BAM both name it |
| `pod5` | yes | file or directory holding that read's signal |
| `anchor_end` | yes | signal index the window ends at (exclusive) |
| `target` | yes | the **resolved** CRF target sequence |
| `label` | no | class name, for evaluation and reporting |
| `group` | no | reporting/balancing bucket (defaults to `label`) |
| `batch` | no | acquisition batch, for leave-one-batch-out holdout |
| `quality_score` / `quality_margin` | no | label quality, gated at *training* time |
| `split` | no | `train`/`test`, when the producer carved one |

```python
from leech.crf import load_manifest, check_geometry

man = load_manifest("manifest.parquet", require=("batch",))
check_geometry(window=3000, target_len=48, samples_per_base=56.0)
print(len(man), man.batches(), man.quality_coverage())
```

Everything above the manifest is vocabulary — panels, codes, oligos, gates —
and belongs to whatever project defines those. Everything below is signal ML.
There is deliberately **no `keep` boolean**: quality travels as numbers so the
gate stays sweepable without re-cutting the corpus.

## Encoder

::: leech.crf.encoder.CrfEncoder
options:
show_root_heading: true
show_source: true

::: leech.crf.encoder.EncoderConfig
options:
show_root_heading: true

::: leech.crf.encoder.encoder_config_from_toml
options:
show_root_heading: true

::: leech.crf.encoder.load_crf_state_dict
options:
show_root_heading: true

## Loss

::: leech.crf.loss.CtcCrfLoss
options:
show_root_heading: true
show_source: true

::: leech.crf.loss.predecessor_index
options:
show_root_heading: true

## Decode

::: leech.crf.decode.decode_batch
options:
show_root_heading: true
show_source: true

::: leech.crf.decode.best_path
options:
show_root_heading: true

## Manifest

::: leech.crf.manifest.load_manifest
options:
show_root_heading: true

::: leech.crf.manifest.CrfManifest
options:
show_root_heading: true

::: leech.crf.manifest.check_geometry
options:
show_root_heading: true

::: leech.crf.manifest.emitted_target
options:
show_root_heading: true

## Configuration

::: leech.crf.config.load_config
options:
show_root_heading: true

The packaged default is `leech/crf/configs/crf_ctc.toml`. It travels with the
package rather than beside a corpus: an architecture config kept only in scratch
means a purge leaves trained weights nobody can load.

## Acceleration

Two optional fast paths, both gated and both falling back to the PyTorch
reference implementation — which stays the correctness oracle.

| Switch | Effect |
|---|---|
| `LEECH_COMPILE=1` | `torch.compile` the CRF tail and the forward-backward scans. CUDA only; the CPU path stays eager because inductor's CPU `tanh` is not bit-exact. |
| `LEECH_NO_TRITON=1` | Disable the Triton lattice kernels and use the PyTorch scans. |
| `LEECH_NO_COMPILE=1` | Disable compilation of the reference scans in `loss.py`. |

Each also answers to the `ESCAPEPOD_` prefix, which is what escapepod-models'
equivalence checks set.
1 change: 1 addition & 0 deletions docs/api/index.md
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ see the [CLI Reference](../reference/cli.md).
### Models and training

- **[Models](models.md)** -- Neural network architectures (ConvLSTMDwell, TransformerDwell, etc.)
- **[CRF](crf.md)** -- CTC-CRF sequence models: encoder, loss, and Viterbi decode
- **[Training](training.md)** -- Trainer class and training loop
- **[Evaluation](evaluation.md)** -- Model evaluation and metrics
- **[Inference](inference.md)** -- Inference engine and bundle inference
Expand Down
Loading