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
42 changes: 42 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,48 @@
All notable changes to `mirpy-lib` (import `mir`). This project follows semantic versioning; the v3 line is a
greenfield ML/embedding rewrite (the classical v1.x/v2 toolkit is frozen on branch `legacy-v2`).

## 3.4.0 — 2026-07-18

Minor: a command-line interface, a uv-based dev setup, and a documentation overhaul. No public
Python API removed; one optional-dependency group split out.

### Added

- **`mir` command-line interface** (`[project.scripts]`, also `python -m mir.cli`) — the two
embedding scales without writing Python:
- `mir embed clonotypes SAMPLE` → a per-clonotype TCREMP embedding table (`e0…`).
- `mir embed repertoires SAMPLE…` → one repertoire vector `Φ(S)` per sample **per chain** on one
shared basis (`phi0…`), with optional `--mmd` pairwise-distance output.

Inputs are any format `vdjtools.io` reads (AIRR/vdjtools/MiXCR/immunoSEQ/parquet); output is TSV
or Parquet. See `mir embed <cmd> -h`. Tests in `tests/test_cli.py`.
- **`mir.repertoire.correct_batch`** — Harmony-like cluster-aware batch correction on a stacked
sample×feature Φ matrix. Removes the batch offset *per soft cluster* (batch-diversity-penalised),
so a batch confounded with a biological cluster is corrected without erasing that biology; reduces
exactly to `mir.cohort.residualize` at `n_clusters=1` / `theta=0` (`prop:batch`).
- **`[ann]` optional-dependency group** for the approximate-NN density backend (`pynndescent`).

### Changed

