Source-linked AI summary

State-specific protein-ligand complex structure prediction with a multi-scale deep generative model

Zhuoran Qiao, Weili Nie, Arash Vahdat, Thomas F. Miller, Anima Anandkumar

arXiv:2209.15171v2q-bio.QMcs.LGq-bio.BM

TL;DR

Existing algorithms do not systematically predict protein-ligand complex structures together with ligand-associated protein conformational changes. NeuralPLexer uses a multi-scale generative model conditioned on protein sequences and ligand molecular graphs to sample complex ensembles, outperforming established docking and structure-prediction methods across reported benchmarks. The approach also supports rapid exploration of protein-ligand interactions and may accelerate protein and ligand design, while broader data are needed for more challenging systems.

  • Problem

    Existing algorithms cannot systematically predict binding ligand structures together with their regulatory effects on protein folding, particularly for complexes involving major receptor conformational changes.

  • Method

    NeuralPLexer directly samples protein-ligand complex ensembles from protein sequences and ligand molecular graphs using a biophysics-informed, multi-scale generative model.

  • Results

    NeuralPLexer consistently outperformed existing docking and protein structure-prediction methods across blind docking, binding-site recovery, and flexible ligand-binding-protein structure benchmarks.

  • Takeaways & Limitations

    NeuralPLexer captures structural cooperativity between proteins and small molecules and supports rapid exploration of protein-ligand interactions and design applications.

  • Takeaways & Limitations

    Broader and improved experimental and bioinformatic data are needed to extend NeuralPLexer to proteins without experimentally determined homologs, post-translational modifications, and large heteromeric multi-state complexes.

Abstract

from arXiv · show

The binding complexes formed by proteins and small molecule ligands are ubiquitous and critical to life. Despite recent advancements in protein structure prediction, existing algorithms are so far unable to systematically predict the binding ligand structures along with their regulatory effects on protein folding. To address this discrepancy, we present NeuralPLexer, a computational approach that can directly predict protein-ligand complex structures solely using protein sequence and ligand molecular graph inputs. NeuralPLexer adopts a deep generative model to sample the 3D structures of the binding complex and their conformational changes at an atomistic resolution. The model is based on a diffusion process that incorporates essential biophysical constraints and a multi-scale geometric deep learning system to iteratively sample residue-level contact maps and all heavy-atom coordinates in a hierarchical manner. NeuralPLexer achieves state-of-the-art performance compared to all existing methods on benchmarks for both protein-ligand blind docking and flexible binding site structure recovery. Moreover, owing to its specificity in sampling both ligand-free-state and ligand-bound-state ensembles, NeuralPLexer consistently outperforms AlphaFold2 in terms of global protein structure accuracy on both representative structure pairs with large conformational changes (average TM-score=0.93) and recently determined ligand-binding proteins (average TM-score=0.89). Case studies reveal that the predicted conformational variations are consistent with structure determination experiments for important targets, including human KRAS$^\textrm{G12C}$, ketol-acid reductoisomerase, and purine GPCRs. Our study suggests that a data-driven approach can capture the structural cooperativity between proteins and small molecules, showing promise in accelerating the design of enzymes, drug molecules, and beyond.

Main

NeuralPLexer addresses limitations of single-structure and case-specific approaches by combining biophysical inductive biases with autoregressive and diffusion generative modeling. It directly samples protein-ligand complex ensembles and achieves strong performance across docking, binding-site recovery, and conformationally flexible structure prediction.

  • Motivation: Existing methods struggle to systematically predict protein-ligand complexes involving substantial receptor conformational changes, often requiring expert intervention or experimental constraints.Single-structure regression is especially limited for non-endogenous small molecules whose identities cannot be inferred from protein motifs.
  • Approach: NeuralPLexer integrates biophysical inductive biases with autoregressive contact prediction and diffusion-based atomistic structure generation.Its design addresses both global ligand-functional context and energetically favorable sub-nanometer inter-atomic interactions.
  • Approach: NeuralPLexer combines protein language-model and template features with ligand molecular graphs to sample ensembles of binding complex structures.The framework uses a multi-scale architecture that mirrors the hierarchical organization of biomolecular complexes.
  • Results: 78% maximum improvement in ligand pose accuracy was achieved over the best-performing existing method on the PDBBind2020 blind-docking benchmark.This result reflects the method’s performance on blind protein-ligand docking.
  • Results: 46% maximum binding-site structure recovery and at least 59% success-rate improvement over Rosetta were reported using truncated AlphaFold2-generated scaffolds.The authors describe this as an application for which no existing machine-learning approaches were applicable to their knowledge.
  • Results: NeuralPLexer outperformed AlphaFold2 on ligand-binding proteins with large structural plasticity, achieving an average TM-score of 0.906 and 11–13% higher accuracy on substantially changing domains.The results support selective prediction of structures subject to induced-fit binding or conformational selection.

