Source-linked AI summary

Harmonic Torsional Diffusion for Protein-Ligand Flexible Docking

Maksim Zhdanov, Pavel Strashnov, Vladislav Kurenkov

arXiv:2608.20366v1q-bio.BMcs.LG

TL;DR

Flexible docking needs to model ligand and protein conformational changes, yet generic torsional heads ignore the periodic geometry of angular variables. Harmony uses harmonic torsional diffusion to represent these variables explicitly, and it improves pose recovery, pocket reconstruction, and generated-complex validity across reported benchmarks. The approach is limited by its reliance on variance-exploding diffusion with a Gaussian noise prior.

  • Problem

    Flexible docking requires joint reasoning about ligand pose and protein flexibility, while existing diffusion models use machinery that does not exploit the geometry of torsional variables.

  • Method

    Harmony parameterizes ligand and side-chain torsional scores as derivatives of learned harmonic potentials in Fourier space within variance-exploding diffusion on the torus.

  • Results

    Harmony improves ligand pose recovery and all-atom pocket reconstruction on PDBBind, increases ligand RMSD < 2 Å accuracy from 39.7% to 64.2%, and improves physical validity on PoseBusters.

  • Takeaways & Limitations

    The results support aligning torsional score parameterization with the intrinsic manifold structure of angular motion as a useful design choice for flexible docking.

  • Takeaways & Limitations

    Harmony relies on variance-exploding diffusion with a Gaussian noise prior, limiting direct adaptation to arbitrary endpoint distributions such as flow-based interpolants.

Abstract

from arXiv · show

Molecular docking requires reasoning jointly about ligand pose and protein flexibility. Most diffusion-based docking models predict torsional updates with generic Euclidean heads that ignore the periodic geometry of angular variables. This mismatch is especially limiting in flexible docking, where ligand conformations and pocket side chains co-adapt to form the bound complex. Here, we introduce Harmony, a harmonic torsional diffusion framework for flexible protein-ligand docking. Harmony parameterizes ligand and side-chain torsional score fields as derivatives of learned harmonic potentials on the circle, whose noise-level dependence is supplied analytically by the heat semigroup of variance-exploding diffusion on the torus. This construction makes periodicity explicit and gives the model a frequency-aware inductive bias over rotameric motion. On the PDBBind benchmark, Harmony improves ligand pose accuracy and pocket all-atom reconstruction over recent flexible docking methods. On PoseBusters, it improves the physical validity of generated complexes. Case studies on EBNA1 and KRAS G12D illustrate the method's behavior on a polar and a shallow binding site, respectively. Together, these results indicate that aligning the score parameterization with the geometry of the diffusion process is a simple and effective lever for improving flexible docking.

1. Introduction

Flexible docking must account for protein conformational changes, but existing approaches often use fixed proteins or geometry-mismatched torsional modeling. Harmony addresses this with harmonic torsional diffusion and reports improved flexible-docking performance.

  • Many docking methods assume protein structures remain fixed, despite proteins adopting different conformations upon ligand association.
  • Flexible docking remains difficult because protein flexibility expands the search space and common generative machinery does not exploit torsional geometry.
  • Harmony parameterizes ligand and protein side-chain torsions in a Fourier basis, with an analytical score formulation and frequency-suppressing noising structure.
  • Harmony improves recent generative flexible-docking methods on PDBBind by approximately 24.5 percentage points and improves chemical validity on PoseBusters.
  • The framework is presented as a generative approach that naturally represents degrees of freedom central to molecular binding.

2. Related Work

Related flexible-docking methods extend generative modeling to ligand and pocket flexibility, while Harmony specifically targets torsional motion with a geometry-aware framework. The present work focuses on side-chain flexibility and supports several docking settings and protein representations.

  • Recent methods extend diffusion, bridge, and flow-based frameworks to model flexibility in ligands and protein binding pockets.
  • Existing flexible-docking approaches differ in scope, with some modeling side chains only and others attempting backbone motion.
  • Harmony introduces a generative framework designed to model torsional motion in both ligands and proteins.
  • Harmony can be used for pocket-based or blind docking and with rigid or flexible protein representations, although this work restricts flexibility to side chains.

