From a94d350046eb965686eb347acd15036b20f0f699 Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Thu, 25 Jun 2020 16:10:43 +0100
Subject: [PATCH 01/16] Adapting the autocorrelation and intermittency
 functions in the new hydrogen bond implementation. Fixes #2547 and we
 depracate the waterdynamic.HydrogenBondLifetime because of #2247.

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 .../analysis/hydrogenbonds/hbond_analysis.py  | 69 ++++++++++++++++++-
 package/MDAnalysis/analysis/waterdynamics.py  |  6 ++
 package/MDAnalysis/lib/correlations.py        | 28 ++++----
 .../analysis/test_hydrogenbonds_analysis.py   | 16 ++++-
 4 files changed, 103 insertions(+), 16 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index d662bfd1eca..c8bac367d4e 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -166,11 +166,13 @@
 .. autoclass:: HydrogenBondAnalysis
    :members:
 """
-import warnings
+import logging
+
 import numpy as np
 
 from .. import base
 from MDAnalysis.lib.distances import capped_distance, calc_angles
+from MDAnalysis.lib.correlations import autocorrelation, correct_intermittency
 from MDAnalysis.exceptions import NoDataError
 
 from ...due import due, Doi
@@ -496,6 +498,71 @@ def _conclude(self):
 
         self.hbonds = np.asarray(self.hbonds).T
 
+    def autocorrelation(self, tau_max=20, window_step=1, intermittency=0):
+        """
+        Computes and returns the time-autocorrelation (HydrogenBondLifetimes) of hydrogen bonds.
+
+        Before calling this method, the hydrogen bonds must first be computed with the `run()` function.
+        The same `start`, `stop` and `step` parameters used in finding hydrogen bonds will be used here
+        for calculating hydrogen bond lifetimes. That is, the same frames will be used in the analysis.
+
+        Unique hydrogen bonds are identified using hydrogen-acceptor pairs. This means an acceptor
+        switching to a different hydrogen atom - with the same donor - from one frame to the next
+        is considered a different hydrogen bond.
+
+        Please see :func:`MDAnalysis.lib.correlations.autocorrelation` and
+        :func:`MDAnalysis.lib.correlations.intermittency` functions for more details.
+
+
+        Parameters
+        ----------
+        window_step : int, optional
+            The number of frames between each t(0).
+        tau_max : int, optional
+            Hydrogen bond lifetime is calculated for frames in the range 1 <= `tau` <= `tau_max`
+        intermittency : int, optional
+            The maximum number of consecutive frames for which a bond can disappear but be counted as present if it
+            returns at the next frame. An intermittency of `0` is equivalent to a continuous autocorrelation, which does
+            not allow for hydrogen bond disappearance. For example, for `intermittency=2`, any given hydrogen bond may
+            disappear for up to two consecutive frames yet be treated as being present at all frames.
+            The default is continuous (intermittency=0).
+
+        Returns
+        -------
+        acf_tau_timeseries : list
+            tau from 1 to `tau_max`. Saved in the field acf_tau_timeseries.
+        acf_timeseries : list
+            autcorrelation value for each value of `tau`. Saved in the field acf_timeseries.
+        acf_timeseries_data: list
+            raw datapoints, of hydrogen and acceptor indices, from which the autocorrelation is calculated.
+             Saved in the field acf_timeseries_data.
+        """
+
+        if not hasattr(self, 'hbonds'):
+            logging.error("Autocorrelation analysis of hydrogen bonds cannot be done before the hydrogen bonds are found")
+            logging.error("Autocorrelation: Please use the .run() before calling this function")
+            return
+
+        if self.step != 1:
+            logging.warning("Autocorrelation: ideally autocorrelation function would be carried out on consecutive frames. ")
+            logging.warning("Autocorrelation: if you would like to allow bonds to break and reform, please use 'intermittency'")
+
+        # Extract the hydrogen bonds IDs only in the format [set(superset(x1,x2), superset(x3,x4)), ..]
+        found_hydrogen_bonds = [set() for _ in self.frames]
+        for frame_index, frame in enumerate(self.frames):
+            for hbond in self.hbonds[self.hbonds[:, 0] == frame]:
+                found_hydrogen_bonds[frame_index].add(frozenset(hbond[2:4]))
+
+        intermittent_hbonds = correct_intermittency(found_hydrogen_bonds, intermittency=intermittency)
+        tau_timeseries, timeseries, timeseries_data = autocorrelation(intermittent_hbonds, tau_max,
+                                                                      window_step=window_step)
+        self.acf_tau_timeseries = tau_timeseries
+        self.acf_timeseries = timeseries
+        # user can investigate the distribution and sample size
+        self.acf_timeseries_data = timeseries_data
+
+        return self
+
     def count_by_time(self):
         """Counts the number of hydrogen bonds per timestep.
 
diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index c5d32264b9b..4776a0471a9 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -513,6 +513,12 @@ class HydrogenBondLifetimes(object):
 
     def __init__(self, universe, selection1, selection2, t0, tf, dtmax,
                  nproc=1):
+        warnings.warn(
+            "This class is depracated. "
+            "Please use MDAnalysis.analysis.hydrogenbonds instead.",
+            category=DeprecationWarning
+        )
+
         self.universe = universe
         self.selection1 = selection1
         self.selection2 = selection2
diff --git a/package/MDAnalysis/lib/correlations.py b/package/MDAnalysis/lib/correlations.py
index bf52dd61707..b0c91aca8c5 100644
--- a/package/MDAnalysis/lib/correlations.py
+++ b/package/MDAnalysis/lib/correlations.py
@@ -223,36 +223,36 @@ def correct_intermittency(list_of_sets, intermittency):
 
     list_of_sets = deepcopy(list_of_sets)
 
-    for i, ids in enumerate(list_of_sets):
+    for i, elements in enumerate(list_of_sets):
         # initially update each frame as seen 0 ago (now)
-        seen_frames_ago = {i: 0 for i in ids}
+        seen_frames_ago = {i: 0 for i in elements}
         for j in range(1, intermittency + 2):
-            for atomid in seen_frames_ago.keys():
+            for element in seen_frames_ago.keys():
                 # no more frames
                 if i + j >= len(list_of_sets):
                     continue
 
-                # if the atom is absent now
-                if not atomid in list_of_sets[i + j]:
+                # if the element is absent now
+                if not element in list_of_sets[i + j]:
                     # increase its absence counter
-                    seen_frames_ago[atomid] += 1
+                    seen_frames_ago[element] += 1
                     continue
 
-                # the atom is found
-                if seen_frames_ago[atomid] == 0:
-                    # the atom was present in the last frame
+                # the element is found
+                if seen_frames_ago[element] == 0:
+                    # the element was present in the last frame
                     continue
 
                 # it was absent more times than allowed
-                if seen_frames_ago[atomid] > intermittency:
+                if seen_frames_ago[element] > intermittency:
                     continue
 
-                # the atom was absent but returned (within <= intermittency_value)
+                # the element was absent but returned (within <= intermittency_value)
                 # add it to the frames where it was absent.
                 # Introduce the corrections.
-                for k in range(seen_frames_ago[atomid], 0, -1):
-                    list_of_sets[i + j - k].add(atomid)
+                for k in range(seen_frames_ago[element], 0, -1):
+                    list_of_sets[i + j - k].add(element)
 
-                seen_frames_ago[atomid] = 0
+                seen_frames_ago[element] = 0
     return list_of_sets
 
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index 80356bcd7a6..1c9077bfd29 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -207,6 +207,20 @@ def test_count_by_time(self, universe):
         assert_array_almost_equal(h.times, ref_times)
         assert_array_equal(counts, ref_counts)
 