Results

NeuralPLexer combines autoregressive contact prediction with equivariant diffusion to generate protein-ligand structures and conformational ensembles from sequence and molecular-graph inputs. Across docking, binding-site recovery, and flexible ligand-binding protein benchmarks, it improves structural accuracy and captures ligand-dependent conformational changes.

  • NeuralPLexer: NeuralPLexer generates protein-ligand structural ensembles from protein sequences and ligand molecular graphs, using auxiliary language-model and template features.The framework jointly samples protein and ligand heavy-atom coordinates through an end-to-end generative model.
  • NeuralPLexer: The model combines hybrid atom-and-frame molecular representations with autoregressive contact-map sampling and equivariant diffusion-based coordinate denoising.CPM predicts multiscale proximity structure, while ESDM generates coordinates through a learned stochastic reverse process.
  • Accurate ligand and binding pocket structure modeling: For human KRASG12C, NeuralPLexer predicted an open-like Switch-II pocket with ligand RMSD 2.08 Å, despite AF2 producing a closed pocket with a steric clash.The binding-site lDDT-BS score decreased from 0.76 for AF2 to 0.69, motivating a combined assessment using ligand RMSD and clash rate.
  • End-to-end structure prediction for ligand-binding proteins: NeuralPLexer achieved average weighted Q-factor 0.608 versus 0.501 for the ligand-free baseline, with 16% of predictions improving over AF2 by at least 0.1.The results indicate selective sampling of ligand-free and ligand-bound states on conformationally variable targets.
  • End-to-end structure prediction for ligand-binding proteins: On recently resolved complexes, average TM-score increased from 0.877 without ligand to 0.893 with ligand, matching AF2 at 0.891; the high-confidence subset reached 0.929 versus AF2 at 0.913.More than 70% of predictions meeting ligand RMSD <2.5 Å and lDDT-BS >0.8 were identified by ligand pLDDT >0.8.
  • End-to-end structure prediction for ligand-binding proteins: NeuralPLexer captured the ligand-associated closure of ketol-acid reductoisomerase’s N-subdomain observed in cryo-EM experiments.The comparison used predicted ligand-bound and ligand-free ensembles for target PDB:6KPE.

Discussions

NeuralPLexer enables rapid, differentiable exploration of protein–ligand structure and conformational space, while remaining limited by the breadth and quality of available training and auxiliary data.

  • Discussions: NeuralPLexer generates protein–ligand structures within seconds on a standard GPU after auxiliary features are gathered, supporting proteome-scale exploration.The authors position this capability for sequence and chemical spaces represented by computational protein and chemical databases.
  • Discussions: Its end-to-end differentiability makes NeuralPLexer suitable for integrating protein, ligand, and structure-based generative models in design workflows.The authors specifically discuss sequence and compound prioritization over conventional screening-based strategies.
  • Discussions: Broader and higher-quality experimental and bioinformatic data could extend analysis to proteins without experimentally determined homologs, post-translational modifications, and large heteromeric complexes.The authors also identify NMR, molecular dynamics, binding affinities, and high-throughput mass spectrometry as promising auxiliary data sources.

Methods

