MobleyLab / alchemical-analysis

An open tool implementing some recommended practices for analyzing alchemical free energy calculations
MIT License
120 stars 60 forks source link

Inconsistent results between TI, EXP, BAR, MBAR #126

Open qinghualiao opened 3 years ago

qinghualiao commented 3 years ago

Hello,

I tried to run the ligand transformation in water using pmemd.cuda in the tutorial: http://ambermd.org/tutorials/advanced/tutorial9/ In the tutorial, there are three steps of transformation, which are called decharge, vdw_bonded and recharge, and here is the results from the tutorial (each step has 11 lambdas (bin size 0.1), and each window was sampled for 200 ps ): ligands/decharge: -6.44585036422 ligands/vdw_bonded: 1.29181321994 ligands/recharge: -33.8336526019

dG sum for ligands = -38.98768974618

To improve the sampling, I used 25 lambdas (bin size 0.05 with extra two windows at each end points), and each window was sampled for 2 ns. The results I got with the python script from the tutorial (TI) is: ligands/decharge: -6.363338637432731 ligands/vdw_bonded: 1.3493867576590386 ligands/recharge: -32.39912095583138

dG sum for ligands = -37.4130728356050724

The results I got are consistent with the ones from the tutorial. However, if I use the alchemical_analysis.py/mbar to do the analysis, then I saw the inconsistenc: ligands/decharge: TI: -6.40023 +- 0.03027 TI-CUBIC: -6.33526 +- 0.04188 DEXP: -15.25098 +- 0.38241 IEXP: 7.29922 +- 0.39631 BAR: -4.37929 +- 0.23774 MBAR: -6.61408 +- 0.66711

ligands/vdw_bonded: TI: 1.37004 +- 0.07038 TI-CUBIC: 1.34016 +- 0.09244 DEXP: -12.25424 +- 0.36173 IEXP: 9.20770 +- 0.34181 BAR: -1.49565 +- 0.20352 MBAR: 1.60530 +- 0.66557

ligands/recharge: TI: -32.26338 +- 0.06157 TI-CUBIC: -32.29185 +- 0.09051 DEXP: -30.09208 +- 0.39278 IEXP: -6.23544 +- 0.39539 BAR: -18.20007 +- 0.23451 MBAR: -20.96509 +- 0.69606

As you can see the results from the alchemical_analysis.py/mbar, the TI results are consistent with the tutorial. But clearly, there are big discrepancy between EXP with TI, BAR and MBAR. Could someone comment on the inconsistence among the methods to help me understand the results? Thanks a lot!

PS: Input and output files could be provided if someone want to check them!

All the best, Qinghua

davidlmobley commented 3 years ago

Please use alchemlyb instead; alchemical-analysis is no longer maintained and is quite out of date.

qinghualiao commented 3 years ago

Please use alchemlyb instead; alchemical-analysis is no longer maintained and is quite out of date.

@davidlmobley Thanks, I will then try with alchemlyb!

tanshy17 commented 2 years ago

hello qinghua, I meet the same problem,which the results between BAR and MBAR are not consistent. Have you try alchemlyb and will the alchemlyb fix this problems?

qinghualiao commented 2 years ago

Hello,

No, I did not fix the issue, then I just stick to TI estimation.

All the best, Qinghua

On 4/12/22 11:00, tanshy17 wrote:

hello qinghua, I meet the same problem,which the results between BAR and MBAR are not consistent. Have you try alchemlyb and will the alchemlyb fix this problems?

— Reply to this email directly, view it on GitHub https://github.com/MobleyLab/alchemical-analysis/issues/126#issuecomment-1096400175, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFROUSBLKMYNQOXXIDC7DUDVEU3SZANCNFSM5GLNOOAA. You are receiving this because you authored the thread.Message ID: @.***>

tanshy17 commented 2 years ago

Hello,

No, I did not fix the issue, then I just stick to TI estimation.

All the best, Qinghua

On 4/12/22 11:00, tanshy17 wrote:

hello qinghua, I meet the same problem,which the results between BAR and MBAR are not consistent. Have you try alchemlyb and will the alchemlyb fix this problems?

— Reply to this email directly, view it on GitHub <#126 (comment)>, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFROUSBLKMYNQOXXIDC7DUDVEU3SZANCNFSM5GLNOOAA. You are receiving this because you authored the thread.Message ID: @.***>

hello Qinghua, thanks for your quick reply! I have another question about the result.txt from amber software. Could you give me the some suggestions? Q1: Is that normal that there is not Coulomb and vdWaals items in the result.txt .The result.txt below are the result from gromacs samples. But when I run the analysis using amber output, there is not Coulomb and vdWaals items. Do you know how can I get the Coulomb and vdWaals informations from amber output? Or what should I set in the amber input script to get the Coulomb and vdWaals informations?

Command line was: ../../alchemical_analysis/alchemical_analysis.py -d data/3-methylindole-11steps -q xvg -p dhdl -u kJ


States TI (kJ/mol) TI-CUBIC (kJ/mol) DEXP (kJ/mol) IEXP (kJ/mol) BAR (kJ/mol) MBAR (kJ/mol)


0 -- 1 13.233 +- 0.113 13.172 +- 0.116 13.215 +- 0.223 12.845 +- 0.168 13.187 +- 0.117 13.311 +- 0.113 1 -- 2 13.813 +- 0.142 13.546 +- 0.165 13.212 +- 0.322 14.586 +- 0.574 13.579 +- 0.150 13.687 +- 0.144 2 -- 3 12.562 +- 0.225 12.277 +- 0.285 13.122 +- 0.304 26.267 +- 2.467 12.105 +- 0.221 12.144 +- 0.196 3 -- 4 -9.557 +- 0.281 -7.168 +- 0.324 -3.288 +- 0.780 -6.079 +- 0.376 -5.332 +- 0.260 -5.241 +- 0.230 4 -- 5 2.817 +- 0.179 3.353 +- 0.244 5.508 +- 0.250 3.563 +- 0.335 3.400 +- 0.191 3.490 +- 0.187 5 -- 6 -0.169 +- 0.238 0.597 +- 0.315 0.351 +- 0.649 0.965 +- 0.442 0.576 +- 0.221 0.700 +- 0.217 6 -- 7 -4.250 +- 0.259 -4.258 +- 0.297 -5.837 +- 1.277 -3.048 +- 0.969 -4.184 +- 0.298 -4.055 +- 0.291 7 -- 8 -7.411 +- 0.238 -8.107 +- 0.269 -7.912 +- 0.650 -5.389 +- 1.857 -8.622 +- 0.268 -8.744 +- 0.239 8 -- 9 -3.203 +- 0.055 -3.251 +- 0.059 -3.169 +- 0.073 -3.214 +- 0.146 -3.239 +- 0.051 -3.218 +- 0.042 9 -- 10 -2.199 +- 0.062 -2.144 +- 0.101 -2.082 +- 0.077 -1.777 +- 0.282 -2.079 +- 0.051 -2.054 +- 0.045


Coulomb: 39.609 +- 0.356 38.994 +- 0.378 39.549 +- 0.496 53.698 +- 2.538 38.870 +- 0.292 39.142 +- 0.332 vdWaals: -23.972 +- 0.725 -20.977 +- 0.711 -16.429 +- 1.776 -14.979 +- 2.222 -19.480 +- 0.564 -19.121 +- 0.647 TOTAL: 15.637 +- 0.808 18.017 +- 0.805 23.120 +- 1.844 38.719 +- 3.374 19.390 +- 0.635 20.021 +- 0.727

####################################################################### thanks again for your reply! Best wishes!