GenAIRR

Synthetic Adaptive Immune Receptor Repertoire Generator

High-performance BCR and TCR sequence simulation with full ground-truth annotations.
Rust kernel · 23 species · constraint-aware sampling · cross-platform wheels

📖 Documentation

Installation

pip install GenAIRR

GenAIRR ships as a single wheel that bundles both the Python API and the Rust simulation kernel — no extra packages, no compiler needed. Pre-built wheels are published for Linux (x86_64, aarch64), macOS (Intel + Apple Silicon), and Windows (x64), supporting Python 3.9+.

Building from source needs a stable Rust toolchain (rustup install stable) — see CONTRIBUTING.md.

Quick Start

import GenAIRR as ga

# Generate 1,000 productive human heavy-chain sequences. Every sequence
# comes back with the full AIRR-format annotation block — gene calls,
# junction, productive flag, identity, mutation counts.
result = (
    ga.Experiment.on("human_igh")
      .recombine()
      .productive_only().run_records(n=1000, seed=42)
)

# `result` is a SimulationResult — list-like over AIRR record dicts.
# Each dict has the 50+ standard AIRR fields per row.
len(result)                 # 1000
rec = result[0]

rec["sequence"]             # 'caggtgcagctggtgcagtctggggctgaggtg...' (nucleotide)
rec["sequence_aa"]          # 'QVQLVQSGAEVKKPGASVKVSCKASGYTFTSYY...' (translated)
rec["locus"]                # 'IGH'
rec["v_call"]               # 'IGHVF6-G23*01'    (generating allele)
rec["recalled_v_call"]      # post-mutation inference; comma-separated when tied
rec["d_call"]               # 'IGHD1-20*01'
rec["j_call"]               # 'IGHJ1*01'
rec["junction_aa"]          # 'CARCRTRLTRAYHITGTAHFQHW'
rec["productive"]           # True
rec["v_identity"]           # 1.0  (matches/total over the V segment)
rec["n_mutations"]          # 0
rec["determinism_version"]  # 4    (which sampling model drew this record)

# Export in any of the standard formats. TSV/FASTA/FASTQ are dependency-free;
# to_dataframe() needs pandas (pip install pandas).
result.to_tsv("repertoire.tsv")        # AIRR-spec TSV (50+ columns)
result.to_fasta("sequences.fasta")     # FASTA with v_call/j_call in the headers
result.to_fastq("sequences.fastq")     # FASTQ with illumina-shaped quality scores
df = result.to_dataframe()             # one row per record, AIRR columns

Python record coordinates are zero-based and half-open. to_tsv() converts them to AIRR's one-based closed intervals, adds CIGAR placement offsets and writes booleans as T/F. Use to_tsv(..., airr_strict=False) for the legacy Python-style TSV. With rev_comp=True, biological annotations index the coding-oriented reverse complement of sequence; paired-read windows index the emitted sequence directly.

For paired-end sequencing pipelines, add .paired_end(...) to the experiment and export to two synchronized FASTQ files:

result = (
    ga.Experiment.on("human_igh")
      .recombine()
      .paired_end(r1_length=150, insert_size=300)
      .run_records(n=100, seed=1)
)
result.to_paired_fastq("reads_R1.fastq", "reads_R2.fastq")

Headers use the AIRR record's sequence_id with the universally-portable /1 and /2 suffix (@seq0/1 and @seq0/2) — no |-pipe metadata that some aligners (STAR < 2.7) split on. R1 / R2 sequences are written verbatim from the AIRR r1_sequence / r2_sequence fields (R2 is already reverse-complemented at projection time; the writer does NOT apply a second flip). Default overwrite=False refuses to clobber existing files; pass overwrite=True to replace. Quality strings use the same pluggable models as to_fastq — quality="illumina" (default trapezoid shape, applied per read) or quality="constant" with q=.... The requested insert_size has to fit inside the molecule: on human IGH a record is roughly 350–440 nt, so an insert of 450 leaves the R2 window off the end of the sequence and validate_records reports ReadWindowOutOfBounds. See audit-docs/fastq_export_design.md.

Experiment.on(...) accepts a config-name string (e.g. "human_igh", "mouse_tcrb"), a DataConfig loaded from the bundled species pickles, or a RefDataConfig for custom reference data. .productive_only() is the constraint-aware bundle — covered in the next section. Drop it to sample without the productivity constraint; records then report true, false, or unresolved productivity according to their observed evidence.

See the full walkthrough in the docs: Quick Start · Your first AIRR record · Interpreting Results

A realistic pipeline — everything in one place

