STRC h01 Phase 9c — Umbrella Sampling PMF (SKELETON)

Placeholder. Computes a converged, equilibrium potential of mean force along the ligand-dissociation coordinate as the orthogonal counterpart to the forced-pulling tauRAMD residence times in Phase 5m / 5q. No production calculation yet.

Question

Along a chosen reaction coordinate (typically COM-COM distance between ligand and pocket Cα atoms):

  1. What is the equilibrium free-energy barrier ΔG‡(unbind) for mutant vs WT?
  2. Does the barrier follow the same ranking as the tauRAMD residence times τ?
  3. Is the bound-state well depth ΔG_bind(PMF) consistent with FEP (9a) and MM-PBSA (9b)?

Method

  1. Define reaction coordinate: ligand-COM to pocket-Cα-COM distance, ξ ∈ [3, 25] Å (or path-CV if linear COM is anatomically wrong for this pocket — decide after smoke).
  2. Generate starting configurations along ξ via slow steered MD pulled apart, sample frames at every 0.5 Å.
  3. Set up umbrella windows: ~40 windows × 0.5 Å spacing × harmonic bias k ≈ 1000-2000 kJ/mol/nm².
  4. Equilibrate each window 2-5 ns, production 10-20 ns.
  5. WHAM (gmx wham, or Grossfield WHAM) → PMF curve W(ξ).
  6. Histogram-overlap diagnostic between adjacent windows must show overlap (no gaps).

Tooling: GROMACS pull code for window generation + gmx wham for analysis.

Why orthogonal

  • vs tauRAMD: tauRAMD applies a stochastic force to pull the ligand out fast; bias the rate via force-constant choice. Umbrella sampling is equilibrium — converged windows give a true free-energy barrier independent of pulling rate.
  • vs FEP (9a): FEP gives ΔG_bind (well depth); umbrella gives ΔG_bind and ΔG‡ (the barrier). Both should agree on well depth.
  • vs MM-PBSA (9b): MM-PBSA is end-point implicit-solvent; PMF is path-based explicit solvent.
  • Caveat: choice of reaction coordinate is decisive. If the unbinding path is non-radial (e.g., cryptic-pocket-mediated, see Phase 5c), 1D COM-distance gives wrong barrier. Path-CV or RAM-based CV may be needed.

Inputs needed

  • Equilibrated holo complex from Phase 9b production (or Phase 5d) for WT and mutant.
  • Force field: same as 9b (AMBER ff19SB / ff14SB + OpenFF/GAFF2 + OPC/TIP3P).
  • Hardware: GPU. ~40 windows × 15 ns × 2 systems × N ligands; embarrassingly parallel across windows.

Smoke test (1-day, theoretical)

Single ligand × WT only × 5 windows × 1 ns each:

  1. Generate 5 starting configs at ξ ∈ {6, 9, 12, 15, 18} Å via 1 ns SMD.
  2. Run 1 ns umbrella per window with bias k=1000 kJ/mol/nm².
  3. Plot histogram of ξ for each window.

Smoke pass: adjacent histograms overlap (≥ 30% bin-overlap in each pair) and no window drifts off its target ξ. → reaction coordinate is sane on this system.

Smoke fail: histograms gap → bias too weak or coordinate doesn’t span; redesign before committing.

Production protocol (theoretical)

  • 40 windows × 0.5 Å × 15 ns each × 2 variants × 1-3 ligands.
  • Discard first 5 ns per window; WHAM on remaining 10 ns.
  • Bootstrap (gmx wham -bs-method b-hist) for error bars on PMF.
  • Convergence diagnostic: PMF stable when 2nd half of each window matches whole-window PMF within 0.5 kcal/mol.

Pass criteria

  • PRIMARY: ΔG‡(unbind, mutant) > ΔG‡(unbind, WT) by ≥ 1 kcal/mol — qualitatively matches tauRAMD ordering.
  • SECONDARY: ΔG_bind(PMF) within 2 kcal/mol of ΔG_bind(MM-PBSA from 9b) — cross-method consistency.
  • TERTIARY: histogram overlap > 30% in all adjacent window pairs; bootstrap error < 1 kcal/mol on barrier height.
  • FAIL state: PMF non-monotonic with no clear barrier, or barrier ordering opposite to tauRAMD → mechanism re-derivation; tauRAMD may have been pulling-rate-biased.

Known artifacts and risks

  • Reaction coordinate choice: 1D COM is the failure mode. If smoke shows the ligand exits via a different geometry than radial unbinding, switch to path-CV (Branduardi 2007) or RAM-CV.
  • Bias-force calibration: too weak → poor histogram; too strong → ligand frozen in window with no overlap. Tune in smoke.
  • Sampling around the barrier top: extra-fine windows (0.25 Å) recommended around the apparent transition state.
  • Slow degrees of freedom: protein side-chain rearrangement on ligand exit may not equilibrate in 15 ns — extend or use replica-exchange umbrella.

References (canonical)

  • Torrie, Valleau 1977. J Comput Phys 23:187 — original umbrella sampling.
  • Kästner 2011. WIREs Comput Mol Sci 1:932 — umbrella sampling review.
  • Hub, de Groot, van der Spoel 2010. J Chem Theory Comput 6:3713 — gmx wham.
  • Roux 1995. Comput Phys Commun 91:275 — WHAM.
  • Branduardi, Gervasio, Parrinello 2007. J Chem Phys 126:054103 — path-CV.

Status

  • 2026-04-27 — skeleton created. Reaction-coordinate decision pending. No runs performed.

Ranking delta

  • No change. Skeleton only.

Connections