Metadata-Version: 2.5
Name: vdjmatch
Version: 1.1.0
Summary: Fast, control-calibrated annotation of T-cell receptor antigen specificity
Project-URL: Homepage, https://github.com/antigenomics/vdjmatch
Project-URL: Repository, https://github.com/antigenomics/vdjmatch
Project-URL: Documentation, https://antigenomics.github.io/vdjmatch/
Author-email: ISALGO laboratory <mikhail.shugay@gmail.com>
License: ISALGO LAB Academic and Non-Commercial Research License
        Version 1.0
        
        Copyright (c) 2026 ISALGO LAB. All rights reserved.
        
        0. DEFINITIONS
        
           "Software" means the software, source code, scripts, models, trained
           parameters, and associated documentation files distributed under this
           license, together with any Modification thereof.
        
           "Modification" means any change to, adaptation of, translation of, or work
           derived from the Software, in whole or in part.
        
           "Output" means any data produced by running the Software, including
           predictions, scores, embeddings, annotations, labels, rankings, alignments,
           simulated or synthetic data, and any file derived from them.
        
           "Model" means any machine-learning or statistical system whose parameters
           are fitted from data, including neural networks, gradient-boosted trees,
           regression models, and any set of trained parameters, whether or not
           distributed together with code.
        
           "You" means the individual or legal entity exercising the rights granted
           below.
        
           "Non-Commercial Use" means use for academic, educational, or scientific
           research purposes from which no commercial advantage or monetary
           compensation is sought or obtained, whether directly or indirectly.
        
           "Commercial Use" means any use that is not Non-Commercial Use, including
           without limitation: (a) use in or to develop a product, service, or
           pipeline that is sold, licensed, rented, or otherwise monetized;
           (b) use by a for-profit entity in connection with its business operations,
           including internal evaluation, benchmarking, or research and development;
           (c) use to provide services to third parties for a fee, including
           consulting, contract research, or hosted or cloud-based provision of the
           Software's functionality; (d) use in research performed under contract for,
           or funded directly by, a for-profit entity for that entity's proprietary
           benefit; and (e) any use described in Section 3.
        
        1. GRANT OF RIGHTS
        
           Permission is hereby granted, free of charge, to any person obtaining a
           copy of the Software, to use, copy, modify, merge, publish, and distribute
           the Software for Non-Commercial Use only, subject to all conditions below.
           No rights are granted other than those expressly stated here; all rights
           not expressly granted are reserved by the copyright holder.
        
        2. PROHIBITED USE
        
           Commercial Use of the Software is prohibited without a separate written
           license signed by the copyright holder.
        
        3. MACHINE-LEARNING AND MODEL TRAINING
        
           3.1 The following are permitted under this license only when the resulting
               Model and its Outputs are themselves put to Non-Commercial Use:
        
               (a) training, fine-tuning, pre-training, continued training, or
                   retraining any Model on the Software, on its source code, or on
                   its Output;
               (b) using the Output as training data, labels, targets, weak
                   supervision, reward or preference signal, validation data, or test
                   data for any Model;
               (c) using the Software or its Output to filter, clean, augment, score,
                   rank, deduplicate, or otherwise curate a data set used to train
                   any Model;
               (d) distilling, pruning, quantizing, compressing, or otherwise
                   deriving a second Model from any Model distributed as part of the
                   Software, or from that Model's Output;
               (e) using the Software or its Output to benchmark, evaluate, validate,
                   tune hyperparameters of, or select among Models;
               (f) embedding the Software, or a Model trained in any of the ways
                   above, in an automated or agent-driven pipeline.
        
           3.2 Each activity listed in Section 3.1 is Commercial Use, and therefore
               prohibited without a separate written license, whenever the resulting
               Model, its Outputs, its parameters, or any system containing it is:
               sold, licensed, rented, hosted, or otherwise monetized; used by or on
               behalf of a for-profit entity in connection with its business
               operations; used to provide services to third parties for a fee; or
               used to support an application for, prosecution of, or exploitation of
               a patent or other proprietary right held by a for-profit entity.
        
           3.3 A Model is subject to this Section regardless of how much of the
               Software or its Output contributed to that Model, and regardless of
               whether the Software or its Output can be recovered from the trained
               parameters. Intermediate steps do not remove a Model from the scope of
               this Section: a Model trained on a data set that was itself produced,
               curated, or filtered using the Software remains subject to it.
        
           3.4 Where a Model is trained partly on the Software or its Output and
               partly on other material, Section 3.2 applies to that Model in full.
        
           3.5 Publication of research results, including trained parameters released
               for Non-Commercial Use under terms no less restrictive than this
               license, is permitted and encouraged.
        
        4. REDISTRIBUTION
        
           You may distribute the Software, with or without Modification, provided
           that You:
        
           (a) retain this license in full, including the copyright notice of
               Section 0 and the disclaimer of Section 9, in all copies and in all
               Modifications;
           (b) do not sublicense the Software or impose any additional or different
               terms, so that every recipient receives exactly the rights and
               restrictions granted here and no more;
           (c) mark each Modification with a prominent notice stating that the file
               has been changed, by whom, and on what date;
           (d) charge no fee for the Software itself, beyond the actual cost of
               physical media or transmission.
        
        5. ATTRIBUTION
        
           Any publication, presentation, preprint, poster, or derived work that makes
           use of the Software must cite this repository, or the corresponding
           publication where one is specified in the repository. Any released Model
           covered by Section 3 must state in its accompanying documentation that the
           Software was used in its construction.
        
        6. NO TRADEMARK OR PATENT RIGHTS
        
           This license grants no right to use the names, trademarks, or logos of
           ISALGO LAB, the copyright holder, or any contributor, whether to endorse or
           promote any work or otherwise, without prior written permission. This
           license grants no patent rights, express or implied, under any patent owned
           or controlled by the copyright holder.
        
        7. TERM AND TERMINATION
        
           The rights granted here apply for as long as You comply with this license
           and terminate automatically, without notice, upon any breach. Rights may be
           reinstated in writing by the copyright holder. Sections 2, 3, 6, 8, and 9
           survive termination. Upon termination You must cease all use of the
           Software and destroy all copies in Your possession, except copies retained
           solely for archival purposes required by law or institutional policy.
        
        8. COMMERCIAL AND TRIAL LICENSING
        
           Parties seeking Commercial Use, use of the Software or its Output in
           training a Model for Commercial Use, an evaluation or trial license, or any
           rights beyond those granted here should contact ISALGO LAB. Such licenses are
           granted only in writing.
        
        9. DISCLAIMER OF WARRANTY AND LIMITATION OF LIABILITY
        
           THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
           IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
           FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
           THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
           LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
           FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
           DEALINGS IN THE SOFTWARE.
        
           The Software is a research tool. It is not validated for clinical,
           diagnostic, therapeutic, or any other use on which health, safety, or
           financial decisions depend.
        
        10. SEVERABILITY
        
           If any provision of this license is held unenforceable, that provision
           shall be modified to the minimum extent necessary to make it enforceable,
           and the remaining provisions shall remain in full effect.