The Experiment DSL is a fluent builder. Each step appends to the pipeline; the same Experiment is returned so calls chain. The example below uses every major feature GenAIRR offers — recombination, BCR clonal lineage trees with affinity maturation, per-division somatic hypermutation, primer-trimming, structural indels, PCR errors, N-base injection, custom metadata, and the productive constraint:

import GenAIRR as ga

result = (
    ga.Experiment.on("human_igh")
      # 1. V(D)J recombination — sample alleles, trim, fill NP1/NP2, assemble.
      .recombine()
      # 2. Clonal lineage trees — grow 50 real B-cell lineages from the
      #    recombined founder: generation-synchronous birth–death over
      #    5 generations, per-division S5F somatic hypermutation at 0.75%
      #    per base, affinity selection, then sample 20 cells per clone.
      #    Passes BEFORE this point shape the founder rearrangement;
      #    passes AFTER apply per sampled cell (library prep). Every clone
      #    shares the founder's V(D)J but carries its own SHM + errors,
      #    and you get a ground-truth tree per clone (see lineage_trees).
      .clonal_lineage(n_clones=50, max_generations=5, n_sample=20,
                      rate=0.0075, selection_strength=10.0)
      # 3. Sequencing artefacts per sampled cell: primer trimming, structural
      #    indels, PCR substitution errors, quality-driven N injection.
      .primer_trim_5prime(length=(0, 8))
      .primer_trim_3prime(length=(0, 4))
      .polymerase_indels(count=(0, 2), insertion_prob=0.5)
      .pcr_amplify(count=(0, 3))
      .ambiguous_base_calls(count=(0, 2))
      # 4. Stamp arbitrary metadata onto every record.
      .with_metadata(experiment_id="exp001", tissue="peripheral_blood")
      # The constraint is carried, not just applied to the founder: each
      # lineage division and each corruption pass draws only substitutions,
      # indels and masks that keep the cell productive, so every emitted
      # record is productive by construction.
      .productive_only().run_records(seed=42)
)

len(result)                                  # 440   (sampled live cells across 50 trees)
sum(1 for r in result if r["productive"])    # 440   (all of them — see the note above)
len(result.lineage_trees)                    # 50    (one ground-truth tree per clone)

# Same clone, different cells — shared founder V(D)J,
# independent SHM + errors:
result[0]["clone_id"], result[1]["clone_id"]              # (0, 0)
result[0]["v_call"],   result[1]["v_call"]                # both 'IGHVF6-G21*11'
result[0]["n_mutations"], result[1]["n_mutations"]        # (11, 17) — independent SHM
result[0]["n_pcr_errors"], result[1]["n_pcr_errors"]      # (1, 0)   — independent errors

# Custom metadata propagated:
result[0]["experiment_id"], result[0]["tissue"]           # ('exp001', 'peripheral_blood')

result.to_tsv("repertoire.tsv")
newick = result.lineage_trees[0].to_newick()  # ground-truth lineage for clone 0

Clonal lineages & repertoires

GenAIRR models clonal structure three ways. Pick by what you need:

import GenAIRR as ga

# BCR clonal lineage trees — affinity maturation, with ground truth
bcr = (ga.Experiment.on("human_igh").recombine()
       .clonal_lineage(n_clones=50, max_generations=6, n_sample=30,
                       rate=0.01, selection_strength=10.0)   # 0 = neutral
       .sequencing_errors(rate=0.001)
       .run_records(seed=0))
bcr.records              # per-cell AIRR records: clone_id, lineage_node_id, lineage_affinity, ...
bcr.lineage_trees[0].to_newick()   # ground-truth tree (branch length = per-edge mutations)

# TCR clone-size repertoire — proliferated rearrangements, no SHM, abundance
tcr = (ga.Experiment.on("human_tcrb").allow_curatable_refdata().recombine()
       .clonal_repertoire(n_clones=200, size_distribution="power_law",
                          exponent=2.0, max_size=500, unexpanded_fraction=0.5)
       .sequencing_errors(rate=0.005)
       .run_records(seed=0))
tcr.records              # AIRR records: clone_id (membership) + duplicate_count (abundance)

clonal_lineage is BCR-only (it applies S5F SHM and rejects TCR loci); clonal_repertoire covers TCR and flat BCR. See the guides: Clonal lineage trees · Clonal repertoires (TCR & abundance).

Other feature flags worth knowing:

Per-individual genotypes

A Genotype is one subject's diploid germline complement — which V/D/J alleles they carry, on which chromosome, with which genes deleted or duplicated. Attach one and recombination becomes haplotype-phased: the V, D and J of a rearrangement are drawn from a single chromosome, so the linkage that haplotype-inference methods exploit is actually there in the data.

import GenAIRR as ga
import GenAIRR.data as gdata
from GenAIRR.genotype import Genotype