+    def test_hydrogen_bond_lifetime(self, universe):
+        h = HydrogenBondAnalysis(universe, **self.kwargs)
+        h.run()
+
+        h.autocorrelation(tau_max=2)
+        assert_array_equal(h.acf_timeseries, [1, 0, 0])
+
+    def test_hydrogen_bond_lifetime_intermittency(self, universe):
+        h = HydrogenBondAnalysis(universe, **self.kwargs)
+        h.run()
+
+        h.autocorrelation(tau_max=2, intermittency=1)
+        assert_array_equal(h.acf_timeseries, 1)
+
 
 class TestHydrogenBondAnalysisTIP3P_GuessAcceptors_GuessHydrogens_UseTopology_(TestHydrogenBondAnalysisTIP3P):
     """Uses the same distance and cutoff hydrogen bond criteria as :class:`TestHydrogenBondAnalysisTIP3P`, so the
@@ -358,4 +372,4 @@ def test_count_by_type(self, h):
         ref_count = 15
 
         counts = h.count_by_type()
-        assert int(counts[0, 2]) == ref_count
\ No newline at end of file
+        assert int(counts[0, 2]) == ref_count

From 7b7badc9a05eeb02cf93fca720d529a62e27e437 Mon Sep 17 00:00:00 2001
From: Mateusz Bieniek <bieniekmateusz@users.noreply.github.com>
Date: Fri, 26 Jun 2020 13:34:03 +0100
Subject: [PATCH 02/16] Update package/MDAnalysis/analysis/waterdynamics.py

Co-authored-by: Oliver Beckstein <orbeckst@gmail.com>
---
 package/MDAnalysis/analysis/waterdynamics.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index 4776a0471a9..a2d0dd7f0d7 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -514,7 +514,7 @@ class HydrogenBondLifetimes(object):
     def __init__(self, universe, selection1, selection2, t0, tf, dtmax,
                  nproc=1):
         warnings.warn(
-            "This class is depracated. "
+            "This class is deprecated. "
             "Please use MDAnalysis.analysis.hydrogenbonds instead.",
             category=DeprecationWarning
         )

From 0b9a3997ca3d95acbbc111d9a1608f4ba738775d Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Sat, 27 Jun 2020 19:44:10 +0100
Subject: [PATCH 03/16] Coverage tests and refactoring for PR #2791

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 package/CHANGELOG                                          | 1 +
 .../MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py    | 2 +-
 package/MDAnalysis/analysis/waterdynamics.py               | 2 +-
 package/MDAnalysis/lib/correlations.py                     | 7 +++++--
 4 files changed, 8 insertions(+), 4 deletions(-)

diff --git a/package/CHANGELOG b/package/CHANGELOG
index a0aea381aac..1bb374bfe76 100644
--- a/package/CHANGELOG
+++ b/package/CHANGELOG
@@ -32,6 +32,7 @@ Enhancements
     token (Issue #2468, PR #2707)
   * Added the `from_smiles` classmethod to the Universe (Issue #2468, PR #2707)
   * Added computation of Mean Squared Displacements (#2438, PR #2619)
+  * Added Hydrogen Bond Lifetime via existing autocorrelation features (PR #2791)
 
 Changes
   * Changes development status from Beta to Mature (Issue #2773)
diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index c8bac367d4e..a103766401a 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -541,7 +541,7 @@ def autocorrelation(self, tau_max=20, window_step=1, intermittency=0):
         if not hasattr(self, 'hbonds'):
             logging.error("Autocorrelation analysis of hydrogen bonds cannot be done before the hydrogen bonds are found")
             logging.error("Autocorrelation: Please use the .run() before calling this function")
-            return
+            raise NoDataError("No .hbonds: use the .run() first")
 
         if self.step != 1:
             logging.warning("Autocorrelation: ideally autocorrelation function would be carried out on consecutive frames. ")
diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index a2d0dd7f0d7..6bec1d5a1f7 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -515,7 +515,7 @@ def __init__(self, universe, selection1, selection2, t0, tf, dtmax,
                  nproc=1):
         warnings.warn(
             "This class is deprecated. "
-            "Please use MDAnalysis.analysis.hydrogenbonds instead.",
+            "Please use MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.autocorellation() instead.",
             category=DeprecationWarning
         )
 
diff --git a/package/MDAnalysis/lib/correlations.py b/package/MDAnalysis/lib/correlations.py
index b0c91aca8c5..5a97f857bd4 100644
--- a/package/MDAnalysis/lib/correlations.py
+++ b/package/MDAnalysis/lib/correlations.py
@@ -223,6 +223,9 @@ def correct_intermittency(list_of_sets, intermittency):
 
     list_of_sets = deepcopy(list_of_sets)
 
+    # an element (a superset) takes the form of:
+    # - an atom pair when computing hydrogen bonds lifetime,
+    # - atom ID in the case of water survival probability,
     for i, elements in enumerate(list_of_sets):
         # initially update each frame as seen 0 ago (now)
         seen_frames_ago = {i: 0 for i in elements}
@@ -238,12 +241,12 @@ def correct_intermittency(list_of_sets, intermittency):
                     seen_frames_ago[element] += 1
                     continue
 
-                # the element is found
+                # the element is present
                 if seen_frames_ago[element] == 0:
                     # the element was present in the last frame
                     continue
 
-                # it was absent more times than allowed
+                # element was absent more times than allowed
                 if seen_frames_ago[element] > intermittency:
                     continue
 

From 42b9ad310d41d3b10d4661379fe6ba03f48e27e6 Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Sat, 27 Jun 2020 19:44:47 +0100
Subject: [PATCH 04/16] Coverage tests (PR #2791)

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 .../analysis/test_hydrogenbonds_analysis.py   | 25 ++++++++++++++++---
 1 file changed, 21 insertions(+), 4 deletions(-)

diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index 1c9077bfd29..7ef09635367 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -21,14 +21,16 @@
 # J. Comput. Chem. 32 (2011), 2319--2327, doi:10.1002/jcc.21787
 #
 import numpy as np
-import MDAnalysis
-from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis
-from MDAnalysis.exceptions import NoDataError
-
+import logging
 import pytest
+
 from numpy.testing import (assert_allclose, assert_equal,
                            assert_array_almost_equal, assert_array_equal,
                            assert_almost_equal)
+
+import MDAnalysis
+from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis
+from MDAnalysis.exceptions import NoDataError
 from MDAnalysisTests.datafiles import waterPSF, waterDCD
 
 
@@ -221,6 +223,21 @@ def test_hydrogen_bond_lifetime_intermittency(self, universe):
         h.autocorrelation(tau_max=2, intermittency=1)
         assert_array_equal(h.acf_timeseries, 1)
 
+    def test_no_attr_hbonds(self, universe):
+        h = HydrogenBondAnalysis(universe, **self.kwargs)
+
+        with pytest.raises(NoDataError, match="No .hbonds"):
+            h.autocorrelation(tau_max=2, intermittency=1)
+
+    def test_logging_step_not_1(self, universe, caplog):
+        h = HydrogenBondAnalysis(universe, **self.kwargs)
+        h.run(step=2)
+
+        caplog.set_level(logging.WARNING)
+        h.autocorrelation(tau_max=2, intermittency=1)
+        assert any("ideally autocorrelation function would be carried out on consecutive frames" in
+                   rec.getMessage() for rec in caplog.records)
+
 
 class TestHydrogenBondAnalysisTIP3P_GuessAcceptors_GuessHydrogens_UseTopology_(TestHydrogenBondAnalysisTIP3P):
     """Uses the same distance and cutoff hydrogen bond criteria as :class:`TestHydrogenBondAnalysisTIP3P`, so the

From 11ed00c66c84f656da5ac31c6f4b81cebdc47dd5 Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Sun, 28 Jun 2020 15:53:54 +0100
Subject: [PATCH 05/16] Parameter "between" for hydrogen binding.

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 .../analysis/hydrogenbonds/hbond_analysis.py  | 97 ++++++++++++++++---
 .../analysis/test_hydrogenbonds_analysis.py   | 26 ++---
 2 files changed, 96 insertions(+), 27 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index a103766401a..7c85dac4641 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -191,7 +191,8 @@ class HydrogenBondAnalysis(base.AnalysisBase):
     """
 
     def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=None,
-                 d_h_cutoff=1.2, d_a_cutoff=3.0, d_h_a_angle_cutoff=150, update_selections=True):
+                 between=None, d_h_cutoff=1.2, d_a_cutoff=3.0, d_h_a_angle_cutoff=150,
+                 update_selections=True):
         """Set up atom selections and geometric criteria for finding hydrogen bonds in a Universe.
 
         Parameters
@@ -208,6 +209,14 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
         acceptors_sel : str
             Selection string for the hydrogen bond acceptor atoms. Leave as `None` to guess which atoms to use in the
             analysis using :attr:`guess_acceptors`
+        between : List (optional)
+            Specify two selection strings for non-updating atom groups between which hydrogen bonds will be calculated.
+            For example, if the donor and acceptor selections include both protein and water, it is possible to
+            find only protein-water hydrogen bonds - and not protein-protein or water-water - by specifying
+            :attr`between=["protein", "SOL"]`. If a two-dimensional list is passed, hydrogen bonds between each pair
+            will be found. For example, :attr`between=[["protein", "SOL"], ["protein", "protein"]]` will calculate
+            all protein-water and protein-protein hydrogen bonds but not water-water hydrogen bonds. If None, hydrogen
+            bonds between all donors and acceptors will be calculated.
         d_h_cutoff : float (optional)
             Distance cutoff used for finding donor-hydrogen pairs [1.2]. Only used to find donor-hydrogen pairs if the
             universe topology does not contain bonding information
@@ -221,7 +230,7 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
         Note
         ----
 
-        It is highly recommended that a universe topology with bonding information is used, as this is the only way
+        It is highly recommended that a universe topology with bond information is used, as this is the only way
         that guarantees the correct identification of donor-hydrogen pairs.
         """
 
@@ -230,6 +239,29 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
         self.donors_sel = donors_sel
         self.hydrogens_sel = hydrogens_sel
         self.acceptors_sel = acceptors_sel
+
+        # If hydrogen bonding groups are selected, then generate corresponding atom groups
+        if between is not None:
+
+            if type(between) is not list or len(between)==0:
+                raise ValueError("between must be a non-empty list")
+            if type(between[0]) == str:
+                between = [between]
+
+            between_ags = []
+            for group1, group2 in between:
+                between_ags.append(
+                    [
+                        self.u.select_atoms(group1, update_selections=False),
+                        self.u.select_atoms(group2, update_selections=False)
+                    ]
+                )
+
+            self.between_ags = between_ags
+        else:
+            self.between_ags = None
+
+
         self.d_h_cutoff = d_h_cutoff
         self.d_a_cutoff = d_a_cutoff
         self.d_h_a_angle = d_h_a_angle_cutoff
@@ -430,6 +462,44 @@ def _get_dh_pairs(self):
 
         return donors, hydrogens
 
+    def _filter(self):
+        """Filter donor, hydrogen and acceptors atoms to consider only hydrogen bonds between two or more
+        specified groups.
+
+           Returns
+           -------
+           None
+        """
+
+        if self.between_ags is None:
+            return
+
+        mask = np.full(self._donors.n_atoms, fill_value=False)
+        for group1, group2 in self.between_ags:
+
+            # Find donors in G1 and acceptors in G2
+            mask[
+                    np.logical_and(
+                        np.in1d(self._donors.indices, group1.indices),
+                        np.in1d(self._accetpors.indices, group2.indices)
+                    )
+            ] = True
+
+            # Find acceptors in G1 and donors in G2
+            mask[
+                np.logical_and(
+                    np.in1d(self._acceptors.indices, group1.indices),
+                    np.in1d(self._donors.indices, group2.indices)
+                )
+            ] = True
+
+        self._donors = self._donors[mask]
+        self._hydrogens = self._hydrogens[mask]
+        self._acceptors = self._acceptors[mask]
+
+        return None
+
+
     def _prepare(self):
         self.hbonds = [[], [], [], [], [], []]
 
@@ -444,6 +514,9 @@ def _prepare(self):
                                               updating=self.update_selections)
         self._donors, self._hydrogens = self._get_dh_pairs()
 
+        # If necessary, filter the donors, hydrogens and acceptors
+        self._filter()
+
     def _single_frame(self):
 
         box = self._ts.dimensions
@@ -452,6 +525,9 @@ def _single_frame(self):
         if self.update_selections:
             self._donors, self._hydrogens = self._get_dh_pairs()
 
+            # If necessary, filter the donors, hydrogens and acceptors
+            self._filter()
+
         # find D and A within cutoff distance of one another
         # min_cutoff = 1.0 as an atom cannot form a hydrogen bond with itself
         d_a_indices, d_a_distances = capped_distance(
@@ -529,13 +605,10 @@ def autocorrelation(self, tau_max=20, window_step=1, intermittency=0):
 
         Returns
         -------
-        acf_tau_timeseries : list
-            tau from 1 to `tau_max`. Saved in the field acf_tau_timeseries.
-        acf_timeseries : list
-            autcorrelation value for each value of `tau`. Saved in the field acf_timeseries.
-        acf_timeseries_data: list
-            raw datapoints, of hydrogen and acceptor indices, from which the autocorrelation is calculated.
-             Saved in the field acf_timeseries_data.
+        tau_timeseries : np.array
+            tau from 1 to `tau_max`
+        timeseries : np.array
+            autcorrelation value for each value of `tau`
         """
 
         if not hasattr(self, 'hbonds'):