3. Preliminaries

Harmony represents docking as joint prediction of ligand pose and pocket side-chain arrangement on a product manifold combining Euclidean, rotational, and toroidal variables. Its input is a heterogeneous, multi-scale graph of ligand and receptor structure.

  • The complex is represented as a heterogeneous graph containing ligand atoms, receptor backbone residues, receptor atoms, and multiple interaction edge types.
  • The pocket-centered graph retains the ligand and nearby receptor residues and atoms while fixing the protein backbone and allowing pocket side-chain motion.
  • The docking state combines ligand translation, ligand rotation, ligand torsions, and flexible pocket side-chain torsions on a product manifold.
  • The network jointly predicts ligand pose and pocket side-chain arrangement from the pocket-centered graph.

4. Harmonic Torsional Diffusion

Harmony builds flexible docking around VE diffusion on ligand rigid-body and periodic torsional variables, using Fourier modes to represent heat-flow dynamics on the circle. The model analytically parameterizes torsional scores, jointly denoises ligand and side-chain motion, and evaluates generated complexes without post hoc relaxation.

  • Harmonic parametrization: Because VE density evolution follows heat flow, Fourier modes provide the natural basis for periodic torsional variables on the circle.The Laplace–Beltrami eigenfunctions are cos(kϑ) and sin(kϑ), and the heat semigroup acts diagonally on finite Fourier expansions.
  • Diffusion formulation: Harmony models ligand rigid-body motion, ligand torsions, and pocket side-chain torsions with independently perturbed components in a VE diffusion process.The protein backbone remains fixed while pocket side chains are allowed to move.
  • Sampling: The network predicts scores for each forward-noising factor, and reverse-time sampling jointly recovers the ligand pose and pocket side-chain arrangement.Ligand translation uses Euclidean Gaussian noise, rotation uses isotropic SO(3) noise, and torsions use wrapped Gaussian noise on the torus.
  • Harmonic parametrization: Harmony predicts torsional scores by differentiating heat-smoothed harmonic functions for every ligand and side-chain rotatable bond.Local bond representations produce harmonic coefficients, which define latent functions on the circle and their torsional score components.
  • Frequency damping: The explicit damping factor exp(−1/2k2σ2_u(λ)) encodes noise-level dependence, preserves periodicity, and captures multimodal torsional landscapes.High-frequency angular structure is suppressed at larger noise levels, while finer modes re-emerge as λ →0.
  • Confidence and relaxation: Without post hoc relaxation, Harmony evaluates the generated complexes directly, while a confidence head predicts pose confidence supervised by interface lDDT.The confidence loss upweights low-noise samples, and the model uses an SE(3)-equivariant graph neural network with standard diffusion score matching.

5. Experiments

Harmony is evaluated on ligand placement, pocket reconstruction, physical validity, and targeted ablations and case studies. It improves flexible-docking accuracy over recent methods, with evidence that harmonic damping and joint ligand–side-chain modeling contribute to performance.

  • Experimental setup: PDBBind evaluation measures top-1 ligand placement and binding-pocket geometry without oracle information, using RMSD-based accuracy metrics.Ligand performance includes median RMSD and the fraction below 2 Å; pocket accuracy uses the fraction below 1 Å.
  • Experimental setup: Post-processing can improve geometric plausibility but may obscure deficiencies in the underlying generative model by correcting molecular inconsistencies.Harmony generates final complexes without auxiliary refinement in the reported comparison.
  • Main results: Harmony improves ligand placement from 39.7% to 64.2% for RMSD < 2 Å and pocket reconstruction from 41.7% to 47.3% for AA-RMSD < 1 Å.The comparison uses the same architectural design and training protocol as DiffDock-Pocket and FlexDock.
  • Main results: Harmony outperforms direct pocket-based flexible competitors DiffDock-Pocket and FlexDock on PoseBusters.Additional comparisons exclude relaxation and post-processing to assess the generative models directly.
  • Ablation studies: Harmony performs best with K = 4, degrades when too many harmonic modes are used, and outperforms FlexDock across all tested K values.Performance drops manageably at minimum K = 2, supporting a preference for coarse torsional dynamics.
  • Ablation studies: Removing frequency damping degrades performance, while applying harmonic learning to one torsional group improves modeling of the other.The result supports noise-dependent suppression of high-frequency components and a joint ligand–side-chain representation.
  • Case studies: In EBNA1 and KRAS G12D case studies, Harmony produces more accurate ligand poses and improved local all-atom reconstruction than the baseline.Harmony preserves EBNA1’s enclosing loop and recovers a more accurate KRAS G12D conformation in a shallow, exposed site.

