The BindEnergy command calculates the binding energy of the selected object with respect to the rest of the soup according to the current force field .

The binding energy is obtained by calculating the energy at infinte separation of the selected object (unbound state) and subtracting the energy of the soup (assumed to be the bound state). The more positive the binding energy, the more favorable the interaction in the context of the chosen force field. All involved structures should be minimized before calculating binding energies.

Results are returned in units of either kJ/mol or kcal/mol, depending on the current EnergyUnit.

If the two binding partners belong the the same object , e.g. a ligand and an enzyme, they have to be separated first following these steps:

• Check that there is a split point (a single vertical line in the sequence selector representing a 'TER' in the PDB file) separating the two partners. If not, add one.
• Split the object into multiple childs.
• If the protein is an oligomer and contained additional split points, you will have to join the peptide chains again.

The BindEnergy command does not include solvation effects. These can be added by performing the calculation steps individually, for example with this macro:

# The scene MyComplex must contain receptor and ligand as objects 1 and 2
# Also the structures must have been energy minimized before
EnergyUnit kcal/mol
LoadSce MyComplex
# Make sure that no water is present, we are using implicit solvent
DelRes HOH
ForceField Yamber2
# Choose a large cell, so that periodic boundaries have minimum impact on energies
Cell Auto,Extension=50
Boundary periodic
# Calculate potential and solvation energy of the complex
epotcomp=Energy
esolcomp=SolvEnergy
# Add Van der Waals solvation energy (empirical formula from
# 'Computational Medicinal Chemistry for Drug Discovery', page 271)
surfacc=Surf Accessible
esolcomp=esolcomp+surfacc*0.005+0.85
# Calculate potential and solvation energy for each of the two
for i=1 to 2
  RemoveObj not (i)
  epot(i)=Energy
  esol(i)=SolvEnergy
  surfacc=Surf Accessible
  esol(i)=esol(i)+surfacc*0.005+0.85
  AddObj all
# Calculate result: energy of separated compounds - energy of complex
ebind=(epot1+epot2+esol1+esol2)-(epotcomp+esolcomp)
print 'Binding energy is (ebind)'

This macro considers the potential energy obtained with the current force field, as well as the electrostatic solvation energy obtained with the SolvEnergy command. It ignores all entropic components, which makes it impossible to obtain values that can be quantitatively compared to experimental results. Energies are also shifted by an unknown amount that depends on the protein, it is thus not possible to say that energies > X (e.g. >0) indicate binding while those < X don't. Nevertheless, the obtained energies are often a helpful qualitative estimate to compare affinities of different ligands and the effect of point mutations on binding, as shown in the example below.

This macro can be extended to predict the effect of point mutations on binding. In the following example, residue Lys 60 of the PQQC enzyme (PDB ID 1OTW, the final catalyst in the production of pyrroloquinoline quinone) is mutated, and the mutant binding energies are calculated:

Clear
# Load the structure
LoadPDB 1otw
# Work only with chain A, delete residues 501, 600, 601, the waters plus molecule (=chain) B
DelRes 501 600 601 HOH or Mol B
# To calculate binding energies, the ligand must be in a separate object
# The 'Split' command splits the object at 'TER' entries in the PDB file, ligand is now object 2
Split
# Set the wild type residue (wtres) to mutate
wtres='Lys 60'
# The mutants, wt must appear again as the first one
mutlist='Lys','Arg','Val','Ala','Glu'
# Set the force field parameters
ForceField Amber99
Cutoff 7.86
Longrange Coulomb
# Add a simulation cell for energy minimization
Cell Auto,Extension=20
Boundary Periodic
# Prepare for simulation, add missing hydrogens etc.
Clean
Style Ribbon,Stick
FixAtom Backbone
# Energy minimize the complex
Experiment Minimization
Experiment On
Wait ExpEnd
# Keep this minimized state for reuse by saving the scene
SaveSce Complex
# First, calculate potential (epot) and solvation energy (esol) of the ligand alone
RemoveObj 1
epotlig=Energy
esollig=SolvEnergy
# Try all mutations
for i=1 to count mutlist
  # Load the original scene saved above
  LoadSce Complex
  mutres=mutlist(i)
  ShowMessage 'Analyzing binding energy for mutant (wtres) -> (mutres)'
  # Fix everything not close by
  FixAtom all with distance>15 from Res (wtres)
  # Make the mutation
  SwapRes (wtres),(mutres)
  # Energy minimize the environment of the residue, so that the mutant can adapt
  Experiment Minimization
  Experiment On
  Wait ExpEnd
  # Calculate potential and solvation energy of the complex
  ChargeObj 1,0
  epotcmp(i)=Energy
  esolcmp(i)=SolvEnergy
  # Energy minimize the environment without the ligand
  RemoveObj 2
  Experiment Minimization
  Experiment On
  Wait ExpEnd
  # Calculate potential and solvation energy of the receptor alone
  ChargeObj 1,0
  epotrec(i)=Energy
  esolrec(i)=SolvEnergy
# Print the results
Console Off
print 'Results of the mutation analysis, values higher than the wildtype indicate better binding:'
print '[Quantitative free energies will be supported in YASARA Structure]'
for i=1 to count mutlist
  print 'Mutation (wtres) -> (mutlist(i)):'
  print '  Receptor: Potential energy = (0.00+epotrec(i)), Solvation energy = (0.00+esolrec(i))'
  print '  Ligand:   Potential energy = (0.00+epotlig), Solvation energy = (0.00+esollig)'
  print '  Complex:  Potential energy = (0.00+epotcmp(i)), Solvation energy = (0.00+esolcmp(i))'
  ebind=epotrec(i)+esolrec(i)+epotlig+esollig-epotcmp(i)-esolcmp(i)
  print '  Resulting binding energy   = (0.00+ebind)'
Console On

A much more sophisticated quantitative analysis of binding energies will become available in YASARA Structure.

BindEnergyObj 3

Calculate binding energy of object 3 to the rest of the soup.

nrg = BindEnergyObj 4

Assign binding energy of object 4 to variable 'nrg'.