Hbond analysis.

[Nawel Mele]

Hello,

Following the results of molecular dynamics, I want to study the evolution of hydrogen bonds during the trajectory.
I therefore uses the Mdanalysis hbond packages.
However, some points are still unclear. First of all here is the command I use to extract my Hbond:


import MDAnalysis.analysis.hbonds

u = MDAnalysis.Universe("file.gro","file.xtc", permissive=True)           
h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, selection1="bynum 5051:6060" , selection2="bynum 1:1010" ,distance=3.0, angle=120.0,start=10000,stop=20001,step=1,update_selection1=True, update_selection2=True)
h.run()
frames=h.timeserie


My questions are:

1- What is the purpose of this option ( red bold ):  h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, selection1=mon1 , selection2=mon2 ,distance=3.0, angle=120.0,start=10000,stop=20001,step=1,update_selection1=True, update_selection2=True)
In fact, when I do not use this option, the results are inconsistent.

2- Why when I did this command : h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, selection1="bynum 1:1010" , selection2="bynum 5051:6060" ,distance=3.0, angle=120.0,start=10000,stop=20001,step=1,update_selection1=True, update_selection2=True)
I get a results.

And when I reverse my selections:  h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, selection1="bynum 5051:6060" , selection2="bynum 1:1010" ,distance=3.0, angle=120.0,start=10000,stop=20001,step=1,update_selection1=True, update_selection2=True)
I receive an error message?

Thanks a lot for your answers.

Nawel

[Oliver Beckstein]

Hi Newal,

It's difficult to help you if you just say "the results are inconsistent" or "there was a error message". You need to tell us how the results are inconsistent or what the error message is or what the difference between outputs is. Otherwise there's not enough context. 

Remember that we don't have your input files in front of us and for many of us the code seems to work for what we're doing. Having enough information is crucial to helping you. 

Thanks,

Oliver

--

Oliver Beckstein

email: orbe...@gmail.com

skype: orbeckst

--
You received this message because you are subscribed to the Google Groups "MDnalysis discussion" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mdnalysis-discus...@googlegroups.com.
To post to this group, send email to mdnalysis-...@googlegroups.com.
Visit this group at http://groups.google.com/group/mdnalysis-discussion?hl=en.
For more options, visit https://groups.google.com/groups/opt_out.
 
 

[Christopher Mielack]

Hi,

I also made the observation that the h-bonding analysis yields inconsistent results with the default parameters. In my case I performed h-bonding analysis between a particular residue and bulk water. No matter what part of my trajectory I choose for this analysis, the number of hydrogen bonds will always start out high (~12) and then decay to zero - even though hydration is obviously pretty constant throughout the trajectory. 

Could it be the case that you perform an optimization in which you only select potential hydrogen bonding partners within some distance of the initial configuration, and that, without updating, these molecules in our cases diffuse away over time so that no bonds get detected? This is what it looks like to me.

Btw great job with the MDAnalysis package, it has been a great help to me throughout my current masters thesis project.

Christopher

To unsubscribe from this group and stop receiving emails from it, send an email to mdnalysis-discussion+unsub...@googlegroups.com.

[Oliver Beckstein]

Hi,

You probably need to set update_selection1=True and/or update_selection2=True.

Oliver

http://pythonhosted.org/MDAnalysis/documentation_pages/analysis/hbonds.html#MDAnalysis.analysis.hbonds.HydrogenBondAnalysis

>> To unsubscribe from this group and stop receiving emails from it, send an email to mdnalysis-discus...@googlegroups.com.

>> To post to this group, send email to mdnalysis-...@googlegroups.com.
>> Visit this group at http://groups.google.com/group/mdnalysis-discussion?hl=en.
>> For more options, visit https://groups.google.com/groups/opt_out.
>>
>>
>

> --
> You received this message because you are subscribed to the Google Groups "MDnalysis discussion" group.

> To unsubscribe from this group and stop receiving emails from it, send an email to mdnalysis-discus...@googlegroups.com.

> To post to this group, send email to mdnalysis-...@googlegroups.com.
> Visit this group at http://groups.google.com/group/mdnalysis-discussion?hl=en.
> For more options, visit https://groups.google.com/groups/opt_out.
>
>

--
Oliver Beckstein * orbe...@gmx.net
skype: orbeckst * orbe...@gmail.com

[Christopher Mielack]

With that it works, thanks! 

I think that update_selection1/2=False may be a troublesome default for many users because it makes hidden assumptions about the analyzed system. As a VMD user, I know I have to enable "update selection" whenever I use distance-dependent selections with solvent because intermolecular distances change all the time. However in this case I did not even suspect that the analysis code depends on the initial geometry of my system, I just told the module to perform hydrogen bonding analysis, and my selections did not depend on distances. In a sense, the code does not respect my selection string if I tell it to consider all water molecules but it only selects those that are closest to selection1 in the first frame and then does not update that selection each frame.

I understand that this is a very efficient optimisation, but I would value safe default behaviors over efficient ones. Why not change the default to True? Or reframe the flag by naming it something that communicates this behavior more clearly, e.g. "assume_rigid_molecule" (you can probably come up with something more fitting than me).

 

Christopher

[Oliver Beckstein]

Hi Christopher,

I think you have a valid point there. Can you file an issue report, please? All you have to is give it a title and pretty much copy & paste what you said.

The reports are our way to keep track of all suggestions because everyone is busy and might not react to a request on the mailing list right away.

If anyone else has suggestions e.g. as to how options should be named etc please chime in.

Many thanks,
Oliver