Theory¶
Components¶
Objective Functions¶
At the core of the objective functions are pure virtual interfaces OEFunc0, OEFunc1 and OEFunc2. A new objective function can be defined by providing appropriate expressions for these interface methods. The objective functions can be used along with the optimizers to find the minima of any function.
Molecule Objective Function¶
The molecule objective functions extend the objective function interfaces to define molecular interactions through OEMolFunc. Similar to the objective functions, a new molecule objective function can also be defined by providing the appropriate expressions for the interface methods. All of the force field components are defined as molecule objective functions to allow easy minimization of molecule systems.
Adaptors¶
Adaptors wrap molecule objective functions to allow working in different coordinate systems, instead of the default all atom Cartesian coordinates. Adaptors can be powerful tool for selective energy evaluation and minimization of specific subdomains. Methods to convert between the adaptor coordinates and the Cartesian coordinates are also provided through these interfaces.
Optimizers¶
A number of general purpose optimizers are provided that can be used to minimize any basic or extended objective functions. Optimizer interfaces are defined through OEOptimizer1 and OEOptimizer2. Pros and cons of various optimizers are established in the literature, and one needs to use their understanding of the problem in choosing the appropriate optimizer. Optimizers that can solve objective functions with known gradients but no Hessian (OEOptimizer1) requires an additional line minimizer (OELineMinimize).
Force Field¶
MMFF¶
The MMFF potential expression is:
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:
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:
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:
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:
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:
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:
where \(R_{ij}\) and \(\epsilon_{ij}\) are defined as follows:
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:
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.
FF14SB¶
FF14SB is Amber primary force field for protein modeling. Its details are described in the paper of Maier et al. ([Maier-2015]) and in http://ambermd.org/AmberModels.php. The functional form of FF14SB is described under the chapter “AMBER force field terms”.
SMIRNOFF¶
SMIRNOFF (SMIRKS Natural Open Force Field) is a kind of small-molecule force field
springing from the Open Force Field Initiative (OFF Initiative or OFFI)
(https://openforcefield.org/). The key feature of this force field is its chemical
representation of parameters by using extended SMARTS patterns. “Extended”
here means that the SMARTS patterns also use atom assignment syntax from SMIRKS. For example,
[*:2]~[*:1]~[#1] specifies three atoms, where the first is assigned to be atom 2,
the second atom is assigned to be atom 1, and the third atom is not assigned. This format
frees the force field from using atom types, instead it uses the direct chemical perception
of the extended SMARTS patterns ([Mobley-2018]).
Compared to traditional atom typing approaches, this representation greatly reduces the
number of parameters needed to specify a general small-molecule force field and makes
each parameter more independent.
While it currently has a traditional “Amber” type functional form as specified below, future versions expect to modify the functional form. The OFF Initiative will be releasing periodic revisions of the parameters which will include changes to the chemical specifications (i.e. the extended SMARTS patterns) as well as the numerical parameters. The revision used by the API points will be given with the results.
Currently SMIRNOFF force field has 6 similar to those in AMBER force field and are given below.
AMBER force field terms¶
Amber force field functional form has 6 terms:
describing respectively bond stretching (b), angle bending (a), out-of-plane displacement (o), torsion (t), Van der Waals (v) and electrostatic (c) interactions.
Bond stretching¶
For a bond b the stretching potential has a simple harmonic form:
where \(k_b\) is the force constant, \(r\), and \(r_0\) are actual and reference bond length, respectively.
Angle bending¶
For bond angle a, the angle bending potential form is:
where \(k_a\) is the force constant, \(\alpha\), and \(\alpha_0\) are actual and reference bond angle values, respectively.
Torsion interactions¶
Every 4 bonded atoms defined a torsion whose energy is represented with 1 to 6 terms of a Fourier series:
where \(V_n\), \(d_n\), \(p_n\) and \(\gamma_n\) are barrier height, degeneracy, periodicity and phase for the n component of a torsion whose dihedral angle is \(\theta\)
Out-of-plane bending¶
It is defined as potential of a displacement of the trigonal center atom a bonded to 3 other out-of-plane atoms, b, c and d. Its functional form is similar to a regular torsion with a single term. Sometimes it is also called an improper torsion.
While in FF14SB this term is used to preserved planarity of the sp2 configuration, any improper torsion moving a central trigonal atom out of plane can be used, in SMIRNOFF this term has been redefined to account for the fact that one can define 3 improper torsions of moving a trigonal atom a out of plane defined by atoms b, c and d. The energy of such a displacement in SMIRNOFF is calculated as 1/3 of the sum of those torsions.
Van der Waals interactions¶
For a pair of nonbonded atoms i and j, SMIRNOFF force field adopts Lennard-Jones potential:
where \(\epsilon_{ij}\) is well depth, \(R_{ij}\) and \(r_{ij}\) are the equilibrium and actual distance between atoms i and j, respectively.
Electrostatic interactions¶
A simple Coulomb term:
represents Coulomb interactions, where \(q_i\) and \(q_j\) are partial atomic charges, \(\epsilon_0\) is the vacuum permittivity, and \(r_{ij}\) is the distance between atoms i and j.
MMFFAmber¶
Intermolecular Amber potential¶
Intermolecular interaction can be described by the Amber force field instead of MMFF potential:
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 work by Wang and coworkers ([Wang-J-2000]).
MMFFIEFF¶
Intermolecular IEFF potential¶
Intermolecular interaction can be described by the IEFF (Intermolecular Force Field) instead of MMFF potential.
The IEFF is an ongoing development effort internally at OpenEye to obtain high accuracy Intermolecular Force Field. The initial development of IEFF was done during 2010. The goal of the initial development was to prove internally at OpenEye as well as to external collaborators that point charge monopole approximation in the Coulomb terms are inadequate to capture intermolecular Coulomb interactions within chemical accuracy (within about 1 kcal/mol) regardless of the charge models used. The higher accuracy in the Coulomb terms is especially important and vital in optimizing and ranking organic crystal structures. The initial IEFF development has had some successful crystal predictions. It has also clearly shown, for the formamide dimer case, that there is a 1D energy curve along a given torsion angle where all tested point charge models fails even qualitatively showing repulsive interactions while accurate ab-initio calculations and IEFF with multipole based Coulomb interactions show a minimum along the torsion with great agreement to each other.
The IEFF has not been improved significantly since those initial academic studies and cannot be considered as a general purpose intermolecular force field due to lack of sufficient coverage. However, major development effort is underway to make IEFF an industry quality general intermolecular force field.
The currently available IEFF supports only a few selected elements that were the most important ones for organic crystal predictions, including H, C, N, O, S, F, and Cl. The atomic multipoles are obtained by an internal analysis program that obtains the SCF converged density matrices from ab-initio (PSI4) calculations. The current IEFF does not use any atom types, all VdW parameters are associated to each element types. This obviously introduces a very significant weakness especially for hydrogens where both polar and apolar hydrogens share the same FF parameters but it is also important deficiency for other elements as well. The new version of the IEFF that is under development will try to address all of these deficiencies.
Sheffield Model¶
Small molecules¶
The “Sheffield Solvation Model” [Grant-2007] is a simple two parameter model for calculating the electrostatic component of molecular solvation energy. The Sheffield model is robust and works well for calculating solvation energies of small drug-like molecules with good accuracy, when compared against the more accurate Poisson–Boltzmann calculations, and requires significantly low computational resource. The model is represented by the following simple expression:
where \(f_{\epsilon}\) is
and \(\epsilon_{in}\), \(\epsilon_{out}\) are he relative permittivities (dielectric constants) in solute and solvent, respectively. \(Q_i\) and \(Q_j\) are charges on atoms i and j, \(\sigma_i\) and \(\sigma_j\) are corresponding atomic radii, and \(R_{ij}\) is an interatomic distance between atoms i and j. \(A\) and \(B\) are parameters to the model.
This model works well in the situation where atoms in molecules are well exposed to the solvent. In macromolecules like proteins, this condition is not met.
Proteins¶
In order to apply Sheffield model for proteins and complexes a modification is required to correct the atomic radii according to their exposure to solvent. The correction to atomic radii is applied by introduction of a third parameter to the original two parameter model.
where \(S_i\) and \(S_{0,i}\) are atomic contributions to the solvent accessible surface of atom \(i\) in the protein and in a residue to which it belongs which was isolated from the rest of the protein. Initial atomic radii \(\sigma_{0,i}\) was assumed to be of ZAP9 [Nicholls-2010] type.