Positive Binding Energy in MMPBSA

136 views
Skip to first unread message

Kushal Singh

unread,
Feb 6, 2023, 1:38:30 PM2/6/23
to gmx_m...@googlegroups.com
Dear Sir

I am new to this field. For a particular MD Simulation of a protein-ligand complex using Gromacs and the system prepared with CHARMM ff, I am trying to do the MMPBSA calculation using gmx_MMPBSA (version 1.5.7). Everytime, I am getting positive value for the binding energy (both with or without entropy correction). I have tried out various reasonable combinations of frame range and frame intervals. But in all the cases, I am getting positive values for the binding energy which I guess indicates unfavourable interaction. 

However, all other parameters indicate towards a stable protein-ligand binding (highly negative docking binding energy, low RMSD of the protein, high number of protein-ligand hydrogen bonding, good results in the wetlab etc).

Is there a way to deal with such a problem? Should I send my files which may help you to analyze my issue?

I will be extremely grateful if you kindly help me with this problem.

Warm regards

Kushal Singh
Research Scholar
West Bengal State University
India

Mario Sergio Valdes

unread,
Feb 7, 2023, 12:47:37 AM2/7/23
to gmx_MMPBSA
Giving information without knowing the system is always a risk, so if you send us the files, we will have a better perspective (you can do it entirely confidentially). I suspect this is the well-known EDISPER term overestimation problem, but I won't know until I check it.
If you can, send us the following files:
- tpr
- top + itps + prm + charmm/toppar folder (if applicable)
- ndx
- 10 frames xtc
- groups
- COMPACT_MMXSA_RESULTS.mmxsa

Mario S.

Mario Sergio Valdes

unread,
Feb 7, 2023, 7:27:25 PM2/7/23
to gmx_MMPBSA
The file you shared with me is corrupted. After downloading 6.1GB (it is quite large), I cannot unzip it. Please share the file again and check that it works. 

Mario Sergio Valdes

unread,
Feb 8, 2023, 3:46:45 PM2/8/23
to gmx_MMPBSA
In this case, you should set the variable inp = 1. This is because you are using a CHARMM topology and the method for computing the non-polar component of the solvation is not optimized for the PBRadii you are using, and overestimates the EDISPER term. According to my test (with a few frames) the energy is negative.
Some remarks:
- Although it has no effect on the calculation, it is advisable to remove rotations and translations from the trajectory. Check this section https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/gmx_MMPBSA_running/#before-running-gmx_mmpbsa
- Your ligand has a rather unstable portion, which could be one of the reasons that you get a relatively low affinity in my test. Check this image
rmsd.png 
- For future doubts, please share a complete folder. This will make your job much easier
- Also, try to write non-confidential content here or on GitHub

Let me know any other questions
Mario S.

Kushal Singh

unread,
Feb 9, 2023, 1:28:59 PM2/9/23
to Mario Sergio Valdes, gmx_MMPBSA
Dear Sir

Thank you so much for your valuable guidance. I have a few doubts which I will be grateful if you may kindly clarify.

1. inp = 1    # Nonpolar solvation method
When to keep it as 1 and when to keep this value as 2 ?

2. radiopt = 1   # Use optimized radii?
When to keep it as 1  ? I have seen users keeping this value as 0. Is there any ideal/standard value here and what is the rationale? Does the value depend on which force field has been used to prepare the system?

3. PBRadii = 3  # Define PBRadii to build amber topology from GROMACS files
Sir, in your previous mail you mentioned this PBRadii. Is there any ideal/standard value here and again what is the rationale? Does the value depend on which force field has been used to prepare the system?

4. For reducing the ligand rmsd and increasing the binding affinity, what do you suggest ? Is repeating the MD with refined ligand topology the only option here?


Warm regards

Kushal Singh
Research Scholar
West Bengal State University
India



--
You received this message because you are subscribed to the Google Groups "gmx_MMPBSA" group.
To unsubscribe from this group and stop receiving emails from it, send an email to gmx_mmpbsa+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/gmx_mmpbsa/f0d32d7e-f013-4285-ae8a-70e201c18da6n%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Mario Sergio Valdes

unread,
Feb 9, 2023, 2:00:07 PM2/9/23
to gmx_MMPBSA
1. inp = 1    # Nonpolar solvation method
When to keep it as 1 and when to keep this value as 2 ?
In most cases, you can use inp = 2. However, in this case, the method used to estimate the nonpolar component of the solvation is optimized for the radii reported by Tan, Yang & Luo (protein and nucleic acids only), which is applied when radiopt = 1. When topology conversion is used (with CHARMM, for example), radiopt = 0 to avoid pre-calculation errors. We are currently working to implement the Tan, Yang & Luo radii so that it can be applied in almost all cases where it is not possible now.

