• Docs »
• SZYBKI Theory

# SZYBKI Theory¶

## Force Field¶

The MMFF potential expression is:

$V_{MMFF} = \sum_{b}V_{b} + \sum_{a}V_{a} + \sum_{s}V_{s} + \sum_{o}V_{o} + \sum_{t}V_{t} + \sum_{v}V_{v} + \sum_{c}V_{c}$

where the seven terms respectively describe bond stretching (b), angle bending (a), stretch-bending (s), out-of-plane bending (o), torsion (t), Van der Waals (v) and electrostatic (c) interactions. Their functional forms are given below.

### Bond stretching¶

For a bond b between atoms i and j the stretching potential is:

$V_b = 143.9325 \frac{k_{ij}}{2} \Delta r_{ij}^2(1 + c_s\Delta r_{ij} + 7/12c_s^2 \Delta r_{ij}^2)$

where $$k_{ij}$$ is the force constant, $$\Delta r_{ij}$$ is the difference between actual and reference bond lengths, and $$c_s$$ is so called “cubic-stretch” constant for which the value is -2 Å$$^{-1}$$.

### Angle bending¶

The bending potential of a bond angle a made by the bonds between atoms i, j and atoms j, k is given by:

$V_a = 0.043844\frac{k_{ijk}}{2} \Delta \vartheta_{ijk}^{2}(1+c_b \Delta \vartheta_{ijk})$

where $$k_{ijk}$$ is the force constant, $$\Delta \vartheta_{ijk}$$ is the difference between actual and reference bond angles, and $$c_b$$ is the so called “cubic-bend” constant for which the value is $$-0.4 rad^{-1}$$.

### Stretch-bend interaction¶

The coupling between the stretching potential of two bonds forming a bond angle and bending that angle is described by:

$V_s = 2.5121(k_{ijk} \Delta r_{ij} + k_{kji} \Delta r_{kj})\Delta \vartheta_{ijk}$

where $$k_{ijk}$$ and $$k_{kji}$$ are force constants which couple stretches of i-j and k-j to the i-j-k bend respectively. $$\Delta r$$ and $$\Delta \vartheta$$ are defined above.

### Out-of-plane bending¶

For a trigonal center j, the potential of displacement for an atom l bonded to atom j out of plane i-j-k is:

$V_o = 0.043844\frac{k_{ijkl}}{2} \chi_{ijkl}^2$

where $$k_{ijkl}$$ is the force constant and $$\chi_{ijkl}$$ is an angle formed by the bond j-l and the plane i-j-k.

### Torsion interactions¶

For every four bonded atoms i-j-k-l the torsion interaction is described by the term:

