STRC Research

Recipe — ABF Adaptive Biasing Force Algorithm

4 min read

Recipe — ABF Adaptive Biasing Force Algorithm

P1 recipe — Algorithms 1 and 2 from 2007-chipot-free-energy-calculations-book §4.6.4 (Darve), pp. 137–139. ABF computes the potential of mean force (PMF) along a chosen reaction coordinate ξ by accumulating instantaneous mean-force samples in bins of ξ and applying the running negative as an external bias. Once the bias balances the underlying free-energy gradient, the system diffuses freely along ξ and barriers are erased. Use this whenever an STRC question reduces to “what is the free-energy profile A(ξ) along this 1D or 2D coordinate?”

Outer loop (verbatim Algorithm 1)

Loop over time steps i = 1, …, n: a1 ← M⁻¹ F(r0) /* Force at time t. M = diagonal mass matrix. F = −∇U. / a1 ← ABF(a1, r0, v) / Call to ABF subroutine; adds the adaptive bias. / v ← v + dt · a1 / Velocity advanced from t − dt/2 to t + dt/2. / r0 ← r0 + dt · v / Position advanced from t to t + dt. */ End of loop

— [Chipot & Pohorille 2007, p. 138, §4.6.4 Algorithm 1]

ABF subroutine essence (Algorithm 2 condensed)

For each step:

  1. Identify bin k corresponding to current ξ(r0).
  2. Compute ramp R(n^k(N)) — a damping factor (typically 0 → 1 as bin sample count rises) so the bias is not applied before statistics accrue. Eq. (4.46) in book.
  3. la ← R(n^k) · F(k) / n(k) — the candidate bias, which is the running mean force in bin k.
  4. a1 ← a1 + la · M⁻¹ ∇ξ(r0) — apply the bias along the coordinate gradient.
  5. Accumulate the instantaneous mean-force sample F(k) ← F(k) − pxi0 where pxi0 estimates −d(m_ξ ξ̇)/dt at the current step (centered finite difference using the saved acceleration history).
  6. Increment bin counter n(k) ← n(k) + 1.

The free-energy derivative is recovered as dA/dξ(k) = −F(k)/n(k), integrated to give A(ξ).

Practical knobs (book §4.6.5)

  • Bin width — small enough to resolve the PMF features, large enough to keep n(k) ≥ a few thousand for σ/√n convergence. The book’s polyalanine example: σ ≈ 13 vs. mean ≈ −1 in a single bin → ~14,000 samples needed for 10% error in dA/dξ.
  • Ramp R(n^k) — turn the bias on gradually. If applied before statistics accumulate, ABF can drive the system into unphysical regions.
  • Partial bias — multiply the bias by a fraction (0 < f < 1) to flatten but not erase barriers; minima/maxima stay in place but reduced. Useful when full-bias ABF over-samples uninteresting regions.
  • Optimal sampling correction — U_opt = U − 2 k_B T ln(σ(ξ)/σ₀) (Eq. 4.47) adds extra sampling where the force std-dev is large. Two-pass: run ABF, fit σ(ξ), restart with the σ-derived correction.
  • Force partitioning — the derivative dA/dξ can be decomposed into contributions from U_bonds, U_elec, U_LJ etc. (Eq. 4.48). Useful for diagnosing which term governs a barrier (electrostatic vs. steric).
  • 2D ABF — extends naturally to ξ = (ξ₁, …, ξ_p) but reconstructing A from the gradient field requires path-independence; book recommends bilinear Q₁ spline fit (§4.6.5 Eq. 4.51) with 4 control points per cell.

Comparison vs other PMF methods (book §4.7)

  • SMD (steered MD): applies F = k(vt − x) external pulling force; gives non-equilibrium work. ABF is preferable when convergence to equilibrium PMF is the goal.
  • Metadynamics (Laio–Parrinello): deposits Gaussian hills at visited points. ABF and metadynamics are mathematically related; ABF tends to converge with smaller-magnitude bias accumulated as direct force averages, metadynamics typically faster for highly multidimensional ξ but with characteristic over-shoot artifacts.

How to use this in STRC

  • h09 hydrogel self-assembly — define ξ as the distance between a peptide monomer and the growing β-sheet edge, or as the number of intermonomer hydrogen bonds. ABF gives A(ξ) for monomer addition; combined with replica exchange (Ch. 8) for sidechain sampling. Tool to add to inventory: any future hydrogel_abf_assembly.py.
  • h26 dimer dissociation PMF — define ξ as the distance between centers of mass of the two stereocilin C-terminal domains; ABF along ξ from contact (~0 nm) to fully separated (~2 nm) gives the binding free energy via integration. Combine with the engineered cysteine designs to compute ΔΔA_dimer.
  • h01 alternative scoring — for scoring a candidate pharmacochaperone, the binding free energy can also be obtained as the PMF along a ligand-receptor pulling coordinate (instead of alchemical λ). ABF here is an independent method, useful as a cross-check.

Connections

Graph View

Backlinks