The methods combine molecular-graph encoding, multi-scale protein–ligand graph representations, contact prediction, and equivariant diffusion to generate atomistic complex structures under hierarchical geometric constraints.

  • Methods: Molecular Heat Transformer encodes ligand and amino-acid graphs with learned edge representations designed to infer inter-atomic geometrical correlations.The encoder uses attention blocks with edge embeddings and is trained against molecular structure and bioactivity information.
  • Methods: The graph representation combines sparse geometry-based edges with dense pair representations over protein and ligand anchor nodes.Protein features derive from language-model embeddings and geometric encodings, while ligand anchors are sampled uniformly from ligand frame nodes.
  • Methods: The equivariant structure denoising module represents nodes with scalar and Cartesian vector features to model three-dimensional coordinate displacements.Its architecture is inspired by point-cloud, rigid-body, and biopolymer representation-learning methods.
  • Methods: Ground-truth contact maps use residue centroid and ligand frame center-atom coordinates, with D0 = 8.0Å and ε = 10^-6.The resulting contact-map refinement uses distance histograms between protein residues and ligand frames.
  • Methods: The contact map is represented with 32 distance bins spanning 2.0Å to 22.0Å and then coarse-grained into protein patches for block-adjacency sampling.The patch-level representation is recycled by the contact prediction network.
  • Methods: A multivariate Ornstein–Uhlenbeck formulation introduces data-dependent drift and covariance, while anisotropic coefficients separate coordinate length scales.The model uses λ values of 6.0 for alpha-carbon atoms and 37.5 for other internal coordinates, with σ = 12.25Å.
  • Methods: The diffusion state concatenates protein alpha-carbon coordinates, remaining protein heavy-atom coordinates, and ligand heavy-atom coordinates.The transformation is designed to remove local chemical-group detail before global protein packing and ligand-interface information.
  • Methods: Sampling uses Langevin simulated annealing, Kabsch alignment correction, and a semi-analytic integrator to control temperature, rotational consistency, and discretization error.The inverse temperature increases from 1.0 to βmax = 10.0, and the integrator uses a harmonic approximation of the score function.

Supplementary Information

The supplementary materials specify the inference interface and implementation conventions used to sample multiple structures and optionally compute confidence estimates.

  • Supplementary Information: Algorithm S1 takes protein sequences, ligand molecular graphs, a requested number of conformations, 40 diffusion steps, template usage, and pLDDT computation as inputs.The algorithm summarizes NeuralPLexer inference for ensemble structure generation.
  • Supplementary Information: Implementation conventions define LinearNoBias as a trainable linear map and MLP as a three-layer GELU network with output layer normalization.These conventions standardize the supplementary architectural descriptions.

B.1 Preliminaries

NeuralPLexer inference encodes inputs, samples contact assignments autoregressively, denoises coordinates through a diffusion schedule, and can return protein- and ligand-level confidence estimates while preserving SE(3) symmetry.

  • B.1 Preliminaries: The reverse-time SDE starts from coordinates sampled from the prior distribution at t = T = 1.0 and iteratively generates the target structure distribution.This is the diffusion basis of the inference procedure.
  • B.1 Preliminaries: Inference computes ESM-2-650M protein features and optionally retrieves template features before generating conformations from prior-sampled protein coordinates.Ligand molecular graphs are incorporated through molecular graph embeddings when present.
  • B.1 Preliminaries: The contact prediction module converts 32-bin distance predictions into contact maps, aggregates them into patches, and samples protein–ligand assignments.The sampled assignments are used as additional signals in subsequent iterations.
  • B.1 Preliminaries: The denoising loop regenerates residue- and atomic-scale graphs, predicts clean coordinates, applies an inverse-temperature schedule, and advances the integrator for each diffusion step.After sampling, a confidence head can produce per-residue and per-ligand-atom pLDDT values.
  • B.1 Preliminaries: The forward diffusion process is constrained to be SE(3)-equivariant so global translations and rotations transform consistently with the molecular coordinates.The reverse process preserves SE(3) equivariance but not parity inversion because the data distribution is physically chiral.

B.3 Derivation of the semi-analytic integrator (Eq. 12)

The section describes a semi-analytic reverse-time SDE integrator and its connection to standard variance-preserving diffusion and DDIM formulations. It also outlines molecular graph representations based on atoms, bond-defined frames, stereochemistry, and short-range graph structure.

  • Semi-analytic integrator: The integrator approximates the score function as linear in the attraction term over each integration interval.The interval runs backward from s = t + ∆t to t.
  • Semi-analytic integrator: Analytic integration and matching of conditional Gaussian transition moments recover Equation 12 in a rotation-corrected multivariate setting.The DDIM integrator is recovered by removing noise and setting β ≡ 0, σ ≡ 1, and θ ≡ 1.
  • Diffusion variance: The variance term is adapted to the forward SDE, corresponding to the standard variance-preserving SDE.The approximation is Var(z_t) ≈ σ².
  • Molecular representations: Molecular inputs combine atom and bond-defined frame representations with stereochemistry encodings, graph-distance features, and frame–atom pair features.Frame representations encode adjacent-bond geometry, while graph features include paths of up to three chemical bonds.

