Skip to content

Move results attributes of helix_analysis to Results #3279

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
orbeckst merged 2 commits into from
May 7, 2021
Merged

Move results attributes of helix_analysis to Results #3279

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
4e3243d
add results
lilyminium May 6, 2021
d63d389
Apply suggestions from code review
lilyminium May 7, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions package/CHANGELOG
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down
58 changes: 30 additions & 28 deletions package/MDAnalysis/analysis/helix_analysis.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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::
Expand Down Expand Up @@ -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.
"""
Expand Down Expand Up @@ -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):
"""
Expand All @@ -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')

Expand All @@ -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
36 changes: 18 additions & 18 deletions testsuite/MDAnalysisTests/analysis/test_helix_analysis.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand All @@ -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,
Expand All @@ -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()
Expand Down Expand Up @@ -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)


Expand Down