cfg = gdata.HUMAN_IGH_OGRDB
g = (
    Genotype.from_dataconfig(cfg)
      .complete_from_reference("homozygous_first_reference")
      .heterozygous("IGHVF1-G1", "IGHVF1-G1*01", "IGHVF1-G1*02")
      .delete_gene("IGHVF3-G7", haplotype="both")
      .with_subject("DONOR01")
)

result = (
    ga.Experiment.on(cfg)
      .with_genotype(g)              # recombination is now haplotype-phased
      .recombine()
      .run_records(n=1000, seed=7)
)

result[0]["subject_id"]              # 'DONOR01' — on every record
result[0]["haplotype"]               # 0 or 1 — which chromosome this read used

# The genotype leaves the library with the records it explains:
result.to_genotype_tsv("genotype_truth.tsv")
# subject_id  locus  segment  gene  allele  chromosome  novel  gene_state

to_genotype_tsv() writes one row per carried allele, so a dataset you hand to someone else arrives with the truth it is meant to be scored against. A deleted gene keeps its row with the literal allele (none); a run with no genotype, or a genotype that leaves reference genes unspecified, raises rather than writing a file that looks complete and is not. Whether a given read drew a novel allele is a join of truth_v_call against that table on (subject_id, allele) — run_records(..., expose_provenance=True) adds the truth_*_call columns.

run_cohort(genotypes, n_per_subject=...) runs many subjects in one call and namespaces their records; CohortResult.to_genotype_tsv(path) writes the whole cohort's truth table.

Docs: Genotypes · Genotype sampling & population priors · Genotype cohorts · Benchmarking genotype inference

Constraint-aware sampling

GenAIRR's signature feature is constraint-aware sampling: contracts that prune the candidate distribution at sample time, not retries after the fact. The canonical bundle is productive() (preserved coding frame, no stops across the observed coding sequence, and intact V/J anchors):

import GenAIRR as ga

# Every sequence is productive by construction. No retry loops, no
# post-hoc filtering — the engine only ever picks NP lengths, NP bases,
# and mutation substitutions that satisfy the bundle.
result = (
    ga.Experiment.on("human_igh")
      .recombine()
      .productive_only().run_records(n=1000, seed=42)
)
assert all(rec["productive"] for rec in result)

The named policy is canonical_observed_vdj_v1, scoped to the observed coding span. It requires V Cys and J Phe/Trp original landmarks, consistent frame, and definitely stop-free complete codons. Missing evidence or possible hidden stops remain nullable in unconstrained records. Reference admission exceptions never change this policy. Both execution modes require an accepted final assessment; this is a sequence-level prediction, not a claim of receptor function.

.productive_only() supports fixed and mixed nonzero P-nucleotide lengths. The engine keeps a candidate only when the remaining allele, trim, N/P, and D-orientation choices can complete a productive sequence. It renormalizes the original weights among those candidates at each draw. Fixed P lengths are preserved; impossible support raises an error even in permissive mode instead of dropping requested additions. This conditions each local draw, rather than sampling from a globally conditioned joint probability law.

Docs: Sample only productive sequences

Strict vs permissive mode

Initial allele/genotype assignment and P-aware productive completion fail explicitly when no admissible support exists, including in permissive mode. Other passes can use their documented no-op or sentinel behavior; strict=True surfaces those sampling failures as exceptions:

import GenAIRR as ga

try:
    ga.Experiment.on("human_igh").recombine().productive_only().run_records(n=10, seed=42, strict=True)
except ga.StrictSamplingError as e:
    pass_name, address, reason = e.args
    # pass_name e.g. "generate_np.np1", address e.g. "np.np1.length",
    # reason in {"empty_admissible_support", "support_unavailable", ...}
    print(f"{pass_name} could not satisfy the contract at {address}: {reason}")

Reproducibility

The promise is narrow and worth stating exactly: one seed and one determinism version give one answer. Two runs at the same GenAIRR.DETERMINISM_VERSION agree byte for byte; two runs at different versions promise each other nothing — not "mostly the same", not "statistically equivalent". Every exported record carries the number as a determinism_version column, and so does every trace file, so a dataset says which sampling model drew it.

import GenAIRR as ga

ga.DETERMINISM_VERSION          # 4 — the installed engine's sampling model

# Same seed + same determinism version → byte-identical records.
a = ga.Experiment.on("human_igh").recombine().run_records(n=100, seed=42)
b = ga.Experiment.on("human_igh").recombine().run_records(n=100, seed=42)
assert a[0]["sequence"] == b[0]["sequence"]
assert a[0]["determinism_version"] == ga.DETERMINISM_VERSION

