{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "561216fa",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "# Explore a TTDPSFLGRY interface\n",
    "\n",
    "Inspect observed atom/residue contacts, chemical types, per-residue potential contributions,\n",
    "structural alanine contrasts and frozen local contact enrichment. All scientific quantities\n",
    "come from TCRen CLI commands or public APIs; the notebook selects and displays their tables.\n",
    "\n",
    "**Input:** supply one alpha/beta TCR–pMHC PDB/mmCIF model whose observed peptide is\n",
    "`TTDPSFLGRY`. This notebook does not distribute that structure or assert that a public\n",
    "TTD coordinate download is available. Use only coordinates you are authorized to access.\n",
    "No experimental binding label is inferred from the input.\n",
    "\n",
    "Install `tcren[notebooks]`, choose that environment's kernel, and prepare references with\n",
    "`tcren build-mhc-ref` and ARDA's documented setup. Reference acquisition/index construction\n",
    "is setup work; no query-cohort fitting or control generation occurs here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f7061426",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Report the kernel and scientific package versions before importing numerical APIs.\n",
    "import os, sys, json, tempfile, subprocess, time\n",
    "from pathlib import Path\n",
    "from importlib.metadata import version\n",
    "from importlib.resources import files\n",
    "import tcren\n",
    "import polars as pl\n",
    "from IPython.display import display\n",
    "print({\"python\": sys.version.split()[0], **{name: version(name) for name in (\"tcren\", \"polars\", \"numpy\", \"biopython\")}})"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "701538ab",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Parameters and a fresh run directory\n",
    "\n",
    "Set `TCREN_TTD_STRUCTURE` or edit `STRUCTURE` below. `THREADS=0` uses allocated CPUs;\n",
    "use a positive number for an explicit bound. Optional paths are existing recorded `screen`\n",
    "output and supplied generator confidence. The single-model representative score is **off**\n",
    "unless explicitly enabled; it is separate from default cohort scoring.\n",
    "\n",
    "Every run creates a new temporary output directory. Requested tables and their provenance\n",
    "remain there for inspection; a later run rebuilds them instead of reusing feature rows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d3fe1b81",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Configure one supplied structure; optional inputs never trigger a cohort benchmark.\n",
    "STRUCTURE = Path(os.environ.get(\"TCREN_TTD_STRUCTURE\", \"model.pdb\")).expanduser().resolve()\n",
    "THREADS = int(os.environ.get(\"TCREN_TTD_THREADS\", \"0\"))\n",
    "COHORT_SCORES = os.environ.get(\"TCREN_TTD_COHORT_SCORES\", \"\")\n",
    "SINGLE_CONFIDENCE = os.environ.get(\"TCREN_TTD_CONFIDENCE\", \"\")\n",
    "RUN_REPRESENTATIVE = os.environ.get(\"TCREN_TTD_REPRESENTATIVE\", \"0\") == \"1\"\n",
    "assert STRUCTURE.is_file(), \"Supply TCREN_TTD_STRUCTURE or edit STRUCTURE to your authorized TTD model.\"\n",
    "assert THREADS >= 0, \"THREADS must be zero or positive.\"\n",
    "work = Path(tempfile.mkdtemp(prefix=\"tcren-ttd-\"))\n",
    "stages = []\n",
    "def cli(label, *arguments):\n",
    "    started = time.perf_counter()\n",
    "    with (work / (label + \".log\")).open(\"w\") as log:\n",
    "        subprocess.run([sys.executable, \"-m\", \"tcren.cli\", *map(str, arguments)],\n",
    "                       check=True, stdout=log, stderr=subprocess.STDOUT, timeout=180)\n",
    "    stages.append({\"stage\": label, \"wall_s\": time.perf_counter() - started})\n",
    "print(\"Fresh local output directory:\", work)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dd17a225",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Export the observed interface\n",
    "\n",
    "Each command receives the complete query set (one model here). The installed CLI owns\n",
    "batched annotation, heavy-atom selection and contact definitions. Atomic contacts retain\n",
    "all six CDR loops; the inventory also retains observed noncontacting atoms. The fixed\n",
    "5 Å cutoff matches the bundled enrichment reference. No threshold sweep is performed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "25459e85",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Export coordinates-derived tables with explicit potential and interface choices.\n",
    "for quantity, filename in ((\"atom-inventory\", \"atoms.parquet\"), (\"contact-atoms\", \"contacts.parquet\"),\n",
    "                           (\"residue\", \"residue.tsv\"), (\"position\", \"position.tsv\")):\n",
    "    arguments = [\"export\", \"-s\", STRUCTURE, \"--quantity\", quantity, \"--interface\", \"tcr_peptide\",\n",
    "                 \"--cutoff\", \"5\", \"-t\", THREADS, \"-o\", work / filename]\n",
    "    if quantity in (\"residue\", \"position\"):\n",
    "        arguments += [\"--potential\", \"tcren2\"]\n",
    "    cli(quantity, *arguments)\n",
    "    ledger = pl.read_csv(str(work / filename) + \".accounting.tsv\", separator=\"\\t\")\n",
    "    assert ledger.height == 1 and ledger[\"status\"].to_list() == [\"ok\"], \"Inspect the local accounting and log before proceeding.\"\n",
    "atoms = pl.read_parquet(work / \"atoms.parquet\")\n",
    "contacts = pl.read_parquet(work / \"contacts.parquet\")\n",
    "identity = atoms.filter(pl.col(\"chain.type\") == \"PEPTIDE\").select(\n",
    "    \"chain.id\", \"residue.index\", \"residue.aa\").unique().sort(\"residue.index\")\n",
    "assert identity[\"chain.id\"].n_unique() == 1\n",
    "assert \"\".join(identity[\"residue.aa\"]) == \"TTDPSFLGRY\", \"The supplied observed peptide is not TTDPSFLGRY.\"\n",
    "MODEL_ID = atoms[\"pdb.id\"].unique().item()\n",
    "display(identity)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fd7eaf2d",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Atom pairs, residue pairs and chemistry\n",
    "\n",
    "These are different observation counts. The CLI summary separates atom and residue pairs,\n",
    "with exact structure support. Chemical flags are multi-label distance/type heuristics;\n",
    "a hydrogen-bond-compatible flag is not explicit-hydrogen/angle validation. The closest\n",
    "atom pair alone need not reveal every chemical interaction between two residues.\n",
    "Tables below can be filtered or sorted in the notebook without recalculating the score."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7ca305a6",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Read the CLI's existing atom/residue chemical summary and inspect physical atom witnesses.\n",
    "summary = pl.read_csv(str(work / \"contacts.parquet\") + \".summary.tsv\", separator=\"\\t\")\n",
    "display(summary)\n",
    "flags = [column for column in contacts.columns if column.startswith(\"is_\")]\n",
    "display(contacts.select(\"chain.id.from\", \"region.type.from\", \"residue.index.from\", \"residue.aa.from\", \"atom.from\",\n",
    "                        \"chain.id.to\", \"residue.index.to\", \"residue.aa.to\", \"atom.to\", \"dist\", *flags).sort(\"dist\").head(40))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "301499c3",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Per-residue potential contributions\n",
    "\n",
    "`phi` is the contact-potential attribution for each contacting residue; `n_contacts` is\n",
    "its residue-pair contact count. More negative `phi` is more favourable on this potential's\n",
    "scale. A contact is attributed to **both** endpoints, so summing both sides doubles the\n",
    "interface total. Attribution is not a causal residue effect or measured binding energy.\n",
    "The position export supplies the peptide-side decomposition and physical residue keys."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c215dfe9",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Display the CLI's potential decomposition, without recomputing or renormalizing it.\n",
    "residue = pl.read_csv(work / \"residue.tsv\", separator=\"\\t\")\n",
    "position = pl.read_csv(work / \"position.tsv\", separator=\"\\t\")\n",
    "display(residue.sort(\"phi\"))\n",
    "display(position.sort(\"peptide.pos\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "33db443b",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Structural single-site alanine contrasts\n",
    "\n",
    "The existing `ddg` CLI truncates one selected side chain at a time and rebuilds the\n",
    "contact map. This is distinct from fixed-map virtual substitution and from the attribution\n",
    "above. The exported `ddG` is **native minus alanine-mutant** statistical-potential value:\n",
    "positive means the alanine mutant has the lower value. It is dimensionless, not kcal/mol.\n",
    "Peptide `pos` is zero-based; receptor rows also retain chain, loop and residue index.\n",
    "These model-based contrasts do not establish an experimental or causal binding effect."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "66c37df5",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Use the library's structural alanine scan for peptide and contacted CDR residues.\n",
    "cli(\"alanine\", \"ddg\", \"-s\", STRUCTURE, \"--native\", \"TTDPSFLGRY\", \"--alanine-scan\", \"--side\", \"both\",\n",
    "    \"--regions\", \"cdr\", \"--interface\", \"tcr_peptide\", \"--potential\", \"tcren2\", \"--cutoff\", \"5\",\n",
    "    \"-t\", THREADS, \"-o\", work / \"alanine.csv\")\n",
    "alanine = pl.read_csv(work / \"alanine.csv\")\n",
    "display(alanine.sort(\"side\", \"pos\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f0b49fcf",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Frozen local enrichment across all six CDR loops\n",
    "\n",
    "`residue_counts` retains zero-contact CDR1/2/3 residues in both chains plus peptide\n",
    "residues. CDR3 uses the library's anchor-excluded annotation. Completeness and observed\n",
    "heavy-atom availability remain explicit; missing atoms are not reconstructed.\n",
    "\n",
    "`score_residues` reads the packaged mock histograms and reports local cycle counts,\n",
    "reference-residue support and inclusive empirical tails. `bundled_enrichment` aggregates\n",
    "through the existing Fisher implementation and calibrates against the packaged whole mock\n",
    "graphs, returning `enrichment`, `null_models` and local-test support. Local tails and whole\n",
    "graph enrichment are descriptive structural evidence, not binding probabilities. This is\n",
    "not a new control fit, a sequence repertoire motif search or a single-loop substitute."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "39c0c2c4",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Call the public frozen-reference APIs; do not reproduce their formulas in the notebook.\n",
    "from tcren.contacts.residue_surprise import residue_counts, score_residues\n",
    "from tcren.contacts.motif import bundled_enrichment\n",
    "local = residue_counts(atoms, contacts, atom_scope=\"all\")\n",
    "reference = json.loads(files(\"tcren.data\").joinpath(\"motif_residue_reference.json\").read_text())\n",
    "local_support = score_residues(local, reference, \"cycles\")\n",
    "enrichment = bundled_enrichment(local)\n",
    "local_support.write_parquet(work / \"local-enrichment.parquet\")\n",
    "enrichment.write_csv(work / \"enrichment.tsv\", separator=\"\\t\")\n",
    "display(local_support.select(\"chain\", \"region\", \"residue_index\", \"aa\", \"available\", \"complete\",\n",
    "                              \"occupied\", \"pairs\", \"partners\", \"cycles\", \"reference_residues\",\n",
    "                              \"reference_mass\", \"tail_probability\", \"bits\").sort(\"chain\", \"region\", \"residue_index\"))\n",
    "display(enrichment)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4b06a1a9",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Read a recorded cohort score, or explicitly select a singleton representative\n",
    "\n",
    "The **default** TCR score uses the complete unlabelled fixed-pMHC cohort's mean pLDDT,\n",
    "with a 30-scorable-model floor. To calculate that cohort elsewhere, supply its complete\n",
    "manifest and generator confidence (unique `complex.id`, `iptm` on [0,1], `plddt` on [0,100]):\n",
    "\n",
    "```sh\n",
    "tcren screen -s ttd-cohort.txt --task tcr --confidence confidence.tsv -t 0 -o cohort-scores.tsv\n",
    "```\n",
    "\n",
    "The cell below only reads an optional existing `screen` table. It does not run a cohort\n",
    "benchmark. Input IDs must match this model; inspect `score_unavailable`, no-contact flags,\n",
    "components and weights alongside the rank score. Binding labels do not enter this inspection.\n",
    "\n",
    "For a singleton, `nearest-plddt` is an **explicit alternative**, enabled only by\n",
    "`TCREN_TTD_REPRESENTATIVE=1` plus `TCREN_TTD_CONFIDENCE`. It selects a frozen representative\n",
    "context, not the independently measured true cohort mean. It is not the default cohort\n",
    "score or a calibrated binding probability. Missing pLDDT withholds this context."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3ea46585",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Optional score context: read supplied output; run a singleton alternative only on explicit opt-in.\n",
    "if COHORT_SCORES:\n",
    "    cohort = pl.read_csv(Path(COHORT_SCORES).expanduser(), separator=\"\\t\")\n",
    "    selected = cohort.filter(pl.col(\"complex.id\") == MODEL_ID)\n",
    "    assert selected.height == 1, \"Expected exactly one matching model in the recorded cohort output.\"\n",
    "    display(selected)\n",
    "if RUN_REPRESENTATIVE:\n",
    "    assert SINGLE_CONFIDENCE, \"Supply generator confidence explicitly for the singleton alternative.\"\n",
    "    cli(\"representative\", \"screen\", \"-s\", STRUCTURE, \"--task\", \"tcr\", \"--confidence-context\", \"nearest-plddt\",\n",
    "        \"--min-cohort\", \"1\", \"--confidence\", Path(SINGLE_CONFIDENCE).expanduser(), \"-t\", THREADS,\n",
    "        \"-o\", work / \"representative.tsv\")\n",
    "    display(pl.read_csv(work / \"representative.tsv\", separator=\"\\t\"))\n",
    "print(\"No cohort scoring is launched by this notebook.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9149ba72",
   "metadata": {
    "language": "markdown"
   },
   "source": [
    "## Run diagnostics and reuse\n",
    "\n",
    "The tables describe this supplied pose, not whether its receptor binds experimentally.\n",
    "Save requested final exports together with their provenance sidecars, accounting and\n",
    "reference hashes. Provenance can contain local paths or private identifiers; review it\n",
    "before sharing. No input coordinates, reference source maps or executed outputs are\n",
    "embedded in the distributed notebook. The documentation renders it without execution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "02d1c4ba",
   "metadata": {
    "language": "python"
   },
   "outputs": [],
   "source": [
    "# Summarize command runtimes and requested final table sizes, without evaluating a benchmark.\n",
    "display(pl.DataFrame(stages))\n",
    "display(pl.DataFrame({\"table\": [\"observed atoms\", \"atom contacts\", \"residue contributions\", \"alanine contrasts\", \"local residue evidence\"],\n",
    "                      \"rows\": [atoms.height, contacts.height, residue.height, alanine.height, local_support.height]}))\n",
    "print(\"Inspect the fresh local output directory above for complete exports and provenance.\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python (TCRen)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}