tcren — structure-based prediction of TCR–epitope recognition

TCRen predicts which epitopes a T-cell receptor recognises from a single TCR–peptide–MHC structure (experimental or modelled). It extracts the TCR–peptide contact map and scores every candidate peptide with a residue-level statistical potential derived from contact preferences in TCR:pMHC crystal structures — answering not "what fancy complex can a model draw?" but "is this binding physically plausible?".

This is a documented, tested, CLI-driven Python library. TCR chains are annotated with the sibling arda; MHC chains are mapped and the groove partitioned against a curated reference; structures are oriented into one canonical frame; and the original contact maps, potential, and scores are reproduced numerically (validated against committed oracles to floating-point precision).

While the original tcren focused on TCR:peptide contacts, the new version brings in features to score TCR:MHC and peptide:MHC interactions, required to get full picture of TCR:pMHC binding mechanics and estimate ddG values.

What it does

From one TCR–peptide–MHC structure (crystal or model), each task is one command or one call:

Scope — ranking, not affinity. TCRen ranks peptide/TCR specificity for a given receptor (and the ddg matrix is a fast triage, not a free energy). It is not an affinity model: on the ATLAS SPR benchmark neither the raw contact energy nor its poly-alanine difference predicts Kd/ΔG/koff/kon (|ρ|≤0.3). The one affinity-adjacent quantity a structure predicts is the off-rate koff, via interface mechanics (tcren mechanics) — not the contact sum.

Install

pip install tcren          # from PyPI — binary wheels ship the C++ extension; pulls in arda-mapper

For development (a repo-local .venv via uv, an editable install, and the reference data fetched into data/):

bash setup.sh                    # uv venv + editable install + arda + fetch data/ (no conda)
source .venv/bin/activate

setup.sh needs only uv and a C++ compiler (macOS: xcode-select --install); it never touches conda. Pass --tests to run the fast suite after install.

tcren ships five small pybind11/C++ extensions, built on install by scikit-build-core (which fetches cmake+ninja automatically): tcren._align (MHC-pseudosequence fitting alignment; a Biopython fallback runs if unbuilt), tcren._refine (DOPE atom-level Monte-Carlo peptide refinement), tcren._relax (DOPE interface energy for tcren energy / ΔΔG), tcren._fold (CCD loop closure) and tcren._geom (interface geometry for tcren binder). TCR annotation is provided by arda, a runtime dependency published to PyPI as arda-mapper (it imports as arda); uv/setup.sh pull it automatically, and from arda-mapper >= 2.5.7 it auto-fetches both its own reference and a static mmseqs2 binary on first use — so no conda/bioconda and no ARDA_HOME to set (override the binary with $ARDA_MMSEQS). setup.sh also runs tcren fetch-data to populate data/ with the reference structure sets (Native2026, Canonical2026) used by orient/superimpose (set TCREN_NO_FETCH=1 to skip).

Command line

# Score structures: the three interface contact energies (TCRen for TCR↔peptide, MJ for
# TCR↔MHC and peptide↔MHC) and their total Φ. One row per structure.
tcren scoring -s complex.pdb.gz -o scores.csv

# Inputs: a file, a directory, a .tar.gz, a quoted glob, a .txt manifest (one path per line),
# a comma-separated list, or a repeated -s. Mix freely.
tcren scoring -s a.pdb.gz -s b.pdb.gz -o scores.csv
tcren scoring -s 'models/*.pdb.gz' -o scores.csv
tcren scoring -s models/ --delta --geometry -t 8 -o scores.csv   # a directory, 8 workers
tcren scoring -s models.txt -o scores.csv

# --delta adds the poly-alanine reference ΔΦ per interface (ΔΦ_TCR:MHC is identically 0).
# Use ΔΦ, not Φ, when each candidate carries its OWN generated pose: raw Φ then partly reads
# the pose the predictor chose rather than the peptide.
tcren scoring -s 'models/*.pdb.gz' --delta -o scores.csv