# `n` runs use seeds [seed, seed+1, ..., seed+n-1], so a seed names a
# *record*, not a run: consecutive batches stitch together by offsetting
# the starting seed by at least `n`.
batch_a = ga.Experiment.on("human_igh").recombine().run_records(n=100, seed=0)
batch_b = ga.Experiment.on("human_igh").recombine().run_records(n=100, seed=100)
# batch_a[50] is byte-equal to a one-off run at seed=50.
# Offsetting by less than `n` overlaps: seeds 0, 1, 2 are not three
# independent replicates, they are one dataset read from three
# starting points.

The reproducibility goldens run in CI on Linux, macOS and Windows across Python 3.10–3.13 and assert that the same seed reproduces the same bytes on each. Linux aarch64 ships a wheel but is not covered by those goldens.

Docs: Trace, replay, reproducibility

Validation posture

GenAIRR's release readiness rests on seven guardrails, each backed by an audit doc and a golden-test file. Every PR runs the full suite; every release additionally runs the golden trace compatibility fixtures so on-disk formats stay byte-stable.

• Productive validity — junction-frame, no-stop, anchor-preserved bundle, with strict / permissive failure modes pinned per pass.
• Provenance correctness — indels, end-loss, allele-call ambiguity, junction fields, and per-pass event ledger each have isolated-scenario tests.
• Reproducibility — every audit slice includes trace-replay round-trips that reproduce sequence, AIRR coordinates, and per-pass event counts.
• Distribution invariants — Monte-Carlo tests (±5σ) prove the constrained samplers draw from natural_weight × admissibility, with explicit negative controls against renormalization bugs.
• Performance budgets — wall-time regression guards on nine representative workloads catch ~10× slowdowns before they ship.
• Postcondition validator — result.validate_records(refdata) runs the engine's own truth oracle over every projected AIRR record (re-derives counters, junction, allele tie-set, structural coords from outcome state). Returns a ValidationReport you can assert on as a one-line CI guard. This is the public output-correctness contract.
• Live-call cache parity — outcome.check_live_call_cache_parity(refdata) compares the cached SegmentLiveCall on the final simulation against a from-scratch recompute, per V/D/J. This is the internal cache-correctness check on the state that feeds projection.

Two independent layers, each one a single call:

result = exp.run_records(n=1000, seed=0)

# Layer 1 — public AIRR output correctness.
report = result.validate_records(refdata)
assert report, report.summary()

# Layer 2 — internal live-call cache correctness.
for outcome in result.outcomes:
    for p in outcome.check_live_call_cache_parity(refdata):
        assert p["tie_set_matches"], p

Troubleshooting rule. If a CI run has both layers failing on the same batch, fix the parity divergence first: a stale cache leaks into projection and produces spurious validator failures downstream. Once parity is green, rerun the validator — remaining failures then point at a genuine projection-layer bug rather than cache fallout.

The navigable index — guarantees → audit docs → test files — lives in audit-docs/validation_matrix.md. Use it to locate the right test for a change you're making, or to scope a follow-up that touches an open drift item in any audit's §6.

make validate-fast      # correctness only (~60s)
make validate-full      # + performance budgets (regression guard)
make validate-release   # + golden trace compat + wheel build

Compile once, run many times

For a hot loop, compile() once and reuse the plan. Contracts (respect=) are baked into the compiled plan, so they only need to be passed once:

import GenAIRR as ga

compiled = (
    ga.Experiment.on("human_igk")
      .recombine()
      .productive_only().compile()
)

# Run 10 batches of 100, seeded so they don't overlap.
for batch in range(10):
    result = compiled.run_records(n=100, seed=batch * 100)
    result.to_tsv(f"batch_{batch:02d}.tsv")

What you get back

.run_records(...) returns a SimulationResult — a list-like wrapper around a batch of AIRR record dicts:

Each record dict has 50+ AIRR fields. The most commonly used:

The full schema (plus the *_sequence_start/end, *_alignment_start/end, *_germline_start/end coordinate fields, vj_in_frame, stop_codon, rev_comp, and others) is documented at The AIRR record.

Advanced: full pipeline state via Outcome

When you need step-by-step IR history or the raw trace of every random draw — debugging an engine bug, building a custom alignment tool, replaying a specific seed — use .run() instead of .run_records(). It returns a list of Outcome objects that carry the full pipeline state:

Each Simulation exposes len(sim) (pool length), sim.bases() → bytes, sim.regions() → list[Region], sim.germline_position(i), sim.v_allele_id() / .d_allele_id() / .j_allele_id(). Each Region carries segment ("V"/"D"/"J"/"NP1"/"NP2"), start/end/len(), frame_phase, and amino_acids() → bytes (codon-rail translation, including stop markers and ambiguous codons).