2. radiopt = 1   # Use optimized radii?
When to keep it as 1  ? I have seen users keeping this value as 0. Is there any ideal/standard value here and what is the rationale? Does the value depend on which force field has been used to prepare the system?
As I mentioned, currently, we have to set radiopt=0 when doing topology conversion, because otherwise, sander fails trying to assign it. The only way to use it conveniently is when the system is entirely protean and the -cp option is not defined. We will be solving this limitation very soon.

3. PBRadii = 3  # Define PBRadii to build amber topology from GROMACS files
Sir, in your previous mail you mentioned this PBRadii. Is there any ideal/standard value here and again what is the rationale? Does the value depend on which force field has been used to prepare the system?
As I mentioned, currently, we have to set radiopt=0 when doing topology conversion, because otherwise, sander fails trying to assign it. The only way to use it conveniently is when the system is only protein and the -cp option is not defined. We will be solving this limitation very soon. You can check radiopt description in this section (https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/input_file/#options-to-define-the-physical-constants)

4. For reducing the ligand rmsd and increasing the binding affinity, what do you suggest ? Is repeating the MD with refined ligand topology the only option here?
The problem is not with the MD or the energy calculation. I think that in this case, the ligand structure does not fit the cavity well. In any case, it is just my observation. 

Sincerely
Mario S.

Kushal Singh

unread,
Feb 9, 2023, 2:14:34 PM2/9/23
to Mario Sergio Valdes, gmx_MMPBSA
Dear Sir

Thank you so much for your detailed explanation. However sir you clarification for question no. 3 (related to PBRadii) seems to be apparently same as that for question no. 2. 

It would be very helpful if you may just help me to understand your response related to the ideal value for PBRadii (question no. 3). 

Warm regards,

Kushal Singh
Research Scholar
West Bengal State University
India
--
You received this message because you are subscribed to the Google Groups "gmx_MMPBSA" group.
To unsubscribe from this group and stop receiving emails from it, send an email to gmx_mmpbsa+...@googlegroups.com.

Mario Sergio Valdes

unread,
Feb 9, 2023, 4:34:38 PM2/9/23
to gmx_MMPBSA
My bad. As far as I know, PB will use the radii only to compute the nonpolar component of the solvation, while GB does have much more weight since it is used to compute the effective Born radius of atoms. Commonly, mbondi is used when inp=1 or Tan&Luo when inp=2. The latter is limited in the current version of the program

Kushal Singh

unread,
Feb 10, 2023, 10:34:34 AM2/10/23
to Mario Sergio Valdes, gmx_MMPBSA
Dear Sir

After going through https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/input_file/, I am under the impression that mbondi can be equal to 5 as well as 6. Also  for MMPBSA calculations with systems prepared with CHARMM ff, we should set PBRadii = 7 in the mmpbsa.in file. I am getting a bit confused.

Can you please guide me that for systems prepared with CHARMM ff and inp=1, what should be the appropriate value assigned to PBRadii? My system doesn't have any halogen atoms it.


And sir, I have one more question. What to do if we get the following warning message related to the entropy correction?
"WARNING: THE INTERACTION ENERGY STANDARD DEVIATION [ σ(Int. Energy) ]
IS GREATER THAN 3.6 kcal/mol (~15 kJ/mol). THUS, THE INTERACTION ENTROPY VALUES ARE
NOT RELIABLE".
How to eliminate this warning?

Thank you


Warm regards,

Kushal Singh
Research Scholar
West Bengal State University
India


Mario Sergio Valdes

unread,
Feb 10, 2023, 11:18:15 AM2/10/23
to gmx_MMPBSA
After going through https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/input_file/, I am under the impression that mbondi can be equal to 5 as well as 6. Also  for MMPBSA calculations with systems prepared with CHARMM ff, we should set PBRadii = 7 in the mmpbsa.in file. I am getting a bit confused.
 The fact that there is a wide variety of radii provides users with options. Some are more recommended than others, depending on the nature of the system you are analyzing. mbondi_pb2 and mbondi_pb3 are the same as mbondi, only they include the radii for halogens. This means, that any of the 3 in non-halogen systems will give the same results. As for charmm_radii, it is simply another option. It is so named because it is the one reported for use with CHARMM. Its use in Amber is not sufficiently tested. However, some publications show that the radii inclusion is relatively small in PB since it modifies only the nonpolar component of the solvation, whereas in GB it may be more significant.

Can you please guide me that for systems prepared with CHARMM ff and inp=1, what should be the appropriate value assigned to PBRadii? My system doesn't have any halogen atoms it.
You can use mbondi or search other references in the literature
 
And sir, I have one more question. What to do if we get the following warning message related to the entropy correction?
"WARNING: THE INTERACTION ENERGY STANDARD DEVIATION [ σ(Int. Energy) ]
IS GREATER THAN 3.6 kcal/mol (~15 kJ/mol). THUS, THE INTERACTION ENTROPY VALUES ARE
NOT RELIABLE".
How to eliminate this warning?
Please read the article mentioned in this warning. In order to perform the Interaction Entropy calculation it is necessary to follow certain guidelines, which are listed and evaluated in this article. 

Sarthak Trivedi SVNIT

unread,
Feb 13, 2023, 1:37:00 AM2/13/23
to Mario Sergio Valdes, gmx_MMPBSA
Hello Sir,

I am facing the same issue regarding positive binding energy, I have completed PBSA calculations for two protein-ligand complexes. But unfortunately, one of them is giving the positive binding energy. I have gone through the discussion of this positive binding energy issue. As you mentioned that this positive value is due to overestimation of EDISPER term. My system also has the same configurations (forcefield). As you mentioned from above discussion about the pulling out the rotational and translational components from the trajectory. Being a beginner to this field I am not getting theoretical knowledge about this as well as not aware of getting rid off these two components from the trajectory. Please guide me in this regard. I can share the files with you, please let me know which files you require to analyze my problem more in detail. I am also going to do an extensive literature survey for theoretical studies about PB/GBSA models for more clarification from the discussion. 

Thank you in advance.

Regards,


--
You received this message because you are subscribed to the Google Groups "gmx_MMPBSA" group.
To unsubscribe from this group and stop receiving emails from it, send an email to gmx_mmpbsa+...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.


--
Sarthak Trivedi
Ph.D. (Research Fellow- D21PH003)
Materials and Biophysics Group,
Department of Physics,
Sardar Vallabhbhai National Institute of Technology, Surat, Gujarat, India,
395007.

Mario Sergio Valdes

unread,
Feb 13, 2023, 2:16:04 AM2/13/23
to gmx_MMPBSA
Please perform both calculations by setting the variable inp = 1. In the next version, this will be the default setting, as inp = 2 should only be used with the Tan, Yang & Luo radii which is not fully implemented.

Mario Sergio Valdes

unread,
Feb 13, 2023, 11:32:41 AM2/13/23
to gmx_MMPBSA
Thank you so much for your response. But i don't know that where i should put the inp=1, is it in mmpbsa.in input file ?
Yes, in the &pb namelist. Please, check this section (https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/input_file/#basic-input-options_2

Sarthak Trivedi SVNIT

unread,
Feb 14, 2023, 10:14:42 AM2/14/23
to Mario Sergio Valdes, gmx_MMPBSA
Hello sir,

Thank you so much for your clarification. Sincerely apologizing for asking silly queries, I needed to look into the manual information provided on the website. I need some valuable suggestions from your expertise, i am a PhD student and I am from a physics background, i found gmxMMPBSA, a vast software to perform different calculations, i have started reading manual fundamentally. My research area is computer aided drug design for neurodegenerative diseases, i have studied pathogenesis of Alzheimer's disease and Parkinson's disease so far, and am trying to screen some series of respective inhibitors. In gmxMMPBSA calculations, I need to study about my system and according to that system configuration, I would need to prepare an input file for gmxMMPBSA calculations. I want to have some guidance so I can make myself better at theoretical understanding, so can you recommend some general articles or books or anything else to focus my work for the particular applications? 

Any comments or opinions are welcomed. 

Thank you in advance 


Mr. Sarthak Trivedi
Research Fellow,
Department of Physics,
Sardar Vallabhbhai National Institute of Technology, Surat.

On Mon, Feb 13, 2023 at 10:02 PM Mario Sergio Valdes <mariosergi...@gmail.com> wrote:
Thank you so much for your response. But i don't know that where i should put the inp=1, is it in mmpbsa.in input file ?
Yes, in the &pb namelist. Please, check this section (https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/input_file/#basic-input-options_2

--
You received this message because you are subscribed to the Google Groups "gmx_MMPBSA" group.
To unsubscribe from this group and stop receiving emails from it, send an email to gmx_mmpbsa+...@googlegroups.com.

gmx_MMPBSA

unread,
Feb 14, 2023, 10:23:55 AM2/14/23
to gmx_MMPBSA
No problem Sarthak, we are happy to help everyone and make our community feel supported.
Regarding studies for the design of drugs against molecular targets in Alzheimer's disease, I don't have much background. A search in PubMed for the most recent review papers and state-of-the-art studies may give you a better perspective. As for the theory behind gmx_MMPBSA, you can check this section (https://valdes-tresanco-ms.github.io/gmx_MMPBSA/dev/introduction/). For parameter modification, you can look for references to other studies based on MMPB/GBSA.
If you have any questions or problems, please contact us.
Best of wishes
gmx_MMPBSA team

Sarthak Trivedi SVNIT

unread,
Feb 14, 2023, 11:43:38 AM2/14/23
to gmx_MMPBSA
Hello sir,

Thank you for immediate response. 

Regards,


Mr. Sarthak Trivedi
Research Fellow,
Department of Physics,
Sardar Vallabhbhai National Institute of Technology, Surat.
Reply all
Reply to author
Forward
0 new messages