C.2 The MHT network architecture

The MHT architecture propagates information across atom and frame nodes using heat-kernel-style pair updates, graph attention with edge bias, and residual node and edge updates. It then returns updated molecular representations, optionally including subsampled ligand-pair features.

  • MHT architecture: MHT propagates both node and edge representations before uniformly sampling ligand frame nodes as anchor nodes for subsequent processing.The sampled anchors are represented by an incidence matrix over all frame nodes.
  • Output representations: Residual note updates refine node and edge features, after which MHT returns atom, frame, atom–atom, and frame–atom representations.The inference algorithm specifies Nblocks = 8 and K = 8.
  • Pair updates: The pair-update block computes normalized nearest-neighbor frame affinities and expands them with an approximate matrix exponential heat kernel.The resulting kernels update frame–atom pair representations.
  • Graph attention: Graph attention concatenates atom and frame node features, combines multiple edge representations, and applies multi-head attention with edge bias.Edge bias enters attention as a relative positional encoding term.

C.3 MHT model pretraining

MHT pretraining combines 3D molecular geometry objectives with chemical-checker regression, masking classification, and masked-token prediction. The encoder uses atom, frame, stereochemical, graph-distance, and protein-residue representations to support molecular structure learning.

  • Pretraining objectives: The MHT pretraining objective combines marginal 3D, denoising score-matching, chemical-checker regression, masked-entry classification, and masked-language-model losses.The combined objective weights the classification loss by 0.01 and the masked-language-model loss by 0.1.
  • 3D geometry learning: A mixture-density head aligns learned pair representations with intramolecular 3D coordinate marginals using frame-aligned coordinates and power-spherical mixture components.Frames are constructed from three points with numerical stabilization in vector-norm calculations.
  • 3D geometry learning: The 3D geometry head denoises coordinates sampled from a variance-preserving SDE using an SE(3)-invariant FAPE-based score-matching loss.The geometry head uses an equivariant graph transformer with receptor nodes removed.
  • Training procedure: MHT training masks atom, bond, edge, and stereochemistry encodings at a 50% ratio with dropout = 0.1 and uses batch size 32 for 1.5×10^6 iterations.These settings are reported for the chemical encoder pretraining procedure.
  • Protein representations: Protein residue features combine amino-acid identities, ESM-2-650M embeddings, and optional template coordinates, while residue edges encode sequence and noisy or template backbone geometry.Residue-edge initialization uses node-feature outer sums, sequence-position encodings, and relative geometric encodings.

E The CPM network architecture

The CPM architecture is a six-block residue-scale network that incorporates diffusion-time and sampled ligand-contact information. It first processes protein-only edges, then integrates ligand interactions across the full residue-scale graph.

  • CPM depth: All reported models use NCPM = 6 stacked contact prediction module blocks.The blocks constitute CPMForward.
  • Conditioning inputs: The network injects a 64-bit random Fourier encoding of diffusion time into frame features and adds sampled block-adjacency encodings to protein–ligand edge representations.The adjacency encodings connect patch-wise protein nodes with selected ligand nodes.
  • Block operations: Figure S1 represents information flow through a CPM block, including triangular gated self-attention and element-wise tensor summation.The local residue–residue neighbor size is kNN = 32 in all reported models.
  • Graph schedule: The first 2 CPM blocks operate on the protein subgraph, while the remaining 4 operate on the complete residue-scale graph containing protein and ligand edge types.The full graph includes BB, BS, BF, SS, SF, and FF edges.

F The ESDM network architecture

The ESDM generates denoised three-dimensional complex structures from noisy coordinates and graph representations through stacked equivariant blocks. Its local updates preserve rigid-motion structure while retaining chirality-sensitive behavior.

  • ESDM overview: ESDM predicts denoised structures from noisy coordinates and graph representations, using a stack of four ESDM blocks.The final denoised coordinates are used to infer the score function.
  • ESDM overview: Diffusion-time Fourier encodings are concatenated with atomic node representations before the ESDM forward pass.The encoded diffusion time is embedded by a standard MLP.
  • Attention and local updates: Point-set attention combines scalar, vector, edge, and coordinate features to compute graph attention weights.Its attention-weight expression is adapted from Invariant Point Attention.
  • Attention and local updates: Local updates transform vector features into reference frames, process concatenated scalar and vector features with an MLP, then rotate the outputs back.The procedure uses reference rotations and channel-wise vector gating.
  • Equivariance and chirality: The reference-frame representation is SE(3)-invariant but generally not parity-invariant, enabling predicted coordinates to capture chiral symmetry-breaking behavior.This property follows from using a pseudovector in the rotation matrix.

