diff --git a/package/CHANGELOG b/package/CHANGELOG index abaf73f9f2a..80b613196df 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -175,6 +175,8 @@ Enhancements checking if it can be used in parallel analysis. (Issue #2996, PR #2950) Changes + * `helix_analysis.HELANAL` now uses the `analysis.base.Results` class to + store results attributes (Issue #3261, #3267) * `analysis.align.AlignTraj` and `analysis.align.AverageStructure` now store their result attributes using the `analysis.base.Results` class (Issues #3278 and #3261) diff --git a/package/MDAnalysis/analysis/helix_analysis.py b/package/MDAnalysis/analysis/helix_analysis.py index 5e8c516ab7e..7518ee6c376 100644 --- a/package/MDAnalysis/analysis/helix_analysis.py +++ b/package/MDAnalysis/analysis/helix_analysis.py @@ -29,7 +29,7 @@ :Year: 2020 :Copyright: GNU Public License v3 -.. versionadded:: 1.0.0 +.. versionadded:: 2.0.0 This module contains code to analyse protein helices using the HELANAL_ algorithm @@ -63,7 +63,9 @@ helanal = hel.HELANAL(u, select='name CA and resnum 161-187') helanal.run() - print(helanal.summary) +All computed properties are available in ``.results``:: + + print(helanal.results.summary) Alternatively, you can analyse several helices at once by passing in multiple selection strings:: @@ -309,42 +311,42 @@ class HELANAL(AnalysisBase): Attributes ---------- - local_twists : array or list of arrays + results.local_twists : array or list of arrays The local twist angle from atom i+1 to i+2. Each array has shape (n_frames, n_residues-3) - local_nres_per_turn : array or list of arrays + results.local_nres_per_turn : array or list of arrays Number of residues per turn, based on local_twist. Each array has shape (n_frames, n_residues-3) - local_axes : array or list of arrays + results.local_axes : array or list of arrays The length-wise helix axis of the local window. Each array has shape (n_frames, n_residues-3, 3) - local_heights : array or list of arrays + results.local_heights : array or list of arrays The rise of each local helix. Each array has shape (n_frames, n_residues-3) - local_helix_directions : array or list of arrays + results.local_helix_directions : array or list of arrays The unit vector from each local origin to atom i+1. Each array has shape (n_frames, n_residues-2, 3) - local_origins :array or list of arrays + results.local_origins :array or list of arrays The projected origin for each helix. Each array has shape (n_frames, n_residues-2, 3) - local_screw_angles : array or list of arrays + results.local_screw_angles : array or list of arrays The local screw angle for each helix. Each array has shape (n_frames, n_residues-2) - local_bends : array or list of arrays + results.local_bends : array or list of arrays The angles between local helix axes, 3 windows apart. Each array has shape (n_frames, n_residues-6) - all_bends : array or list of arrays + results.all_bends : array or list of arrays The angles between local helix axes. Each array has shape (n_frames, n_residues-3, n_residues-3) - global_axis : array or list of arrays + results.global_axis : array or list of arrays The length-wise axis for the overall helix. This points at the first helix window in the helix, so it runs opposite to the direction of the residue numbers. Each array has shape (n_frames, 3) - global_tilts : array or list of arrays + results.global_tilts : array or list of arrays The angle between the global axis and the reference axis. Each array has shape (n_frames,) - summary : dict or list of dicts + results.summary : dict or list of dicts Summary of stats for each property: the mean, the sample standard deviation, and the mean absolute deviation. """ @@ -418,49 +420,49 @@ def _prepare(self): for key, dims in self.attr_shapes.items(): empty = [self._zeros_per_frame( dims, n_positions=n) for n in n_res] - setattr(self, key, empty) + self.results[key] = empty - self.global_axis = [self._zeros_per_frame((3,)) for n in n_res] - self.all_bends = [self._zeros_per_frame((n-3, n-3)) for n in n_res] + self.results.global_axis = [self._zeros_per_frame((3,)) for n in n_res] + self.results.all_bends = [self._zeros_per_frame((n-3, n-3)) for n in n_res] def _single_frame(self): _f = self._frame_index for i, ag in enumerate(self.atomgroups): results = helix_analysis(ag.positions, ref_axis=self.ref_axis) for key, value in results.items(): - attr = getattr(self, key) + attr = self.results[key] attr[i][_f] = value def _conclude(self): # compute tilt of global axes - self.global_tilts = [] + self.results.global_tilts = tilts = [] norm_ref = (self.ref_axis**2).sum() ** 0.5 - for axes in self.global_axis: + for axes in self.results.global_axis: cos = np.matmul(self.ref_axis, axes.T) / \ (mdamath.pnorm(axes)*norm_ref) cos = np.clip(cos, -1.0, 1.0) - self.global_tilts.append(np.rad2deg(np.arccos(cos))) + tilts.append(np.rad2deg(np.arccos(cos))) global_attrs = ['global_axis', 'global_tilts', 'all_bends'] attrnames = list(self.attr_shapes.keys()) + global_attrs # summarise - self.summary = [] + self.results.summary = [] for i in range(len(self.atomgroups)): stats = {} for name in attrnames: - attr = getattr(self, name) + attr = self.results[name] mean = attr[i].mean(axis=0) dev = np.abs(attr[i]-mean) stats[name] = {'mean': mean, 'sample_sd': attr[i].std(axis=0, ddof=1), 'abs_dev': dev.mean(axis=0)} - self.summary.append(stats) + self.results.summary.append(stats) # flatten? if len(self.atomgroups) == 1 and self._flatten: for name in attrnames + ['summary']: - attr = getattr(self, name) - setattr(self, name, attr[0]) + attr = self.results[name] + self.results[name] = attr[0] def universe_from_origins(self): """ @@ -471,7 +473,7 @@ def universe_from_origins(self): Universe or list of Universes """ try: - origins = self.local_origins + origins = self.results.local_origins except AttributeError: raise ValueError('Call run() before universe_from_origins') @@ -485,6 +487,6 @@ def universe_from_origins(self): atom_resindex=np.arange(n_res), trajectory=True).load_new(xyz) universe.append(u) - if not isinstance(self.local_origins, list): + if not isinstance(self.results.local_origins, list): universe = universe[0] return universe diff --git a/testsuite/MDAnalysisTests/analysis/test_helix_analysis.py b/testsuite/MDAnalysisTests/analysis/test_helix_analysis.py index f96691d5292..1e39147be81 100644 --- a/testsuite/MDAnalysisTests/analysis/test_helix_analysis.py +++ b/testsuite/MDAnalysisTests/analysis/test_helix_analysis.py @@ -351,7 +351,7 @@ def helanal(self, psf_ca): return ha.run(start=10, stop=80) def test_regression_summary(self, helanal): - bends = helanal.summary['all_bends'] + bends = helanal.results.summary['all_bends'] old_helanal = read_bending_matrix(HELANAL_BENDING_MATRIX_SUBSET) assert_almost_equal(np.triu(bends['mean'], 1), old_helanal['Mean'], decimal=1) @@ -368,11 +368,11 @@ def test_regression_values(self): for key, value in HELANAL_SINGLE_DATA.items(): if 'summary' in key: - data = getattr(ha, key.split()[0]) + data = ha.results[key.split()[0]] calculated = [data.mean(), data.std(ddof=1), np.fabs(data-data.mean()).mean()] else: - calculated = getattr(ha, key)[0] + calculated = ha.results[key][0] msg = 'Mismatch between calculated and reference {}' assert_almost_equal(calculated, value, @@ -385,10 +385,10 @@ def test_multiple_selections(self, psf_ca): ha.run() n_frames = len(psf_ca.universe.trajectory) assert len(ha.atomgroups) == 2 - assert len(ha.summary) == 2 - assert len(ha.all_bends) == 2 - assert ha.all_bends[0].shape == (n_frames, 8, 8) - assert ha.all_bends[1].shape == (n_frames, 18, 18) + assert len(ha.results.summary) == 2 + assert len(ha.results.all_bends) == 2 + assert ha.results.all_bends[0].shape == (n_frames, 8, 8) + assert ha.results.all_bends[1].shape == (n_frames, 18, 18) def test_universe_from_origins(self, helanal): u = helanal.universe_from_origins() @@ -464,22 +464,22 @@ def test_len_groups_short(self, psf_ca): def test_helanal_zigzag(self, zigzag, ref_axis, screw_angles): ha = hel.HELANAL(zigzag, select="all", ref_axis=ref_axis, flatten_single_helix=True).run() - assert_almost_equal(ha.local_twists, 180, decimal=4) - assert_almost_equal(ha.local_nres_per_turn, 2, decimal=4) - assert_almost_equal(ha.global_axis, [[0, 0, -1]], decimal=4) + assert_almost_equal(ha.results.local_twists, 180, decimal=4) + assert_almost_equal(ha.results.local_nres_per_turn, 2, decimal=4) + assert_almost_equal(ha.results.global_axis, [[0, 0, -1]], decimal=4) # all 0 vectors - assert_almost_equal(ha.local_axes, 0, decimal=4) - assert_almost_equal(ha.local_bends, 0, decimal=4) - assert_almost_equal(ha.all_bends, 0, decimal=4) - assert_almost_equal(ha.local_heights, 0, decimal=4) - assert_almost_equal(ha.local_helix_directions[0][0::2], + assert_almost_equal(ha.results.local_axes, 0, decimal=4) + assert_almost_equal(ha.results.local_bends, 0, decimal=4) + assert_almost_equal(ha.results.all_bends, 0, decimal=4) + assert_almost_equal(ha.results.local_heights, 0, decimal=4) + assert_almost_equal(ha.results.local_helix_directions[0][0::2], [[-1, 0, 0]]*49, decimal=4) - assert_almost_equal(ha.local_helix_directions[0][1::2], + assert_almost_equal(ha.results.local_helix_directions[0][1::2], [[1, 0, 0]]*49, decimal=4) origins = zigzag.atoms.positions[1:-1].copy() origins[:, 0] = 0 - assert_almost_equal(ha.local_origins[0], origins, decimal=4) - assert_almost_equal(ha.local_screw_angles[0], + assert_almost_equal(ha.results.local_origins[0], origins, decimal=4) + assert_almost_equal(ha.results.local_screw_angles[0], screw_angles*49, decimal=4)