License-File: LICENSE
Keywords: AIRR,E-value,TCR,VDJdb,antigen-specificity,immunology,repertoire
Classifier: Development Status :: 3 - Alpha
Classifier: Intended Audience :: Science/Research
Classifier: Programming Language :: Python :: 3
Classifier: Topic :: Scientific/Engineering :: Bio-Informatics
Requires-Python: >=3.10
Requires-Dist: polars<2,>=1.41.2
Requires-Dist: scipy>=1.9
Requires-Dist: seqtree>=1.3.2
Requires-Dist: threadpoolctl
Requires-Dist: vdjtools>=5.6.0
Provides-Extra: bench
Requires-Dist: huggingface-hub; extra == 'bench'
Requires-Dist: matplotlib; extra == 'bench'
Requires-Dist: psutil; extra == 'bench'
Requires-Dist: rapidfuzz; extra == 'bench'
Requires-Dist: scikit-learn; extra == 'bench'
Requires-Dist: seaborn; extra == 'bench'
Provides-Extra: cli
Requires-Dist: tqdm; extra == 'cli'
Provides-Extra: control
Requires-Dist: huggingface-hub; extra == 'control'
Requires-Dist: seqtree[control]; extra == 'control'
Provides-Extra: diagnostics
Requires-Dist: leidenalg<0.13,>=0.11; extra == 'diagnostics'
Provides-Extra: ndn
Requires-Dist: arda-mapper>=2.38.2; extra == 'ndn'
Provides-Extra: notebook
Requires-Dist: huggingface-hub; extra == 'notebook'
Requires-Dist: marimo; extra == 'notebook'
Requires-Dist: matplotlib; extra == 'notebook'
Provides-Extra: ods
Requires-Dist: ezodf; extra == 'ods'
Provides-Extra: precursor
Requires-Dist: numpy; extra == 'precursor'
Provides-Extra: test
Requires-Dist: biopython>=1.84; extra == 'test'
Requires-Dist: hypothesis; extra == 'test'
Requires-Dist: pytest>=7.4; extra == 'test'
Provides-Extra: vdjam
Requires-Dist: biopython>=1.84; extra == 'vdjam'
Description-Content-Type: text/markdown

<p align="center">
  <picture>
    <source media="(prefers-color-scheme: dark)" srcset="assets/vdjmatch_dark.svg">
    <source media="(prefers-color-scheme: light)" srcset="assets/vdjmatch_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.
         Use the light (dark-on-transparent) raster so it reads on PyPI's white page. -->
    <img alt="vdjmatch" src="https://raw.githubusercontent.com/antigenomics/vdjmatch/master/assets/vdjmatch_light.png" width="340">
  </picture>
</p>

<h1 align="center">vdjmatch — TCR peptide–MHC annotation</h1>

<p align="center">
  <a href="https://pypi.org/project/vdjmatch/"><img alt="PyPI" src="https://img.shields.io/pypi/v/vdjmatch"></a>
  <a href="https://github.com/antigenomics/vdjmatch/actions/workflows/tests.yml"><img alt="tests" src="https://github.com/antigenomics/vdjmatch/actions/workflows/tests.yml/badge.svg"></a>
  <a href="https://docs.isalgo.dev/vdjmatch/"><img alt="docs" src="https://github.com/antigenomics/vdjmatch/actions/workflows/docs.yml/badge.svg"></a>
  <img alt="python" src="https://img.shields.io/badge/python-3.10%2B-blue">
  <a href="LICENSE"><img alt="license" src="https://img.shields.io/badge/license-Academic%20%26%20Non--Commercial-green"></a>
</p>

