Fix HydrogenBondAnalysis when no atom pair is with 3A #1325

xiki-tempula · 2017-04-28T15:20:53Z

When there is no atom pair between the two selections of HydrogenBondAnalysis that are within 3A, HydrogenBondAnalysis throws an exception instead of giving a empty list.
Code to reproduce the error

import MDAnalysis as mda
import MDAnalysis.analysis.hbonds

text = '''TITLE     Two atoms far away
CRYST1   52.763   52.763   52.763  90.00  90.00  90.00 P 1           1
MODEL         1
ATOM      1  N   LEU     1      32.310  13.778  14.372  1.00  0.00      SYST N 0
ATOM      2  OW  SOL     2       3.024   4.456   4.147  1.00  0.00      SYST H 0'''
with open('temp.pdb', 'w') as f:
    f.write(text)
u = mda.Universe('temp.pdb')
h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, 'resname SOL', 'protein')
h.run()

Gives AttributeError: 'list' object has no attribute 'select_atoms'

After analysing the code, it seems that there is no check after filtering to make sure that a selection of atom can be done.

This PR fixed this problem.

PR Checklist

Tests?
Docs?
CHANGELOG updated?
Issue raised/referenced?

orbeckst · 2017-04-28T18:30:40Z

@xiki-tempula thanks for report and fix.

You probably need to rebase your branch against develop and force push.

I also don't understand why travis did not run; there should be a proper unit test run. Hopefully that should show up after you push again.

As a minor thing: If you find a bug, raise an issue and describe the bug (as you did in the first part of your report). Then open a PR that references the issue. The reason is that sometimes PRs get thrown away or changed. Then you are also loosing a good reference point. It's generally considered good practice to separate the issue report from the PR (or PRs). No big deal, but for the future.

orbeckst

Looking good, especially that it comes with a test case. Revisions:

I think we need the same check for _s1: can you comment and if necessary implement, too?
suggestions for code clean-up

Btw, instead of bundling a mini-file you can probably use the new "stream reading" capability

import MDAnalysis as mda
import StringIO

pdb = '''TITLE     Two atoms far away
CRYST1   52.763   52.763   52.763  90.00  90.00  90.00 P 1           1
MODEL         1
ATOM      1  N   LEU     1      32.310  13.778  14.372  1.00  0.00      SYST N 0
ATOM      2  OW  SOL     2       3.024   4.456   4.147  1.00  0.00      SYST H 0'''

u = mda.Universe(StringIO.StringIO(pdb), format="pdb")

and then you don't need to store a tiny file.

orbeckst · 2017-04-28T18:34:19Z

package/CHANGELOG


  * 0.16.1

 Enhancements
  * Universe now works with StringIO objects for topologies and trajectories.
+  * Residues with the same residue ids are not merged by default now 


There are lots of pieces that should not show up in your PR, they come from elsewhere. Something went wrong when you merged. Rebase against develop.

orbeckst · 2017-04-28T18:41:07Z

package/MDAnalysis/analysis/hbonds/hbond_analysis.py

-                    self._s2_donors_h[i] = tmp
-            self.logger_debug("Selection 2 donors: {0:d}".format(len(self._s2_donors)))
-            self.logger_debug("Selection 2 donor hydrogens: {0:d}".format(len(self._s2_donors_h)))
+        if self._s2:


Don't you need a similar check for _s1?

You can probably re-use your test case by just switching selections.

In my opinion, it would be cleaner to just return here instead of doing a big if block. If you re-order the code so that all data structures are initialized and then use the if not self._s2 where the logger is already reporting this problem and return from there.

self._s2_donors = {} self._s2_donors_h = {} self._s2_acceptors = {} if not self._s2: # nothing to do here, just keep empty data structures logger.warn('Selection 2 "{0}" did not select any atoms.'.format( str(self.selection2)[:80])) return if self.selection1_type in ('donor', 'both'): self._s2_acceptors = self._s2.select_atoms(...) ...

This is reasonably clean (and the diff will also be smaller).

xiki-tempula · 2017-04-29T10:05:58Z

@orbeckst Thank you for your advice.
I have checked and it seems that this bug occurred because there was no check after filtering.
There are checks at the beginning which dealt with an empty selection. Then the script compares the pairwise distance between s1 and s2 and then get rid of the atoms in s2 which are too far away from s1. Thus, s2 was changed while s1 was not. So if the s2 becomes empty after filtering, it escaped the check and caused error.
This can be check using
h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, 'protein and not protein', 'resname SOL')

orbeckst

Thanks for the clarification regarding testing s1 and s2.

Minor requests:

clean up (see inline)
rebase against develop (do you know how to do this?)