- **Development now uses a repo-local `.venv` via [uv](https://docs.astral.sh/uv/)** instead of
conda. `setup.sh` is rewritten (bash/zsh portable; `--dev-parents`, `--docs`, `--tests`); the
conda `environment.yml` is removed. Runtime is unchanged — still a pure-Python `py3-none-any` wheel.
- **`pynndescent` moved from `[bench]` to the new `[ann]` extra.** `[bench]` is now all pure-wheel
(no numba/llvmlite), so `pip install "mirpy-lib[bench]"` resolves cleanly on any Python. Users of
`density.neighbor_enrichment(backend="ann")` should install `"mirpy-lib[ann]"`.
- **`vdjtools>=3.0.0`** (was `>=2.3.0`).
- Documentation overhaul: a use-case-driven user guide, the two CLI commands documented, an
examples/notebooks page, `mir.cohort` and `mir.bench.eval` added to the API reference, a logo, and
the sample-embedding schematic + real depth-robustness figure. Zero-warning Sphinx build.

### Fixed

- `DEFAULT_GAP_POSITIONS = (3, 4, -4, -3)` had three independent definitions
(`distances/junction`, `embedding/tcremp`, `ml/bundle`); now defined once in `distances.junction`
and imported, so the coordinate constant cannot drift.
- `cohort.cluster_samples` docstring described itself; `AntigenMetric` and `mir.bench.eval` gained
the docstrings/module-map entries they were missing.

## 3.3.0 — 2026-07-17

Minor: one new public parameter, nothing removed or changed.
Expand Down
39 changes: 35 additions & 4 deletions CLAUDE.md
Original file line number Diff line number Diff line change
Expand Up @@ -88,6 +88,20 @@ fix were added to `vdjtools` under the owner's direction (this is that owner's e
Opt-in `fit_repertoire_space(n_eigs=r)` swaps the second-moment block's full `D₂(D₂+1)/2` upper-triangle
for its top-`r` eigenvalues (rotation-invariant spectrum; default `None` = upper-tri, unchanged — but
lossy for the *directional* HLA imprint, so the full triangle stays the recommended block; `benchmark_repertoire_spectral.py`).
Also `fit_repertoire_spaces` (one basis **per locus** — the multi-chain fit), `centroid_atypicality`
(per-sample cosine distance to a group centroid — a Φ-geometry op feeding the digital-donor atypicality channel),
and `correct_batch` (Harmony-like cluster-aware batch correction on a stacked Φ matrix; reduces to
`cohort.residualize` at `n_clusters=1`/`theta=0` — `prop:batch`; test in `test_repertoire.py`).
- `cohort.py` — **the digital donor** (T.7): `fit_donor_embeddings`/`DonorCohort` fuse per-chain identity
(kernel-mean, cross-sample PCA-reduced) ‖ diversity ‖ coverage across loci through **one `ChannelBuilder`**,
with an `extra_channels` hook for the analysis's own tissue/clinical blocks; comparability bites twice
(per-locus `prototype_hash` **and** the stored identity PCA) so `save`/`load` verify every hash and
`transform` is the only held-out path. Plus `residualize` (batch cookbook), `cluster_samples`
(MMD→precomputed-metric cluster), `incidence_biomarkers` (subject-incidence Fisher, delegating to
`vdjtools.biomarker.fisher`). Generalizes the analysis-repo `_tcga_embedding.build_embedding` glue.
- `bench/eval.py` — the scorers `channel_report` consumes (kept out of `explain.py` so it stays scorer-free):
`cv_auc` / `held_out_auc` (classification), `cv_cindex` / `km_logrank` (Cox survival, `[bench]`→lifelines),
`kmer_matrix` (baseline).
- `explain.py` — **explainable readouts over any feature matrix** (T7). `ChannelSpec`/`ChannelBuilder`/
`stack_embeddings` attach the name→column map `Φ.vector` does not carry (per-chain blocks merge by
name); `channel_report(X, spec, scorer, base=, mode="in"|"out"|"both", n_permutations=)` ablates each
Expand All @@ -99,14 +113,22 @@ fix were added to `vdjtools` under the owner's direction (this is that owner's e
has no clonotype pre-image and it raises. Depends one-way on `repertoire.py`; nothing there changed.
- `ml/` — Part 2 (torch), neural codecs + `set_encoder.py` (learned repertoire track: Set-Transformer/DeepRC
attention pooling + `SetEncoderBundle`).
- `cli.py` — the `mir` console script (argparse, stdlib): `embed clonotypes` (TCREmp table) /
`embed repertoires` (per-locus `fit_repertoire_space`→`sample_embedding` Φ(S), optional `--mmd`).
Reads via `vdjtools.io.read`; writes TSV/parquet.
- `resources/` — `prototypes/` (TSVs + manifest), `gene_library/` (region_annotations.txt),
`germline_dist/` (baked `.npz`, from `build_germline_dist.py`).

## Build / test / run
- Conda env **`mirpy`** (Python 3.12; do NOT use `.venv` here). `pip install -e .`
(pure-Python hatchling; no C build). Extras: `[bench] [annotate] [build] [ml] [docs] [dev]`.
- Tests: `python -m pytest tests/ -q` (78 pass; `-m "not integration"` for the ~5s fast tier —
the pynndescent ANN parity test carries a one-time JIT cost). All self-contained on bundled resources.
- **Repo-local `.venv` via uv** (Python 3.12; conda retired 2026-07-18). `bash setup.sh`
(bash/zsh; `--dev-parents` editable-installs `../seqtree ../vdjtools ../vdjmatch`, `--docs`,
`--tests`), or `uv pip install -e ".[dev,bench]"`. Pure-Python hatchling; no C build for `mir`
itself. Extras: `[bench] [ann] [annotate] [build] [ml] [docs] [dev] [examples]` (`[ann]` =
pynndescent, split out of `[bench]` so `[bench]` stays numba-free / clean to resolve).
- CLI: `mir embed clonotypes SAMPLE` (→ clonotype embedding table) / `mir embed repertoires
SAMPLE… ` (→ per-sample-per-chain Φ(S)); `mir/cli.py`, `tests/test_cli.py`.
- Tests: `python -m pytest tests/ -q` (`-m "not integration and not benchmark"` for the fast tier;
torch/pynndescent tests skip unless `[ml]`/`[ann]` installed). All self-contained on bundled resources.
- Experiments: `python benchmarks/reproduce_supplementary.py` (theory S1–S3),
`python benchmarks/benchmark_vdjdb.py` (Table S1). Analyses: `analyze_prototype_counts.py`
(geometry saturates by K≈100 — T.1/S4), `analyze_pc_decomposition.py` (V/J η² ≈0.44/0.49,
Expand All @@ -126,6 +148,15 @@ fix were added to `vdjtools` under the owner's direction (this is that owner's e
GPU only in `mir.ml`: `pick_device()` = **CUDA → MPS → CPU** auto, override `device=`/`MIR_DEVICE`.

## Open loops / next steps
- **Embedding-tier roadmap** (`ROADMAP.md`, 2026-07-17) — the "vdjtools at the embedding level" audit +
plan (three verbs: make / measure / generate-decode). **Phase 0** (robustness + optimization quick wins)
and **Phase 1** (cohort tier: `bench/eval.py`, `repertoire.{fit_repertoire_spaces,centroid_atypicality}`,
`cohort.py` digital donor) are **DONE**. Next: **Phase 2** generative loop (`generate.py` `DescriptorDensity`
+ `evolve`/`sample`, `CodecBundle.from_unified/from_decoder`), then multimodal encoders + embedding
trajectory. **Analysis-repo follow-up:** refactor `_tcga_embedding.build_embedding` onto
`cohort.fit_donor_embeddings` (+ `extra_channels` for isotype/composition/atypicality) and re-verify the
pan-cancer ΔC numbers. NB Phase 0 flipped the density default to `backend="kdtree"` — re-verify any
recorded balloon-mode baselines (±1 boundary counts).
- **v3.0 remaining**: 10X paired benchmark; docs (Sphinx theory section + notebooks); CI; publish
`py3-none-any` wheel; regenerate `generate_prototypes.py` via `vdjtools.model.generate`.
- **Bench tuning**: raw kneedle eps over-merges; `cluster(eps_factor=0.4)` recovers the paper
Expand Down
62 changes: 51 additions & 11 deletions README.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,13 @@
<p align="center">
<picture>
<source media="(prefers-color-scheme: dark)" srcset="assets/mirpy_dark.svg">
<source media="(prefers-color-scheme: light)" srcset="assets/mirpy_light.svg">
<!-- Absolute PNG fallback: PyPI strips <picture>/<source> and cannot render a relative or
raw-served SVG, so the logo must be an absolute-URL raster here. GitHub uses the SVG sources. -->
<img alt="mirpy" src="https://raw.githubusercontent.com/antigenomics/mirpy/master/assets/mirpy_light.png" width="360">
</picture>
</p>

<h1 align="center">mirpy — ML embeddings for immune repertoires</h1>

[![PyPI](https://img.shields.io/pypi/v/mirpy-lib.svg)](https://pypi.org/project/mirpy-lib/)
Expand Down Expand Up @@ -50,6 +60,24 @@ labels = cluster(pca_denoise(X, n_components=50))
Paired chains concatenate per-chain embeddings via `PairedTCREmp`. Input/output are AIRR polars
frames keyed by `vdjtools.io.schema` column names.

## Command line

`pip install mirpy-lib` also installs a `mir` command for the two embedding scales — no Python
needed. Inputs are any format `vdjtools.io` reads (AIRR TSV, vdjtools, MiXCR, immunoSEQ, parquet).

```bash
# one repertoire -> per-clonotype embedding table (e0…), the input to clustering / ML
mir embed clonotypes sample.tsv --pca 50 -o clonotypes.parquet

# a dataset of repertoires -> one fingerprint Φ(S) per sample, per chain (phi0…), on one
# shared basis so the rows are mutually comparable; --mmd also writes the pairwise MMD matrix
mir embed repertoires cohort/*.tsv.gz -o phi.tsv --mmd mmd.tsv
```

`mir embed clonotypes -h` / `mir embed repertoires -h` list every flag (species, locus,
prototype count, weight, Φ blocks, …). Sample id defaults to the filename stem; the locus is
inferred per file (or restrict with `--locus`).

## Recommended presets

`TCREmp.from_defaults(species, locus)` uses the per-chain preset when `n_prototypes` is
Expand Down Expand Up @@ -87,9 +115,11 @@ Xr = pca_denoise(X, n_components=p.n_components_recon) # codec reconstruction
| `mir.embedding.tcremp` | `TCREmp` / `PairedTCREmp` — the prototype embedding |
| `mir.embedding.pca` | PCA denoising of embeddings |
| `mir.distances` | junction distance (`seqtree.gapblock`; `metric`/`matrix`/`alignment` options) + baked germline distances |
| `mir.bench` | VDJdb loader, clustering (`cluster(method=…)`: DBSCAN/HDBSCAN/OPTICS) + F1/retention, theory experiments (incl. `codec_losslessness`) |
| `mir.bench` | VDJdb loader, clustering (`cluster(method=…)`: DBSCAN/HDBSCAN/OPTICS) + F1/retention, theory experiments (incl. `codec_losslessness`), cohort scorers (`bench.eval`: `cv_auc`/`cv_cindex`/`km_logrank`) |
| `mir.density` | continuous-density TCRNET/ALICE — enrichment (+ clonal-abundance channel, `backend=` exact/kdtree/ann) + noise-filtering (Theory T6) |
| `mir.repertoire` | sample-level (repertoire) embedding — RFF kernel mean ‖ Hill diversity ‖ second moment; MMD / HLA-stratified distance; motif witness (Theory §T.7) |
| `mir.repertoire` | sample-level (repertoire) embedding — RFF kernel mean ‖ Hill diversity ‖ second moment; MMD / HLA-stratified distance; motif witness; `centroid_atypicality`, multi-locus `fit_repertoire_spaces` (Theory §T.7) |
| `mir.explain` | named-channel fusion (`ChannelBuilder`) + scorer-agnostic ablation (`channel_report`/`channel_drivers`) — which part of Φ carries the signal (§T.7) |
| `mir.cohort` | the **digital donor** — multi-chain `fit_donor_embeddings`/`DonorCohort` (+ `transform`/`save`/`load`) + `residualize` / `cluster_samples` / `incidence_biomarkers` (§T.7) |
| `mir.ml` | neural codecs (forward/inverse/Pgen/unified) + learned repertoire `set_encoder` (Set-Transformer/DeepRC) — Part 2, experimental; `[ml]` extra |

## Background subtraction & clustering (`mir.density`)
Expand All @@ -116,11 +146,12 @@ Use a **biological control** as the background when you have one (e.g. pre- vs p
patient vs healthy) — differential enrichment cancels generic public convergence and isolates the
antigen-specific response. With no control, `generate_background(locus, n)` samples the vdjtools
P_gen model (the ALICE regime); the "water level" of a naive repertoire is handled by the
empirical-null calibration. See `experiments/benchmark_density_{yfv,ankspond,tcrnet}.py`.
empirical-null calibration. The density benchmarks (YFV, ankylosing-spondylitis B27, TCRNET)
live in the companion [`2026-mirpy-analysis`](https://github.com/antigenomics) repo.

At whole-repertoire scale, pass `neighbor_enrichment(..., backend="kdtree")` (exact scipy cKDTree,
5–9× faster than the default BallTree) or `backend="ann"` (approximate pynndescent, ~30× faster
past ~10⁵ clones, trading a small conservative undercount); see `experiments/benchmark_ann.py`.
past ~10⁵ clones, trading a small conservative undercount; `pip install "mirpy-lib[ann]"`).

## Sample-level (repertoire) embedding (`mir.repertoire`)

Expand Down Expand Up @@ -155,16 +186,23 @@ a batch offset is first-order and cancels, while a batch-orthogonal signal (e.g.
empirical rule of thumb — **diversity for how-even, the embedding for which-clones**: clone-size
phenotypes (age, CMV) are a diversity summary's turf, while clonotype identity (HLA — strongest in
TRA and class II) lives in the second moment / witness. A learned co-equal set encoder
(Set-Transformer / DeepRC) is in `mir.ml.set_encoder` (`[ml]` extra). See
`experiments/benchmark_repertoire_*.py`, `BENCHMARKS.md`, and `THEORY.md` T7.
(Set-Transformer / DeepRC) is in `mir.ml.set_encoder` (`[ml]` extra). See [`BENCHMARKS.md`](BENCHMARKS.md)
and [`THEORY.md`](THEORY.md) T7 (the benchmark scripts themselves live in the analysis repo).

## Reproduce the paper

The self-contained theory notebooks run on bundled data:

```bash
python experiments/reproduce_supplementary.py # supplementary S1–S3
python experiments/benchmark_vdjdb.py # Table S1 (needs a VDJdb dump)
pip install "mirpy-lib[examples]"
marimo edit notebooks/theory.py # supplementary S1–S3 (distance laws, D↔d, prototype robustness)
marimo edit notebooks/quickstart.py # embed + cluster VDJdb antigens
```

The full benchmark suite (VDJdb Table S1, density, repertoire/TCGA) and result docs live in the
companion analysis repo [`2026-mirpy-analysis`](https://github.com/antigenomics) — this repo is the
library + CI tests only.

Method: Kremlyakova *et al.*, *TCREMP: a bioinformatic pipeline for efficient embedding of
T-cell receptor sequences*, **J Mol Biol** 437 (2025) 169205.

Expand All @@ -178,12 +216,14 @@ mirpy is CPU-parallel by default and uses the GPU for the neural codecs. Knobs,
| Density kNN / balloon | `neighbor_enrichment(..., backend=…)` | `"exact"` (BallTree, **1 core**) | `backend="kdtree"` = exact scipy cKDTree, **all cores** (`workers=-1`), 5–9× faster; `backend="ann"` = pynndescent, auto all-core, ~30× at ≥1e5. Prefer `kdtree` for multicore exact. |
| Clustering | `cluster(..., n_jobs=-1)` | sklearn default (1) | forwarded to DBSCAN/OPTICS/HDBSCAN via `**kwargs`; parallelizes the neighbour search. |
| BLAS (PCA, RFF, matmul) | `OMP_NUM_THREADS` / `OPENBLAS_NUM_THREADS` env | all cores | numpy/sklearn use the platform BLAS; cap via env if oversubscribed. |
| Neural codecs (`mir.ml`) | `pick_device()` / `device=` / `MIR_DEVICE` env | **CUDA → MPS → CPU**, auto | every `train_*` / codec / bundle takes `device=`; e.g. `MIR_DEVICE=cuda:1 python experiments/train_forward_encoder.py`. Torch-free paths (`density`, `repertoire`) never touch the GPU. |
| Neural codecs (`mir.ml`) | `pick_device()` / `device=` / `MIR_DEVICE` env | **CUDA → MPS → CPU**, auto | every `train_*` / codec / bundle takes `device=`; e.g. `MIR_DEVICE=cuda:1` pins the second GPU. Torch-free paths (`density`, `repertoire`) never touch the GPU. |

Rule of thumb: leave `threads=0` (all cores) for embedding; switch density to `backend="kdtree"`
for exact multicore or `"ann"` at whole-repertoire scale; the GPU is used only by `mir.ml`.

## Development

Conda env `mirpy` (Python 3.12): `bash setup.sh`. Tests: `python -m pytest tests/ -q`.
See [`CLAUDE.md`](CLAUDE.md) for the architecture and reuse map.
Repo-local `.venv` via [uv](https://docs.astral.sh/uv/) (bash/zsh): `bash setup.sh` — add
`--dev-parents` to editable-install the sibling `seqtree` / `vdjtools` / `vdjmatch` checkouts and
`--tests` to run the fast suite. Tests: `python -m pytest tests/ -q`. See [`CLAUDE.md`](CLAUDE.md)
for the architecture and reuse map.
Loading
Loading