outcome.trace() supports find(address), prefix_query(prefix), and prefix_count(prefix) — every random draw is keyed by a hierarchical address ("sample_allele.v", "np.np1.length", "np.np1.bases[3]", …). This is the same trace the engine uses internally for replay determinism.

.run_records(...) also exposes these via result.outcomes[i] — so you can have both the AIRR records and the deep introspection from a single call.

Docs: The Experiment builder · Trace, replay, reproducibility · The AIRR record

Supported Species & Chains

GenAIRR ships with 106 built-in configurations covering 23 species (sourced from IMGT and OGRDB).

import GenAIRR as ga
print(ga.list_configs())  # all available configs

import GenAIRR as ga

ga.Experiment.on("rabbit_tcrb").recombine().run_records(n=500)

# 80 of the 106 bundled catalogues carry unresolved anchor or
# pseudogene/ORF findings, so strict compilation quarantines them. Opt in
# explicitly to simulate from one — the records are then drawn from a
# catalogue whose landmarks the project has not been able to confirm.
ga.Experiment.on("mouse_igh").allow_curatable_refdata().recombine().run_records(n=500)
ga.Experiment.on("rhesus_igk").allow_curatable_refdata().recombine().run_records(n=500)

The 26 catalogues that compile strictly today are the three human OGRDB cartridges (human_igh, human_igk, human_igl), HUMAN_IGH_EXTENDED, and 22 IMGT catalogues whose anchors all resolved. Everything else — including every mouse catalogue and human_tcrb — needs allow_curatable_refdata(). The compile error lists the findings allele by allele, so you can also filter the catalogue or curate it with Experiment.curate_refdata(...) instead.

Docs: Choose your path · Reference evidence and quarantine · API Reference

Custom reference data

For non-builtin alleles (custom IMGT pulls, in-house references, etc.) you can build a RefDataConfig directly and pass it to Experiment.on(...):

import GenAIRR as ga

# Real bytes, so the example runs as written. Substitute your own.
igk = ga.HUMAN_IGK_OGRDB
v1, v2 = [a for g in igk.v_alleles.values() for a in g][:2]
j1 = next(a for g in igk.j_alleles.values() for a in g)

cfg = ga.RefDataConfig.vj()
cfg.add_v_allele("v_custom*01", "v_custom", v1.ungapped_seq.upper().encode(), anchor=v1.anchor)
cfg.add_v_allele("v_custom*02", "v_custom", v2.ungapped_seq.upper().encode(), anchor=v2.anchor)
cfg.add_j_allele("j_custom*01", "j_custom", j1.ungapped_seq.upper().encode(), anchor=j1.anchor)

result = (
    ga.Experiment.on(cfg)
      .recombine(np1_lengths=[(n, 1.0) for n in range(7)])  # no empirical NP model
      .trim(enabled=False)                                  # no trim distributions
      .run_records(n=100, seed=42)
)
result[0]["v_call"], result[0]["junction_aa"]   # ('v_custom*01', 'CLQHNSYPPWRFTF')

RefDataConfig.vdj() builds a heavy-chain-shaped refdata (with a D pool); add_d_allele(...) populates it. Anchors are 0-based offsets of the V Cys / J W or F codon's first base, used to keep the junction frame-aligned during recombination — the engine checks that the codon at the anchor really translates to Cys (V) or Trp/Phe (J) and refuses the reference otherwise. Sequence bytes must be A/C/G/T/N; a raw RefDataConfig carries no empirical NP-length or trim distributions, so pass them explicitly (as above) or accept the two UserWarnings that say the engine fell back to a uniform NP length and a no-op trim.

For non-human or otherwise non-standard references — custom J anchor amino acids, extended sequence alphabets, hand-authored NP-length or trim distributions, pseudogene-aware allele filtering — configure a reference cartridge instead of relying on the built-in locus defaults. The cartridge has four typed planes (identity, catalogue, rules, empirical models) plus an explicit curation policy. See audit-docs/reference_cartridge.md for the model, the authoring surface (ReferenceRulesSpec, ReferenceEmpiricalModels, Experiment.curate_refdata), and end-to-end examples.

