Predicting Protein-Ligand Binding Explained

Binding prediction covers two different questions that get merged: where does the ligand sit, and how tightly does it hold on? The second is far harder, and the arithmetic below shows how little room a scoring function has to be wrong before its ranking is noise.

What is being predicted

The experimental quantities are a dissociation constant Kd, an inhibition constant Ki, or a half-maximal inhibitory concentration IC50. The first two are equilibrium constants; the third is assay-dependent and shifts with substrate concentration, which is why mixing IC50 values from different assays into one training set introduces label noise that no model can remove.

Values are usually modelled on a log scale: pKd = -log10(Kd) with Kd in molar. A Kd of 1 µM is pKd 6; 1 nM is pKd 9. Regression targets are in pKd or pKi units and errors are reported the same way, so “RMSE 1.4” means the typical prediction is off by more than a factor of twenty in concentration.

Worked: the error budget in kcal/mol

Binding affinity is a free energy, and the conversion sets the whole accuracy discussion.

dG = -RT ln Ka = RT ln Kd

at 298 K:  R = 0.001987 kcal/(mol K),  RT = 0.592 kcal/mol
           ln(10) * RT = 2.303 * 0.592 = 1.364 kcal/mol

so   dG (kcal/mol) = -1.364 * pKd

  Kd = 1 mM   -> pKd 3 -> dG = -4.1 kcal/mol
  Kd = 1 uM   -> pKd 6 -> dG = -8.2 kcal/mol
  Kd = 1 nM   -> pKd 9 -> dG = -12.3 kcal/mol

one log unit of potency = 1.36 kcal/mol
one methyl group in a hit-to-lead campaign is often worth
about that much, or less.

A scoring function that wants to rank two analogues differing by a single log unit therefore needs to be accurate to well under 1.4 kcal/mol on the difference. Empirical docking scores are not, and the published comparative assessments — the CASF series built on PDBbind is the standard one — report scoring-power correlations and errors that leave the ranking of close analogues unresolved. Consult the current CASF tables for the exact figures rather than a remembered number; they move as the benchmark and the entrants change.

This is the mechanical reason a virtual screen uses the score to select a set and not to order it. Separating micromolar from millimolar is a 4 kcal/mol question and is within reach; separating 10 nM from 100 nM is not.

What the model sees in a pocket

Classical empirical scoring functions decompose the interaction into weighted terms. AutoDock Vina, from Oleg Trott and Arthur Olson in the Journal of Computational Chemistry (2010), uses two Gaussian steric terms, a repulsion term, a hydrophobic term and a hydrogen-bonding term, all functions of the surface distance between heavy-atom pairs, plus a penalty proportional to the number of active rotatable bonds in the ligand — a crude stand-in for the entropy lost on binding. The published weights are in the paper.

Machine-learned scoring functions replace the fixed functional form. Three representations are common:

• Interaction counts. Count atom-type pairs within distance shells — carbon-to-carbon within 4 Å, nitrogen-to-oxygen within 3.5 Å — and regress on the counts. This is the RF-Score family, and it is a strong baseline for how simple it is.
• Voxel grids. Rasterise the pocket into a 3D grid with channels for atom types and run a convolutional network. Captures shape complementarity; needs pose-aligned inputs and is not rotation invariant unless you augment.
• Interaction graphs. Build a graph whose nodes are ligand and pocket atoms and whose edges include non-covalent contacts under a distance cutoff. Naturally handles variable pocket size and makes the cutoff an explicit hyperparameter.

What none of these represent well is water. Displacing an ordered water molecule from a pocket can dominate the binding energy, and structural waters are usually stripped before featurisation. Entropy is similarly absent beyond the rotatable-bond penalty, and protonation states — which histidine tautomer, whether a carboxylate is charged at assay pH — are assigned by preparation software and are frequently wrong.

The shortcuts models learn instead

The most important result in this area is negative. Mikhail Volkov and colleagues, writing in the Journal of Medicinal Chemistry (2022), dissected what deep affinity models actually use by retraining them on ligand-only and protein-only inputs. Performance held up disturbingly well without the interaction — the models were largely predicting from ligand identity and protein identity, exploiting the fact that PDBbind contains many complexes per protein family and that affinity correlates with simple ligand properties such as size.

That is shortcut learning, and it is invisible to any evaluation that does not include the control. A model that scores well on a random split of PDBbind may be recognising the protein family and recalling the typical affinity of its ligands. It will do nothing useful on a new target, which is the only case anyone cares about.

A good score with an impossible pose

Pose prediction has its own failure. Learned docking and co-folding methods can place a ligand in the right region and produce geometry that is not chemically possible — distorted bond lengths, impossible angles, atoms clashing with protein atoms, stereochemistry inverted from the input.

Martin Buttenschoen, Garrett Morris and Charlotte Deane made this concrete with PoseBusters in Chemical Science (2024), a suite of physical and chemical checks — bond geometry, planarity, internal clashes, protein-ligand distances, stereochemistry retention — applied to predicted poses. Applying it to several deep-learning docking methods showed that a substantial share of poses within the usual 2 Å RMSD success criterion failed basic physical validity. RMSD alone was hiding it.

The operational lesson is to report validity checks alongside RMSD, and to run a force-field minimisation of any predicted pose before believing a score computed on it.

How to validate one honestly

1. Split by target, not at random. Every complex for a given protein goes entirely into train or entirely into test, and better, split by sequence-identity cluster so close homologues do not straddle the boundary.
2. Run the ligand-only control. Train the identical pipeline with the protein input removed. If it performs comparably, your model is not learning the interaction, whatever its headline number.
3. Report error in the units of the label — RMSE in pKd, and the equivalent in kcal/mol — not only a correlation coefficient, which depends on the affinity range of your test set.
4. Check pose validity, not just RMSD, for anything that predicts geometry.
5. State the applicability domain: which target families are represented, which ligand size range, and what the model does outside it.

Held to that standard, the honest use of an affinity model is coarse classification rather than ranking: separating plausible binders from implausible ones at a scale no assay could cover, and handing the survivors to people who can measure them. That is a smaller claim than the literature often makes and it is the one the arithmetic at the top of this page supports.