MDAnalysis / mdanalysis

MDAnalysis is a Python library to analyze molecular dynamics simulations.
https://mdanalysis.org
Other
1.29k stars 647 forks source link

trouble with HBL water dynamics analysis #1632

Closed grauph closed 7 years ago

grauph commented 7 years ago

Expected behaviour

To get the hydrogen bond lifetimes of water in an enzyme active site from a trajectory XTC file I used waterdynamics.

Actual behaviour

After run the analysis I get a warning (below), after that the process seem correct. Finally I print the results but all are zero. I tested whit different selections but its the same. How can I manage this issue? I have to specify a donor and acceptor indices?

HBL_analysis.run()
/anaconda/envs/mdaenv/lib/python2.7/site-packages/MDAnalysis/analysis/hbonds/hbond_analysis.py:594: DeprecationWarning: The donor and acceptor indices being 1-based is deprecated in favor of a zero-based index. These can be accessed by 'donor_index' or 'acceptor_index', removal of the 1-based indices is targeted for version 0.17.0
  " version 0.17.0", category=DeprecationWarning)
/anaconda/envs/mdaenv/lib/python2.7/site-packages/MDAnalysis/analysis/hbonds/hbond_analysis.py:718: SelectionWarning: No donors found in selection 1. You might have to specify a custom 'donors' keyword. Selection will update so continuing with fingers crossed.
  warnings.warn(errmsg, category=SelectionWarning)
/anaconda/envs/mdaenv/lib/python2.7/site-packages/MDAnalysis/analysis/hbonds/hbond_analysis.py:718: SelectionWarning: No acceptors found in selection 1. You might have to specify a custom 'acceptors' keyword. Selection will update so continuing with fingers crossed.
  warnings.warn(errmsg, category=SelectionWarning)
/anaconda/envs/mdaenv/lib/python2.7/site-packages/MDAnalysis/analysis/hbonds/hbond_analysis.py:718: SelectionWarning: No acceptors found in selection 2. You might have to specify a custom 'acceptors' keyword. Selection will update so continuing with fingers crossed.
  warnings.warn(errmsg, category=SelectionWarning)
/anaconda/envs/mdaenv/lib/python2.7/site-packages/MDAnalysis/analysis/hbonds/hbond_analysis.py:718: SelectionWarning: No donors found in selection 2. You might have to specify a custom 'donors' keyword. Selection will update so continuing with fingers crossed.
  warnings.warn(errmsg, category=SelectionWarning)
HBonds frame  5000:  5001/5001 [100.0%]

Code to reproduce the behaviour

import numpy as np
import MDAnalysis
from MDAnalysis.analysis.waterdynamics import HydrogenBondLifetimes as HBL
u = MDAnalysis.Universe('ref_2ypi.pdb', '2ypi_center.xtc')

water = "byres name SOL and sphzone 6.0 protein"
activesite = "resid 12 171 210 212 230 232 233 234"

HBL_analysis = HBL(u, water, activesite, 0, 5000, 30)
HBL_analysis.run()
time = 0

for HBLc, HBLi in HBL_analysis.timeseries:
    print("{time} {HBLc} {time} {HBLi}".format(time=time, HBLc=HBLc, HBLi=HBLi))
    time += 1
0 1.0 0 0.0
1 0.0 1 0.0
2 0.0 2 0.0
3 0.0 3 0.0
4 0.0 4 0.0
5 0.0 5 0.0
6 0.0 6 0.0
7 0.0 7 0.0
8 0.0 8 0.0
9 0.0 9 0.0
10 0.0 10 0.0
11 0.0 11 0.0
12 0.0 12 0.0
13 0.0 13 0.0
14 0.0 14 0.0
15 0.0 15 0.0
16 0.0 16 0.0
17 0.0 17 0.0
18 0.0 18 0.0
19 0.0 19 0.0
20 0.0 20 0.0
21 0.0 21 0.0
22 0.0 22 0.0
23 0.0 23 0.0
24 0.0 24 0.0
25 0.0 25 0.0
26 0.0 26 0.0
27 0.0 27 0.0
28 0.0 28 0.0
29 0.0 29 0.0

Currently version of MDAnalysis:

0.16.2

orbeckst commented 7 years ago

Maybe @alejob can help here?

(Possibly see also #1603).

alejob commented 7 years ago

I'm on it. Anyway, the second selection should be:

activesite = "resid 12 or resid 171 or resid 210 or resid 212 or resid 230
 or resid 232 or resid 233 or resid 234"

Does it changes the 4th and 5th warning?

orbeckst commented 7 years ago

The selection "resid 12 171 210 212 230 232 233 234" should be equivalent to ""resid 12 or resid 171 or resid 210 or resid 212 or resid 230 or resid 232 or resid 233 or resid 234" because we have an implicit "or" in the selections – admittedly, that isn't documented in the docs for Boolean selections, only in the tutorial... needs to be fixed (edit: ... in #1643).

alejob commented 7 years ago

Ah, ok, I searched in docs but didn't find it.

grauph commented 7 years ago

I fixed the selection as you say, but the warnings and the results are the same.

alejob commented 7 years ago

@grauph could you send us part of the pdb?, where we can see the solvent SOL and its atoms.

orbeckst commented 7 years ago

We would need more input to help, closing for now.

We will re-open if we have more information.