6. Discussion

Harmony’s main scope boundary is its reliance on variance-exploding diffusion with a Gaussian noise prior. The framework nevertheless consistently improves pose recovery and all-atom reconstruction across the reported benchmarks and case studies.

  • Harmony requires a Gaussian noise prior because it is based on a variance-exploding diffusion process.
  • Adapting the harmonic framework to arbitrary endpoint distributions, including flow-based interpolants, is identified as future work for fully all-atom docking with backbone motion.
  • Harmony makes angular periodicity explicit and improves modeling of flexible binding sites through Fourier-space torsional scores and variance-exploding diffusion on the torus.
  • Across standard benchmarks and case studies, Harmony consistently improves pose recovery and all-atom reconstruction.

A. Proofs of Propositions

The propositions establish the heat-equation and harmonic-mode behavior underlying Harmony’s torsional diffusion. The PoseBusters breakdown reports improved physically meaningful validity checks for Harmony over the baseline.

  • Propositions: The driftless diffusion’s forward density evolves according to a heat equation with time-dependent diffusivity ν(λ) = g(λ)^2/2.
  • Propositions: On the circle, cosine and sine modes are Laplacian eigenfunctions with eigenvalue −k^2, so the heat semigroup acts on each harmonic mode independently.
  • Propositions: Before wrapping, the angular perturbation is Gaussian; quotienting by R/2πZ produces the wrapped normal kernel.
  • Propositions: The heat kernel on the torus is a periodized Gaussian with variance 2t, obtained by setting t = σ^2(λ)/2.
  • PoseBusters: Each PoseBusters bar reports a criterion-specific pass rate, with Harmony improving checks involving internal steric validity, ligand strain, and protein–ligand separation.

C. Training and Inference

The training and inference procedures include implementation details for Harmony, including data preparation, pocket selection, conformer matching, frame alignment, and flexibility formalization.

  • The paper provides additional details on training, inference, data preparation, pocket selection, conformer matching, protein-frame alignment, and flexibility formalization.

C.1. Algorithms for Training and Inference of Harmony

Harmony’s algorithms prepare pocket-centered apo graphs, add torsional and rigid-body noise, predict harmonic torsional scores, and iteratively reverse the variance-exploding diffusion process. Data preparation aligns apo and holo structures and constructs rotatable-bond fragments for ligands and side chains.

  • Data preparation: Pocket selection retains nearby residues and enlarges the crop with a 10 Å buffer to improve robustness to pocket choice.
  • Training algorithm: Training samples noise levels and translation, rotational, ligand-torsional, and side-chain-torsional perturbations before predicting rigid-body scores and harmonic coefficients.
  • Training algorithm: The training procedure evaluates harmonic torsional scores, computes score targets from the forward kernel, and optimizes matching losses with Adam.
  • Data preparation: Conformer matching optimizes ligand and receptor side-chain torsions so generated apo inputs better correspond to reference crystal structures.
  • Flexibility representation: Flexibility is represented with directed torsional edges, rotatable-bond masks, affected-atom fragments, and dependency-based torsion application order.
  • Inference algorithm: Inference begins from a pocket-centered apo graph and iteratively predicts scores, samples reverse-time noise, applies reverse VE updates, and returns ligand and side-chain coordinates.
  • Data preparation: Apo–holo alignment removes global rigid-body offsets and aligns side chains into local backbone frames; inference uses the apo structure as the working frame.