$V_t = 0.5(V_1(1 + \cos\Phi) + V_2(1 - \cos2\Phi) +V_3(1+ \cos3\Phi)$

where $$V_1$$, $$V_2$$ and $$V_3$$ are constants depending on atoms i, j, k, l, and $$\Phi$$ is the dihedral angle formed by bonds i-j and k-l.

### Van der Waals interactions¶

For a pair of atoms i and j separated by three or more bonds, where the distance between them is $$r_{ij}$$, MMFF adopts the following Van der Waals potential:

$V_v = \epsilon_{ij} \left(\frac{1.07R_{ij}}{r_{ij}+0.07R_{ij}} \right) ^7 \left( \frac{1.12R_{ij}^7}{r_{ij}^7+0.12R_{ij}} -2 \right)$

where $$R_{ij}$$ and $$\epsilon_{ij}$$ are defined as follows:

$R_{i} = A_i\alpha_i^{0.25}$
$R_{ij} = 0.5(R_i + R_j)(1+B(1-\exp(-12\gamma_{ij}^2)))$
$\gamma_{ij} = (R_i - R_j)/(R_i + R_j)$
$\epsilon_{ij} = \frac{181.16G_iG_j \alpha_i \alpha_j}{(\alpha_i/N_i)^{0.5} + (\alpha_j/N_j)^{0.5}} \frac{1}{R_{ij}^6}$

where $$\alpha_i$$ is atomic polarizability of atom i, B is 0.2 or 0.0 if one of the atoms is polar hydrogen, $$N_i$$ and $$N_j$$ are the Slater-Kirkwood effective numbers of valence electrons, $$G_i$$, $$G_j$$ and $$A_i$$ are scale factors.

### Electrostatic interactions¶

The electrostatic interaction between two charged atoms i and j separated by at least three bonds is calculated from the standard Coulombic expression:

$V_c = f\frac{332.0716q_iq_j}{D(r_{ij} + \delta)}$

where D is the dielectric constant for which the default value is 1, $$q_i$$ and $$q_j$$ are the MMFF partial charges on atoms i and j, $$r_{ij}$$ is the interatomic distance, $$\delta$$ is the “electrostatic buffering” constant of 0.05 Å. Scaling factor f is 0.75 for 1,4 interactions, and 1.0 otherwise.

### Protein-ligand Amber potential¶

Protein-ligand interaction can be described by the Amber force field instead of MMFF potential:

$V_{Amber} = \sum_{i,j} \left\{ \epsilon_{ij}\left[\left(\frac{R_{ij}}{r_{ij}}\right)^{12} -2\left(\frac{R_{ij}}{r_{ij}}\right)^6\right] + \frac{q_iq_j}{4\pi\epsilon_0r_{ij}}\right\}$

where the summation is over protein-ligand atom pairs. $$r_{ij}$$ is the interatomic distance, $$R_{ij}$$ is the Vdw distance for a pair of atoms, $$q_i$$ and $$q_j$$ are the Amber partial charges on atoms i and j, and $$\epsilon_0$$ is the vacuum permittivity. Vdw parameters $$R_{ij}$$ and $$\epsilon_{ij}$$ and partial charges are taken from ([Wang-2000]).

## Entropy evaluation¶

Ligand entropy is evaluated as a sum of configurational entropy ($$S_c$$) and solvation entropy ($$\Delta S_s$$):

$S = S_c + \Delta S_s$

### Configurational entropy¶

Configurational entropy is calculated as:

$S_c = kN \Bigg[ 1 + \ln \Bigg(\frac{q}{N} \Bigg) + \frac{T}{q} \frac{\partial q}{\partial T} \Bigg]$

where q is the conformation dependent partition function:

$q = q_t \sum_{i=1}^{n_c} e^{-\frac{\epsilon_i}{kT}} q_{iv} q_{ir}$

Here $$q_t$$, $$q_ir$$ and $$q_iv$$ are the translational, rotational and vibrational partition functions respectively, $$n_c$$ is the number of unique conformations in the ensemble. All 3 partition functions are calculated from the classical statistical mechanics expressions which could be found in [McQuarrie-1976]. Vibrational frequencies for each conformation, needed for evaluation of $$q_{ir}$$ are derived from diagonalization of a Hessian matrix obtained from Quasi-Newton optimization when convergence is achieved. Eigenvalues $$\lambda_i$$ of the mass-weighted Hessian:

$\mathbf{H}^{m} = \mathbf{M^{-1/2}HM^{-1/2}}$

are converted into wavenumbers $$\tilde{\nu}_i$$ according to:

$\tilde{\nu}_i = \frac{1}{2\pi c} \sqrt{\lambda_i}$

### Solvation entropy¶

Solvation entropy is split into electrostatic and hydrophobic parts:

$\Delta S_s = \Delta S_{s,elec} + \Delta S_{s,hyd}$

The electrostatic part of solvation entropy is divided in to the bulk component and tight electrostatic polar solute - water interactions (hydrogen bonds). The bulk contribution is estimated from the temperature dependence of the solvent dielectric constant as:

$\Delta S_{s,elec\_bulk} = -\bigg(\frac{\partial \Delta G_s}{\partial \epsilon_{solv}}\bigg) \bigg(\frac{\partial \epsilon_{solv}}{\partial T}\bigg)$

The second term of the electrostatic solvation entropy is estimated as a constant of 28 J/(mol K).

The hydrophobic term, $$\Delta S_{s,hyd}$$, is evaluated as:

$\Delta S_{s,hyd} = -\bigg(\frac{\partial \Delta G_{s,hyd}}{\partial T}\bigg)$

where $$\Delta G_{s,hyd}$$ consists of 3 components:

$\Delta G_{s,hyd} = \Delta G_{cav} + \Delta G_{VdW} + \Delta G_{Ind}$

describing the free energy of cavitation, solute-solvent van der Waals and inductive terms respectively. The cavity formation term is calculated from Scaled Particle Theory [Pierotti-1976]. Analytical expressions for $$\Delta G_{VdW}$$ and $$\Delta G_{Ind}$$ terms are also taken from the 1976 Pierotti review.

### Protein-bound ligand entropy¶

Configurational entropy of a protein bound ligand is calculated totally as vibrational entropy for 3N modes, assuming that 3 rotational and 3 translational degrees of freedom of a solution ligand are transformed into low-vibrational degrees of freedom for the bound ligand.

Solvation entropy for a ligand in the active site is assumed to be a sum of its fractional value in solution determined by the percentage of the ligand surface exposed to the solvent, $$f\Delta S_s$$, and a partial desolvation entropy of the protein active site, $$\Delta S_{des}$$

$S_{protein} = S_v + f\Delta S_s + \Delta S_{des}$

where f is the fraction of ligand surface exposed to the solvent. It is important to notice that $$S_{protein}$$ is not an experimentally measurable value, and that only the difference between $$S_{protein}$$ and $$S = S_c + \Delta S_s$$ might be compared with experimental binding entropy.