# --geometry adds the interface descriptors and Q, the directional decorrelated
# interface-quality score (native-crystal calibrated, so it is defined for a single structure).
tcren scoring -s complex.pdb.gz --delta --geometry -o scores.csv

# Configurable per-interface potential: swap a bundled name (tcren|mj|keskin), a CSV, or
# None for any interface; default reproduces the built-in per-interface families exactly.
tcren scoring -s complex.pdb -o scores.csv --tcr-mhc-potential keskin

# Opt-in TCR framework regions: --regions {all,cdr,cdr+fr} chooses which TCR regions
# contribute on the TCR side (cdr = CDR1-3 only; cdr+fr adds FR1-3; all = unfiltered, default).
tcren score -s complex.pdb -c candidates.txt -o ranked.csv --regions cdr+fr

# Percentile-rank the native (or candidate) peptide's TCRen energy against a random pMHC
# background — small rank_pct = the peptide scores among the best binders.
tcren rank -s complex.pdb -o rank.csv

# Fast ΔΔG of peptide point mutations (virtual-matrix path: no atoms move, no re-docking).
# Requires --native (the peptide) and exactly one mode: --alanine-scan or --mutant.
# ddG = E(native) - E(mutant), and lower energy binds better, so POSITIVE = stabilising.
tcren ddg -s complex.pdb --native EPITOPE --alanine-scan -o ddg.csv

# Predict a combinatorial-peptide-library (CPL) response matrix from ONE template TCR:pMHC
# structure: every peptide position x all 20 residues, threaded on the template's own contact map.
# Every cell sums BOTH peptide-bearing interfaces (TCRen over TCR:peptide + Miyazawa-Jernigan over
# peptide:MHC), because the assay reads activation, which needs presentation as well as engagement.
tcren cpl -s complex.pdb -o cpl_matrix.csv
# Two reference states, both emitted, and a cell means nothing except against one of them:
#   effect_equimolar  vs the 1/20 mixture  -> the CPL background; compare against a measured matrix
#   effect_wild_type  vs the template residue -> the mutation-scan / neoantigen question
# Positive is favourable on both. Three narrower questions off the same matrix:
tcren cpl -s complex.pdb --position 5                  # every substitution at position 5, best first
tcren cpl -s complex.pdb --position 5 --mutation W     # just that one cell
tcren cpl -s complex.pdb --position 5 --to-mixture     # cost of giving position 5 up to the mixture

# Binder vs non-binder from AF-orthogonal interface geometry + the CDR1/2-vs-CDR3a TCRen term —
# ranks candidate TCRs against a fixed pMHC, on par with AlphaFold/TCRmodel2 confidence with no
# external tool (raw-label macro AUC ~0.80 vs AF ipTM 0.79). PREFER the fit-free Q = tcren.cohort.
# q_score, which matches this and generalises across cohorts where the fitted p_bind does not; with
# ipTM, z(ipTM)+z(Q) is the fit-free synergy (macro 0.83 vs 0.79). `tcren binder` emits the fitted
# p_bind (retained for reproducibility); `tcren recognize --scores` adds q_bind + s_strain.
tcren binder -s complex.pdb -o binder.csv

# One TSV per structure: every interface descriptor (geometry + energies) + joint P(real).
tcren recognize -s my_pdbs/ -o recognize.tsv          # descriptors + p_real + p_real_bn, one row/PDB

# End-to-end candidate-epitope scoring from a structure
tcren score -s complex.pdb -c candidates.txt -o ranked.csv

# Wrong-TCR decoys: keep each ORIENTED complex's pMHC, graft on 10 other complexes' TCRs (within
# MHC class, no real pairing). Real-vs-decoy trains a label-free TCR-recognition classifier.
tcren orient -s natives/ -o oriented/          # inputs must share the canonical MHC frame
tcren shuffle -s oriented/ -o shuffled/ --n 10

# Substitute a peptide and refine its pose (knowledge-based MC scored by the DOPE atom-level
# statistical potential — independent of the TCRen/MJ scoring potentials, restrained to the input).
# Not physics relaxation — use Rosetta FlexPepDock for that.
tcren refine -s complex.pdb -o refined/ --substitute KQWLVWLFL

