Versions Compared

Key

  • This line was added.
  • This line was removed.
  • Formatting was changed.

...

To generate 2-D TorsionDrives of protein backbone dihedrals in dipeptides (2021-11-18-OpenFF-Protein-Dipeptide-2D-TorsionDrive), the heavy-atom dihedral angles of the amino acid sidechains (which were not driven in the TorsionDrive) were constrained to values corresponding to highly populated sidechain rotamers. These constraints prevent large sidechain rearrangements that lead to discontinuous jumps in the QM energy profile (2021-12-16 Meeting notes). There is concern that using CONstraints instead of REstraints may introduce artificially high energy barriers during the torsion scan.

...

The categories of atom subsets are all atoms, driven atoms (atoms in phi or psi), non-driven atoms (atoms in chi1 or chi2), non-driven and not driven atoms (atoms in chi1 or chi2 but not n phi or psi), and non-constrained atoms (atoms not in phi, psi, chi1 or chi2). The units for RMS force are

Mathinline
body--uriencoded--\mathsf%7Bm%7D E_h a_0%5e%7B-1%7D
(millihartree per bohr).

Atom category

N Atoms

RMS Force, no chi constraints

RMS Force, with chi constraints

All

36

2.422

2.870

Driven

5

6.497

7.108

Non-driven

5

4.982

6.259

Non-driven and not driven

3

0.083

3.822

Non-constrained

28

0.043

0.041

The RMS force on non-constrained atoms is similar with and without constraints on sidechain dihedrals, and the RMS force on non-driven and not driven constrained atoms is only 25 % 50 times higher with constraints on sidechain dihedrals than without.

Energy differences between TorsionDrives

We can also compare the QM energy profiles with and without constraints on sidechain dihedrals. First, we shift each energy profile by its own minimum energy so that the minimum energy is zero for both profiles. The minimum energy occurs at the same grid point for both profiles. Then, we compute an RMS error over grid points.

Mathblock
\mathsf{RMS\ Error} = \left( \frac {1} {N_{\mathsf{grid}}} \sum_{\mathsf{grid}} \left( E_1 - E_2 \right)^2 \right)

The RMS error is 2.97 kcal/mol−1. However, in a ForceBalance fit, grid points are weighted as a function of their shifted QM energies in kcal mol−1.

Mathblock
w(E) = \begin{cases} 1 & E < 1.0 \\ \left( 1 + (E - 1)^2 \right)^{-1/2} & 1.0 \leq E < 5.0 \\ 0 & 5.0 \leq E \end{cases}

We can thus compute a weighted RMS error.

Mathblock
\mathsf{Weighted\ RMS\ Error} = \left( \frac {1} {N_{\mathsf{grid}}} \sum_{\mathsf{grid}} \left( w(E_1) \left( E_1 - E_2 \right) \right)^2 \right)

With the weights defined using the shifted QM energies of the TorsionDrive without constraints on sidechain dihedrals, the weighted RMS error is 0.0779 kcal mol−1. Thus, the energy differences occur almost entirely in high energy regions which do not contribute to the ForceBalance objective function during torsion fitting.

Code to reproduce

The results in the table above can be reproduced with the following python code. Arguments are

...