{ "cells": [ { "cell_type": "markdown", "id": "892511b5", "metadata": {}, "source": [ "# Case Study 10 — Ligand Pose Sensitivity and Local Rescue\n", "\n", "[](https://colab.research.google.com/github/uw-ipd/tmol/blob/master/docs/tutorial/10_ligand_pose_sensitivity.ipynb)\n", "\n", "A bound ligand can look plausible while making poor local contacts. This case study creates controlled rigid-body decoys of one crystallographic ligand, scores every pose in one batch, and locally minimizes three diagnostic states.\n", "\n", "## Biological question\n", "\n", "Does the one-complex TMol interaction score distinguish the deposited ligand placement from controlled displacements, and can short local minimization rescue selected poses?\n", "\n", "## Learning objectives\n", "\n", "- Reuse authoritative ligand chemistry through one build context.\n", "- Build and score a matched ligand-decoy batch.\n", "- Relate interaction score to ligand heavy-atom displacement.\n", "- Minimize a local ligand/pocket shell in one batched call.\n", "- Separate pose sensitivity from docking and binding-affinity claims.\n", "\n", "## Before you begin\n", "\n", "Complete [07 — Ligands and Parameter Files](07_ligand_and_params.ipynb) first. This case study reuses its pinned ADA/LG1 fixture and chemistry but keeps the experimental question separate from parameter-file mechanics.\n" ] }, { "cell_type": "markdown", "id": "7ff5f3cf", "metadata": {}, "source": [ "## Setup\n", "\n", "The checked-in CIF and matching `.tmol` file are the authoritative coordinate and chemical inputs. The notebook runs on CPU or CUDA, fixes every random seed, and performs no live structure download.\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "dbec29b1", "metadata": {}, "outputs": [], "source": [ "try:\n", " import google.colab # noqa: F401\n", "except ImportError:\n", " IN_COLAB = False\n", "else:\n", " IN_COLAB = True\n", "\n", "if IN_COLAB:\n", " from urllib.request import urlopen\n", "\n", " exec(\n", " urlopen(\n", " \"https://raw.githubusercontent.com/uw-ipd/tmol/\"\n", " \"master/docs/tutorial/colab_setup.py\"\n", " ).read(),\n", " globals(),\n", " )\n", " setup_colab(\n", " [\n", " \"tmol/tests/data/protein_ligand_test/ada.tmol.nomin.cif\",\n", " \"tmol/tests/data/protein_ligand_test/ada.xtal-lig.mmff94.tmol\",\n", " ]\n", " )" ] }, { "cell_type": "code", "execution_count": 2, "id": "cc968d37", "metadata": {}, "outputs": [], "source": [ "from contextlib import redirect_stderr, redirect_stdout\n", "from io import StringIO\n", "from pathlib import Path\n", "import warnings\n", "\n", "import biotite.structure.io\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import pandas as pd\n", "import torch\n", "from IPython.display import display\n", "\n", "import tmol\n", "from tmol.database import ParameterDatabase\n", "from tmol.io import build_context_from_biotite, pose_stack_from_biotite\n", "from tmol.ligand import inject_params_file\n", "from tmol.ops import build_sidechain_coord_mask, res_mask_to_coord_mask\n", "from tmol.optimization import run_cart_min\n", "from tmol.pose import PoseStackBuilder\n", "from tmol.score import beta2016_score_function\n", "\n", "SEED = 20260810\n", "np.random.seed(SEED)\n", "torch.manual_seed(SEED)\n", "if torch.cuda.is_available():\n", " torch.cuda.manual_seed_all(SEED)\n", "warnings.filterwarnings(\n", " \"ignore\", message=r\"Sparse invariant checks are implicitly disabled.*\"\n", ")\n", "\n", "device = (\n", " torch.device(\"cuda\", torch.cuda.current_device())\n", " if torch.cuda.is_available()\n", " else torch.device(\"cpu\")\n", ")\n", "LIGAND_NAME = \"LG1\"\n", "\n", "\n", "def show_table(frame):\n", " \"\"\"Display a sortable table when available, with a pandas fallback.\"\"\"\n", " try:\n", " from itables import show\n", " except ImportError:\n", " return display(frame)\n", " return show(frame)\n", "\n", "\n", "def block_mask_for_name3(pose, name3):\n", " \"\"\"Select blocks with one residue name across a PoseStack.\"\"\"\n", " mask = torch.zeros_like(pose.block_type_ind, dtype=torch.bool)\n", " for pose_index in range(pose.n_poses):\n", " for block_index in range(pose.max_n_blocks):\n", " if int(pose.block_type_ind64[pose_index, block_index]) < 0:\n", " continue\n", " mask[pose_index, block_index] = (\n", " pose.block_type(pose_index, block_index).name3 == name3\n", " )\n", " return mask\n", "\n", "\n", "def heavy_atom_mask_for_blocks(pose, block_mask):\n", " \"\"\"Expand selected blocks to their non-hydrogen coordinate atoms.\"\"\"\n", " mask = torch.zeros_like(pose.real_atoms)\n", " for pose_index, block_index in torch.nonzero(block_mask, as_tuple=False).tolist():\n", " block_type_index = int(pose.block_type_ind64[pose_index, block_index].item())\n", " block_type = pose.block_type(pose_index, block_index)\n", " offset = int(pose.block_coord_offset64[pose_index, block_index])\n", " n_atoms = len(block_type.atoms)\n", " is_hydrogen = pose.packed_block_types.atom_is_hydrogen[\n", " block_type_index, :n_atoms\n", " ].bool()\n", " mask[pose_index, offset : offset + n_atoms] = ~is_hydrogen\n", " return mask & pose.real_atoms\n", "\n", "\n", "def ligand_protein_interactions(pose, score_function):\n", " \"\"\"Return both-orientation weighted ligand–protein scores per pose.\"\"\"\n", " ligand = block_mask_for_name3(pose, LIGAND_NAME)\n", " protein = (pose.block_type_ind64 >= 0) & ~ligand\n", " pair_mask = (ligand[:, :, None] & protein[:, None, :]) | (\n", " protein[:, :, None] & ligand[:, None, :]\n", " )\n", " scorer = score_function.render_block_pair_scoring_module(pose)\n", " with torch.no_grad():\n", " weighted = scorer(\n", " pose.coords, sum_terms=False, apply_weights=True\n", " ).sum(dim=0)\n", " return (weighted * pair_mask).sum(dim=(1, 2))\n", "\n", "\n", "def ligand_rmsd_from_native(pose, native_heavy_coords):\n", " \"\"\"Measure ligand heavy-atom RMSD in the fixed protein coordinate frame.\"\"\"\n", " ligand = block_mask_for_name3(pose, LIGAND_NAME)\n", " heavy = heavy_atom_mask_for_blocks(pose, ligand)\n", " values = []\n", " for pose_index in range(pose.n_poses):\n", " delta = pose.coords[pose_index, heavy[pose_index]] - native_heavy_coords\n", " values.append(torch.sqrt(torch.mean(torch.sum(delta * delta, dim=-1))))\n", " return torch.stack(values)" ] }, { "cell_type": "markdown", "id": "824d24e2", "metadata": {}, "source": [ "## Build one ligand-aware context\n", "\n", "The pose and score function use the same extended parameter database. Hydrogen optimization is disabled because the experiment isolates controlled rigid-body ligand placement and short local Cartesian refinement; every compared state starts from the same prepared coordinates and chemistry.\n" ] }, { "cell_type": "code", "execution_count": 3, "id": "b42f443f", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "device: cpu\n", "coordinate input: ada.tmol.nomin.cif\n", "authoritative chemistry: ada.xtal-lig.mmff94.tmol\n", "ligand heavy atoms: 19\n" ] } ], "source": [ "repo_root = Path.cwd()\n", "data_dir = repo_root / \"tmol\" / \"tests\" / \"data\" / \"protein_ligand_test\"\n", "if not (data_dir / \"ada.tmol.nomin.cif\").exists():\n", " repo_root = Path(tmol.__file__).resolve().parents[1]\n", " data_dir = repo_root / \"tmol\" / \"tests\" / \"data\" / \"protein_ligand_test\"\n", "\n", "complex_path = data_dir / \"ada.tmol.nomin.cif\"\n", "params_path = data_dir / \"ada.xtal-lig.mmff94.tmol\"\n", "parameter_database = inject_params_file(ParameterDatabase.get_default(), params_path)\n", "atom_array = biotite.structure.io.load_structure(\n", " str(complex_path), model=1, include_bonds=True\n", ")\n", "\n", "diagnostics = StringIO()\n", "try:\n", " with redirect_stdout(diagnostics), redirect_stderr(diagnostics):\n", " build_context = build_context_from_biotite(\n", " atom_array,\n", " device,\n", " param_db=parameter_database,\n", " prepare_ligands=False,\n", " )\n", " native_pose = pose_stack_from_biotite(\n", " atom_array,\n", " device,\n", " context=build_context,\n", " no_optH=True,\n", " )\n", "except Exception:\n", " print(diagnostics.getvalue())\n", " raise\n", "\n", "score_function = beta2016_score_function(\n", " device, param_db=build_context.parameter_database\n", ")\n", "native_ligand = block_mask_for_name3(native_pose, LIGAND_NAME)\n", "if int(native_ligand.sum().item()) != 1:\n", " raise RuntimeError(\"Expected exactly one LG1 ligand block\")\n", "native_ligand_coords = res_mask_to_coord_mask(native_pose, native_ligand)\n", "native_heavy_mask = heavy_atom_mask_for_blocks(native_pose, native_ligand)\n", "native_heavy_coords = native_pose.coords[native_heavy_mask]\n", "\n", "print(f\"device: {device}\")\n", "print(f\"coordinate input: {complex_path.name}\")\n", "print(f\"authoritative chemistry: {params_path.name}\")\n", "print(\"ligand heavy atoms:\", int(native_heavy_mask.sum().item()))" ] }, { "cell_type": "markdown", "id": "601cd2e7", "metadata": {}, "source": [ "## Generate a matched decoy series\n", "\n", "Every decoy preserves the protein and the ligand's internal geometry. Only the ligand receives a declared rotation about its centroid and translation in the protein frame. These seven states are sensitivity probes, not samples from a docking search distribution.\n" ] }, { "cell_type": "code", "execution_count": 4, "id": "1a9239ed", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "
| \n", "\n", " Loading ITables v2.9.1 from the internet...\n", " (need help?)\n", " | \n", "
| ⓘpose_index | \n", "state | \n", "rotation_degrees | \n", "translation_A | \n", "ligand_heavy_atom_RMSD_A | \n", "ligand_protein_interaction_score | \n", "whole_pose_score | \n", "
|---|---|---|---|---|---|---|
| 0 | deposited | 0.0 | 0.000000 | 0.000000 | 12.211265 | 954.375000 |
| 1 | small rotation | 15.0 | 0.250000 | 0.832808 | 607.751770 | 1549.915527 |
| 2 | small shift | -20.0 | 0.500000 | 1.193737 | 415.838959 | 1358.002686 |
| 3 | mixed 1 | 30.0 | 0.790569 | 1.747893 | 1741.305176 | 2683.468994 |
| 4 | mixed 2 | -45.0 | 1.145644 | 2.672955 | 1377.492310 | 2319.656982 |
| 5 | large shift | 60.0 | 1.500000 | 3.402580 | 2660.247314 | 3602.409912 |
| 6 | far decoy | 90.0 | 2.061553 | 4.788392 | 3613.181641 | 4555.344238 |
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
| \n", "\n", " Loading ITables v2.9.1 from the internet...\n", " (need help?)\n", " | \n", "
| ⓘstate | \n", "movable_atoms | \n", "interaction_before | \n", "interaction_after | \n", "interaction_change | \n", "whole_pose_before | \n", "whole_pose_after | \n", "whole_pose_change | \n", "ligand_RMSD_before_A | \n", "ligand_RMSD_after_A | \n", "optimizer_budget | \n", "
|---|---|---|---|---|---|---|---|---|---|---|
| deposited | 199 | 12.211265 | -20.625813 | -32.837078 | 954.374939 | 903.116699 | -51.258240 | 0.000000 | 0.198972 | 15-iteration smoke test; convergence not assessed |
| small shift | 198 | 415.838959 | -7.526464 | -423.365417 | 1358.002319 | 924.598877 | -433.403442 | 1.193737 | 0.900713 | 15-iteration smoke test; convergence not assessed |
| far decoy | 213 | 3613.181641 | 59.484619 | -3553.697021 | 4555.344238 | 1346.806152 | -3208.538086 | 4.788392 | 4.342620 | 15-iteration smoke test; convergence not assessed |
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.
3Dmol.js failed to load for some reason. Please check your browser console for error messages.