Skip to content

Commit

Permalink
Accept AtomGroup alongwith Universe in the DistanceMatrix (#3541)
Browse files Browse the repository at this point in the history
* Accept AtomGroup instead of Universe in the DistanceMatrix

* Modified Tests

* Added a test to assert same values of distmatrix

* Accept both AG and Universe

* Resolve PEP8

* Final changes

* updated AUTHOR and CHANGELOG

* corrected and updated docstrings

* modified and displaced import statement

* updated changelog and added self to authors

* updated changelog and authors

* updated docstrings

* fixed typo and edited docstrings

* fixed merge issues with authors and changelog

* update docstring

* remove github id from previous release

* tests on an updating AG

* upstream changes merged

* revoked the redundant testcase and user warning if updating atomgroup is parsed

* minor docstring edit

* warning test for updating atomgroup

* modified warning message

* added test where ag and select string is parsed
  • Loading branch information
AnirG authored Apr 11, 2022
1 parent 8eda52f commit 0eb42ca
Show file tree
Hide file tree
Showing 4 changed files with 50 additions and 11 deletions.
1 change: 1 addition & 0 deletions package/AUTHORS
Original file line number Diff line number Diff line change
Expand Up @@ -184,6 +184,7 @@ Chronological list of authors
- Tamandeep Singh
- Mingyi Xue
- Meghan Osato
- Anirvinya G

External code
-------------
Expand Down
4 changes: 3 additions & 1 deletion package/CHANGELOG
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@ The rules for this file:

------------------------------------------------------------------------------
??/??/?? IAlibay, BFedder, inomag, Agorfa, aya9aladdin, shudipto-amin, cbouy,
HenokB, umak1106, tamandeeps, Mrqeoqqt, megosato
HenokB, umak1106, tamandeeps, Mrqeoqqt, megosato, AnirG

* 2.2.0

Expand All @@ -36,6 +36,8 @@ Fixes
* Fixed FrameIteratorIndices doesn't rewind. (Issue #3416)

Enhancements
* `DistanceMatrix` class in MDAnalysis.analysis.diffusionmap now also
accepts `AtomGroup` as a valid argument (Issue #3486, PR #3541)
* Improves the RDKitConverter's accuracy (PR #3044): AtomGroups containing
monatomic ion charges or edge cases with nitrogen, sulfur, phosphorus and
conjugated systems should now have correctly assigned bond orders and
Expand Down
27 changes: 20 additions & 7 deletions package/MDAnalysis/analysis/diffusionmap.py
Original file line number Diff line number Diff line change
Expand Up @@ -156,6 +156,7 @@
import numpy as np

from MDAnalysis.core.universe import Universe
from MDAnalysis.core.groups import AtomGroup, UpdatingAtomGroup
from .rms import rmsd
from .base import AnalysisBase

Expand All @@ -171,7 +172,7 @@ class DistanceMatrix(AnalysisBase):
Parameters
----------
universe : `~MDAnalysis.core.universe.Universe`
universe : `~MDAnalysis.core.universe.Universe` or `~MDAnalysis.core.groups.AtomGroup`
The MD Trajectory for dimension reduction, remember that
computational cost of eigenvalue decomposition
scales at O(N^3) where N is the number of frames.
Expand Down Expand Up @@ -243,12 +244,20 @@ class DistanceMatrix(AnalysisBase):
.. versionchanged:: 2.0.0
:attr:`dist_matrix` is now stored in a
:class:`MDAnalysis.analysis.base.Results` instance.
.. versionchanged:: 2.2.0
:class:`DistanceMatrix` now also accepts `AtomGroup`.
"""
def __init__(self, universe, select='all', metric=rmsd, cutoff=1E0-5,
weights=None, **kwargs):
# remember that this must be called before referencing self.n_frames
super(DistanceMatrix, self).__init__(universe.trajectory, **kwargs)
super(DistanceMatrix, self).__init__(universe.universe.trajectory,
**kwargs)

if isinstance(universe, UpdatingAtomGroup):
wmsg = ("U must be a static AtomGroup. Parsing an updating AtomGroup "
"will result in a static AtomGroup with the current frame "
"atom selection.")
warnings.warn(wmsg)

self.atoms = universe.select_atoms(select)
self._metric = metric
Expand Down Expand Up @@ -309,14 +318,18 @@ class DiffusionMap(object):
transform(n_eigenvectors, time)
Perform an embedding of a frame into the eigenvectors representing
the collective coordinates.
.. versionchanged:: 2.2.0
:class:`DiffusionMap` now also accepts `AtomGroup`.
"""

def __init__(self, u, epsilon=1, **kwargs):
"""
Parameters
-------------
u : MDAnalysis Universe or DistanceMatrix object
Can be a Universe, in which case one must supply kwargs for the
u : MDAnalysis Universe or AtomGroup or DistanceMatrix object.
Can be a Universe or AtomGroup, in which case one must supply kwargs for the
initialization of a DistanceMatrix. Otherwise, this can be a
DistanceMatrix already initialized. Either way, this will be made
into a diffusion kernel.
Expand All @@ -328,12 +341,12 @@ def __init__(self, u, epsilon=1, **kwargs):
Parameters to be passed for the initialization of a
:class:`DistanceMatrix`.
"""
if isinstance(u, Universe):
if isinstance(u, AtomGroup) or isinstance(u, Universe):
self._dist_matrix = DistanceMatrix(u, **kwargs)
elif isinstance(u, DistanceMatrix):
self._dist_matrix = u
else:
raise ValueError("U is not a Universe or DistanceMatrix and"
raise ValueError("U is not a Universe or AtomGroup or DistanceMatrix and"
" so the DiffusionMap has no data to work with.")
self._epsilon = epsilon

Expand Down
29 changes: 26 additions & 3 deletions testsuite/MDAnalysisTests/analysis/test_diffusionmap.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@
import numpy as np
import pytest
from MDAnalysisTests.datafiles import PDB, XTC
from numpy.testing import assert_array_almost_equal
from numpy.testing import assert_array_almost_equal, assert_allclose


@pytest.fixture(scope='module')
Expand Down Expand Up @@ -69,6 +69,22 @@ def test_dist_weights(u):
[.707, -.707, 0, 0]]), 2)