@@ -556,12 +629,8 @@ def autocorrelation(self, tau_max=20, window_step=1, intermittency=0):
         intermittent_hbonds = correct_intermittency(found_hydrogen_bonds, intermittency=intermittency)
         tau_timeseries, timeseries, timeseries_data = autocorrelation(intermittent_hbonds, tau_max,
                                                                       window_step=window_step)
-        self.acf_tau_timeseries = tau_timeseries
-        self.acf_timeseries = timeseries
-        # user can investigate the distribution and sample size
-        self.acf_timeseries_data = timeseries_data
 
-        return self
+        return np.vstack([tau_timeseries, timeseries])
 
     def count_by_time(self):
         """Counts the number of hydrogen bonds per timestep.
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index 7ef09635367..6f1d1970b2e 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -210,31 +210,31 @@ def test_count_by_time(self, universe):
         assert_array_equal(counts, ref_counts)
 
     def test_hydrogen_bond_lifetime(self, universe):
-        h = HydrogenBondAnalysis(universe, **self.kwargs)
-        h.run()
+        hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
+        hbonds.run()
 
-        h.autocorrelation(tau_max=2)
-        assert_array_equal(h.acf_timeseries, [1, 0, 0])
+        tau_timeseries, timeseries = hbonds.lifetime(tau_max=2)
+        assert_array_equal(timeseries, [1, 0, 0])
 
     def test_hydrogen_bond_lifetime_intermittency(self, universe):
-        h = HydrogenBondAnalysis(universe, **self.kwargs)
-        h.run()
+        hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
+        hbonds.run()
 
-        h.autocorrelation(tau_max=2, intermittency=1)
-        assert_array_equal(h.acf_timeseries, 1)
+        tau_timeseries, timeseries = hbonds.lifetime(tau_max=2, intermittency=1)
+        assert_array_equal(timeseries, 1)
 
     def test_no_attr_hbonds(self, universe):
-        h = HydrogenBondAnalysis(universe, **self.kwargs)
+        hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
 
         with pytest.raises(NoDataError, match="No .hbonds"):
-            h.autocorrelation(tau_max=2, intermittency=1)
+            hbonds.lifetime(tau_max=2, intermittency=1)
 
     def test_logging_step_not_1(self, universe, caplog):
-        h = HydrogenBondAnalysis(universe, **self.kwargs)
-        h.run(step=2)
+        hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
+        hbonds.run(step=2)
 
         caplog.set_level(logging.WARNING)
-        h.autocorrelation(tau_max=2, intermittency=1)
+        hbonds.lifetime(tau_max=2, intermittency=1)
         assert any("ideally autocorrelation function would be carried out on consecutive frames" in
                    rec.getMessage() for rec in caplog.records)
 

From d2b947a055b06efe3626f3556058cca39e88576c Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Mon, 29 Jun 2020 16:24:26 +0100
Subject: [PATCH 06/16] Missing function rename

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index 7c85dac4641..24db963ec24 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -574,7 +574,7 @@ def _conclude(self):
 
         self.hbonds = np.asarray(self.hbonds).T
 
-    def autocorrelation(self, tau_max=20, window_step=1, intermittency=0):
+    def lifetime(self, tau_max=20, window_step=1, intermittency=0):
         """
         Computes and returns the time-autocorrelation (HydrogenBondLifetimes) of hydrogen bonds.
 

From 865f35feb3bbda0c1d729c49edfefc7253421896 Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Mon, 29 Jun 2020 17:50:23 +0100
Subject: [PATCH 07/16] Added ability to select groups between which to find
 hydrogen bonds.

The `between` keyword can be used to specify pairs of groups between
which hydrogen bonds will be calculated. Hydrogen bonds found
other pairs of atom groups will be discarded.

Added test cases for filtering hydrogen bonds.

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 .../analysis/hydrogenbonds/hbond_analysis.py  |  44 ++++----
 .../analysis/test_hydrogenbonds_analysis.py   | 105 ++++++++++++++++++
 2 files changed, 125 insertions(+), 24 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index 24db963ec24..178da0686df 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -252,8 +252,8 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
             for group1, group2 in between:
                 between_ags.append(
                     [
-                        self.u.select_atoms(group1, update_selections=False),
-                        self.u.select_atoms(group2, update_selections=False)
+                        self.u.select_atoms(group1, updating=False),
+                        self.u.select_atoms(group2, updating=False)
                     ]
                 )
 
@@ -462,42 +462,37 @@ def _get_dh_pairs(self):
 
         return donors, hydrogens
 
-    def _filter(self):
-        """Filter donor, hydrogen and acceptors atoms to consider only hydrogen bonds between two or more
+    def _filter_atoms(self, donors, hydrogens, acceptors):
+        """Filter donor, hydrogen and acceptor atoms to consider only hydrogen bonds between two or more
         specified groups.
 
+        Groups are specified with the `between` keyword when creating the HydrogenBondAnalysis object.
+
            Returns
            -------
-           None
+           donors, hydrogens, acceptors: Filtered AtomGroups
         """
 
-        if self.between_ags is None:
-            return
-
-        mask = np.full(self._donors.n_atoms, fill_value=False)
+        mask = np.full(donors.n_atoms, fill_value=False)
         for group1, group2 in self.between_ags:
 
             # Find donors in G1 and acceptors in G2
             mask[
                     np.logical_and(
-                        np.in1d(self._donors.indices, group1.indices),
-                        np.in1d(self._accetpors.indices, group2.indices)
+                        np.in1d(donors.indices, group1.indices),
+                        np.in1d(acceptors.indices, group2.indices)
                     )
             ] = True
 
             # Find acceptors in G1 and donors in G2
             mask[
                 np.logical_and(
-                    np.in1d(self._acceptors.indices, group1.indices),
-                    np.in1d(self._donors.indices, group2.indices)
+                    np.in1d(acceptors.indices, group1.indices),
+                    np.in1d(donors.indices, group2.indices)
                 )
             ] = True
 
-        self._donors = self._donors[mask]
-        self._hydrogens = self._hydrogens[mask]
-        self._acceptors = self._acceptors[mask]
-
-        return None
+        return donors[mask], hydrogens[mask], acceptors[mask]
 
 
     def _prepare(self):
@@ -514,9 +509,6 @@ def _prepare(self):
                                               updating=self.update_selections)
         self._donors, self._hydrogens = self._get_dh_pairs()
 
-        # If necessary, filter the donors, hydrogens and acceptors
-        self._filter()
-
     def _single_frame(self):
 
         box = self._ts.dimensions
@@ -525,9 +517,6 @@ def _single_frame(self):
         if self.update_selections:
             self._donors, self._hydrogens = self._get_dh_pairs()
 
-            # If necessary, filter the donors, hydrogens and acceptors
-            self._filter()
-
         # find D and A within cutoff distance of one another
         # min_cutoff = 1.0 as an atom cannot form a hydrogen bond with itself
         d_a_indices, d_a_distances = capped_distance(
@@ -544,6 +533,10 @@ def _single_frame(self):
         tmp_hydrogens = self._hydrogens[d_a_indices.T[0]]
         tmp_acceptors = self._acceptors[d_a_indices.T[1]]
 
+        # Remove donor-acceptor pairs between pairs of AtomGroups we are not interested in
+        if self.between_ags is not None:
+            tmp_donors, tmp_hydrogens, tmp_acceptors = self._filter_atoms(tmp_donors, tmp_hydrogens, tmp_acceptors)
+
         # Find D-H-A angles greater than d_h_a_angle_cutoff
         d_h_a_angles = np.rad2deg(
             calc_angles(
@@ -592,6 +585,9 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
 
         Parameters
         ----------
+        of : list, optional
+            A tuple ("A", "B") such as ("Protein", "Ligand"), or a list of tuples, e.g [("A", "B"), ("A", "A and around 5 B")]
+            for which the lifetime should be calculated.
         window_step : int, optional
             The number of frames between each t(0).
         tau_max : int, optional
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index 6f1d1970b2e..84cfa9d1c7e 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -239,6 +239,111 @@ def test_logging_step_not_1(self, universe, caplog):
                    rec.getMessage() for rec in caplog.records)
 
 
+class TestHydrogenBondAnalysisBetween(object):
+
+    kwargs = {
+        'donors_sel': 'name O P',
+        'hydrogens_sel': 'name H1 H2 PH',
+        'acceptors_sel': 'name O P',
+        'd_h_cutoff': 1.2,
+        'd_a_cutoff': 3.0,
+        'd_h_a_angle_cutoff': 120.0
+    }
+
+    @staticmethod
+    @pytest.fixture(scope='class')
+    def universe():
+        # create two water molecules and two "protein" molecules
+        # P1-PH1 are the two atoms that comprise the toy protein PROT1
+        """
+                       H4
+                        \
+            O1-H2 .... O2-H3 ... P1-PH1 ... P2-PH2
+           /
+          H1
+        """
+        n_residues = 4
+        n_sol_residues = 2
+        n_prot_residues = 2
+        u = MDAnalysis.Universe.empty(
+            n_atoms=n_sol_residues*3 + n_prot_residues*2,
+            n_residues=n_residues,
+            atom_resindex=[0, 0, 0, 1, 1, 1, 2, 2, 3, 3],
+            residue_segindex=[0, 0, 1, 1],
+            trajectory=True,  # necessary for adding coordinates
+            )
+
+        u.add_TopologyAttr('name', ['O', 'H1', 'H2'] * n_sol_residues + ['P', 'PH'] * n_prot_residues)
+        u.add_TopologyAttr('type', ['O', 'H', 'H'] * n_sol_residues + ['P', 'PH'] * n_prot_residues)
+        u.add_TopologyAttr('resname', ['SOL'] * n_sol_residues + ['PROT'] * n_prot_residues)
+        u.add_TopologyAttr('resid', list(range(1, n_residues + 1)))
+        u.add_TopologyAttr('id', list(range(1, (n_sol_residues * 3 + n_prot_residues * 2) + 1)))
+
+        # Atomic coordinates with hydrogen bonds between:
+        #     O1-H2---O2
+        #     O2-H3---P1
+        #     P1-PH1---P2
+        pos = np.array([[0, 0, 0],             # O1
+                        [-0.249, -0.968, 0],    # H1
+                        [1, 0, 0],              # H2
+                        [2.5, 0, 0],            # O2
+                        [3., 0, 0],             # H3
+                        [2.250, 0.968, 0],      # H4
+                        [5.5, 0, 0],             # P1
+                        [6.5, 0, 0],             # PH1
+                        [8.5, 0, 0],             # P2
+                        [9.5, 0, 0],             # PH2
+                        ])
+
+        coordinates = np.empty((1,  # number of frames
+                                u.atoms.n_atoms,
+                                3))
+        coordinates[0] = pos
+        u.load_new(coordinates, order='fac')
+
+        return u
+
+    def test_between_all(self, universe):
+        # don't specify groups between which to find hydrogen bonds
+        hbonds = HydrogenBondAnalysis(universe, between=None, **self.kwargs)
+        hbonds.run()
+
+        # indices of [donor, hydrogen, acceptor] for each hydrogen bond
+        expected_hbond_indices = [
+            [0, 2, 3], # water-water
+            [3, 4, 6], # protein-water
+            [6, 7, 8]  # protein-protein
+        ]
+        assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices)
+
+    def test_between_PW(self, universe):
+        # Find only protein-water hydrogen bonds
+        hbonds = HydrogenBondAnalysis(universe, between=["resname PROT", "resname SOL"], **self.kwargs)
+        hbonds.run()
+
+        # indices of [donor, hydrogen, acceptor] for each hydrogen bond
+        expected_hbond_indices = [
+            [3, 4, 6] # protein-water
+        ]
+        assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices)
+
+    def test_between_PW_PP(self, universe):
+        # Find protein-water and protein-protein hydrogen bonds (not water-water)
+        hbonds = HydrogenBondAnalysis(
+            universe,
+            between=[["resname PROT", "resname SOL"], ["resname PROT", "resname PROT"]],
+            **self.kwargs
+        )
+        hbonds.run()
+
+        # indices of [donor, hydrogen, acceptor] for each hydrogen bond
+        expected_hbond_indices = [
+            [3, 4, 6], # protein-water
+            [6, 7, 8] # protein-protein
+        ]
+        assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices)
+
+
 class TestHydrogenBondAnalysisTIP3P_GuessAcceptors_GuessHydrogens_UseTopology_(TestHydrogenBondAnalysisTIP3P):
     """Uses the same distance and cutoff hydrogen bond criteria as :class:`TestHydrogenBondAnalysisTIP3P`, so the
     results are identical, but the hydrogens and acceptors are guessed whilst the donor-hydrogen pairs are determined