Cartridges can also carry V-region substructure annotations — per-V-allele IMGT FWR1 / CDR1 / FWR2 / CDR2 / FWR3 intervals derived from the IMGT-gapped reference sequence. The bundled OGRDB cartridges populate them at load time for all human IGH / IGK / IGL V alleles; user cartridges can override the dict explicitly. Subregion intervals are inspectable on the bridged RefDataConfig, surfaced in DataConfig.cartridge_manifest()["models"]["shm"]["v_subregion_support"], and folded into refdata_content_hash so two cartridges that differ only in subregion boundaries hash differently. They also drive SHM targeting via Experiment.mutate(v_subregion_rates={…}) — see the targeted-SHM example in the pipeline-steps table above — and they bucket the per-region AIRR mutation counters (n_cdr1_mutations …) and their recoverable counterparts (n_cdr1_observable_mutations …), both described in the pipeline-steps table. See audit-docs/v_region_substructure_audit.md and audit-docs/v_subregion_shm_rate_design.md for the audits and per-slice scope.

Typed NP base models (uniform / empirical / Markov)

Cartridges can also author per-NP-region base distributions via the typed ReferenceEmpiricalModels.np_bases plane. Three kind values are supported end-to-end:

• "uniform" — equivalent to the engine default (4-way A/C/G/T, byte-identical to the pre-typed-model baseline).
• "empirical_first_base" — position-independent weighted categorical. Every NP slot samples from the same first_base distribution.
• "markov" — 1-step previous-base-conditional sampling. Position 0 uses first_base; positions 1+ select a transition row keyed by the previous emitted base.

from GenAIRR.reference_models import (
    NpBaseModelSpec,
    ReferenceEmpiricalModels,
)

cfg.reference_models = ReferenceEmpiricalModels(
    np_bases={
        "NP1": NpBaseModelSpec(
            kind="markov",
            first_base={"G": 1.0},
            transitions={
                "A": {"T": 1.0},
                "C": {"G": 1.0},
                "G": {"A": 1.0},
                "T": {"C": 1.0},
            },
        )
    }
)

Both rows must cover every canonical from-base — the Python validator rejects partial matrices at construction time. Plan signatures fold the full Markov payload (first-base row + 4 transition rows), so replay against a different matrix fails the signature gate before any choice is consumed; same-cartridge + same-seed replay is byte-identical. Productive-only sampling composes through the existing admit-mask intersection.

Legacy NP_transitions / NP_first_bases are NOT auto-lifted yet. The bundled cartridges still carry those dicts as documented orphan fields, and the manifest reports legacy_fallback=False. Authors who want Markov today populate ReferenceEmpiricalModels.np_bases explicitly; an opt-in auto-lift slice is a separate cartridge-migration decision. See audit-docs/np_markov_base_generator_design.md and audit-docs/junction_n_addition_audit.md for the full design.

Templated P-nucleotide (palindromic) additions

For an opt-in model tied to the original coding tips, configure .process_coding_ends(v_3=[(-1, 1.0)], reference_ends_are_coding_ends=True) after .recombine(). Positive signed outcomes delete reference bases, negative outcomes retain original-tip P extensions, and omitted ends are fixed at zero. This replaces automatic trim/P settings and cannot be combined with an explicit .trim() call. The boundary confirmation requires complete junction-facing coding ends in the reference. These are caller-supplied net outcome probabilities, not calibrated hairpin-opening probabilities. See original coding-end processing for VDJ orientation, provenance, and productive-sampling semantics.

The existing cartridge P model below remains available when process_coding_ends() is not configured.

Cartridges can also author per-end P-nucleotide length distributions via ReferenceEmpiricalModels.p_nucleotide_lengths, keyed by junction side label ("V_3", "D_5", "D_3", "J_5"). P-bases are templated palindromic complements of the source allele's post-trim coding flank — only the per-end length is sampled; the bytes themselves derive deterministically from (allele, trim, orientation, length) via complement_base. The four ends interleave with NP1 / NP2 in pool order between V/D and D/J coding bases:

from GenAIRR.reference_models import (
    EmpiricalDistributionSpec,
    ReferenceEmpiricalModels,
)

cfg.reference_models = ReferenceEmpiricalModels(
    p_nucleotide_lengths={
        "V_3": EmpiricalDistributionSpec([(0, 0.8), (1, 0.2)]),
        "J_5": EmpiricalDistributionSpec([(0, 0.9), (1, 0.1)]),
    }
)

Empty dict (the bundled-cartridge default) means the pipeline omits every P-addition pass — byte-identical to the pre-slice baseline. VJ cartridges reject D-end keys at the validation layer (no D segment to extend); VDJ cartridges accept all four ends. Each end's length distribution folds into the plan signature, so replay against a different P-length distribution fails the signature gate before any choice is consumed.

P-nucleotide provenance is surfaced on every AIRR record via four int fields:

p_v_3_length, p_d_5_length, p_d_3_length, p_j_5_length

…each counting the number of templated P-bytes emitted at that V(D)J coding-end junction side. Zero on records from cartridges without a typed P-plane. The validator's PLengthMismatch issue kind catches downstream record tampering; cache parity is unaffected (P-bytes don't move structural region boundaries).