# Structures: any of .pdb / .cif / .pdb.gz / .cif.gz, a directory, or a .tar.gz batch
tcren contacts -s batch.tar.gz -o contacts.csv --interface tcr_peptide

# Per-residue markup: TCR (CDR/FR) + MHC groove (helix/floor) + peptide in one table.
# --regions all|tcr|mhc|peptide filters; --pseudo also marks NetMHCpan groove residues (MPS).
tcren annotate -s complex.cif.gz -o markup.csv --regions mhc --pseudo

# Superimpose structure(s) onto the canonical frame, by MHC, against the canonical database
# (data/Canonical2026, fetched at install). Detects MHC class + species and averages the
# superposition over every database structure of that class/species. Chains -> A=Vα B=Vβ
# C=peptide D=MHCα E=MHCβ/β2m. -s takes a file / directory / .tar.gz / glob; -o is a directory,
# or a single structure file (one input) whose extension must match --mmCIF/--compress; -t threads.
tcren superimpose -s complex.pdb -o oriented.pdb           # single file
tcren superimpose -s 'data/*.pdb' -o oriented/ -t 8        # glob -> directory, threaded

# Build a canonical database from native complexes (how Canonical2026 is produced). Annotation
# is one batched mmseqs call; -t threads only the structural alignment + write.
tcren orient -s data/Native2026 -o data/Canonical2026 -t 8

# Structure outputs are plain .pdb by default; add --mmCIF for .cif and --compress for .gz.
tcren superimpose -s complex.pdb -o oriented/ --mmCIF --compress   # -> oriented/<id>.cif.gz

# Fetch recent TCR-pMHC structures from RCSB -> data/pdb_recent (mmCIF .cif.gz, 5-chain validated)
tcren fetch-recent --discover --after 2024-01-01

# Build the MHC reference once (IMGT/HLA + mouse H-2; cached, not committed)
tcren build-mhc-ref

tcren info
tcren --install-completion        # shell tab-completion (bash/zsh)

tcren orient and tcren superimpose need the reference sets in data/ (Native2026, Canonical2026); setup.sh fetches them at install via tcren fetch-data (re-run it any time).

One table per structure: descriptors, energies & the joint recognizer

Give tcren recognize a list of complexes (a file, directory, .tar.gz, or glob) and it writes one TSV row per structure with the full interface descriptor set and the joint recognition probability P(real):

tcren recognize -s my_pdbs/ -o recognize.tsv               # 35 descriptors + p_real + p_real_bn
tcren recognize -s my_pdbs/ -o scored.tsv --scores         # + q_bind, s_strain (recommended) + p_bind, p_forced
tcren recognize -s my_pdbs/ -o feats.tsv --features-only   # descriptors only, skip the models

Where the joint model lives. p_real is the frozen recognizer we derive from real crystals vs wrong-TCR shuffled decoys: code in tcren.recognition (recognition_features → real_probability), coefficients shipped in src/tcren/data/shuffle_logistic.json.gz, and the full derivation (PyMC fit, encoding, ROC/PR, posterior forest) in the appendix appendix/logistic_stan/. Decoys come from tcren shuffle; the Gaussian-BN companion is appendix/shuffle_bn/.

(c) physics of the interaction. The koff proxies fold into the same table with --mechanics; only the mutation scan, which is per-residue rather than per-structure, needs its own command:

tcren recognize -s models/ --scores --mechanics -t 0 -o out.tsv   # every per-structure descriptor, one table
tcren ddg       -s complex.pdb -o ddg.csv     # per-residue alanine / neoantigen ΔΔF (fast virtual matrix)

--mechanics is how to ask for the stiffness tensor, steered rupture and coupling residues on a cohort. tcren mechanics still exists and gives the same numbers, but as a second command it repeats the parse and both mmseqs searches to return a second table — CSV, keyed pdb.id rather than complex.id — that then has to be joined. Inside recognize the structures are already annotated, so the flag costs only the mechanics arithmetic (12 crystals: 19.0 s → 19.5 s, against 22.5 s for the two commands).

