diff --git a/package/CHANGELOG b/package/CHANGELOG index 35569ec940e..c009b86c35a 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -174,6 +174,8 @@ Enhancements checking if it can be used in parallel analysis. (Issue #2996, PR #2950) Changes + * `analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalyis` now stores + `hbonds` data using the `analysis.base.Results class (Issues #3271, #3261) * `analysis.rms.RMSD` and `analysis.rms.RMSF` classes now store `rmsd` and `rmsf` data using the `analysis.base.Results` class (Issues #3274 #3261) * `analysis.dihedrals` classes now store angle data using the @@ -235,6 +237,9 @@ Changes * Added OpenMM coordinate and topology converters (Issue #2863, PR #2917) Deprecations + * The `hbonds` attribute of + `hydrogenbonds.hbond_analysis.HydrogenBondAnalysis.hbonds` is now + deprecated in favour of `results.hbonds` (Issues #3271, #3261) * The `analysis.rms.RMSD.rmsd` and `analysis.rms.RMSF.rmsf` attributes are now deprecated in favour of `analysis.rms.RMSD.results.rmsd` and `analysis.rms.RMSF.results.rmsf` respectively. They will be removed in diff --git a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py index 020871d99f0..17d907aa309 100644 --- a/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py +++ b/package/MDAnalysis/analysis/hydrogenbonds/hbond_analysis.py @@ -67,7 +67,7 @@ - *angle* (degrees): angle of the hydrogen bond Hydrogen bond data are returned in a :class:`numpy.ndarray` on a "one line, one observation" basis -and can be accessed via :attr:`HydrogenBondAnalysis.hbonds`:: +and can be accessed via :attr:`HydrogenBondAnalysis.results.hbonds`:: results = [ [ @@ -209,13 +209,29 @@ .. autoclass:: HydrogenBondAnalysis :members: + + .. attribute:: results.hbonds + + A :class:`numpy.ndarray` which contains a list of all observed hydrogen + bond interactions. See `Output`_ for more information. + + .. versionadded:: 2.0.0 + + .. attribute:: hbonds + + Alias to the :attr:`results.hbonds` attribute. + + .. deprecated:: 2.0.0 + Will be removed in MDAnalysis 3.0.0. Please use + :attr:`results.hbonds` instead. """ import logging +import warnings from collections.abc import Iterable import numpy as np -from .. import base +from ..base import AnalysisBase, Results from MDAnalysis.lib.distances import capped_distance, calc_angles from MDAnalysis.lib.correlations import autocorrelation, correct_intermittency from MDAnalysis.exceptions import NoDataError @@ -231,7 +247,7 @@ del Doi -class HydrogenBondAnalysis(base.AnalysisBase): +class HydrogenBondAnalysis(AnalysisBase): """ Perform an analysis of hydrogen bonds in a Universe. """ @@ -330,7 +346,8 @@ 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 + self.results = Results() + self.results.hbonds = None def guess_hydrogens(self, select='all', @@ -563,7 +580,7 @@ def _filter_atoms(self, donors, hydrogens, acceptors): def _prepare(self): - self.hbonds = [[], [], [], [], [], []] + self.results.hbonds = [[], [], [], [], [], []] # Set atom selections if they have not been provided if not self.acceptors_sel: @@ -625,16 +642,25 @@ def _single_frame(self): hbond_angles = d_h_a_angles[hbond_indices] # Store data on hydrogen bonds found at this frame - self.hbonds[0].extend(np.full_like(hbond_donors, self._ts.frame)) - self.hbonds[1].extend(hbond_donors.indices) - self.hbonds[2].extend(hbond_hydrogens.indices) - self.hbonds[3].extend(hbond_acceptors.indices) - self.hbonds[4].extend(hbond_distances) - self.hbonds[5].extend(hbond_angles) + self.results.hbonds[0].extend(np.full_like(hbond_donors, + self._ts.frame)) + self.results.hbonds[1].extend(hbond_donors.indices) + self.results.hbonds[2].extend(hbond_hydrogens.indices) + self.results.hbonds[3].extend(hbond_acceptors.indices) + self.results.hbonds[4].extend(hbond_distances) + self.results.hbonds[5].extend(hbond_angles) def _conclude(self): - self.hbonds = np.asarray(self.hbonds).T + self.results.hbonds = np.asarray(self.results.hbonds).T + + @property + def hbonds(self): + wmsg = ("The `hbonds` attribute was deprecated in MDAnalysis 2.0.0 " + "and will be removed in MDAnalysis 3.0.0. Please use " + "`results.hbonds` instead.") + warnings.warn(wmsg, DeprecationWarning) + return self.results.hbonds def lifetime(self, tau_max=20, window_step=1, intermittency=0): """Computes and returns the time-autocorrelation @@ -681,7 +707,7 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0): autcorrelation value for each value of `tau` """ - if self.hbonds is None: + if self.results.hbonds is None: logging.error( "Autocorrelation analysis of hydrogen bonds cannot be done" "before the hydrogen bonds are found" @@ -709,7 +735,7 @@ def lifetime(self, tau_max=20, window_step=1, intermittency=0): # [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]: + for hbond in self.results.hbonds[self.results.hbonds[:, 0] == frame]: found_hydrogen_bonds[frame_index].add(frozenset(hbond[2:4])) intermittent_hbonds = correct_intermittency( @@ -735,7 +761,7 @@ def count_by_time(self): the number of hydrogen bonds over time. """ - indices, tmp_counts = np.unique(self.hbonds[:, 0], axis=0, + indices, tmp_counts = np.unique(self.results.hbonds[:, 0], axis=0, return_counts=True) indices -= self.start @@ -760,8 +786,8 @@ def count_by_type(self): acceptor atoms in a hydrogen bond. """ - d = self.u.atoms[self.hbonds[:, 1].astype(int)] - a = self.u.atoms[self.hbonds[:, 3].astype(int)] + d = self.u.atoms[self.results.hbonds[:, 1].astype(int)] + a = self.u.atoms[self.results.hbonds[:, 3].astype(int)] tmp_hbonds = np.array([d.resnames, d.types, a.resnames, a.types], dtype=str).T @@ -789,9 +815,9 @@ def count_by_ids(self): in a hydrogen bond. """ - d = self.u.atoms[self.hbonds[:, 1].astype(int)] - h = self.u.atoms[self.hbonds[:, 2].astype(int)] - a = self.u.atoms[self.hbonds[:, 3].astype(int)] + d = self.u.atoms[self.results.hbonds[:, 1].astype(int)] + h = self.u.atoms[self.results.hbonds[:, 2].astype(int)] + a = self.u.atoms[self.results.hbonds[:, 3].astype(int)] tmp_hbonds = np.array([d.ids, h.ids, a.ids]).T hbond_ids, ids_counts = np.unique(tmp_hbonds, axis=0, diff --git a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py index d90a2ce0488..0d27429bcd7 100644 --- a/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py +++ b/testsuite/MDAnalysisTests/analysis/test_hydrogenbonds_analysis.py @@ -59,18 +59,22 @@ def h(self, universe): def test_hbond_analysis(self, h): - assert len(np.unique(h.hbonds[:, 0])) == 10 - assert len(h.hbonds) == 32 + assert len(np.unique(h.results.hbonds[:, 0])) == 10 + assert len(h.results.hbonds) == 32 reference = { 'distance': {'mean': 2.7627309, 'std': 0.0905052}, 'angle': {'mean': 158.9038039, 'std': 12.0362826}, } - assert_allclose(np.mean(h.hbonds[:, 4]), reference['distance']['mean']) - assert_allclose(np.std(h.hbonds[:, 4]), reference['distance']['std']) - assert_allclose(np.mean(h.hbonds[:, 5]), reference['angle']['mean']) - assert_allclose(np.std(h.hbonds[:, 5]), reference['angle']['std']) + assert_allclose(np.mean(h.results.hbonds[:, 4]), + reference['distance']['mean']) + assert_allclose(np.std(h.results.hbonds[:, 4]), + reference['distance']['std']) + assert_allclose(np.mean(h.results.hbonds[:, 5]), + reference['angle']['mean']) + assert_allclose(np.std(h.results.hbonds[:, 5]), + reference['angle']['std']) def test_count_by_time(self, h): @@ -103,6 +107,11 @@ def test_count_by_ids(self, h, universe): assert_allclose(counts, ref_counts) + def test_hbonds_deprecated_attr(self, h): + wmsg = "The `hbonds` attribute was deprecated in MDAnalysis 2.0.0" + with pytest.warns(DeprecationWarning, match=wmsg): + assert_equal(h.hbonds, h.results.hbonds) + class TestHydrogenBondAnalysisIdeal(object): @@ -196,10 +205,10 @@ def test_no_bond_info_exception(self, universe): h._get_dh_pairs() def test_first_hbond(self, hydrogen_bonds): - assert len(hydrogen_bonds.hbonds) == 2 + assert len(hydrogen_bonds.results.hbonds) == 2 frame_no, donor_index, hydrogen_index, acceptor_index, da_dst, angle =\ - hydrogen_bonds.hbonds[0] + hydrogen_bonds.results.hbonds[0] assert_equal(donor_index, 0) assert_equal(hydrogen_index, 2) assert_equal(acceptor_index, 3) @@ -333,7 +342,8 @@ def test_between_all(self, universe): [3, 4, 6], # protein-water [6, 7, 8] # protein-protein ] - assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices) + assert_array_equal(hbonds.results.hbonds[:, 1:4], + expected_hbond_indices) def test_between_PW(self, universe): # Find only protein-water hydrogen bonds @@ -348,7 +358,8 @@ def test_between_PW(self, universe): expected_hbond_indices = [ [3, 4, 6] # protein-water ] - assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices) + assert_array_equal(hbonds.results.hbonds[:, 1:4], + expected_hbond_indices) def test_between_PW_PP(self, universe): # Find protein-water and protein-protein hydrogen bonds (not @@ -368,7 +379,8 @@ def test_between_PW_PP(self, universe): [3, 4, 6], # protein-water [6, 7, 8] # protein-protein ] - assert_array_equal(hbonds.hbonds[:, 1:4], expected_hbond_indices) + assert_array_equal(hbonds.results.hbonds[:, 1:4], + expected_hbond_indices) class TestHydrogenBondAnalysisTIP3P_GuessAcceptors_GuessHydrogens_UseTopology_(TestHydrogenBondAnalysisTIP3P): @@ -396,7 +408,7 @@ def test_no_hydrogens(self, universe): assert h._hydrogens.n_atoms == 0 assert h._donors.n_atoms == 0 - assert h.hbonds.size == 0 + assert h.results.hbonds.size == 0 class TestHydrogenBondAnalysisTIP3P_GuessDonors_NoTopology(object): """Guess the donor atoms involved in hydrogen bonds using the partial charges of the atoms. @@ -507,18 +519,22 @@ def h(self, universe): def test_hbond_analysis(self, h): - assert len(np.unique(h.hbonds[:, 0])) == 5 - assert len(h.hbonds) == 15 + assert len(np.unique(h.results.hbonds[:, 0])) == 5 + assert len(h.results.hbonds) == 15 reference = { 'distance': {'mean': 2.73942464, 'std': 0.05867924}, 'angle': {'mean': 157.07768079, 'std': 9.72636682}, } - assert_allclose(np.mean(h.hbonds[:, 4]), reference['distance']['mean']) - assert_allclose(np.std(h.hbonds[:, 4]), reference['distance']['std']) - assert_allclose(np.mean(h.hbonds[:, 5]), reference['angle']['mean']) - assert_allclose(np.std(h.hbonds[:, 5]), reference['angle']['std']) + assert_allclose(np.mean(h.results.hbonds[:, 4]), + reference['distance']['mean']) + assert_allclose(np.std(h.results.hbonds[:, 4]), + reference['distance']['std']) + assert_allclose(np.mean(h.results.hbonds[:, 5]), + reference['angle']['mean']) + assert_allclose(np.std(h.results.hbonds[:, 5]), + reference['angle']['std']) def test_count_by_time(self, h):