From b0d341b9dbd9bcd960b8f50e6ba162bbdabe1cb4 Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Sun, 5 Jul 2020 17:10:51 +0100
Subject: [PATCH 08/16] Addressing feedback in PR #2791

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 package/CHANGELOG                             |  1 +
 .../analysis/hydrogenbonds/hbond_analysis.py  | 52 +++++++++++++++----
 package/MDAnalysis/analysis/waterdynamics.py  |  5 +-
 .../analysis/test_hydrogenbonds_analysis.py   | 45 +++++++---------
 4 files changed, 66 insertions(+), 37 deletions(-)

diff --git a/package/CHANGELOG b/package/CHANGELOG
index 1bb374bfe76..229e8275445 100644
--- a/package/CHANGELOG
+++ b/package/CHANGELOG
@@ -33,6 +33,7 @@ Enhancements
   * Added the `from_smiles` classmethod to the Universe (Issue #2468, PR #2707)
   * Added computation of Mean Squared Displacements (#2438, PR #2619)
   * Added Hydrogen Bond Lifetime via existing autocorrelation features (PR #2791)
+  * Added Hydrogen Bond Lifetime keyword "between" (PR #2791)
 
 Changes
   * Changes development status from Beta to Mature (Issue #2773)
diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index 178da0686df..a3fcab1d205 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -141,6 +141,35 @@
   hbonds.acceptors_sel = f"({protein_acceptors_sel}) or ({water_acceptors_sel} and around 10 not resname TIP3})"
   hbonds.run()
 
+To calculate the hydrogen bonds between different groups, for example a protein and a ligand, one can use the
+:attr:`between` keyword:
+
+  import MDAnalysis
+  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
+
+  u = MDAnalysis.Universe(psf, trajectory)
+
+  hbonds = HBA(universe=u, hydrogens_sel='resname TIP3 and name H1 H2', acceptors_sel='resname TIP3 and name OH2',
+                between=['resname LIG', 'protein'])
+  hbonds.run()
+
+It is further possible to compute hydrogen bonds between several groups with with use of :attr:`between`
+
+  import MDAnalysis
+  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
+
+  u = MDAnalysis.Universe(psf, trajectory)
+
+  hbonds = HBA(universe=u, hydrogens_sel='resname TIP3 and name H1 H2', acceptors_sel='resname TIP3 and name OH2',
+                between=[['resname LIG', 'protein'], ['resname LIG', 'resname WAT']])
+  hbonds.run()
+
+In order to compute the hydrogen bond lifetime, after computing one can use the :attr:`lifetime` function
+
+    ...
+    hbonds.run()
+    tau_timeseries, timeseries = hbonds.lifetime()
+
 It is highly recommended that a topology with bonding information is used to generate the universe, e.g `PSF`, `TPR`, or
 `PRMTOP` files. This is the only method by which it can be guaranteed that donor-hydrogen pairs are correctly identified.
 However, if, for example, a `PDB` file is used instead, a :attr:`donors_sel` may be provided along with a
@@ -167,6 +196,7 @@
    :members:
 """
 import logging
+from collections.abc import Iterable
 
 import numpy as np
 
@@ -209,14 +239,14 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
         acceptors_sel : str
             Selection string for the hydrogen bond acceptor atoms. Leave as `None` to guess which atoms to use in the
             analysis using :attr:`guess_acceptors`
-        between : List (optional)
+        between : List (optional),
             Specify two selection strings for non-updating atom groups between which hydrogen bonds will be calculated.
             For example, if the donor and acceptor selections include both protein and water, it is possible to
             find only protein-water hydrogen bonds - and not protein-protein or water-water - by specifying
             :attr`between=["protein", "SOL"]`. If a two-dimensional list is passed, hydrogen bonds between each pair
             will be found. For example, :attr`between=[["protein", "SOL"], ["protein", "protein"]]` will calculate
-            all protein-water and protein-protein hydrogen bonds but not water-water hydrogen bonds. If None, hydrogen
-            bonds between all donors and acceptors will be calculated.
+            all protein-water and protein-protein hydrogen bonds but not water-water hydrogen bonds. If `None`, hydrogen
+            bonds between all donors and acceptors will be calculated. See examples
         d_h_cutoff : float (optional)
             Distance cutoff used for finding donor-hydrogen pairs [1.2]. Only used to find donor-hydrogen pairs if the
             universe topology does not contain bonding information
@@ -232,6 +262,9 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
 
         It is highly recommended that a universe topology with bond information is used, as this is the only way
         that guarantees the correct identification of donor-hydrogen pairs.
+
+        .. versionadded:: 2.0.0
+            Added `between` keyword
         """
 
         self.u = universe
@@ -242,10 +275,9 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
 
         # If hydrogen bonding groups are selected, then generate corresponding atom groups
         if between is not None:
-
-            if type(between) is not list or len(between)==0:
-                raise ValueError("between must be a non-empty list")
-            if type(between[0]) == str:
+            if not isinstance(between, Iterable) or len(between) == 0:
+                raise ValueError("between must be a non-empty list/iterable")
+            if isinstance(between[0], str):
                 between = [between]
 
             between_ags = []
@@ -585,9 +617,6 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
 
         Parameters
         ----------
-        of : list, optional
-            A tuple ("A", "B") such as ("Protein", "Ligand"), or a list of tuples, e.g [("A", "B"), ("A", "A and around 5 B")]
-            for which the lifetime should be calculated.
         window_step : int, optional
             The number of frames between each t(0).
         tau_max : int, optional
@@ -613,7 +642,8 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
             raise NoDataError("No .hbonds: use the .run() first")
 
         if self.step != 1:
-            logging.warning("Autocorrelation: ideally autocorrelation function would be carried out on consecutive frames. ")
+            logging.warning("Autocorrelation: Hydrogen bonds were computed with step > 1.  ")
+            logging.warning("Autocorrelation: We recommend recomputing hydrogen bonds with step = 1. ")
             logging.warning("Autocorrelation: if you would like to allow bonds to break and reform, please use 'intermittency'")
 
         # Extract the hydrogen bonds IDs only in the format [set(superset(x1,x2), superset(x3,x4)), ..]
diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index 6bec1d5a1f7..6ad23b8613a 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -508,6 +508,9 @@ class HydrogenBondLifetimes(object):
     .. versionchanged:: 1.0.0
        The ``nproc`` keyword was removed as it linked to a portion of code that
        may have failed in some cases.
+    .. deprecated:: 2.0.0
+       Please use MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime() instead.
+
     """
 
 
@@ -515,7 +518,7 @@ def __init__(self, universe, selection1, selection2, t0, tf, dtmax,
                  nproc=1):
         warnings.warn(
             "This class is deprecated. "
-            "Please use MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.autocorellation() instead.",
+            "Please use MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime() instead.",
             category=DeprecationWarning
         )
 
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index 84cfa9d1c7e..05f56d20d89 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -168,6 +168,14 @@ def universe():
 
         return u
 
+    @staticmethod
+    @pytest.fixture(scope='class')
+    def hydrogen_bonds(universe):
+        h = HydrogenBondAnalysis(universe, **TestHydrogenBondAnalysisMock.kwargs)
+        h.run()
+        return h
+
+
     def test_no_bond_info_exception(self, universe):
 
         kwargs = {
@@ -183,54 +191,41 @@ def test_no_bond_info_exception(self, universe):
             h = HydrogenBondAnalysis(universe, **kwargs)
             h._get_dh_pairs()
 
-    def test_first(self, universe):
+    def test_first(self, hydrogen_bonds):
+        assert len(hydrogen_bonds.hbonds) == 2
 
-        h = HydrogenBondAnalysis(universe, **self.kwargs)
-        h.run()
-
-        assert len(h.hbonds) == 2
-
-        frame_no, donor_index, hydrogen_index, acceptor_index, da_dst, dha_angle = h.hbonds[0]
+        frame_no, donor_index, hydrogen_index, acceptor_index, da_dst, dha_angle = hydrogen_bonds.hbonds[0]
         assert_equal(donor_index, 0)
         assert_equal(hydrogen_index, 2)
         assert_equal(acceptor_index, 3)
         assert_almost_equal(da_dst, 2.5)
         assert_almost_equal(dha_angle, 180)
 
-    def test_count_by_time(self, universe):
-
-        h = HydrogenBondAnalysis(universe, **self.kwargs)
-        h.run()
-
+    def test_count_by_time(self, hydrogen_bonds):
         ref_times = np.array([0, 1, 2]) # u.trajectory.dt is 1
         ref_counts = np.array([1, 0, 1])
 
-        counts = h.count_by_time()
-        assert_array_almost_equal(h.times, ref_times)
+        counts = hydrogen_bonds.count_by_time()
+        assert_array_almost_equal(hydrogen_bonds.times, ref_times)
         assert_array_equal(counts, ref_counts)
 
-    def test_hydrogen_bond_lifetime(self, universe):
-        hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
-        hbonds.run()
-
-        tau_timeseries, timeseries = hbonds.lifetime(tau_max=2)
+    def test_hydrogen_bond_lifetime(self, hydrogen_bonds):
+        tau_timeseries, timeseries = hydrogen_bonds.lifetime(tau_max=2)
         assert_array_equal(timeseries, [1, 0, 0])
 
-    def test_hydrogen_bond_lifetime_intermittency(self, universe):
-        hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
-        hbonds.run()
-
-        tau_timeseries, timeseries = hbonds.lifetime(tau_max=2, intermittency=1)
+    def test_hydrogen_bond_lifetime_intermittency(self, hydrogen_bonds):
+        tau_timeseries, timeseries = hydrogen_bonds.lifetime(tau_max=2, intermittency=1)
         assert_array_equal(timeseries, 1)
 
     def test_no_attr_hbonds(self, universe):
         hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
-
+        # hydrogen bonds are not computed
         with pytest.raises(NoDataError, match="No .hbonds"):
             hbonds.lifetime(tau_max=2, intermittency=1)
 
     def test_logging_step_not_1(self, universe, caplog):
         hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
+        # using step 2
         hbonds.run(step=2)
 
         caplog.set_level(logging.WARNING)

From edb717bc5fac38fbd512d755ac937b6bbfe03196 Mon Sep 17 00:00:00 2001
From: Paul Smith <paul.smith@kcl.ac.uk>
Date: Mon, 6 Jul 2020 13:26:28 +0100
Subject: [PATCH 09/16] Silence pep8speaks

Fixed the line length and whitespace issues
---
 .../analysis/hydrogenbonds/hbond_analysis.py  | 175 ++++++++++++------
 package/MDAnalysis/analysis/waterdynamics.py  |   6 +-
 package/MDAnalysis/lib/correlations.py        |   7 +-
 .../analysis/test_hydrogenbonds_analysis.py   |  90 ++++++---
 4 files changed, 186 insertions(+), 92 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index a3fcab1d205..eea9f97fc4b 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -141,30 +141,42 @@
   hbonds.acceptors_sel = f"({protein_acceptors_sel}) or ({water_acceptors_sel} and around 10 not resname TIP3})"
   hbonds.run()
 
-To calculate the hydrogen bonds between different groups, for example a protein and a ligand, one can use the
-:attr:`between` keyword:
+To calculate the hydrogen bonds between different groups, for example a
+protein and a ligand, one can use the :attr:`between` keyword:
 
   import MDAnalysis
-  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
+  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
+  HydrogenBondAnalysis as HBA)
 
   u = MDAnalysis.Universe(psf, trajectory)
 
-  hbonds = HBA(universe=u, hydrogens_sel='resname TIP3 and name H1 H2', acceptors_sel='resname TIP3 and name OH2',
-                between=['resname LIG', 'protein'])
+  hbonds = HBA(
+      universe=u,
+      hydrogens_sel='resname TIP3 and name H1 H2',
+      acceptors_sel='resname TIP3 and name OH2',
+      between=['resname LIG', 'protein']
+      )
   hbonds.run()
 
-It is further possible to compute hydrogen bonds between several groups with with use of :attr:`between`
+It is further possible to compute hydrogen bonds between several groups with
+with use of :attr:`between`
 
   import MDAnalysis
-  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
+  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
+  HydrogenBondAnalysis as HBA)
 
   u = MDAnalysis.Universe(psf, trajectory)
 
-  hbonds = HBA(universe=u, hydrogens_sel='resname TIP3 and name H1 H2', acceptors_sel='resname TIP3 and name OH2',
-                between=[['resname LIG', 'protein'], ['resname LIG', 'resname WAT']])
+  hbonds = HBA(
+      universe=u,
+      hydrogens_sel='resname TIP3 and name H1 H2',
+      acceptors_sel='resname TIP3 and name OH2',
+      between=[['resname LIG', 'protein'], ['resname LIG', 'resname WAT']]
+      )
   hbonds.run()
 
-In order to compute the hydrogen bond lifetime, after computing one can use the :attr:`lifetime` function
+In order to compute the hydrogen bond lifetime, after computing one can use
+the :attr:`lifetime` function
 
     ...
     hbonds.run()
@@ -220,33 +232,45 @@ class HydrogenBondAnalysis(base.AnalysisBase):
     Perform an analysis of hydrogen bonds in a Universe.
     """
 
-    def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=None,
-                 between=None, d_h_cutoff=1.2, d_a_cutoff=3.0, d_h_a_angle_cutoff=150,
+    def __init__(self, universe,
+                 donors_sel=None, hydrogens_sel=None, acceptors_sel=None,
+                 between=None, d_h_cutoff=1.2,
+                 d_a_cutoff=3.0, d_h_a_angle_cutoff=150,
                  update_selections=True):
-        """Set up atom selections and geometric criteria for finding hydrogen bonds in a Universe.
+        """Set up atom selections and geometric criteria for finding hydrogen
+        bonds in a Universe.
 
         Parameters
         ----------
         universe : Universe
             MDAnalysis Universe object
         donors_sel : str
-            Selection string for the hydrogen bond donor atoms. If the universe topology contains bonding information,
-            leave :attr:`donors_sel` as `None` so that donor-hydrogen pairs can be correctly identified.
+            Selection string for the hydrogen bond donor atoms. If the
+            universe topology contains bonding information, leave
+            :attr:`donors_sel` as `None` so that donor-hydrogen pairs can be
+            correctly identified.
         hydrogens_sel :  str
-            Selection string for the hydrogen bond hydrogen atoms. Leave as `None` to guess which hydrogens to use in
-            the analysis using :attr:`guess_hydrogens`. If :attr:`hydrogens_sel` is left as `None`, also leave
-            :attr:`donors_sel` as None so that donor-hydrogen pairs can be correctly identified.
+            Selection string for the hydrogen bond hydrogen atoms. Leave as
+            `None` to guess which hydrogens to use in the analysis using
+            :attr:`guess_hydrogens`. If :attr:`hydrogens_sel` is left as
+            `None`, also leave :attr:`donors_sel` as None so that
+            donor-hydrogen pairs can be correctly identified.
         acceptors_sel : str
-            Selection string for the hydrogen bond acceptor atoms. Leave as `None` to guess which atoms to use in the
-            analysis using :attr:`guess_acceptors`
+            Selection string for the hydrogen bond acceptor atoms. Leave as
+            `None` to guess which atoms to use in the analysis using
+            :attr:`guess_acceptors`
         between : List (optional),
-            Specify two selection strings for non-updating atom groups between which hydrogen bonds will be calculated.
-            For example, if the donor and acceptor selections include both protein and water, it is possible to
-            find only protein-water hydrogen bonds - and not protein-protein or water-water - by specifying
-            :attr`between=["protein", "SOL"]`. If a two-dimensional list is passed, hydrogen bonds between each pair
-            will be found. For example, :attr`between=[["protein", "SOL"], ["protein", "protein"]]` will calculate
-            all protein-water and protein-protein hydrogen bonds but not water-water hydrogen bonds. If `None`, hydrogen
-            bonds between all donors and acceptors will be calculated. See examples
+            Specify two selection strings for non-updating atom groups between
+            which hydrogen bonds will be calculated. For example, if the donor
+            and acceptor selections include both protein and water, it is
+            possible to find only protein-water hydrogen bonds - and not
+            protein-protein or water-water - by specifying
+            :attr`between=["protein", "SOL"]`. If a two-dimensional list is
+            passed, hydrogen bonds between each pair will be found. For
+            example, :attr`between=[["protein", "SOL"], ["protein", "protein"]]`
+            will calculate all protein-water and protein-protein hydrogen
+            bonds but not water-water hydrogen bonds. If `None`, hydrogen
+            bonds between all donors and acceptors will be calculated.
         d_h_cutoff : float (optional)
             Distance cutoff used for finding donor-hydrogen pairs [1.2]. Only used to find donor-hydrogen pairs if the
             universe topology does not contain bonding information
@@ -260,8 +284,9 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
         Note
         ----
 
-        It is highly recommended that a universe topology with bond information is used, as this is the only way
-        that guarantees the correct identification of donor-hydrogen pairs.
+        It is highly recommended that a universe topology with bond
+        information is used, as this is the only way that guarantees the
+        correct identification of donor-hydrogen pairs.
 
         .. versionadded:: 2.0.0
             Added `between` keyword
@@ -273,7 +298,8 @@ def __init__(self, universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=
         self.hydrogens_sel = hydrogens_sel
         self.acceptors_sel = acceptors_sel
 
-        # If hydrogen bonding groups are selected, then generate corresponding atom groups
+        # If hydrogen bonding groups are selected, then generate
+        # corresponding atom groups
         if between is not None:
             if not isinstance(between, Iterable) or len(between) == 0:
                 raise ValueError("between must be a non-empty list/iterable")
@@ -495,10 +521,11 @@ def _get_dh_pairs(self):
         return donors, hydrogens
 
     def _filter_atoms(self, donors, hydrogens, acceptors):
-        """Filter donor, hydrogen and acceptor atoms to consider only hydrogen bonds between two or more
-        specified groups.
+        """Filter donor, hydrogen and acceptor atoms to consider only hydrogen
+        bonds between two or more specified groups.
 
-        Groups are specified with the `between` keyword when creating the HydrogenBondAnalysis object.
+        Groups are specified with the `between` keyword when creating the
+        HydrogenBondAnalysis object.
 
            Returns
            -------
@@ -565,9 +592,11 @@ def _single_frame(self):
         tmp_hydrogens = self._hydrogens[d_a_indices.T[0]]
         tmp_acceptors = self._acceptors[d_a_indices.T[1]]
 
-        # Remove donor-acceptor pairs between pairs of AtomGroups we are not interested in
+        # Remove donor-acceptor pairs between pairs of AtomGroups we are not
+        # interested in
         if self.between_ags is not None:
-            tmp_donors, tmp_hydrogens, tmp_acceptors = self._filter_atoms(tmp_donors, tmp_hydrogens, tmp_acceptors)
+            tmp_donors, tmp_hydrogens, tmp_acceptors = \
+                self._filter_atoms(tmp_donors, tmp_hydrogens, tmp_acceptors)
 
         # Find D-H-A angles greater than d_h_a_angle_cutoff
         d_h_a_angles = np.rad2deg(
@@ -601,18 +630,23 @@ def _conclude(self):
 
     def lifetime(self, tau_max=20, window_step=1, intermittency=0):
         """
-        Computes and returns the time-autocorrelation (HydrogenBondLifetimes) of hydrogen bonds.
+        Computes and returns the time-autocorrelation (HydrogenBondLifetimes)
+         of hydrogen bonds.
 
-        Before calling this method, the hydrogen bonds must first be computed with the `run()` function.
-        The same `start`, `stop` and `step` parameters used in finding hydrogen bonds will be used here
-        for calculating hydrogen bond lifetimes. That is, the same frames will be used in the analysis.
+        Before calling this method, the hydrogen bonds must first be computed
+        with the `run()` function. The same `start`, `stop` and `step`
+        parameters used in finding hydrogen bonds will be used here for
+        calculating hydrogen bond lifetimes. That is, the same frames will be
+        used in the analysis.
 
-        Unique hydrogen bonds are identified using hydrogen-acceptor pairs. This means an acceptor
-        switching to a different hydrogen atom - with the same donor - from one frame to the next
-        is considered a different hydrogen bond.
+        Unique hydrogen bonds are identified using hydrogen-acceptor pairs.
+        This means an acceptor switching to a different hydrogen atom - with
+        the same donor - from one frame to the next is considered a different
+        hydrogen bond.
 
         Please see :func:`MDAnalysis.lib.correlations.autocorrelation` and
-        :func:`MDAnalysis.lib.correlations.intermittency` functions for more details.
+        :func:`MDAnalysis.lib.correlations.intermittency` functions for more
+        details.
 
 
         Parameters
@@ -620,13 +654,17 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
         window_step : int, optional
             The number of frames between each t(0).
         tau_max : int, optional
-            Hydrogen bond lifetime is calculated for frames in the range 1 <= `tau` <= `tau_max`
+            Hydrogen bond lifetime is calculated for frames in the range
+            1 <= `tau` <= `tau_max`
         intermittency : int, optional
-            The maximum number of consecutive frames for which a bond can disappear but be counted as present if it
-            returns at the next frame. An intermittency of `0` is equivalent to a continuous autocorrelation, which does
-            not allow for hydrogen bond disappearance. For example, for `intermittency=2`, any given hydrogen bond may
-            disappear for up to two consecutive frames yet be treated as being present at all frames.
-            The default is continuous (intermittency=0).
+            The maximum number of consecutive frames for which a bond can
+            disappear but be counted as present if it returns at the next
+            frame. An intermittency of `0` is equivalent to a continuous
+            autocorrelation, which does not allow for hydrogen bond
+            disappearance. For example, for `intermittency=2`, any given
+            hydrogen bond may disappear for up to two consecutive frames yet
+            be treated as being present at all frames. The default is
+            continuous (intermittency=0).
 
         Returns
         -------
@@ -637,24 +675,45 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
         """
 
         if not hasattr(self, 'hbonds'):
-            logging.error("Autocorrelation analysis of hydrogen bonds cannot be done before the hydrogen bonds are found")
-            logging.error("Autocorrelation: Please use the .run() before calling this function")
+            logging.error(
+                "Autocorrelation analysis of hydrogen bonds cannot be done"
+                "before the hydrogen bonds are found"
+            )
+            logging.error(
+                "Autocorrelation: Please use the .run() before calling this"
+                "function"
+            )
             raise NoDataError("No .hbonds: use the .run() first")
 
         if self.step != 1:
-            logging.warning("Autocorrelation: Hydrogen bonds were computed with step > 1.  ")
-            logging.warning("Autocorrelation: We recommend recomputing hydrogen bonds with step = 1. ")
-            logging.warning("Autocorrelation: if you would like to allow bonds to break and reform, please use 'intermittency'")
+            logging.warning(
+                "Autocorrelation: Hydrogen bonds were computedwith step > 1."
+            )
+            logging.warning(
+                "Autocorrelation: We recommend recomputing hydrogen bonds with"
+                "step = 1."
+            )
+            logging.warning(
+                "Autocorrelation: if you would like to allow bonds to break"
+                "and reform, please use 'intermittency'"
+            )
 
-        # Extract the hydrogen bonds IDs only in the format [set(superset(x1,x2), superset(x3,x4)), ..]
+        # Extract the hydrogen bonds IDs only in the format
+        # [set(superset(x1,x2), superset(x3,x4)), ..]
         found_hydrogen_bonds = [set() for _ in self.frames]
         for frame_index, frame in enumerate(self.frames):
             for hbond in self.hbonds[self.hbonds[:, 0] == frame]:
                 found_hydrogen_bonds[frame_index].add(frozenset(hbond[2:4]))
 
-        intermittent_hbonds = correct_intermittency(found_hydrogen_bonds, intermittency=intermittency)
-        tau_timeseries, timeseries, timeseries_data = autocorrelation(intermittent_hbonds, tau_max,
-                                                                      window_step=window_step)
+        intermittent_hbonds = correct_intermittency(
+            found_hydrogen_bonds,
+            intermittency=intermittency
+        )
+        tau_timeseries, timeseries, timeseries_data = autocorrelation(
+            intermittent_hbonds,
+            tau_max,
+            window_step=window_step
+        )
 
         return np.vstack([tau_timeseries, timeseries])
 
diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index 6ad23b8613a..b5a372069f5 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -509,7 +509,8 @@ class HydrogenBondLifetimes(object):
        The ``nproc`` keyword was removed as it linked to a portion of code that
        may have failed in some cases.
     .. deprecated:: 2.0.0
-       Please use MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime() instead.
+       Instead, please use
+       MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime
 
     """
 
@@ -518,7 +519,8 @@ def __init__(self, universe, selection1, selection2, t0, tf, dtmax,
                  nproc=1):
         warnings.warn(
             "This class is deprecated. "
-            "Please use MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime() instead.",
+            "Instrad, please use"
+            "MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime",
             category=DeprecationWarning
         )
 
diff --git a/package/MDAnalysis/lib/correlations.py b/package/MDAnalysis/lib/correlations.py
index 5a97f857bd4..9b12823b576 100644
--- a/package/MDAnalysis/lib/correlations.py
+++ b/package/MDAnalysis/lib/correlations.py
@@ -236,7 +236,7 @@ def correct_intermittency(list_of_sets, intermittency):
                     continue
 
                 # if the element is absent now
-                if not element in list_of_sets[i + j]:
+                if element not in list_of_sets[i + j]:
                     # increase its absence counter
                     seen_frames_ago[element] += 1
                     continue
@@ -250,8 +250,9 @@ def correct_intermittency(list_of_sets, intermittency):
                 if seen_frames_ago[element] > intermittency:
                     continue
 
-                # the element was absent but returned (within <= intermittency_value)
-                # add it to the frames where it was absent.
+                # The element was absent but returned, i.e.
+                # within <= intermittency_value.
+                # Add it to the frames where it was absent.
                 # Introduce the corrections.
                 for k in range(seen_frames_ago[element], 0, -1):
                     list_of_sets[i + j - k].add(element)
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index 05f56d20d89..d0d565ef2fb 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -29,7 +29,8 @@
                            assert_almost_equal)
 
 import MDAnalysis
-from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis
+from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
+    HydrogenBondAnalysis)
 from MDAnalysis.exceptions import NoDataError
 from MDAnalysisTests.datafiles import waterPSF, waterDCD
 
@@ -171,7 +172,10 @@ def universe():
     @staticmethod
     @pytest.fixture(scope='class')
     def hydrogen_bonds(universe):
-        h = HydrogenBondAnalysis(universe, **TestHydrogenBondAnalysisMock.kwargs)
+        h = HydrogenBondAnalysis(
+            universe,
+            **TestHydrogenBondAnalysisMock.kwargs
+        )
         h.run()
         return h
 
@@ -194,7 +198,8 @@ def test_no_bond_info_exception(self, universe):
     def test_first(self, hydrogen_bonds):
         assert len(hydrogen_bonds.hbonds) == 2
 
-        frame_no, donor_index, hydrogen_index, acceptor_index, da_dst, dha_angle = hydrogen_bonds.hbonds[0]
+        frame_no, donor_index, hydrogen_index, acceptor_index, da_dst, angle =\
+            hydrogen_bonds.hbonds[0]
         assert_equal(donor_index, 0)
         assert_equal(hydrogen_index, 2)
         assert_equal(acceptor_index, 3)
@@ -214,7 +219,9 @@ def test_hydrogen_bond_lifetime(self, hydrogen_bonds):
         assert_array_equal(timeseries, [1, 0, 0])
 
     def test_hydrogen_bond_lifetime_intermittency(self, hydrogen_bonds):
-        tau_timeseries, timeseries = hydrogen_bonds.lifetime(tau_max=2, intermittency=1)
+        tau_timeseries, timeseries = hydrogen_bonds.lifetime(
+            tau_max=2, intermittency=1
+        )
         assert_array_equal(timeseries, 1)
 
     def test_no_attr_hbonds(self, universe):
@@ -230,8 +237,10 @@ def test_logging_step_not_1(self, universe, caplog):
 
         caplog.set_level(logging.WARNING)
         hbonds.lifetime(tau_max=2, intermittency=1)
-        assert any("ideally autocorrelation function would be carried out on consecutive frames" in
-                   rec.getMessage() for rec in caplog.records)
+
+        message = "ideally autocorrelation function would be carried out on" \
+            "consecutive frames"
+        assert any(message in rec.getMessage() for rec in caplog.records)
 
 
 class TestHydrogenBondAnalysisBetween(object):
@@ -268,26 +277,41 @@ def universe():
             trajectory=True,  # necessary for adding coordinates
             )
 
-        u.add_TopologyAttr('name', ['O', 'H1', 'H2'] * n_sol_residues + ['P', 'PH'] * n_prot_residues)
-        u.add_TopologyAttr('type', ['O', 'H', 'H'] * n_sol_residues + ['P', 'PH'] * n_prot_residues)
-        u.add_TopologyAttr('resname', ['SOL'] * n_sol_residues + ['PROT'] * n_prot_residues)
-        u.add_TopologyAttr('resid', list(range(1, n_residues + 1)))
-        u.add_TopologyAttr('id', list(range(1, (n_sol_residues * 3 + n_prot_residues * 2) + 1)))
+        u.add_TopologyAttr(
+            'name',
+            ['O', 'H1', 'H2'] * n_sol_residues + ['P', 'PH'] * n_prot_residues
+        )
+        u.add_TopologyAttr(
+            'type',
+            ['O', 'H', 'H'] * n_sol_residues + ['P', 'PH'] * n_prot_residues
+        )
+        u.add_TopologyAttr(
+            'resname',
+            ['SOL'] * n_sol_residues + ['PROT'] * n_prot_residues
+        )
+        u.add_TopologyAttr(
+            'resid',
+            list(range(1, n_residues + 1))
+        )
+        u.add_TopologyAttr(
+            'id',
+            list(range(1, (n_sol_residues * 3 + n_prot_residues * 2) + 1))
+        )
 
         # Atomic coordinates with hydrogen bonds between:
         #     O1-H2---O2
         #     O2-H3---P1
         #     P1-PH1---P2
         pos = np.array([[0, 0, 0],             # O1
-                        [-0.249, -0.968, 0],    # H1
-                        [1, 0, 0],              # H2
-                        [2.5, 0, 0],            # O2
-                        [3., 0, 0],             # H3
-                        [2.250, 0.968, 0],      # H4
-                        [5.5, 0, 0],             # P1
-                        [6.5, 0, 0],             # PH1
-                        [8.5, 0, 0],             # P2
-                        [9.5, 0, 0],             # PH2
+                        [-0.249, -0.968, 0],   # H1
+                        [1, 0, 0],             # H2
+                        [2.5, 0, 0],           # O2
+                        [3., 0, 0],            # H3
+                        [2.250, 0.968, 0],     # H4
+                        [5.5, 0, 0],           # P1
+                        [6.5, 0, 0],           # PH1
+                        [8.5, 0, 0],           # P2
+                        [9.5, 0, 0],           # PH2
                         ])
 
         coordinates = np.empty((1,  # number of frames
@@ -305,36 +329,44 @@ def test_between_all(self, universe):
 
         # indices of [donor, hydrogen, acceptor] for each hydrogen bond
         expected_hbond_indices = [
-            [0, 2, 3], # water-water
-            [3, 4, 6], # protein-water
-            [6, 7, 8]  # protein-protein
+            [0, 2, 3],  # water-water
+            [3, 4, 6],  # protein-water
+            [6, 7, 8]   # protein-protein
         ]
         assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices)
 
     def test_between_PW(self, universe):
         # Find only protein-water hydrogen bonds
-        hbonds = HydrogenBondAnalysis(universe, between=["resname PROT", "resname SOL"], **self.kwargs)
+        hbonds = HydrogenBondAnalysis(
+            universe,
+            between=["resname PROT", "resname SOL"],
+            **self.kwargs
+        )
         hbonds.run()
 
         # indices of [donor, hydrogen, acceptor] for each hydrogen bond
         expected_hbond_indices = [
-            [3, 4, 6] # protein-water
+            [3, 4, 6]  # protein-water
         ]
         assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices)
 
     def test_between_PW_PP(self, universe):
-        # Find protein-water and protein-protein hydrogen bonds (not water-water)
+        # Find protein-water and protein-protein hydrogen bonds (not
+        # water-water)
         hbonds = HydrogenBondAnalysis(
             universe,
-            between=[["resname PROT", "resname SOL"], ["resname PROT", "resname PROT"]],
+            between=[
+                ["resname PROT", "resname SOL"],
+                ["resname PROT", "resname PROT"]
+            ],
             **self.kwargs
         )
         hbonds.run()
 
         # indices of [donor, hydrogen, acceptor] for each hydrogen bond
         expected_hbond_indices = [
-            [3, 4, 6], # protein-water
-            [6, 7, 8] # protein-protein
+            [3, 4, 6],  # protein-water
+            [6, 7, 8]   # protein-protein
         ]
         assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices)
 

From f6c3e2a7427469c5386807dbc7f3dbce29e83e76 Mon Sep 17 00:00:00 2001
From: Paul Smith <paul.smith@kcl.ac.uk>
Date: Mon, 6 Jul 2020 13:44:47 +0100
Subject: [PATCH 10/16] Fix `between` argument examples.

Previously, a selection involving water and protein was used, but
the example suggested hydrogen bonds between a ligand and protein
would be found.
---
 .../analysis/hydrogenbonds/hbond_analysis.py  | 57 ++++++++++---------
 1 file changed, 30 insertions(+), 27 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index eea9f97fc4b..1eff55e3af3 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -142,53 +142,56 @@
   hbonds.run()
 
 To calculate the hydrogen bonds between different groups, for example a
-protein and a ligand, one can use the :attr:`between` keyword:
+protein and water, one can use the :attr:`between` keyword. The
+following will find protein-water hydrogen bonds but not protein-protein
+or water-water hydrogen bonds:
 
   import MDAnalysis
   from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
-  HydrogenBondAnalysis as HBA)
+    HydrogenBondAnalysis as HBA)
 
   u = MDAnalysis.Universe(psf, trajectory)
 
   hbonds = HBA(
-      universe=u,
-      hydrogens_sel='resname TIP3 and name H1 H2',
-      acceptors_sel='resname TIP3 and name OH2',
-      between=['resname LIG', 'protein']
-      )
-  hbonds.run()
+    universe=u,
+    between=['resname TIP3', 'protein']
+    )
 
-It is further possible to compute hydrogen bonds between several groups with
-with use of :attr:`between`
+  protein_hydrogens_sel = hbonds.guess_hydrogens("protein")
+  protein_acceptors_sel = hbonds.guess_acceptors("protein")
 
-  import MDAnalysis
-  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
-  HydrogenBondAnalysis as HBA)
+  water_hydrogens_sel = "resname TIP3 and name H1 H2"
+  water_acceptors_sel = "resname TIP3 and name OH2"
 
-  u = MDAnalysis.Universe(psf, trajectory)
+  hbonds.hydrogens_sel = f"({protein_hydrogens_sel}) or ({water_hydrogens_sel}"
+  hbonds.acceptors_sel = f"({protein_acceptors_sel}) or ({water_acceptors_sel}"
 
-  hbonds = HBA(
-      universe=u,
-      hydrogens_sel='resname TIP3 and name H1 H2',
-      acceptors_sel='resname TIP3 and name OH2',
-      between=[['resname LIG', 'protein'], ['resname LIG', 'resname WAT']]
-      )
   hbonds.run()
 
-In order to compute the hydrogen bond lifetime, after computing one can use
-the :attr:`lifetime` function
+It is further possible to compute hydrogen bonds between several groups with
+with use of :attr:`between`. If in the above example,
+`between=[['resname TIP3', 'protein'], ['protein', 'protein']]`, all
+protein-water and protein-protein hydrogen bonds will be found, but
+no water-water hydrogen bonds.
+
+In order to compute the hydrogen bond lifetime, after finding hydrogen bonds
+one can use the :attr:`lifetime` function:
 
     ...
     hbonds.run()
     tau_timeseries, timeseries = hbonds.lifetime()
 
-It is highly recommended that a topology with bonding information is used to generate the universe, e.g `PSF`, `TPR`, or
-`PRMTOP` files. This is the only method by which it can be guaranteed that donor-hydrogen pairs are correctly identified.
-However, if, for example, a `PDB` file is used instead, a :attr:`donors_sel` may be provided along with a
-:attr:`hydrogens_sel` and the donor-hydrogen pairs will be identified via a distance cutoff, :attr:`d_h_cutoff`::
+It is **highly recommended** that a topology with bond information is used to
+generate the universe, e.g `PSF`, `TPR`, or `PRMTOP` files. This is the only
+method by which it can be guaranteed that donor-hydrogen pairs are correctly
+identified. However, if, for example, a `PDB` file is used instead, a
+:attr:`donors_sel` may be provided along with a :attr:`hydrogens_sel` and the
+donor-hydrogen pairs will be identified via a distance cutoff,
+:attr:`d_h_cutoff`::
 
   import MDAnalysis
-  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
+  from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
+    HydrogenBondAnalysis as HBA)
 
   u = MDAnalysis.Universe(pdb, trajectory)
 

From 0698e22911326129121eea0be97066191739e614 Mon Sep 17 00:00:00 2001
From: Paul Smith <paul.smith@kcl.ac.uk>
Date: Mon, 6 Jul 2020 14:15:24 +0100
Subject: [PATCH 11/16] Fix docstring formatting for pep8

---
 package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index 1eff55e3af3..ba5da70b06f 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -268,9 +268,9 @@ def __init__(self, universe,
             and acceptor selections include both protein and water, it is
             possible to find only protein-water hydrogen bonds - and not
             protein-protein or water-water - by specifying
-            :attr`between=["protein", "SOL"]`. If a two-dimensional list is
+            between=["protein", "SOL"]`. If a two-dimensional list is
             passed, hydrogen bonds between each pair will be found. For
-            example, :attr`between=[["protein", "SOL"], ["protein", "protein"]]`
+            example, between=[["protein", "SOL"], ["protein", "protein"]]`
             will calculate all protein-water and protein-protein hydrogen
             bonds but not water-water hydrogen bonds. If `None`, hydrogen
             bonds between all donors and acceptors will be calculated.

From d16d227bfe9068822249adb4419c9a29a84db1ab Mon Sep 17 00:00:00 2001
From: Paul Smith <paul.smith@kcl.ac.uk>
Date: Mon, 6 Jul 2020 19:24:35 +0100
Subject: [PATCH 12/16] Fix tests that I broke

---
 .../MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py   | 8 ++++----
 .../analysis/test_hydrogenbonds_analysis.py               | 8 ++++----
 2 files changed, 8 insertions(+), 8 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index ba5da70b06f..005dc963c0f 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -686,19 +686,19 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
                 "Autocorrelation: Please use the .run() before calling this"
                 "function"
             )
-            raise NoDataError("No .hbonds: use the .run() first")
+            raise NoDataError("No .hbonds attribute: use .run() first")
 
         if self.step != 1:
             logging.warning(
-                "Autocorrelation: Hydrogen bonds were computedwith step > 1."
+                "Autocorrelation: Hydrogen bonds were computed with step > 1."
             )
             logging.warning(
                 "Autocorrelation: We recommend recomputing hydrogen bonds with"
-                "step = 1."
+                " step = 1."
             )
             logging.warning(
                 "Autocorrelation: if you would like to allow bonds to break"
-                "and reform, please use 'intermittency'"
+                " and reform, please use 'intermittency'"
             )
 
         # Extract the hydrogen bonds IDs only in the format
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index d0d565ef2fb..eb3ef05263b 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -204,7 +204,7 @@ def test_first(self, hydrogen_bonds):
         assert_equal(hydrogen_index, 2)
         assert_equal(acceptor_index, 3)
         assert_almost_equal(da_dst, 2.5)
-        assert_almost_equal(dha_angle, 180)
+        assert_almost_equal(angle, 180)
 
     def test_count_by_time(self, hydrogen_bonds):
         ref_times = np.array([0, 1, 2]) # u.trajectory.dt is 1
@@ -238,9 +238,9 @@ def test_logging_step_not_1(self, universe, caplog):
         caplog.set_level(logging.WARNING)
         hbonds.lifetime(tau_max=2, intermittency=1)
 
-        message = "ideally autocorrelation function would be carried out on" \
-            "consecutive frames"
-        assert any(message in rec.getMessage() for rec in caplog.records)
+        warning = \
+            "Autocorrelation: Hydrogen bonds were computed with step > 1."
+        assert any(warning in rec.getMessage() for rec in caplog.records)
 
 
 class TestHydrogenBondAnalysisBetween(object):

From 33f260c9788937782654d3ca33488e04bb96fe24 Mon Sep 17 00:00:00 2001
From: Paul Smith <paul.smith@kcl.ac.uk>
Date: Mon, 6 Jul 2020 22:44:36 +0100
Subject: [PATCH 13/16] Fixed example in docs

Sphinx docs were previously not building
---
 package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index 005dc963c0f..05d9fbd0a01 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -144,7 +144,7 @@
 To calculate the hydrogen bonds between different groups, for example a
 protein and water, one can use the :attr:`between` keyword. The
 following will find protein-water hydrogen bonds but not protein-protein
-or water-water hydrogen bonds:
+or water-water hydrogen bonds::
 
   import MDAnalysis
   from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (

From 9b087a464a96da7afa25b988b77932390798b3e9 Mon Sep 17 00:00:00 2001
From: Mateusz Bieniek <bieniekmateusz@users.noreply.github.com>
Date: Tue, 7 Jul 2020 13:58:24 +0100
Subject: [PATCH 14/16] Update package/MDAnalysis/analysis/waterdynamics.py

Co-authored-by: Oliver Beckstein <orbeckst@gmail.com>
---
 package/MDAnalysis/analysis/waterdynamics.py | 5 +++--
 1 file changed, 3 insertions(+), 2 deletions(-)

diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index b5a372069f5..5ef1820c581 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -508,8 +508,9 @@ class HydrogenBondLifetimes(object):
     .. versionchanged:: 1.0.0
        The ``nproc`` keyword was removed as it linked to a portion of code that
        may have failed in some cases.
-    .. deprecated:: 2.0.0
-       Instead, please use
+    .. deprecated:: 1.0.1
+       ``waterdynamics.HydrogenBondLifetimes`` is deprecated and will be removed in 2.0.0.
+       Instead, please use (available in 2.0.0)
        MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime
 
     """

From 904e75dafa8ad9fc3d18b7c1f1cfc46aa3746da2 Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Tue, 7 Jul 2020 14:13:02 +0100
Subject: [PATCH 15/16] Addressing feedback in PR #2791

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py | 5 +++--
 .../MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py | 6 +++---
 2 files changed, 6 insertions(+), 5 deletions(-)

diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
index 05d9fbd0a01..0ce9062ded2 100644
--- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
+++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py
@@ -327,6 +327,7 @@ def __init__(self, universe,
         self.d_a_cutoff = d_a_cutoff
         self.d_h_a_angle = d_h_a_angle_cutoff
         self.update_selections = update_selections
+        self.hbonds = None
 
     def guess_hydrogens(self,
                         select='all',
@@ -677,7 +678,7 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
             autcorrelation value for each value of `tau`
         """
 
-        if not hasattr(self, 'hbonds'):
+        if self.hbonds is None:
             logging.error(
                 "Autocorrelation analysis of hydrogen bonds cannot be done"
                 "before the hydrogen bonds are found"
@@ -686,7 +687,7 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0):
                 "Autocorrelation: Please use the .run() before calling this"
                 "function"
             )
-            raise NoDataError("No .hbonds attribute: use .run() first")
+            raise NoDataError(".hbonds attribute is None: use .run() first")
 
         if self.step != 1:
             logging.warning(
diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
index eb3ef05263b..d107f9932fe 100644
--- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
+++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py
@@ -104,7 +104,7 @@ def test_count_by_ids(self, h, universe):
         assert_allclose(counts, ref_counts)
 
 
-class TestHydrogenBondAnalysisMock(object):
+class TestHydrogenBondAnalysisIdeal(object):
 
     kwargs = {
         'donors_sel': 'name O',
@@ -174,7 +174,7 @@ def universe():
     def hydrogen_bonds(universe):
         h = HydrogenBondAnalysis(
             universe,
-            **TestHydrogenBondAnalysisMock.kwargs
+            **TestHydrogenBondAnalysisIdeal.kwargs
         )
         h.run()
         return h
@@ -227,7 +227,7 @@ def test_hydrogen_bond_lifetime_intermittency(self, hydrogen_bonds):
     def test_no_attr_hbonds(self, universe):
         hbonds = HydrogenBondAnalysis(universe, **self.kwargs)
         # hydrogen bonds are not computed
-        with pytest.raises(NoDataError, match="No .hbonds"):
+        with pytest.raises(NoDataError, match=".hbonds attribute is None"):
             hbonds.lifetime(tau_max=2, intermittency=1)
 
     def test_logging_step_not_1(self, universe, caplog):

From 93b67078641469bbd245cb5574796f33556fd74e Mon Sep 17 00:00:00 2001
From: bieniekmat <bieniekmat@gmail.com>
Date: Tue, 7 Jul 2020 14:14:59 +0100
Subject: [PATCH 16/16] Addressing feedback in PR #2791

co-authored-by: p-j-smith <paul.smith@kcl.ac.uk>
---
 package/MDAnalysis/analysis/waterdynamics.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/package/MDAnalysis/analysis/waterdynamics.py b/package/MDAnalysis/analysis/waterdynamics.py
index 5ef1820c581..fd26da5f021 100644
--- a/package/MDAnalysis/analysis/waterdynamics.py
+++ b/package/MDAnalysis/analysis/waterdynamics.py
@@ -509,8 +509,8 @@ class HydrogenBondLifetimes(object):
        The ``nproc`` keyword was removed as it linked to a portion of code that
        may have failed in some cases.
     .. deprecated:: 1.0.1
-       ``waterdynamics.HydrogenBondLifetimes`` is deprecated and will be removed in 2.0.0.
-       Instead, please use (available in 2.0.0)
+       ``waterdynamics.HydrogenBondLifetimes`` is deprecated and will be
+       removed in 2.0.0. Instead, please use (available in 2.0.0)
        MDAnalysis.analysis.hydrogenbonds.HydrogenBondAnalysis.lifetime
 
     """