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):
- What is the equilibrium free-energy barrier ΔG‡(unbind) for mutant vs WT?
- Does the barrier follow the same ranking as the tauRAMD residence times τ?
- Is the bound-state well depth ΔG_bind(PMF) consistent with FEP (9a) and MM-PBSA (9b)?
Method
- 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).
- Generate starting configurations along ξ via slow steered MD pulled apart, sample frames at every 0.5 Å.
- Set up umbrella windows: ~40 windows × 0.5 Å spacing × harmonic bias k ≈ 1000-2000 kJ/mol/nm².
- Equilibrate each window 2-5 ns, production 10-20 ns.
- WHAM (gmx wham, or Grossfield WHAM) → PMF curve W(ξ).
- 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:
- Generate 5 starting configs at ξ ∈ {6, 9, 12, 15, 18} Å via 1 ns SMD.
- Run 1 ns umbrella per window with bias k=1000 kJ/mol/nm².
- 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
[parent]STRC h01 Phase 9 Orthogonal Cross-Checks Plan 2026-04-27[see-also]h01 hub[see-also]STRC h01 Phase 5m Production tauRAMD Ranking 2026-04-26 — forced-pull counterpart this phase replaces with equilibrium PMF[see-also]STRC h01 Phase 5q v5.3 STRC Within-Target tauRAMD 2026-04-26[about]Misha