(Per the affinity scope caveat above, structures predict the off-rate koff via the mechanics columns, not Kd/ΔG/kon.) From Python:

from tcren.recognition import recognition_features, real_probability
feats = recognition_features("complex.pdb")    # dict of the 35 descriptors (RECOGNITION_FEATURES)
p = real_probability(feats)                     # {"logistic": P(real), "bn": P(real)}

Library

from tcren import run_pipeline, parse_structure, import_structure, ContactMap, score_peptides
from tcren.annotation import classify_chains
from tcren.potential import tcren

# One call: annotate -> superimpose -> contacts -> per-interface energies + total
res = run_pipeline("complex.pdb")              # res.scores, res.markup, res.contacts, res.oriented
res = run_pipeline("complex.pdb", reference_aa="A")  # + delta_* : the poly-alanine ΔΦ per interface

# Oracle facade: one structure -> a bundle of ready-to-tabulate frames for the paper
# notebooks (scores, percentile rank, ΔΔG alanine scan, markup, contacts). Configurable
# per-interface potentials and TCR-region selection are forwarded to every milestone.
from tcren import summarize_structure
bundle = summarize_structure("complex.pdb", alanine=True)   # bundle["scores"], ["rank"], ["ddg"], …

# …or the individual steps:
s = parse_structure("complex.pdb.gz")          # also .cif/.cif.gz; import_structure trims the C-gene
classify_chains(s, organism="human")           # TRA/TRB via arda, peptide, MHC
cm = ContactMap.from_structure(s)              # 5 Å contacts + interface partitioning
ranked = score_peptides(cm, ["KQWLVWLFL", "RLLHPHHPL"], tcren())

CPL response matrices from one template structure

A positional-scanning combinatorial peptide library fixes position i to residue a and leaves every other position an equimolar 1/20 mixture, so a measured cell is an ensemble mean, R[i,a] = E[response | x_i = a]. tcren.cpl predicts that matrix from a single template complex — each of the twenty residues threaded through the template's own contact map, nothing re-docked, nothing fitted to any assay.

from tcren import (ContactMap, parse_structure, response_matrix,
                   mutation_effect, position_scan, equimolar_effect)
from tcren.annotation import classify_chains
from tcren.mhc import annotate_mhc

s = parse_structure("3HG1.pdb", pdb_id="3HG1")
classify_chains(s, organism="human")
annotate_mhc(s)                       # REQUIRED: without it peptide:MHC is empty and anchors zero out
rm = response_matrix(ContactMap.from_structure(s, cutoff=5.0))

rm.to_frame()                         # the whole matrix, one row per (position, amino acid) cell
position_scan(rm, 5)                  # every substitution at position 5, best first
mutation_effect(rm, 5, "W")           # one cell
equimolar_effect(rm, 5)               # cost of giving position 5 up to the 1/20 mixture

Every cell sums both peptide-bearing interfaces — TCRen over TCR:peptide plus Miyazawa–Jernigan over peptide:MHC — because the assay reads activation, which needs the peptide presented as well as the receptor engaged. A position the receptor never touches is an anchor; its TCR term is constant along the row, so the sum degrades to presentation alone rather than to a special case.

Two reference states, and a cell is meaningless except against one of them. A raw Φ carries a large per-position offset that says only how many contacts the position makes:

They differ by a per-position constant — how far the template's residue sits above its column mean. Positive is favourable on both, since lower energy is the better binder. Under "wild_type" the template's own cell is identically zero; under "equimolar" it is an ordinary measurement.

Batch inputs, gzip, archives

from tcren.structure import iter_structures
for pdb_id, structure in iter_structures("batch.tar.gz"):   # file | directory | .tar.gz
    classify_chains(structure, organism="human")
    ...

Canonical orientation, contacts, docking geometry

from tcren.mhc import annotate_mhc
from tcren.orient import canonicalize_structure, superimpose, docking_angles
from tcren.contacts import multi_contacts, ContactDefinition

