Introduction of drugs puts an evolutionary pressure on viruses and makes the viral proteins evolve. On the other hand, the change in the protein structure should be compatible with the protein function. The degree of the protein sequence conservativity is thus a measure of the structural importance of a particular region on the protein surface. Certain binding sites are evolving and mutate often, certain are not. Finding the most conservative fractions of a target protein sequence can help identify important and druggable binding sites for inhibitors search, find drugs with the least potential for drug resistance development.

Below we provide an analysis of H5N1 neuromidase protein surface next to the tamiflu binding region. The Figure shows the level of protein structure conservation, red portions corresponds to the most conservative residues. A small yellow patch of the pocket has evolved since 1990 and gave way to develop tamiflu resistance.

The analysis shows that in spite of a lot of pressure from the drug application, the virus was not able to change the red part of the pocket. Future drugs should target the binding sites with the least mutable residues.

The free binding energy (F) of a small molecule and a protein is a non-additive complex function of individual interatomic interactions. There are two major contributing quantities leading to non-additivity in F: the electrostatic energy, and the entropy.

A common approach to molecular docking is to develop a simple (generally additive) model of intermolecular interactions and train it to reproduce experimental values of binding energy for a range of known inhibitors. The obvious advantage is the calculations speed.

Any more sophisticated approach takes requires more computational resources and may hardly lead to an immense improvement in calculations accuracy. The question thus is: do non-additive force fields have any advantages in docking situtations?

To investigate the issue the following experiment was performed:

1. A simple force Molecular Mechanics (MM) force field, containing a reasonable approximation for (distant-dependent media polarization) electrostatics, van-der-waals and hydrogen bonding, was developed.
2. Same van-der-waals and hydrogen bonds were paired with vacuum electrostatics and a (non-additive) water model to simulate solvation effects.

Both models were transformed into simple linear regressions (to avoid sophisticated thermodynamic integration) and trained to reproduce the same amount of experimental values.

Although the two models show roughly the same level of accuracy in predicting the binding energies, the selectivity is drustically different. The figure above shows the results of the energy calculations for a set of few hundreds decoys. Both graphs represent the relative binding energy (counted from the minimum position) vs. the r.m.s. deviation of the calculated ligand positions from the known native position. The lines are the energy offset values averaged over the conformers (decoys) with similar r.m.s. positions.

The conclusion?
The solvation energy is the major source of non-additivity in ligand binding. Though being one of the most complicated quantities to be acounted properly in a docking run, the solvation energy is one of the major mechanisms of molecular recognition