orbeckst · 2017-04-29T16:30:05Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py

@@ -31,8 +31,12 @@

 import itertools
 import warnings
+try:


Why is this needed?? I think it can be removed, can't it?

Thank you for your comment. There is no StringIO in py3 so I was trying to make it compatible with py3. I still need some time to figure out how this stream reading works.

Oh, I didn't think about this. Then we should probably use six.StringIO: something like

from six import StringIO

which replaces

from StringIO import StringIO

Thank you. It seems that six is more powerful than I thought after spending an hour playing with py3.io.StringIO.

orbeckst · 2017-04-29T16:30:31Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py

@@ -87,6 +91,18 @@ def test_generate_table(self):
        assert_array_equal(h.table.donor_resid, self.values['donor_resid'])
        assert_array_equal(h.table.acceptor_resnm, self.values['acceptor_resnm'])

+    def test_atoms_too_far(self):


Minor thing: Make this a static method

@staticmethod def test_atoms_too_far(): ...

orbeckst · 2017-04-29T16:30:46Z

testsuite/MDAnalysisTests/data/HBond_atoms_far.pdb

@@ -0,0 +1,5 @@
+TITLE     Two atoms far away


remove (you're using the inline test)

orbeckst · 2017-04-29T16:30:52Z

testsuite/MDAnalysisTests/datafiles.py

@@ -144,6 +144,7 @@
    "MMTF", "MMTF_gz",
    "ALIGN_BOUND",  # two component bound system
    "ALIGN_UNBOUND", # two component unbound system
+    "HBond_atoms_far", # Universe with two atoms far away


remove (you're using the inline test)

orbeckst · 2017-04-29T16:30:56Z

testsuite/MDAnalysisTests/datafiles.py

@@ -221,6 +222,7 @@
 XTC_sub_sol = resource_filename(__name__, 'data/cobrotoxin.xtc')
 PDB_sub_sol = resource_filename(__name__, 'data/cobrotoxin.pdb')
 PDB_xlserial = resource_filename(__name__, 'data/xl_serial.pdb')
+HBond_atoms_far = resource_filename(__name__, 'data/HBond_atoms_far.pdb')


remove (you're using the inline test)

orbeckst · 2017-04-29T16:31:37Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py

@@ -88,6 +92,18 @@ def test_generate_table(self):
        assert_array_equal(h.table.acceptor_resnm, self.values['acceptor_resnm'])

    @staticmethod
+    def test_atoms_too_far(self):


remove self from the argument list of the static method

orbeckst · 2017-04-29T16:36:51Z

I resolved the conflict in CHANGELOG using GitHub's friendly Resolve button (I hadn't done this before). Apparently that merged develop into your branch. Normally that is not what I would do but this should be ok. Remember to git pull before committing to this branch.

When it comes to merging I will probably do any squashing manually from the command line, just to make sure that the merge of develop into your branch did not mess up anything.

orbeckst · 2017-04-29T16:54:52Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py

@@ -103,7 +104,6 @@ def test_atoms_too_far(self):
        h.run(verbose=False)
        assert_equal(h.timeseries, [[]])

-    @staticmethod
    def test_true_traj():


This @staticmethod has to stay.

orbeckst

use six.StringIO
fixes as in comments

orbeckst · 2017-04-29T17:03:28Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py


-from MDAnalysisTests.datafiles import PDB_helix, GRO, XTC
+from MDAnalysisTests.datafiles import PDB_helix, GRO, XTC, HBond_atoms_far


remove HBond_atoms_far

orbeckst · 2017-04-29T17:04:05Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py

+        h = MDAnalysis.analysis.hbonds.HydrogenBondAnalysis(u, 'resname SOL', 'protein')
+        h.run(verbose=False)
+        assert_equal(h.timeseries, [[]])
+


There has to be a @staticmethod in front of the test_true_traj() function; Travis fails https://travis-ci.org/MDAnalysis/mdanalysis/jobs/227170225 , too, and QC flags it as a critical issue.

orbeckst · 2017-04-29T17:06:07Z

testsuite/MDAnalysisTests/analysis/test_hbonds.py

@@ -31,8 +31,12 @@

 import itertools
 import warnings
+try:


Use six.StringIO (as discussed in previous comment)

xiki-tempula · 2017-04-29T20:58:14Z

@orbeckst Thank you for your advice. I guess I have learnt about decorators today.

orbeckst · 2017-04-30T06:57:25Z

I am rebasing/squashing/merging this manually.

* add check that s2 selection is empty * add test * add change log entry

orbeckst · 2017-04-30T07:15:59Z

See PR #1327 (I needed to make a new PR for the rebasing/history rewriting).

xiki-tempula · 2017-04-30T07:22:36Z

@orbeckst I have approved it. Thank you for all the help and advice.

orbeckst · 2017-04-30T07:29:33Z

Thanks for fixing a bug!

We'll have to wait for another coredev to approve #1327 – I can't do it because I submitted the PR.

* add check that s2 selection is empty * add test * add change log entry

xiki-tempula added 4 commits April 28, 2017 16:10

Fix hbond cannot deal with atoms too far

53e1d94

add change log

b1f71e1

resolve conflict

dff9ad1

resolve conflict

b6a8a07

orbeckst added Component-Analysis defect labels Apr 28, 2017

orbeckst requested changes Apr 28, 2017

View reviewed changes

xiki-tempula added 4 commits April 29, 2017 10:20

rebase

c1d3eef

Merge remote-tracking branch 'origin/h_bond' into h_bond

0f89262

Fix test and refine hbond analysis

65536ea

Fix minor error and add py3 compatibility

98aba1f

change to static method

cc947a4

orbeckst requested changes Apr 29, 2017

View reviewed changes

Merge branch 'develop' into h_bond

a44bbc5

orbeckst self-assigned this Apr 29, 2017

delete the file and lines

a9b238d

orbeckst reviewed Apr 29, 2017

View reviewed changes

orbeckst requested changes Apr 29, 2017

View reviewed changes

fix StringIO and static method

120f380

orbeckst approved these changes Apr 30, 2017

View reviewed changes

orbeckst pushed a commit that referenced this pull request Apr 30, 2017

Fix hbond cannot deal with atoms too far (PR #1325)

f17247e

* add check that s2 selection is empty * add test * add change log entry

orbeckst mentioned this pull request Apr 30, 2017

Fix HydrogenBondAnalysis cannot deal with atoms too far (PR #1325) #1327

Merged

2 tasks

orbeckst closed this Apr 30, 2017

richardjgowers pushed a commit that referenced this pull request May 1, 2017

Fix hbond cannot deal with atoms too far (PR #1325)

f6414bf

* add check that s2 selection is empty * add test * add change log entry

richardjgowers pushed a commit that referenced this pull request May 1, 2017

Fix hbond cannot deal with atoms too far (PR #1325)

d5bc8f5

* add check that s2 selection is empty * add test * add change log entry

xiki-tempula deleted the h_bond branch May 16, 2017 20:19

orbeckst mentioned this pull request May 21, 2017

Hydrogen bond analysis exits with AttributeError: 'list' object has no attribute 'select_atoms' #1268

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Fix HydrogenBondAnalysis when no atom pair is with 3A #1325

Fix HydrogenBondAnalysis when no atom pair is with 3A #1325

xiki-tempula commented Apr 28, 2017 •

edited

Loading

orbeckst commented Apr 28, 2017

orbeckst left a comment

orbeckst Apr 28, 2017

orbeckst Apr 28, 2017

orbeckst Apr 28, 2017

xiki-tempula commented Apr 29, 2017 •

edited

Loading

orbeckst left a comment

orbeckst Apr 29, 2017

xiki-tempula Apr 29, 2017

orbeckst Apr 29, 2017 •

edited

Loading

xiki-tempula Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst commented Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst left a comment

orbeckst Apr 29, 2017

orbeckst Apr 29, 2017

orbeckst Apr 29, 2017

xiki-tempula commented Apr 29, 2017

orbeckst commented Apr 30, 2017

orbeckst commented Apr 30, 2017

xiki-tempula commented Apr 30, 2017

orbeckst commented Apr 30, 2017


		from MDAnalysisTests.datafiles import PDB_helix, GRO, XTC
		from MDAnalysisTests.datafiles import PDB_helix, GRO, XTC, HBond_atoms_far

Fix HydrogenBondAnalysis when no atom pair is with 3A #1325

Fix HydrogenBondAnalysis when no atom pair is with 3A #1325

Conversation

xiki-tempula commented Apr 28, 2017 • edited Loading

PR Checklist

orbeckst commented Apr 28, 2017

orbeckst left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

xiki-tempula commented Apr 29, 2017 • edited Loading

orbeckst left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

orbeckst Apr 29, 2017 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

orbeckst commented Apr 29, 2017

Choose a reason for hiding this comment

orbeckst left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

xiki-tempula commented Apr 29, 2017

orbeckst commented Apr 30, 2017

orbeckst commented Apr 30, 2017

xiki-tempula commented Apr 30, 2017

orbeckst commented Apr 30, 2017

xiki-tempula commented Apr 28, 2017 •

edited

Loading

xiki-tempula commented Apr 29, 2017 •

edited

Loading

orbeckst Apr 29, 2017 •

edited

Loading