search

Valdes-Tresanco-MS / gmx_MMPBSA

gmx_MMPBSA is a new tool based on AMBER's MMPBSA.py aiming to perform end-state free energy calculations with GROMACS files.
https://valdes-tresanco-ms.github.io/gmx_MMPBSA/
GNU General Public License v3.0
226 stars 66 forks source link

[Bug-gmx_MMPBSA]: mmpbsa produces a different enthalpy when repeated on the same trajectory #380

Closed laurenreid1 closed 1 year ago

laurenreid1 commented 1 year ago

Bug summary

I have run 4 ST mmpbsa calculations on a protein-protein complex using the same trajectory and the same general settings. 1 calculation was run on the wild type (WT) trajectory and the other 3 calculations were run on the same WT trajectory with the same settings but also incorporated an alanine scan calculation in each one. My question relates to the WT calculation produced by each gmx_mmpbsa run. I observed that 2 of the calculations (the original WT-only calc and the WT calc produced by one of the alanine scan runs) produced exactly the same enthalpy, however the other 2 calculations produced quite different numbers for the enthaply. Specifically, the EEL and EPB contributions are different, and the remaining terms are the same.

If an mmpbsa calculation is run on the same trajectory with exactly the same settings, shouldn't the result be the same each time, or is there are a random seed used in the calculation that could be causing this discrepancy? If it's caused by a difference in a random seed, how should I interpret the results?... is it a sign of lack of convergence in my trajectory?

The 4 WT results are as follows:

1) WT-only calculation:

-------------------------------------------------------------------------------
Delta (Complex - Receptor - Ligand):
Energy Component       Average     SD(Prop.)         SD   SEM(Prop.)        SEM
-------------------------------------------------------------------------------
ΔBOND                     0.00          8.89       0.00         0.40       0.00
ΔANGLE                    0.00         11.85       0.00         0.53       0.00
ΔDIHED                    0.00          6.17       0.00         0.28       0.00
ΔUB                       0.00          2.89       0.00         0.13       0.00
ΔIMP                     -0.00          3.31       0.00         0.15       0.00
ΔCMAP                     0.00          3.68       0.00         0.16       0.00
ΔVDWAALS                -61.14          7.59       4.99         0.34       0.22
ΔEEL                   -460.68         43.29      40.77         1.93       1.82
Δ1-4 VDW                  0.00          4.79       0.00         0.21       0.00
Δ1-4 EEL                  0.00         20.99       0.00         0.94       0.00
ΔEPB                    486.02         24.22      35.35         1.08       1.58
ΔENPOLAR                 -8.87          0.26       0.32         0.01       0.01
ΔEDISPER                  0.00          0.00       0.00         0.00       0.00

ΔGGAS                  -521.82         44.32      39.57         1.98       1.77
ΔGSOLV                  477.15         24.22      35.22         1.08       1.57

ΔTOTAL                  -44.67         50.51      12.70         2.25       0.57
-------------------------------------------------------------------------------

2) WT calculation from alanine scan 1:

-------------------------------------------------------------------------------
Delta (Complex - Receptor - Ligand):
Energy Component       Average     SD(Prop.)         SD   SEM(Prop.)        SEM
-------------------------------------------------------------------------------
ΔBOND                     0.00          8.89       0.00         0.40       0.00
ΔANGLE                    0.00         11.85       0.00         0.53       0.00
ΔDIHED                    0.00          6.17       0.00         0.28       0.00
ΔUB                       0.00          2.89       0.00         0.13       0.00
ΔIMP                     -0.00          3.31       0.00         0.15       0.00
ΔCMAP                     0.00          3.68       0.00         0.16       0.00
ΔVDWAALS                -61.14          7.59       4.99         0.34       0.22
ΔEEL                   -460.68         43.29      40.77         1.93       1.82
Δ1-4 VDW                  0.00          4.79       0.00         0.21       0.00
Δ1-4 EEL                  0.00         20.99       0.00         0.94       0.00
ΔEPB                    486.02         24.22      35.35         1.08       1.58
ΔENPOLAR                 -8.87          0.26       0.32         0.01       0.01
ΔEDISPER                  0.00          0.00       0.00         0.00       0.00

ΔGGAS                  -521.82         44.32      39.57         1.98       1.77
ΔGSOLV                  477.15         24.22      35.22         1.08       1.57

ΔTOTAL                  -44.67         50.51      12.70         2.25       0.57
-------------------------------------------------------------------------------

3) WT calculation from alanine scan 2:

-------------------------------------------------------------------------------
Delta (Complex - Receptor - Ligand):
Energy Component       Average     SD(Prop.)         SD   SEM(Prop.)        SEM
-------------------------------------------------------------------------------
ΔBOND                     0.00          8.89       0.00         0.40       0.00
ΔANGLE                    0.00         11.85       0.00         0.53       0.00
ΔDIHED                    0.00          6.17       0.00         0.28       0.00
ΔUB                       0.00          2.89       0.00         0.13       0.00
ΔIMP                     -0.00          3.31       0.00         0.15       0.00
ΔCMAP                     0.00          3.68       0.00         0.16       0.00
ΔVDWAALS                -61.14          7.59       4.99         0.34       0.22
ΔEEL                   -153.56         14.43      13.59         0.64       0.61
Δ1-4 VDW                  0.00          4.79       0.00         0.21       0.00
Δ1-4 EEL                  0.00         20.99       0.00         0.94       0.00
ΔEPB                    157.97          7.89      11.34         0.35       0.51
ΔENPOLAR                 -8.87          0.26       0.32         0.01       0.01
ΔEDISPER                  0.00          0.00       0.00         0.00       0.00

ΔGGAS                  -214.70         17.28      13.00         0.77       0.58
ΔGSOLV                  149.10          7.89      11.20         0.35       0.50

ΔTOTAL                  -65.60         19.00       5.45         0.85       0.24
-------------------------------------------------------------------------------

4) WT calculation from alanine scan 3:

-------------------------------------------------------------------------------
Delta (Complex - Receptor - Ligand):
Energy Component       Average     SD(Prop.)         SD   SEM(Prop.)        SEM
-------------------------------------------------------------------------------
ΔBOND                     0.00          8.89       0.00         0.40       0.00
ΔANGLE                    0.00         11.85       0.00         0.53       0.00
ΔDIHED                    0.00          6.17       0.00         0.28       0.00
ΔUB                       0.00          2.89       0.00         0.13       0.00
ΔIMP                     -0.00          3.31       0.00         0.15       0.00
ΔCMAP                     0.00          3.68       0.00         0.16       0.00
ΔVDWAALS                -61.14          7.59       4.99         0.34       0.22
ΔEEL                    -92.14          8.66       8.15         0.39       0.36
Δ1-4 VDW                  0.00          4.79       0.00         0.21       0.00
Δ1-4 EEL                  0.00         20.99       0.00         0.94       0.00
ΔEPB                     92.38          4.62       6.55         0.21       0.29
ΔENPOLAR                 -8.87          0.26       0.32         0.01       0.01
ΔEDISPER                  0.00          0.00       0.00         0.00       0.00

ΔGGAS                  -153.28         12.86       8.19         0.57       0.37
ΔGSOLV                   83.51          4.63       6.42         0.21       0.29

ΔTOTAL                  -69.77         13.67       4.75         0.61       0.21
-------------------------------------------------------------------------------

Terminal output

N/A

gmx_MMPBSA.log

-

Operating system

Linux

gmx_MMPBSA Version

v1.6.1

Python version

3.9.15

Installation

None

Valdes-Tresanco-MS commented 1 year ago

This is because you set the variable cas_intdiel=1, which assigns internal dielectric constant values depending on the residue being mutated. We discussed related aspects in this issue. Let me know if you have any questions.

laurenreid1 commented 1 year ago

Hi @Valdes-Tresanco-MS,

Thanks for the speedy reply! Having read the documentation and the paper linked in the issue, I'd like to check that my understanding is correct: To compare the deltaH(bind) from multiple alanine scanning simulations, I need to use the default cas_intdiel=0 and PB indi = 1.0 for each system? This will make the different alanine scan mutant calculations comparable but some will have larger errors than others depending on the nature of the residue in question?

Thanks!

Valdes-Tresanco-MS commented 1 year ago

Yes and no. The results you get can be comparable to each other in any instance. However, it is sensible to use the same internal dielectric constant (intdiel). In the article, they observed a better correlation of the mutant DDG vs experimental (for the systems they tested) if the intdiel was varied according to the mutated residue. The remainder is sacrificed for trying to predict the importance of the contribution of that particular residue. In general, they found that the results fit better under these conditions than if they used the same intdiel. Along these same lines, it is important to remember that these methods are good for comparisons with relative values, not absolute ones. However, these may be specific cases, and may not be suitable for your system. In certain cases, as discussed in #280, the variation in the electrostatic component (EEL and EPB) is so large that it creates more problems than solutions. In any case, we implement tools that can be useful in some cases and not in others. It is important to understand how they work and in which cases to apply them. Let me know any doubts

laurenreid1 commented 1 year ago

Thanks a lot for your detailed reply, that has clarified the setting for me!