Legacy p_nucleotide_length_probs is NOT auto-lifted yet. The bundled cartridges still carry the orphan dict; the manifest reports legacy_fallback=False. Per-base P strings (p_v_3, …) and an aggregate n_p_nucleotides are still not exported, but each NP region now reports how its bases arose: np1_n_tdt / np1_n_p / np1_n_inserted (and the np2_* triple) partition np1_length into TdT N additions, templated P bases and bases an indel inserted. See audit-docs/p_nucleotide_design.md for the full design.

Build a cartridge from FASTA

Beyond the 106 bundled cartridges, the 4-path manual construction surfaces, and the typed ReferenceEmpiricalModels / ReferenceRulesSpec authoring layer, GenAIRR ships an audit-trail builder for new cartridges from raw FASTA. ReferenceCartridgeBuilder stages FASTA parsing, identity attribution, V-subregion derivation, and rules/models attachment into a single fluent chain; every stage writes a structured entry to the build report so downstream consumers can audit how the cartridge was constructed.

import GenAIRR as ga

builder = (
    ga.ReferenceCartridgeBuilder
    .from_fasta(
        v_fasta="v.fa",
        d_fasta="d.fa",
        j_fasta="j.fa",
        chain_type="BCR_HEAVY",
    )
    .infer_identity(
        species="human",
        locus="IGH",
        reference_set="custom",
        name="my_igh",
    )
    .infer_v_subregions()
)

cfg = builder.build()
print(cfg.cartridge_manifest())
print(builder.report().to_dict())

.from_fasta(...) accepts file paths, open text files, or raw FASTA strings (string with \n and leading >). Duplicate allele names, empty sequences, and constructor failures land in report.rejected with structured {stage, segment, allele_name, reason} entries rather than corrupting the cartridge silently. .infer_v_subregions() reuses the same compute_v_region_boundaries helper the bundled-cartridge bridge uses, so a builder-produced cartridge with IMGT-gapped V FASTA gets the five canonical subregion intervals for free.

.build() finalises the cartridge: it stamps the canonical checksum onto schema_sha256, attaches the build_report as a typed CartridgeBuildReport dataclass, runs verify_integrity(), and returns a plain DataConfig you can drop straight into ga.Experiment.on(cfg). The build report's .to_dict() round-trips through json.dumps so CI artifacts can capture the full provenance.

Estimate allele usage from observed rearrangements. The first data-derived estimator landed: ReferenceCartridgeBuilder.estimate_allele_usage(rearrangements, *, min_count=1.0, ambiguous="fractional", replace=True) reads AIRR v_call / d_call / j_call columns from a list[dict], AIRR-C TSV path, or open text handle, and writes per-segment allele weights into the typed ReferenceEmpiricalModels.allele_usage plane.

builder = (
    ga.ReferenceCartridgeBuilder
    .from_fasta(v_fasta="v.fa", d_fasta="d.fa", j_fasta="j.fa", chain_type="BCR_HEAVY")
    .infer_identity(species="human", locus="IGH", reference_set="custom", name="my_igh")
    .estimate_allele_usage(rearrangements, ambiguous="fractional")
)
cfg = builder.build()

Recombination on the resulting cartridge uses the estimated weights by default — ga.Experiment.on(cfg).recombine() now samples alleles according to the typed plane unless an explicit recombine(v_allele_weights=..., d_allele_weights=..., j_allele_weights=...) kwarg is passed (kwarg > cartridge plane > uniform). Ambiguity policy is "fractional" (splits credit across tie-set entries), "truth_first" (uses the first call only), or "reject" (drops ambiguous rows into report.rejected). Unknown alleles, missing-D-on-VDJ rows, and below-min_count drops all surface as structured entries on the build report. Per-segment weights are normalised to sum to 1.0; the manifest["models"]["allele_usage"] block surfaces available / nonempty_segments / legacy_gene_use_dict_present / in_plan_signature=False (inherited soft gap from v_allele_weights). See audit-docs/allele_usage_estimation_design.md.

The other estimators have landed since. estimate_trim_distributions, estimate_np_length_distributions, estimate_np_base_model, estimate_p_nucleotide_lengths and estimate_genotype_priors all ship alongside estimate_allele_usage and append their own stage entries to the build report; set_genotype_priors attaches a population model you already have. Only estimate_shm_rates is still deferred. See audit-docs/reference_cartridge_authoring_audit.md for the full design and Estimate models from data for the per-estimator walkthrough.

Key Features