`vdjmatch` annotates selected T-cell receptors and repertoires against
[VDJdb](https://github.com/antigenomics/vdjdb-db), preserving query identities and
reporting peptide–MHC candidates, competing evidence and optional background calibration.
Native sequence search uses [`seqtree`](https://github.com/antigenomics/seqtree);
repertoire ingestion reuses published `vdjtools` readers.
Requires released `seqtree>=1.3.2` and `vdjtools>=5.6.0`, including allocation-aware
native thread validation.

<p align="center">
  <img alt="VDJmatch workflow: shared receptor-table validation and curated VDJdb search, followed by receptor annotation or matched-reference repertoire enrichment" src="https://raw.githubusercontent.com/antigenomics/vdjmatch/master/assets/algorithm-overview.png" width="100%">
</p>

Receptor annotation asks which peptide–MHC targets have sequence support.
Repertoire inference asks whether target support exceeds matched mock-reference
expectations. The [repertoire method guide](https://docs.isalgo.dev/vdjmatch/repertoire-method.html)
connects the matching rule, null population and counting unit to that second question.
Full-data discovery curves retain clone abundance while measuring matched
clonotype richness over whole-library sampling effort; raw counts, reads and
Hill diversity remain separate outputs.

The Python implementation follows the legacy Java/Groovy package, preserved on the
`legacy-java` branch and [Java tag `1.3.1`](https://github.com/antigenomics/vdjmatch/tree/1.3.1).
Existing Java consumers use the separate, pinned JitPack artifact
`com.github.antigenomics:vdjmatch:1.3.1`; Python releases do not replace that artifact.
See the [tutorial](https://docs.isalgo.dev/vdjmatch/tutorial.html),
[workflow guide](https://docs.isalgo.dev/vdjmatch/how-to.html) and
[score explanation](https://docs.isalgo.dev/vdjmatch/explanation.html).

## Pick an entry point

| Task | Entry point | Guide |
| --- | --- | --- |
| Annotate a repertoire or cohort | `vdjmatch annotate` | [Tutorial](https://docs.isalgo.dev/vdjmatch/tutorial.html) |
| Compare retained geometries from raw inputs | `sh examples/compare.sh NEW_DIRECTORY` | [CLI examples](https://docs.isalgo.dev/vdjmatch/examples.html) |
| Calibrate the combined score from one control | `vdjmatch annotate --control CONTROL` | [Calibration guide](https://docs.isalgo.dev/vdjmatch/control-calibration.html) |
| Check literal database presence and unusual resemblance | `vdjmatch database-evidence` | [Score explanation](https://docs.isalgo.dev/vdjmatch/explanation.html) |
| Choose a pinned or offline reference | `--vdjdb LOCAL_REFERENCE` or `vdjmatch update` | [Workflow guide](https://docs.isalgo.dev/vdjmatch/how-to.html) |
| Test full-depth matched-clonotype enrichment | `vdjmatch repertoire-exact --pgen-features FEATURES --abundance-metrics ...` | [Repertoire method](https://docs.isalgo.dev/vdjmatch/repertoire-method.html) |
| Plot completed enrichment against every mock | `vdjmatch repertoire-report COUNTS SUMMARY ...` | [Offline notebook](examples/repertoire_evidence.ipynb) |

## Runnable offline examples

Browse the [examples gallery](https://docs.isalgo.dev/vdjmatch/examples.html) or
the [notebooks gallery](https://docs.isalgo.dev/vdjmatch/notebooks.html) for
rendered previews, source downloads and raw inputs.
The [versioned Hugging Face asset mirror](https://huggingface.co/datasets/isalgo/airr_benchmark/tree/314b72d2955e9578157c6767845741444d89f540/vdjmatch/releases/0.6.3)
also provides the synthetic offline examples and bundled static models with a checksum manifest.

The [examples directory](examples/) contains tiny raw synthetic inputs, a pinned
offline reference, [annotation](examples/annotate.sh) and
[comparison](examples/compare.sh) shell scripts, and output-free
[annotation](examples/annotation.ipynb) and
[comparison](examples/comparison.ipynb) notebooks. Activate your installed
environment and run from the repository root:

```sh
sh examples/annotate.sh annotation-demo
sh examples/compare.sh comparison-demo
```

Both use four native threads, reject existing output directories and fail on CLI
errors. The comparison explicitly runs exact-excluded detailed stem4 and advanced
vcdr12; notebooks analyze only fresh CLI tables. `DEMO_A`/`DEMO_B` labels are
synthetic software examples, not binding truth. See the
[example guide](https://docs.isalgo.dev/vdjmatch/examples.html).

The offline [repertoire evidence notebook](examples/repertoire_evidence.ipynb)
uses synthetic final-table fixtures to demonstrate the report CLI, all 200 mock
endpoints and preservation of recorded P/q values. It performs no biological inference.

## Scientific worked examples

Three longer notebooks connect raw biological inputs to fresh CLI tables and
publication-size figures, with editable PDF/SVG and 300-dpi PNG exports:

| Question | Notebook |
| --- | --- |
| How do common epitopes separate supplied positives and controls? | [ROC/PR, score overlap and TP/FP counts](examples/accuracy.ipynb) |
| Draft: LLW neighbourhood evidence across donors (validation deferred) | [Vaccination dataset walkthrough](examples/vaccination.ipynb) |
| How do model mass, precursor frequency and expected cells differ? | [Precursor estimates](examples/precursor_estimates.ipynb) |

Install the matching checkout from its repository root with
`python -m pip install -r examples/requirements-notebooks.txt`; keep the complete
`examples` directory together. Inputs remain in the local `airr_benchmark` and
`airr_yfv19` mirrors, configured through `AIRR_BENCHMARK` and `AIRR_YFV19`.
The notebooks preserve supplied labels and never reuse predictions as inputs.
Browse the [rendered gallery](https://docs.isalgo.dev/vdjmatch/notebooks.html)
and [measured laptop resource table](https://docs.isalgo.dev/vdjmatch/performance.html).
Four threads sped up the 2,670-query example from 13.71 to 4.69 seconds
(2.92×), with identical evidence and about 0.85 GB peak CLI process-tree RSS.

## Features

- Latest/pinned GitHub releases and offline legacy, rich or AIRR reference inputs.
- Source-aware junction ingestion, stable row identities and explicit ingestion counts.
- Native batched searches by locus, with retained reference observations.
- Targeted stem4 kernel annotation; frozen alignment scores, PSSM-density, fixed-scope unit-edit distance
  and neighbor counts, and other retained scorers remain comparison options.
- All peptide–MHC candidates, NED-v1 ranking and competing-label support.
- Optional finite-background calibration with a stated counting unit and search predicate.
- Paired evidence from linked reference complexes plus unlinked same-target chain marginals.
- CLI hit/candidate/call/summary/ingestion tables and a reproducibility manifest.
- Separate precursor-frequency APIs and their established worked examples.

### Targeted annotation

```console
$ vdjmatch annotate sample.tsv --threads 4 -o annotation
$ vdjmatch annotate airr_directory/ --threads 4 -o cohort
$ vdjmatch annotate --metadata samples.tsv --vdjdb reference.zip --threads 4 -o cohort
```

The reference is fetched automatically; `--vdjdb` selects a pinned offline table
or archive. A metadata TSV has `path` and optional `sample_id` columns; relative
paths resolve beside the metadata file. Directory discovery accepts TSV, CSV,
text, gzip and Parquet tables. Sample names must be unique.

AIRR input uses anchor-inclusive `junction_aa` and V calls. Version 0.6.0 defaults
to `--reference-population curated-benchmark`, matching the manuscript primary
reference view. Use `--reference-population all-valid` to reproduce the accepted
all-valid baseline. Curation excludes KLGGALQAK and records sourced to PMID
40498839; natural observations from other sources remain eligible. Exact hits
remain included unless `--exclude-exact` is supplied. The default scores all
peptide–MHC targets with the stem4 BLOSUM62 kernel ladder. `--locus paired` also
retains independent TRA/TRB evidence from unlinked references. Scores are relative
evidence; whole-selection null calibration is required for false-call recommendations.
The original PSSM-density scorer remains available through `historical-density`;
`trained-score` uses the bundled frozen alignment coefficients and hit informativeness.
Advanced alternatives remain available through `historical-density` and `neighbours`.

| Query format | Required sequence/gene fields | Identity/linkage |
| --- | --- | --- |
| Single AIRR | `junction_aa`, `v_call`, `locus` (or locus inferred from V) | Optional `sequence_id`, `query_id`, `clonotype_id`, `clone_id` |
| Simple repertoire | `junction_aa`, `v`, `locus`; optional `j` | Optional `clonotype_id`, `clone_id` or `cloneid` |
| Linked AIRR | One TRA/TRB row per linked receptor, `junction_aa`, `v_call`, `locus` | Shared `cell_id`, `pair_id`, `clone_id` or `cloneid`; cell/pair linkage takes priority |
| Wide paired repertoire | `junction_aa_alpha`, `junction_aa_beta`, `v_alpha`, `v_beta`; optional `j_alpha`, `j_beta` | Optional `clone_id` or `cloneid`; otherwise row identity links the pair |

`--locus auto` infers single or paired populations; `--locus TRA`, `TRB` or
`paired` selects explicitly. Missing chains and invalid rows retain dispositions;
multiple chains in one clone/locus require resolution before annotation.
`clonotype_id` and `clone_id` are retained as separate identities. Bare AIRR
`cdr3_aa` excludes conserved anchors and cannot replace `junction_aa`.
For amino-acid clonotype inputs, collapse by the supplied V, junction and J
calls, summing counts before annotation. Distinct J calls remain distinct
clonotypes even when their search geometry is shared.

The default writes `.queries.tsv`, `.selection.tsv` and `.manifest.json`, retaining
sample and source identities, competing scores, ties and no-support rows.
`--details` additionally writes candidates and per-target evidence; its memory
grows with queries × targets. `--epitope PEPTIDE` restricts the requested targets.
`--exclude-exact` is for honest benchmarks; ordinary annotation retains exact hits.
V-loop lookup accepts model-backed dual-use gene spellings such as
`TRAV14_DV4*02` and `TRAV14/DV4*02` as the same allele, while retaining
the submitted call in outputs. Unknown alleles are not silently renamed.

Advanced scorer and geometry recipes remain in the [CLI reference](https://docs.isalgo.dev/vdjmatch/cli.html).
The manuscript's curated primary comparison has higher equal-task mean average
precision than each of seven available rivals on that rival's fixed matched task
set. This is a matched-mean result, not a claim of winning every individual task.
The accepted all-valid baseline and its curves and intervals remain separate;
changing the reference population does not transfer calibration between views.
Calibration is specific to the complete score and selection procedure.
`--threads N` budgets native compute threads; searches batch queries by population
and do not create a worker process per receptor or input file.
Mixed single/paired cohorts share marginal chain searches. V-loop tables are
prepared once per invocation, including the selected matrix and V regions.
The library also enforces this ceiling for the Polars table pool. If Polars was
already imported with a larger pool, start a fresh process with
`POLARS_MAX_THREADS=N` set before importing Polars; the API fails explicitly
rather than exceeding the requested budget. Input metadata whose names conflict
with generated output fields is retained under `input_` prefixes.

```python
from vdjmatch import annotate_repertoires

result = annotate_repertoires("airr_directory", reference="reference.zip",
                             output_prefix="cohort", threads=4)
queries, selections = result["queries"], result["selection"]
```

### Complete-score calibration from one control

```sh
vdjmatch annotate sample.tsv --vdjdb reference.zip --control control.tsv --threads 4 -o calibrated
```

With the default `--background vdjdb-other`, supplying one raw control selects the
opt-in combined score and computes both conditional randomization tails on demand.
It writes `.queries.tsv`, `.database.tsv`, `.epitopes.tsv` and `.calibration.tsv`;
no fitted calibration file or JSON input is needed. `--mock-replicates` defaults
to 200 and `--seed` to 0. Reference selectors (`--epitope`, MHC,
`--evidence-type`, `--min-score`) apply before calibration.

Benchmark `--exclude-exact` filters exact hits only; reference models, counts,
priors and graph weights remain fixed across queries.

Your own reference can be a `junction_aa,locus,epitope` table, with `pair_id`
linking long paired rows, or a wide `junction_aa_alpha/beta,v_alpha/beta,j_alpha/beta,epitope`
table. Supplied pairing and metadata are retained; calibration uses that selected table.

The control requires full junctions and unambiguous V/J calls. Productive flags
and valid amino-acid junctions define eligibility; excluded counts are reported.
Rearrangement-event weights, rather than read abundance, define sampling.
For paired queries, the same control table supplies TRA and TRB, drawn independently.
Missing partners retain their raw rows and use eligible single-chain evidence.

Database tails include label and paired-axis selection; epitope tails test a
fixed full peptide–MHC group. `p_both = max(p_reference, p_query)` requires evidence
in both directions. Null match rates are separate `k/B` estimates; none of these
quantities is a binding probability. This candidate recipe derives the original
component prior from the matching chain control (TRB for paired queries), which
differs from historical benchmark model assets. It remains opt-in pending full
nonregression checks. The no-control stem4 default and older real/generative
control routes remain available. See the [calibration guide](https://docs.isalgo.dev/vdjmatch/control-calibration.html).

The original `Annotator` and `annotate` APIs remain available. Agent usage is
documented in the repository [vdjmatch skill](skills/vdjmatch/SKILL.md).
The synthetic library/CLI batch resource gate runs with test/bench extras:
`python -m pytest tests/test_annotation_resources.py -q`.

Offline call-policy experiments remain available through `reference-calls`,
`calibrate-calls` and `apply-call-calibration`; see the
[CLI reference](https://docs.isalgo.dev/vdjmatch/cli.html). Their pinned operating
points belong to the exact scorer, reference, null population and paired selection
event used to fit them. They are not inputs to the one-control workflow above.

### Precursor frequency

How much repertoire mass can see an epitope, and how many cognate clonotypes exist — seen and unseen.
Use the documented precursor setup (`pip install 'vdjmatch[precursor]'`) for these examples:

```console
$ vdjmatch precursor --vdjdb --mhc-class MHCI --min-junctions 10 \
      --n-eff 1e8 --selection auto -o precursor.txt
```

```python
from vdjmatch import precursor as P

model = P.load_model("TRB")
junctions = ["CASSIRSSYEQYF", "CASSLGQAYEQYF", ...]     # junctions, not IMGT CDR3s

P.union_mass(model, junctions)      # exact union of the 1-mm balls, and the overlap a sum invents
P.occupancy(model, junctions, n_eff=1e8, selection=P.SELECTION_BY_CHAIN["TRB"])
# -> {"S": ..., "F": ..., "n_seen": ..., "n_unseen": ..., "seen_fraction": ...}
```

`S` is the effective cognate-set size, `F` the fraction of the naive repertoire, and `n_seen`
saturates in depth — a cognate set concentrated on a few high-`Pgen` junctions is exhausted at
shallow sequencing, a broad one keeps accumulating clonotypes. `SELECTION_BY_CHAIN` carries the
measured per-chain depth factors (TRB 4.62, TRA 1.07); they differ and should not be pooled.

Two things worth knowing before trusting a number:

- **`--source` picks the recombination model, and it matters for a mass over a set.** The default
  `olga` is a bit-faithful import of OLGA's published models — faithful to OLGA's deletion-bin grid
  too, under which **5.7% of human TRA junctions in VDJdb score `Pgen` exactly zero** (up to 15% for
  some epitopes) and vanish from the total without an error. `isalgo` selects current gene-level fits; `learned` and `arda` alias those same
  weights. Mouse TRA/TRB/IGH/IGK/IGL are available; mouse TRG/TRD fail explicitly.
  Historical learned/ARDA results require their original pinned weights. The `n_zero_pgen` output column reports the
  loss per group and the CLI warns above 1%.
- **`--species` selects records, `--organism` selects the model, and they must agree.** A mismatch
  used to run and return a plausible number that was wrong by an epitope-dependent factor; since
  0.3.0 it is refused.

#### Checking it against what laboratories measured

`F(e)` estimates a **naive repertoire mass**, and that quantity has been measured directly — by
tetramer and dextramer enrichment of unexposed donors, cord blood and naive mice. Both sides of that
comparison ship, and are deliberately kept apart:

```python
from vdjmatch.precursor import load_compendium, load_estimates

comp = load_compendium()                                  # experimental: 147 published measurements,
                                                          # 14 studies, 31 originating laboratories
est = load_estimates(species="human", chain="TRB", model="arda")   # derived: what the model predicts
est.select("epitope_seq", (est["union"] * 1e6).alias("predicted"), "measured_per_1e6")
```

Both fetch from `isalgo/airr_benchmark` (`vdjmatch/precursor_freq/`) and cache next to the VDJdb
releases. Merging them into one column would conflate a measurement with a prediction of it.

Beyond a point estimate, the same two inputs — the `Pgen` spectrum and the catalogued junction count
— answer more than one question, and the answers are not interchangeable:

```python
u = P.union_mass(model, junctions, r=1, count_members=True)
u["n_union"] - u["n_seqs"]                     # uncatalogued neighbours: a census, finite and exact
P.observed_mass(model, junctions) / u["union"] # share of the neighbourhood's mass already catalogued
u["overlap"]                                   # what naive per-sequence summing would double-count
```

**Unseen mass converges; unseen count does not.** The Horvitz–Thompson weight `p/π → 1/N` as
`p → 0`, which is what makes `coverage_corrected_mass` work without guessing the shape of the tail.
But `1/π` diverges, so a richness extrapolation on real TCR data returns 10¹¹–10¹⁴ unseen junctions
— not a number. `unseen_junctions` marks `richness_reliable` false when the rarest observed junction
was captured too seldom, and the ball census above is the finite alternative.

An interactive [marimo](https://marimo.io) notebook walks the whole comparison — the scatter against
measurement, the `Pgen` spectrum inside one cognate set, the ball census, and the step from a mass to
a number of cells:

```console
$ pip install 'vdjmatch[precursor,notebook]'
$ marimo edit docs/notebooks/precursor_compendium.py
```

The method, its benchmarks and the paper live outside this repository:
[repseq/2026-precursor-freq](https://github.com/repseq/2026-precursor-freq) (benchmarks, result
tables, the literature compendium) and
[repseq/2026-precursor-freq-ms](https://github.com/repseq/2026-precursor-freq-ms) (the manuscript).
This repository holds only the software, its CLI, tests, docs and the worked examples
(`docs/notebooks/`).

## Install

```console
python -m pip install vdjmatch
```

For development, use an isolated environment:

```console
python -m venv .venv
source .venv/bin/activate
python -m pip install -e ".[test,bench]"
```

## Python API

```python
import polars as pl
from vdjmatch import Annotator

ann = Annotator.latest()  # GitHub latest; .version(tag) pins a release
# ann = Annotator.from_path("reference.zip")  # existing local input, offline
queries = pl.DataFrame({
    "sequence_id": ["query-a", "query-b"],
    "junction_aa": ["CASSIRSSYEQYF", "CWWWWF"],
    "locus": ["TRB", "TRB"],
})
calls = ann.annotate(queries, threads=4)  # every input row remains
candidates = ann.candidates(queries, threads=4)  # all peptide–MHC alternatives
calibrated = ann.annotate(queries, calibrate=True, species="human", threads=4)
```

AIRR `junction_aa` contains the conserved anchors; AIRR `cdr3_aa` excludes them.
True AIRR input without a junction is rejected. Legacy VDJdb `cdr3` can already
contain the junction. Custom ingestion requires an explicit junction convention;
the software never fabricates anchors.

The API uses uncalibrated ranking unless `calibrate=True` or an explicit control
index is supplied. NED ranks candidates; `E` is an expected background match count,
and `p_enrichment` tests the stated reference set. None is a posterior probability
of the winning specificity. Exact ranking ties abstain as `ambiguous`.

`ann.annotate_paired(cell_df, cdr3a="cdr3a", cdr3b="cdr3b")` requires both chains
to match the same reference complex and peptide–MHC label. Paired background
calibration states an independent-chain null explicitly. Agreement between
separate alpha and beta top labels alone does not establish paired-reference evidence.

## Command line

```console
vdjmatch update --asset default --cache reference-inputs
vdjmatch match --input-format airr --threads 4 -o results/run sample.tsv
vdjmatch match --vdjdb reference.zip --no-evalue -o results/offline sample.tsv
vdjmatch match --vdjdb reference.zip --paired --link cell_id -o results/cells cells.tsv
vdjmatch match --vdjdb reference.zip --search-mode ball -o results/balls sample.tsv
```

CLI matching requests calibration by default; `--no-evalue` selects uncalibrated
ranking. The command emits hits, candidates, row-preserving calls, descriptive
summaries and ingestion diagnostics as TSV, plus a JSON manifest.

Use `--input-format custom --sequence-convention junction` for an ambiguous custom
table. `--scope s,i,d,t` sets substitution/insertion/deletion/total budgets. Hard
`--match-v` or `--match-j` requests require `--no-evalue` with sequence-only backgrounds.
Ordinary CLI search defaults to allocated CPUs (`--threads 0`); use `--threads 1` for serial execution. See `vdjmatch match -h` and the
[CLI reference](https://docs.isalgo.dev/vdjmatch/cli.html).

Ranked mode scans the complete eligible reference with a native restricted gap-block scorer,
retaining at most K unique junctions or genuine paired keys before expanding observations.
Its global fixed-K enrichment test is separate from descriptive epitope ranks. Ball mode
reports radii1–5 separately. Neither selects the best K/radius after seeing the results.
These new modes require seqtree1.0.2; see the CLI reference for scoring and null assumptions.

## Scoring: what works (and what doesn't)

For new fuzzy-matching comparisons, start from the established BLOSUM62 and
end-anchored positional-BLOSUM62 implementations. Keep the one-mismatch baseline.
The M4 conclusion is that BLOSUM62 is optimal among the tested substitution
matrices: the substitution patterns reflect V(D)J rearrangement, N-insertion
and positional conservation rather than amino-acid chemistry. The observed
similar accuracy near the BLOSUM62/conservation upper cap supports the expectation
that fully optimized VDJAM would converge to this same generative structure;
this does not claim that optimization is complete. Keep the measured constructions distinct: the
genetic-code substitution-accessibility comparator **VDJAMr** matches BLOSUM62,
while the data-derived **VDJAM** matrix underperformed in the recorded retrieval
comparison. Positional weighting provides the additional improvement below.
Reuse these recorded comparisons before proposing another matrix.

An empirical study on VDJdb (see `appendix/vdjmatch_scoring.tex`; regenerated on the **2026-06-11-ZENODO**
release with composition controls and balanced metrics) settles the scoring question honestly:

- **Hamming distance 1 is the signal:noise optimum** — macro purity (per-epitope mean) falls
  0.49 → 0.07 across edit distance 1–5, a 56× → 2.5× enrichment over chance (reproducing the original
  VDJdb observation, Shugay et al. NAR 2018). The search-ball radius, not the substitution matrix, is
  the dominant lever. (On the dense 2026 release this only shows up once a few 10× mega-studies are
  capped and the random tail is treated as an admixed control — naive pooled purity reads a flat ~0.9.)
- **Central (NDN) substitutions carry the specificity signal** — a mismatch in the CDR3 core most
  often changes specificity (P(same epitope) ≈ 0.31) while near-anchor mismatches are germline noise
  (≈ 0.75); the NDN core is also ~31–34% glycine (insertion/D-gene signature).
- **No amino-acid matrix clearly beats BLOSUM62 — and a genetic-code null ties it.** BLOSUM62 ≈ PAM250 ≈
  structural > Hamming > data-derived VDJAM. Strikingly, **VDJAMr** — a matrix built from the *genetic
  code alone* (how mutationally accessible one AA is from another; `loo_vdjam.codon_dissim`) — matches
  BLOSUM62 (0.581 vs 0.564 @≤2), so TCR CDR3 substitution structure is **generative, not chemical**.
  A published TCR-specific matrix, [tcrBLOSUM](https://doi.org/10.1093/bib/bbae602), does *not* beat
  BLOSUM62 either (0.555 vs 0.567; `bench/tcrblosum_refute.py`) — its same-epitope counts don't transfer.
- **Position-weighting BLOSUM62 does beat it.** Encoding the central-substitution finding as a
  seqtree positional matrix (`PositionalMatrix.from_weights(BLOSUM62, …)`, centre ~2× the borders)
  raises leave-one-out retrieval (balanced PR-AUC) above flat BLOSUM62 in **7/8 held-out epitopes**
  (0.564 → 0.598 @≤2; 6/8 @≤4). The weight is **end-anchored** (offset from each germline anchor,
  Beta-Binomial-smoothed) and BLOSUM severity discriminates exactly where the weight is high (the core).
  For CDR3, *where* a mismatch falls matters more than *which* residue it is; the first-order statistic
  is still the control-calibrated E-value.
- **The V gene is a strong, near-binary prior — partly recovered at near-exact germline identity.**
  Same-V neighbours share the epitope **~45–64%** of the time vs **~6–17%** cross-V (ratio up to ~7×).
  *Loose* CDR1/CDR2 similarity barely predicts it (point-biserial r ≈ 0), but cross-V co-specificity
  **rises monotonically as germline CDR1+CDR2 approach identity** — ~11% at ≥6 mismatches to ~32% at
  edit-0: a real ~3× lift, but only **~half** the same-V level (32% vs 64%), recovered *whole-loop*
  rather than via a sparse pseudosequence (per-position lift flat except CDR2 pos 5; edit-0 bin n≈25,
  noisy). So the germline loops carry part of the V prior at near-exact tolerance; the rest is
  gene-identity-specific (`bench/vregion_decompose.py`, `bench/vpseudo.py`; `vdjmatch.match.vgene`).

## Benchmark

The standing benchmark is the **2026-06-11-ZENODO VDJdb release** (mirrored on the
[`isalgo/airr_benchmark`](https://huggingface.co/datasets/isalgo/airr_benchmark) HF dataset, fetched via
`db.fetch_hf`). For the scoring studies it is **composition-controlled** (`_bench.long_list`: keep
epitopes ≥30 clonotypes, cap mega-epitopes to a random 3000, drop spectratype-anomalous spike studies),
with imbalance-robust metrics (`bench/metrics.py`: ROC-AUC + balanced PR/F1). Its highest-confidence
subset — the **shortlist** of clonotype–epitope pairs in **≥2 independent references** (`db.replicated`,
~3000 TRB + ~2000 TRA) — is the gold standard (as in `mhcmatch`), kept separate from the long-list.
Leave-one-out NN annotation (`bench/shortlist_accuracy.py`, subs 1, exact self excluded) reaches
**~44–47% top-1**, modestly above single-reference controls (the dramatic 3× gap on a sparse older
export was an artefact — the dense release makes controls findable too).

## Repertoire-level enrichment (M3)

Start with the [algorithm and parameter guide](docs/repertoire-method.rst).
It explains the two null models, matching/V/J choices, chain and epitope
aggregation, Fisher evidence and the completed manuscript comparisons.

### Full-depth discovery from raw marginal Pgen

```sh
vdjmatch repertoire-exact samples.tsv --vdjdb reference.zip \
  --control observed-control.tsv --pgen-features marginal-pgen.tsv \
  --abundance-metrics --output-prefix full
vdjmatch repertoire-report full.counts.tsv full.summary.tsv \
  --metric curve_mean --output-prefix figures/richness
```

Use verified junction-only marginal Pgen from the current vdjtools model and
the observed control repertoire. Probabilities enter the existing bins directly,
without Q, V/J-frequency corrections, AA-max reduction or a model prior.
The full source-clonotype abundance vectors supply expected matched richness
over a common whole-library read-depth range. This curve-average endpoint is
distinct from raw counts, reads, Chao1 and matched-subset Hill diversity; all
remain available for comparison. Declare test families with `--fdr-group`
metadata columns. The report displays recorded empirical P and BH q values
alongside every mock endpoint; it does not refit inference.

### Historical matcher comparisons and reproduction

There is no universal best radius: positional BLOSUM62 radius 0.30 improved
bulk day15 clonotype detection (5/7 libraries; 7/7 by reads, 0/12 day0), while
Hamming<=1 retained the broader paired single-cell exposure coverage
(5/6 day14 donors, 6/6 ever post, 0/6 prechallenge). These are different endpoints.
Keep exact Pgen, original correction families and the full reference scope
attached to each result. Chain P-values support inference; paired exposure and
cross-epitope unions are descriptive. See the guide before transferring settings.


`vdjmatch repertoire-exact samples.tsv --vdjdb reference.zip --control probes.tsv --model-vj-counts model-vj.tsv --output-prefix enrichment`
preserves the compiled reference-panel bulk calculation: exact junction/V/J matching,
source clonotype rows and read abundance, joint Pgen with the original AA-max
reducer and within-V J-usage correction. Use `--sampler r --pool-mhc` for the
original single-epitope sampling convention; declare replicate-level correction
families explicitly with repeated `--fdr-group` options. NumPy is the default
sampler. Complete-zero tests and original-notebook statistics are separate columns.

See [CLI definitions](docs/cli.rst) for optional panel matchers and input filtering.

The following flexible-matching mode uses a different null and counting unit;
it is not a reproduction of the compiled reference-panel calculation.


`vdjmatch repertoire-enrichment samples.tsv --null reference --vdjdb reference.zip --control control.tsv.gz --threads 4 --output-prefix enrichment`
holds each repertoire fixed and compares its unique receptor support against fresh
Pgen-matched mock epitope corpora. `--mock-policy pgen-vj` preserves V/J usage;
`pgen-different-v` changes V and junction identities. Paired inputs retain orphan
chains and count a receptor once per target. `--significance nominal|bh|both`
selects explicit inference options; no binary significance call is the default.
All no-hit targets remain in the outputs. See the
[CLI contract](https://docs.isalgo.dev/vdjmatch/cli.html) for matching, Pgen,
replicate resolution, test families and input schemas. The historical
`--null repertoire` mode remains separate.

For paired reference-null comparisons, `--chain-mode tra|trb|union|intersection|all`
preserves receptor identity throughout real and mock matching. The asymmetric
`--tra-matcher stem4 --mismatches 1 --chain-mode intersection` combines strict
TRB support with fuzzy TRA evidence for the same target; chain scores remain
separate from repertoire significance.

The established `unit`, `v-unit`, `weighted` and `v-weighted` families are also
available through these matcher options. They reuse the graded edit-bounded
alignment engine. An explicit chain `--*-radius-factor` bounds the best individual
alignment cost divided by 1 (unit) or 400 (weighted); `--mismatches` sets the
1–3 total-edit ceiling. Count each accepted receptor and its abundance once.
Weak alignments cannot add up to cross this bound. The fourteen-setting bulk
comparison is complete: hard-V exact and hard-V weighted bounds 0.1/0.2 retain 4/7 day15
clonotype and 7/7 read calls, with 0/12 day0 calls. Weighted bound 0.1 is exact
membership by construction. The separate positional BLOSUM62 preset at 0.30
remains the strongest recorded bulk matcher in that historical comparison
(5/7 and 7/7), using pooled J-usage correction within V. These results retain
their producing null and are separate from the current raw-marginal discovery
endpoint. Paired four-family outputs retain their own completed-arm audit and
scope. Summed continuous evidence is deferred.

Both `--tra-matcher` and `--trb-matcher` accept `hamming`, `stem4`, `gapblock`,
`blosum62` or `positional-blosum62`. The direct BLOSUM presets require
`--tra-radius-factor`/`--trb-radius-factor`: a dimensionless budget relative to
the matrix substitution scale, with fixed positional weights and at most three
substitutions, no indels, exact V and optional exact J. The bulk panel command
uses the same names with `--radius-factor` and retains exact V/J. These
pairwise budgets are distinct from kernel scores. `--tra-min-score` and `--trb-min-score` set dimensionless strict
lower bounds for fuzzy support; the same bounds apply to real and every mock
panel. Keep calibrated scorer parameters fixed. The default zero accepts weak
positive support and is a historical diagnostic, not a precision operating point. Use `--trb-matcher stem4 --chain-mode intersection` for the reverse
asymmetric comparison. Cutoffs require calibration for their own kernel and
locus; positive support alone is not calibrated binding confidence. Defaults
remain Hamming on both chains. Full reference models are preserved.

An independent counts-only comparison filters both the sample and each pMHC
reference with the existing TCRnet neighbourhood test, then runs flexible
VDJmatch matching:

```sh
vdjmatch repertoire-tcrnet samples.tsv --vdjdb reference.zip \
  --control control.tsv.gz --tcrnet-fdr 0.1 --scope 3,1,1,3 \
  --threads 4 --output-prefix filtered
```

This mode exports filter evidence, matching hits and raw-denominator counts;
it does not produce repertoire P-values or assume random backgrounds filter to
an empty set. Paired inputs retain orphan chains and count each receptor once.


For inputs too large to retain as one sample, `vdjmatch.io.iter_rearrangement`
yields normalized batches with file-wide positional IDs. The underlying
`vdjtools.io.iter_table_batches` supports projected raw CSV/TSV/gzip/Parquet
batches. Prefer the simpler one-sample-at-a-time path when profiling shows it
fits; within-sample batching does not by itself bound reference/control memory.

`--chain-mode all` computes the four chain rules from the same native searches and mock draws.
Summary, mock and hit tables include a chain_mode column. BH correction is
performed separately for each rule within the selected target or cohort family.

Reference-null repertoire enrichment accepts ``--model-source isalgo`` or
``--model-source olga-legacy``. Fresh probabilities, model identity checks and
prepared feature validation use that same family. Omission retains historical
model selection until the calibrated-default release is finalized.

Reference-null significance reports both unique-receptor and abundance endpoints.
``q_count``, ``significant_nominal_count`` and ``significant_bh_count`` describe
summed abundance (reads in bulk data, supplied cell counts in single-cell data).
Existing unsuffixed call flags describe receptor counts. Each endpoint uses the
selected correction family separately; positive excess is required for a call.

For cohorts that fit in memory, reference-null ``repertoire-enrichment --batch-samples``
shares each seeded mock draw and native search across samples. It preserves
per-sample receptor identities, counts, null statistics and the selected BH family.
The default still loads and processes one sample at a time.

## License

ISALGO LAB Academic and Non-Commercial Research License v1.0. See [`LICENSE`](LICENSE).
