Epitope (pMHC) search#

TCR cross-reactivity is driven by a shared central, TCR-facing motif, not the MHC anchors (Dolton et al., Cell 2023: one HLA-A*02:01 TCR recognizes EAAGIGILTV / LLLGIGILVL / NLSALGIFST via x-x-x-A/G-I/L-G-I-x-x-x). seqtree.pmhc therefore masks anchor positions and searches anchor-masked TCR-facing k-mers over the C++ KmerIndex (seed-and-gather).

MHC class I (global, length 8–11): anchors at P2 and the C-terminus (PΩ) are masked; the central bulge + P1 drive homology across lengths.
MHC class II (local, trim/shift): the 9-mer core register is handled register-agnostically — shared central k-mers match regardless of flanks or register.

Anchor positions are parametrized (presets per class, overridable) and masking is realized by a zero-cost PositionalMatrix (weight 0 = masked/free; up-weight a hotspot for a graded TCR-facing distance).

Homology search#

from seqtree.pmhc import PMHCStore

store = PMHCStore.from_pmhc("pmhc_full.tsv.gz")        # IEDB-style epitope table
# or PMHCStore.from_records([{ "epitope": ..., "mhc": ..., "mhc_class": "MHCI"}, ...])

for h in store.search_homologs("EAAGIGILTV", "mhc1", mhc="HLA-A*02:01", max_subs=2):
    print(h.epitope, h.mhc, h.shared_kmers, h.score)

mhc= restricts to a single allele (cross-allele groove similarity is out of scope). The Dolton trio is a built-in positive control: each epitope recovers the other two as mutual homologs.

Neoantigen mimic discovery + E-values#

Significance is presentation-aware — calibrated per allele against an anchor-masked presented background, so anchor sharing (presentation, not recognition) does not inflate hits:

from seqtree.pmhc import find_mimics

res = find_mimics(neoantigen, self_set=human_presented, bacterial_sets={"ecoli": ecoli_presented})
for source, r in res.items():
    print(source, "E=%.3g" % r["E"], "p=%.3g" % r["p_enrichment"],
          [h.epitope for h in r["hits"][:5]])

self_set / bacterial_sets are presented peptides for the relevant allele (from pmhc_data or any caller-supplied proteome slice). The math — per-allele null, the E=(N/M)·n_control estimator, the cross-reactivity distance, and the limitation that distinct alleles are distinct nulls — is in appendix/evalue.tex §”Epitopes”.

Reverse problem: peptide → allele#

Flipping the weights to score the anchors ranks the likely presenting alleles (a lightweight presentation prior, not a trained predictor):

store.assign_allele("KLEEEEEEV", "mhc1")   # -> [(allele, score, n_match, n_allele), ...]

A neighbour-voting variant (bench/bench_mhc_guess.py) attaches an aggregate E-value / confidence to the guess and shows that length-matched random peptides get a much higher E-value, so they are rejected as noise by a threshold — see Benchmarks. This works on anchor (presentation) features only; TCR-facing homology carries no allele information.

Filtering non-binders#

The per-allele E-value also doubles as a non-binder filter, at two granularities:

Binds no MHC. A peptide whose best per-allele E-value stays high across the whole panel (equivalently, assign_allele returns no allele scoring above the marginal background) is flagged as a non-binder and dropped. The real-vs-random arm of bench/bench_mhc_guess.py measures exactly this: it is the noise-rejection AUROC in Benchmarks.
Does not bind a specific MHC. For a candidate allele a, a per-allele E-value E_a above a threshold α means the peptide’s anchors are atypical of a’s motif, so a is excluded as a restriction element even if the peptide binds something else.

Note

MHC class II is promiscuous. A single class-II peptide is routinely presented by many alleles (shared, degenerate P1/P4/P6/P9 pockets and an open-ended groove), so allele assignment is genuinely multi-label: the “true” set for a peptide is several alleles, not one. The benchmark treats it as such (positives = every observed restricting allele), which is why class-II precision–recall baselines sit higher than class I. A non-binder filter for class II should therefore test “binds none of the panel” rather than “binds exactly one”.

Note

This module is a reference implementation and benchmark for the E-value methodology. The production allele-guessing and non-binder filtering — with the optimized index, tuned thresholds, and ROC/PR calibration — live in the downstream vdjmatch (TCR–antigen specificity) and mhcmatch (peptide–MHC) packages.