• Rust simulation kernel — persistent IR with full revision history, addressed-trace introspection, cargo test-grade unit coverage.
• Constraint-aware sampling — contracts prune candidate distributions at sample time so productive sequences come out of the engine by construction; no retry loops.
• Strict-mode opt-in — surface unsatisfiable plans as StrictSamplingError instead of silently relaxing the bundle.
• Deterministic seeds — one seed and one determinism_version reproduce every byte of the pool and every entry of the trace; the goldens check that on Linux, macOS and Windows across Python 3.10–3.13.
• Full revision history — outcome.revision(i) exposes the IR after each pass for fine-grained debugging.
• Addressed trace — every random draw is keyed by a hierarchical string ("np.np1.bases[3]") and survives end-to-end into the returned Outcome.
• 23 species, 106 configs — built-in IMGT + OGRDB reference pickles ship with the wheel.
• Zero mandatory Python dependencies — one wheel, everything in the box.

Optional Extras

pip install GenAIRR[all]          # numpy, scipy, graphviz, tqdm, fastmcp
pip install GenAIRR[dataconfig]   # numpy + scipy (custom DataConfig analysis)
pip install GenAIRR[viz]          # graphviz
pip install GenAIRR[mcp]          # fastmcp (for the MCP server, see next section)

MCP server — drive GenAIRR from an LLM agent

GenAIRR ships an MCP server that exposes 14 tools an LLM agent (Claude, Cursor, etc.) can call to discover configs, simulate repertoires, validate AIRR records, and replay specific seeds — all without writing Python. Install the extra, then point your MCP client at python -m GenAIRR.mcp_server:

pip install GenAIRR[mcp]

Config snippets

Claude Code — .mcp.json in the project root (or ~/.claude/mcp.json globally):

{
  "mcpServers": {
    "genairr": {
      "type": "stdio",
      "command": "python",
      "args": ["-m", "GenAIRR.mcp_server"]
    }
  }
}

Claude Desktop — ~/Library/Application Support/Claude/claude_desktop_config.json (macOS) or %APPDATA%\Claude\claude_desktop_config.json (Windows):

{
  "mcpServers": {
    "genairr": {
      "command": "/path/to/venv/bin/python",
      "args": ["-m", "GenAIRR.mcp_server"]
    }
  }
}

Cursor — ~/.cursor/mcp.json or per-project:

{
  "mcpServers": {
    "genairr": {
      "command": "python",
      "args": ["-m", "GenAIRR.mcp_server"]
    }
  }
}

Use the full path to the venv's python if GenAIRR isn't installed in the system interpreter — the MCP server inherits the launching process's Python environment.

What you get

After reloading the client, the agent has 14 tools available under the genairr namespace:

Every tool returns a uniform {ok, tool, elapsed_ms, result | error} envelope; failures carry a stable error-code token (config_not_found, allele_not_found, invalid_preset, invalid_parameter, malformed_record, seed_replay_mismatch) the agent can branch on.

A quick smoke-test prompt to verify the install: "List the available GenAIRR configs, then simulate 100 productive human heavy-chain sequences with moderate SHM and summarise the V-gene usage." — the agent should chain list_configs → simulate_repertoire(config="human_igh", n=100, productive_only=true, mutation_model="s5f", mutation_count_min=5, mutation_count_max=15) and read v_usage_top from the result.

Documentation

The full documentation site is at mutejester.github.io/GenAIRR. Useful starting points:

• Getting started — Installation · Quick start · Your first AIRR record · Export the results
• Concepts — AIRR Record · Reference Cartridge
• Guides — The Experiment builder · Productive sampling · Clonal lineage trees · Clonal repertoires · Genotypes · Paired-end reads and FASTQ · Corruption + sequencing artefacts · Targeted SHM rates · Trace, replay, reproducibility
• Validation — validate_records · Validation overview
• Reference — API + Configs + AIRR fields

Citing GenAIRR

If GenAIRR is useful in your research, please cite:

Konstantinovsky T, Peres A, Polak P, Yaari G. An unbiased comparison of immunoglobulin sequence aligners. Briefings in Bioinformatics. 2024;25(6):bbae556. doi:10.1093/bib/bbae556

Contributing

Contributions are welcome. See CONTRIBUTING.md for development setup and guidelines. If you're adding a pass, a contract, or any code that mutates the simulation IR:

• Read first: audit-docs/engine_architecture.md — codifies the invariants (contracts constrain support, trace = choices, events = consequences, live-call refresh follows events) and lists the anti-patterns CI catches.
• Then copy from: audit-docs/adding_a_pass.md — minimal pass template, three required test patterns with passes::test_support helpers, and a crib sheet of which existing pass to model the new one on.

License

GPL-3.0. See LICENSE.