annotate_mhc(s)
oriented, info = canonicalize_structure(s)     # frame: z=MHC→TCR, y=peptide, x=thin; chains A–E
oriented, info = superimpose(s)                # orient onto data/Canonical2026 by MHC (class+species ensemble)
layers = multi_contacts(s, ContactDefinition(d1=5, d2=8, d3=12))   # heavy-atom / Cβ / Cα
d = docking_angles(s)                          # crossing (~20–70° αβ) + incident angle

2D complementarity maps & region-pair contacts

from tcren.project2d import (project_structure, residue_markup_table, contacts_table,
                             region_pair_summary)
from tcren.viz import render_complementarity_map, view_pocket_cdr

proj = project_structure(s)                                   # canonical groove plane
svg  = render_complementarity_map(residue_markup_table(s, proj),
                                  contacts=contacts_table(s, threshold=5.0))
region_pair_summary(s, kind="closest")        # contacts per region pair + bond types (cb/ca too)
view_pocket_cdr(s).show()                      # interactive 3D pocket + CDR overlay (py3Dmol)

Modules

Data

Structures live in the Hugging Face dataset isalgo/tcren_structures, all gzipped:

tcren reads .pdb/.cif/.pdb.gz/.cif.gz and .tar.gz batches; an installed library lazily fetches the canonical reference structures from the Hub when orienting a new complex. The root data/ holds Native2026 (+ Canonical2026, gitignored, fetched on demand), PDB_date.tsv, orient_metadata.json, and TCRen_potential.csv — the current potential derived from the Native2026 set (use it with tcren score -p data/TCRen_potential.csv).

Notebooks

Runnable examples under notebooks/ (rendered in the docs):

• complementarity_map_2d — 2D interface maps, multiple structural + map views of 1ao7
• contact_thresholds_and_bondtypes — region-pair contact counts (closest/Cβ/Cα) + bond types
• canonical_frame_figures — canonical-frame QC across the Native2026 set
• pymol_canonical_figures — ray-traced PyMOL panels (overlay, groove, interface) by class/species
• mhc_pseudosequence_mps — NetMHCpan MHC pseudosequence (MPS) residues vs. peptide contacts
• example_gil_a02_rs_motif — GILGFVFTL/HLA-A*02 and the public CDR3β Arg–Ser motif
• natcompsci2022/ — full reproduction of the Nat Comput Sci 2022 analyses

Performance

Per-stage wall time (best of n) on a TCR-pMHC complex (1ao7), Apple M-series, single thread (RUN_BENCHMARK=1 pytest -k benchmark -s to reproduce the core stages):

Annotation is the only network/compute-heavy step and is always batched (one mmseqs2 search over all chains; mmseqs2 parallelises internally — never per-structure, never Python-threaded). Threads are used only for the embarrassingly-parallel, mmseqs-free stages (structural alignment, write, rendering): tcren orient -t N. Screening a peptide/TCR panel is embarrassingly parallel — references are annotated and oriented once, so the hot loop is just refine + contacts + score per complex.

Tests

pytest -m "not slow"          # unit + fast regression (the CI gate)
pytest                        # add the arda/mmseqs-backed regression tests
RUN_BENCHMARK=1 pytest -k benchmark -s

Methods appendix

The coordinate-level extensions — backbone-preserving peptide substitution and the potential-guided Monte-Carlo refinement kernel (energy function, the restraint-necessity argument, sampler, and citations) — are written up in the technical appendix appendix/tcren.tex (built with make -C appendix → appendix/tcren.pdf).

Citing

TCRen is free for academic and non-commercial use. If you use it, please cite our latest Nature Computational Science 2024 paper:

Karnaukhov VK, Shcherbinin DS, Chugunov AO, Chudakov DM, Efremov RG, Zvyagin IV, Shugay M. Structure-based prediction of T cell receptor recognition of unseen epitopes using TCRen. Nat Comput Sci. 2024 Jul;4(7):510-521. doi: 10.1038/s43588-024-00653-0. Epub 2024 Jul 10. PMID: 38987378.