G Confidence estimation heads

NeuralPLexer estimates confidence for generated protein and ligand structures with a separate pLDDT head. The head predicts local distance-error distributions around protein residues and ligand atoms.

  • Confidence-head architecture: The pLDDT head uses the ESDM architecture with independent parameters after a structure is generated by the main network.It first uses the pretrained contact prediction module to generate residue-scale graph embeddings.
  • Confidence prediction: A 6-layer MLP predicts distance-error histograms for each protein residue centroid and ligand atom.The histograms describe deviations between generated and reference pairwise distances for target points within 15.0 Å.
  • Training setup: Training used six NVIDIA-Tesla-V100-SXM2-32GB GPUs with a total batch size of 6, data parallelism, and automatic mixed precision.The training details also specify cosine annealing and exponential moving average procedures.

H Model training details

NeuralPLexer training combines contact and distance prediction, denoising, geometric, clash, and confidence-estimation objectives. The procedure samples templates, diffusion states, contacts, and perturbed structures to train its multi-scale generative pipeline.

  • Objective functions: FAPE aligns relative point coordinates across protein and ligand frames, while distance-based RMSD is evaluated globally, near binding sites, across protein-ligand pairs, or in template-variable regions.These terms provide complementary structural supervision across molecular scales.
  • Geometric regularization: Distance-geometry and steric-clash losses regularize atom-pair distances and penalize geometrically invalid structures.The definitions use indicator functions and atom-specific van der Waals radii.
  • Model scale and evaluation: The NeuralPLexer main network has 65.8×10^6 trainable parameters, and the pLDDT head adds 5.0×10^6 parameters.The main network includes MHT, CPM, ESDM, and projection layers.
  • Training procedure: The training algorithm supports template-conditioned or template-free inputs, randomly selecting templates for non-rigid receptors.Templates may be retrieved from AlphaFold2 predictions and PL2019-74k structures.
  • Training procedure: Training encodes ligand molecular graphs, samples diffusion states, perturbs reference structures, and constructs residue-scale and atomic-scale graphs.The procedure also derives distance and contact-map targets from the reference structure.
  • Objective functions: The main loss combines distogram, contact-map, FAPE, distance-based RMSD, distance-geometry, and clash terms with time-dependent weighting.Prior training emphasizes distogram and contact-map losses, while later training uses the full weighted objective.
  • Model scale and evaluation: Blind docking and binding-site recovery use 25 diffusion integrator steps, whereas end-to-end structure prediction uses 100 steps.Evaluation also specifies procedures for AlphaFold2, EquiBind, RosettaLigand, and structural metrics.

J Supplementary results

Supplementary evaluations compare NeuralPLexer with AlphaFold2 across apo, holo, and ligand-induced conformational-change targets. NeuralPLexer obtains the strongest average accuracy on both reported datasets, with additional gains from selecting the best sampled structure per target.

  • AlphaFold2 comparison: The AlphaFold2 comparison uses six models with specified MSA sizes, extra-MSA sizes, recycling cycles, and model versions.Model 0 uses AlphaFold-ptm:3, while models 1–5 use AlphaFold-ptm:1 through AlphaFold-ptm:5.
  • Sampling strategy: NeuralPLexer samples eight structures for each apo protein or protein-ligand pair using independent random seeds.The six AlphaFold2 structures serve as template inputs for six independent NeuralPLexer prediction sets.
  • Results: NeuralPLexer achieves the highest average prediction accuracy on both datasets, measured by the highest TM-score and lowest backbone RMSD.This comparison is reported across the apo and ligand-bound evaluation sets.
  • Results: Selecting the highest-TM-score sampled structure for each target further improves performance against AlphaFold2.The authors interpret this as better coverage of experimental structures than randomizing AlphaFold2 predictions.
  • Evaluation datasets: The supplementary end-to-end evaluation reports predictions for the PocketMiner dataset and 118 recent targets with ligand-induced conformational changes.PocketMiner contains 33 holo and 29 apo structures.
Loading 2209.15171v2…