def test_distvalues_ag_universe(u):
dist_universe = diffusionmap.DistanceMatrix(u, select='backbone').run()
ag = u.select_atoms('backbone')
dist_ag = diffusionmap.DistanceMatrix(ag).run()
assert_allclose(dist_universe.results.dist_matrix,
dist_ag.results.dist_matrix)


def test_distvalues_ag_select(u):
dist_universe = diffusionmap.DistanceMatrix(u, select='backbone').run()
ag = u.select_atoms('protein')
dist_ag = diffusionmap.DistanceMatrix(ag, select='backbone').run()
assert_allclose(dist_universe.results.dist_matrix,
dist_ag.results.dist_matrix)


def test_different_steps(u):
dmap = diffusionmap.DiffusionMap(u, select='backbone')
dmap.run(step=3)
Expand All @@ -92,9 +108,16 @@ def test_long_traj(u):
dmap.run()


def test_not_universe_error(u):
def test_updating_atomgroup(u):
with pytest.warns(UserWarning, match='U must be a static AtomGroup'):
resid_select = 'around 5 resname ALA'
ag = u.select_atoms(resid_select, updating=True)
dmap = diffusionmap.DiffusionMap(ag)
dmap.run()

def test_not_universe_atomgroup_error(u):
trj_only = u.trajectory
with pytest.raises(ValueError, match='U is not a Universe'):
with pytest.raises(ValueError, match='U is not a Universe or AtomGroup'):
diffusionmap.DiffusionMap(trj_only)


Expand Down

0 comments on commit 0eb42ca

Please sign in to comment.