D. Model Architecture

Harmony combines heterogeneous protein–ligand graphs with an equivariant backbone and specialized torsional prediction heads. Its preprocessing and confidence module support flexible ligand and side-chain modeling.

  • Graph representation: Harmony represents ligand atoms, receptor residues, and receptor atoms as three coupled node sets with static and dynamically rebuilt geometric graphs.Dynamic neighborhoods are recomputed from current noised coordinates during message passing.
  • Graph representation: The static graphs encode ligand covalent bonds, receptor atom-to-residue assignments, and residue-wise side-chain torsional topology.Ligand rotatable edge masks define ligand torsional degrees of freedom.
  • Equivariant backbone: A heterogeneous SE(3)-equivariant e3nn backbone processes ligand, residue, and receptor-atom graphs with shared tensor-product convolution layers.The architecture uses six layers, scalar width ns = 60, higher-order width nv = 15, and spherical harmonics through lmax = 2.
  • Prediction heads: Three prediction modules estimate rigid-body ligand scores, harmonic ligand torsion coefficients, and noise-conditioned harmonic side-chain torsion coefficients.The predicted coefficients are converted into periodic torsional score fields.
  • Confidence estimation: The confidence head jointly predicts a scalar pose confidence score using pooled complex representations and PLI-lDDT supervision.PLI-lDDT measures preservation of native receptor–ligand contact distances and emphasizes local binding-site quality.

E.1. Data

The evaluation uses a time-based PDBBind split with paired holo and predicted apo receptors, and measures ligand accuracy, all-atom pocket reconstruction, and physical plausibility. These metrics cover both pose-level and full-complex quality.

  • Dataset: PDBBind complexes deposited before 2019 form the training and validation sets, while 363 post-2019 complexes are reserved for testing.The pipeline follows the time-based split used by FlexDock.
  • Dataset: Holo references are repaired with PDBFixer, while ESMFold predicts the corresponding apo receptor structures used for flexible docking.The corrected complex structures define holo references before apo prediction.
  • Physical validity: PoseBusters evaluates physical plausibility through ligand internal validity and protein–ligand interaction checks, including steric-clash detection.This complements geometric docking metrics with chemical and structural validity checks.
  • Docking metrics: Ligand pose quality is measured by L-RMSD after optimal rigid alignment of the pockets and summarized by a ligand success rate.The success rate is the fraction of complexes satisfying the reported ligand criterion.
  • Docking metrics: Full-complex reconstruction is assessed with receptor-pocket AA-RMSD and the fraction of predictions meeting its stricter all-atom criterion.AA-RMSD evaluates retained pocket receptor atoms, including flexible binding-site side chains.

E.3. Similarity Analysis for Case Studies

The case-study similarity analysis compares target complexes with the PDBBind training split using ligand, protein, and contact-environment descriptors. A combined score averages these three similarity components.

  • Analysis design: Each target complex is compared against every complex in the PDBBind training split using a descriptor-based similarity analysis.The analysis is intended to contextualize qualitative case studies.
  • Ligand similarity: Ligand similarity uses radius-2, 1024-bit Morgan fingerprints after hydrogen removal and measures the Dice coefficient.The fingerprint comparison uses sets of active bits for each ligand.
  • Protein similarity: Protein similarity is global sequence identity computed from Needleman–Wunsch alignment of receptor amino-acid sequences.The alignment includes matched residues and gaps in its aligned length.
  • Contact similarity: Contact-interaction similarity records receptor residues within 4.5 Å of the ligand and separates contacts into close and near distance buckets.Close contacts are at most 3.5 Å.
  • Combined score: The combined similarity score averages ligand, protein, and contact-interaction similarity components against the training complexes.This produces a coarse descriptor-level comparison of each target complex.
Loading 2608.20366v1…