0.
...
1. Tracing the origin of the amide issue
...
(1) Addition of new amide-specific torsion parameters
Torsion profile of N-methylacetamide amide bond generated using v1.2.0 fails failed to fit the expected minimum (parabola shape), giving reproduce the QM torsion profile(left plots), forming a small hump at the minimum. Re180 degree, which doesn’t appear in the QM profile. Energy decoupling(bottom left plot) showed that torsion term contributes the most to the formation of the hump, which indicated re-fitting the related torsion parameters might solve the issue.
Before this amide issue was raised, I found that there’s a missing torsion term for dialkyl amide in the current parameter set and to fix the dialkyl amide issue(for more details, please refer to here), we’ve decided to add new amide-specific torsion parameters. The preliminary study showed that the addition of the new torsion terms improved the performance of reproducing QM optimized geometries of dialkyl amides.
To see if the new amide-specific torsion terms would also solve the amide issue, re-fitting to the selected torsion targets which rotate an amide bond couldn’t remove the hump. targets(targets scanning amide bond dihedral) was carried out with the new amide-specific torsion terms, then the test calculation was done to see if the re-fitted parameter set improved the torsion profile of N-methylacetamide amide bond(removal of the hump + reproduction of QM basin depth at the minimum). The test calculation(right plots) showed that while the re-fitting improved the torsion energy term(bottom right plot), it deformed other terms to force the formation the hump at 180 degree.
(2) Investigation of fitting targets
The next step of the tracing was investigating targets.
...
From a closer look at the targets used to train the amide bond parameters(t69s, t70s), in the re-fitting, I found that some targets behave not as expected at around its their planar geometry geometries (having a small peak, instead of parabolic shape near the minimum. right figure), which leads the parameter set fitted to be trained to fail to reproduce the torsion profiles with the targets to make a small hump even when scanning amide bonds in simple molecules whose torsion profiles have parabolic shape near the minimum (like N-methyacetamide case).
So one of the easy fixes One easy fix of the problem will be is using simple targets which behaves exactly as we expect , whose torsion profiles have parabolic shape near the minimum (planar at 180 and small steric hindrance) and here’s what I did:
1. filtering non-planar structure at minimum geometry
explanation of scheme2 and scheme3
scheme 2 (fb-fit2):
scheme 3 (fb-fit3)
fb-fit2 vs fb-fit2-2
...
fb-fit2
...
SMIRKS
...
initial guess
...
final
...
t69a
...
[*:1]-[#7X3:2]-[#6X3$(*=[#8,#16,#7]):3]~[*:4]
...
2.5 (1+cos(2x-180))
...
2.5 (1+cos(2x-180))
...
t69b
...
[*:1]-[#7X3:2]-!@[#6X3$(*=[#8,#16,#7]):3]~[*:4]
...
2.5 (1+cos(2x-180))
...
1.652160092455e+00 (1+cos(2x-180))
...
t70
...
[#1:1]-[#7X3:2]-[#6X3:3]=[#8,#16,#7:4]
...
3.459249459574e+00 (1+cos(2x-180)) + 1.356955617521e+00 (1+cosx)
...
3.459249459574e+00 (1+cos(2x-180)) + 1.356955617521e+00 (1+cosx)
...
t70a
...
[#1:1]-[#7X3:2]-!@[#6X3:3]=[#8,#16,#7:4]
...
2.5 (1+cos(2x-180)) + 2.0 (1+cosx)
...
-2.081991778915e-01 (1+cos(2x-180)) + 1.400460484703e+00 (1+cosx)
...
t70b
...
[*:1]-[#7X3:2]-!@[#6X3:3](=[#8,#16,#7:4])-[#6,#1]
...
2.5 (1+cos(2x-180))
...
4.091207931356e+00 (1+cos(2x-180))
...
t70c
...
[#1:1]-[#7X3:2]-!@[#6X3:3](=[#8,#16,#7:4])-[#6,#1]
...
2.5 (1+cos(2x-180)) + 2.0 (1+cosx)
...
1.301446828323e+00 (1+cos(2x-180)) + 9.476830965983e-01 (1+cosx)
...
t70d
...
[*:1]-[#7X3:2]-!@[#6X3:3](=[#8,#16,#7:4])-[#7X3]
...
2.5 (1+cos(2x-180))
...
1.271374777009e+00 (1+cos(2x-180))
...
fb-fit2-2
...
SMIRKS
...
initial guess
...
final
...
t69a
geometry at the minimum). + Expecting that this experiment to be one proof of the need for using simple-as-possible torsiondrive targets in general torsion parameter fitting.
2. Filtering out non-planar structures at the minimum geometry
Several schemes to filter non-planar molecules:
Scheme 2
Check a planarity of the geometry at the minimum point of a QM torsion profile;
Planar geometry: C-center improper dihedral angle, N-center improper dihedral angle between -5.0, 5.0 degree;
number of targets selected: 51 - 1(N-methylacetamide) = 50
Scheme 3 (rough)
Check if the QM profile has local minima at 0 and 180 degree;
Number of targets selected: 38
Scheme 4
Targets passing either scheme 2 or scheme 3
Number of targets selected: 63 -1(N-methylacetamide) = 62
(Since the 2nd generation torsion training set is small to filter, i pulled Roche torsion set (1st generation) and filtered with various scheme. )
...
For example, the target #540 is filtered out in scheme 3 (hump at 0 degree), but not in scheme 2, since the geometry at the minimum is planar.
3. Comparison of different filtering schemes
Test set 1: (1) the simplest molecule, N-methyacetamide, (2) non-planar molecules filtered out in both scheme 2 and scheme 3;
Test set 2: all amide torsions in the 2nd generation training dataset;
Overall quality of fitting and test calculation:
v1.2.0 | fb-fit2 (scheme 2) | fb-fit3 (scheme3) | fb-fit4 (scheme 4) | Unfiltered | |
---|---|---|---|---|---|
test cal 1 | 144.59 | 91.81 | 82.39 | 88.84 | - |
test cal 2 | 25.37 | 15.62 | 14.75 | 14.15 | 23.21 |
test cal 2-1 simple | 7.10 | 6.21 | 7.49 | 8.29 | 11.04 |
test cal 2-2 complex | 18.26 | 9.40 | 7.26 | 6.87 | 12.17 |
final X2 (fitting) | - | 1.31e+02 → 2.22e+01 | 1.30e+02 → 1.73e+01 | 2.09e+02 → 3.01e+01 | 5.28e+02 → 6.90e+01 |
(1) N-methylacetamide
...
(2) Test set 1 (N-methyacetamide, non-planar molecules in the 1st gen. Roche torsion set)
...
(3) Test set 2 (All amide torsions in the 2nd gen. training set)
...
(4) Check performance on non-planar molecules
...
(5) Direct comparison of final parameters
F(θ) = k1*(1+cos(2θ-180)) + k2*(1+cosθ)
ID, SMIRKS | Initial guess | fb-fit2 | fb-fit3 | fb-fit4 | |||
---|---|---|---|---|---|---|---|
| k1 = 2.5 ( | k1= 1+cos(2x-180)) | 1.679285484776e+00 (1+cos(2x-180)) | t70 | .6793e+0 | k1 = 1.3647e+0 | k1 = 1.7786e+0 |
| k1 = 3. | 459249459574e+00 (1+cos(2x-180)) + 1.356955617521e+00 (1+cosx)3.566161051804e-01 (1+cos(2x-180)) + 1.354248559088e+00 (1+cosx) | t70b | 4592e+0 k2 = 1.3570e+0 (from 1.2.0) | k1 = 3.5662e-1 k2 =1.3542e+0 | k1 = 2.4868e-1 K2 = 9.4280e-1 | k1 = -1.3577e+0 k2 = 1.2663e+0 |
| k1 = 2.5 (1+cos(2x-180)) | 4.169159533591e+00 (1+cos(2x-180)) | t70c | k1 = 4.1692e+0 | k1 = 4.4151e+0 | k1 = 4.1071e+0 | |
| k1 = 2.5 (1+cos(2x-180)) + k2 = 2.0 ( | k1 = 1+cosx) | 1.240609067243e+00 (1+cos(2x-180)) + 9.065552775144e-01 (1+cosx) | t70d | .2406e+0 k2 = 9.0656e-1 | k1 = 1.5283e+0 k2 = 7.3406e-1 | k1 = 1.1159e+0 k2 = 7.6465e-1 |
| k1 = 2.5 (1+cos(2x-180)) 1.303139175986e+00 (1+cos(2x-180)) | k1= 1.3031e+0 | k1 = 1.6149e+0 | k1 = 1.1680e+0 |
* About the counter-intuitive negative k1 value of t70
, [#1:1]-[#7X3:2]-[#6X3:3]=[#8,#16,#7:4]
...
Although the k1 value for t70a seems unphysical, the result below shows the plots got improved after the optimization. checked improvement for all of the four fitting targets.
...
4. Conclusions
+ Selecting simple geometries (no strong electrostatic interaction) may need to be considered in the next torsion training dataset selection procedure.