Skip to content

Commit

Permalink
density_from_PDB: deprecated, test, and bugfix
Browse files Browse the repository at this point in the history
- density_from_PDB is terrible code (I wrote it) and not used: deprecated and to be
  removed for 2.0.0
- added test so that it is covered
- fixed bug that was uncovered by test (us ag.tempfactors)
- added test for Bfactors2RMSF
- test for deprecation warnings
- organized in classes
- shortened test names for density_from_universe()
  • Loading branch information
orbeckst committed Feb 20, 2020
1 parent 9764455 commit 1154532
Show file tree
Hide file tree
Showing 3 changed files with 66 additions and 22 deletions.
3 changes: 2 additions & 1 deletion package/CHANGELOG
Original file line number Diff line number Diff line change
Expand Up @@ -129,8 +129,9 @@ Deprecations
* analysis.hbonds.HydrogenBondAnalysis is deprecated in 1.0 (remove in 2.0)
* analysis.density.density_from_Universe() (remove in 2.0)
* analysis.density.notwithin_coordinates_factory() (remove in 2.0)
* analysis.density.density_from_PDB and BfactorDensityCreator (remove in 2.0)



09/05/19 IAlibay, richardjgowers

* 0.20.1
Expand Down
18 changes: 14 additions & 4 deletions package/MDAnalysis/analysis/density.py
Original file line number Diff line number Diff line change
Expand Up @@ -1268,6 +1268,7 @@ def Bfactor2RMSF(B):
return np.sqrt(3. * B / 8.) / np.pi


@deprecate(release="1.0.0", remove="2.0.0")
def density_from_PDB(pdb, **kwargs):
"""Create a density from a single frame PDB.
Expand Down Expand Up @@ -1313,10 +1314,15 @@ def density_from_PDB(pdb, **kwargs):
--------
:func:`Bfactor2RMSF` and :class:`BfactorDensityCreator`
.. deprecated: 1.0.0
This function is not well tested or optimized and will be removed
in 2.0.0.
"""
return BfactorDensityCreator(pdb, **kwargs).Density()


# REMOVE in 2.0.0
class BfactorDensityCreator(object):
"""Create a density grid from a pdb file using MDAnalysis.
Expand All @@ -1337,6 +1343,10 @@ class BfactorDensityCreator(object):
.. * Using a temporary Creator class with the
.. :meth:`BfactorDensityCreator.Density` helper method is clumsy.
.. deprecated: 1.0.0
This function is not well tested or optimized and will be removed
in 2.0.0.
"""

@deprecate(release="1.0.0", remove="2.0.0")
Expand Down Expand Up @@ -1405,11 +1415,11 @@ def __init__(self, pdb, delta=1.0, select='resname HOH and name O',
if sigma is None:
# histogram individually, and smear out at the same time
# with the appropriate B-factor
if np.any(group.bfactors == 0.0):
if np.any(group.tempfactors == 0.0):
wmsg = "Some B-factors are Zero (will be skipped)."
logger.warning(wmsg)
warnings.warn(wmsg, category=MissingDataWarning)
rmsf = Bfactor2RMSF(group.bfactors)
rmsf = Bfactor2RMSF(group.tempfactors)
grid *= 0.0 # reset grid
self.g = self._smear_rmsf(coord, grid, self.edges, rmsf)
else:
Expand Down Expand Up @@ -1450,7 +1460,7 @@ def _smear_sigma(self, grid, sigma):
for iwat in range(len(pos[0])): # super-ugly loop
p = tuple([wp[iwat] for wp in pos])
g += grid[p] * np.fromfunction(self._gaussian, grid.shape, dtype=np.int, p=p, sigma=sigma)
print("Smearing out atom position {0:4d}/{1:5d} with RMSF {2:4.2f} A\r".format(iwat + 1, len(pos[0]), sigma),)
# print("Smearing out atom position {0:4d}/{1:5d} with RMSF {2:4.2f} A\r".format(iwat + 1, len(pos[0]), sigma),)
return g

def _smear_rmsf(self, coordinates, grid, edges, rmsf):
Expand All @@ -1463,7 +1473,7 @@ def _smear_rmsf(self, coordinates, grid, edges, rmsf):
continue
g += np.fromfunction(self._gaussian_cartesian, grid.shape, dtype=np.int,
c=coord, sigma=rmsf[iwat])
print("Smearing out atom position {0:4d}/{1:5d} with RMSF {2:4.2f} A\r".format(iwat + 1, N, rmsf[iwat]),)
# print("Smearing out atom position {0:4d}/{1:5d} with RMSF {2:4.2f} A\r".format(iwat + 1, N, rmsf[iwat]),)
return g

def _gaussian(self, i, j, k, p, sigma):
Expand Down
67 changes: 50 additions & 17 deletions testsuite/MDAnalysisTests/analysis/test_density.py
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,7 @@
import MDAnalysis as mda
from MDAnalysis.analysis import density

from MDAnalysisTests.datafiles import TPR, XTC, GRO
from MDAnalysisTests.datafiles import TPR, XTC, GRO, PDB_full


class TestDensity(object):
Expand Down Expand Up @@ -282,7 +282,7 @@ def test_density_from_Universe(self, universe, tmpdir):
tmpdir=tmpdir
)

def test_density_from_Universe_sliced(self, universe, tmpdir):
def test_sliced(self, universe, tmpdir):
self.check_density_from_Universe(
self.selections['static'],
self.references['static_sliced']['meandensity'],
Expand All @@ -291,7 +291,7 @@ def test_density_from_Universe_sliced(self, universe, tmpdir):
tmpdir=tmpdir
)

def test_density_from_Universe_update_selection(self, universe, tmpdir):
def test_update_selection(self, universe, tmpdir):
self.check_density_from_Universe(
self.selections['dynamic'],
self.references['dynamic']['meandensity'],
Expand All @@ -300,7 +300,7 @@ def test_density_from_Universe_update_selection(self, universe, tmpdir):
tmpdir=tmpdir
)

def test_density_from_Universe_notwithin(self, universe, tmpdir):
def test_notwithin(self, universe, tmpdir):
self.check_density_from_Universe(
self.selections['static'],
self.references['notwithin']['meandensity'],
Expand All @@ -310,7 +310,7 @@ def test_density_from_Universe_notwithin(self, universe, tmpdir):
tmpdir=tmpdir
)

def test_density_from_Universe_userdefn_eqbox(self, universe, tmpdir):
def test_userdefn_eqbox(self, universe, tmpdir):
self.check_density_from_Universe(
self.selections['static'],
self.references['static_defined']['meandensity'],
Expand All @@ -322,7 +322,7 @@ def test_density_from_Universe_userdefn_eqbox(self, universe, tmpdir):
tmpdir=tmpdir
)

def test_density_from_Universe_userdefn_neqbox(self, universe, tmpdir):
def test_userdefn_neqbox(self, universe, tmpdir):
self.check_density_from_Universe(
self.selections['static'],
self.references['static_defined_unequal']['meandensity'],
Expand All @@ -334,14 +334,14 @@ def test_density_from_Universe_userdefn_neqbox(self, universe, tmpdir):
tmpdir=tmpdir
)

def test_density_from_Universe_userdefn_boxshape(self, universe):
def test_userdefn_boxshape(self, universe):
D = density.density_from_Universe(
universe, select=self.selections['static'],
delta=1.0, xdim=8.0, ydim=12.0, zdim=17.0,
gridcenter=self.gridcenters['static_defined'])
assert D.grid.shape == (8, 12, 17)

def test_density_from_Universe_userdefn_padding(self, universe):
def test_userdefn_padding(self, universe):
regex = ("Box padding \(currently set at 1\.0\) is not used "
"in user defined grids\.")
with pytest.warns(UserWarning, match=regex):
Expand All @@ -350,7 +350,7 @@ def test_density_from_Universe_userdefn_padding(self, universe):
delta=self.delta, xdim=100.0, ydim=100.0, zdim=100.0, padding=1.0,
gridcenter=self.gridcenters['static_defined'])

def test_density_from_Universe_userdefn_selwarning(self, universe):
def test_userdefn_selwarning(self, universe):
regex = ("Atom selection does not fit grid --- "
"you may want to define a larger box")
with pytest.warns(UserWarning, match=regex):
Expand All @@ -359,7 +359,7 @@ def test_density_from_Universe_userdefn_selwarning(self, universe):
delta=self.delta, xdim=1.0, ydim=2.0, zdim=2.0, padding=0.0,
gridcenter=self.gridcenters['static_defined'])

def test_density_from_Universe_userdefn_ValueErrors(self, universe):
def test_userdefn_ValueErrors(self, universe):
# Test len(gridcenter) != 3
with pytest.raises(ValueError):
D = density.density_from_Universe(
Expand All @@ -379,13 +379,14 @@ def test_density_from_Universe_userdefn_ValueErrors(self, universe):
delta=self.delta, xdim="MDAnalysis", ydim=10.0, zdim=10.0,
gridcenter=self.gridcenters['static_defined'])

def test_density_from_Universe_Deprecation_warning(self, universe):
def test_has_DeprecationWarning(self, universe):
with pytest.warns(DeprecationWarning,
match="will be removed in release 2.0.0"):
D = density.density_from_Universe(
universe, select=self.selections['static'],
delta=self.delta)


class TestNotWithin(object):
# tests notwithin_coordinates_factory
# only checks that KDTree and distance_array give same results
Expand Down Expand Up @@ -418,10 +419,42 @@ def test_not_within(self, u):
use_kdtree=False)()
assert_equal(vers1, vers2)

def test_has_DeprecationWarning(self, u):
with pytest.warns(DeprecationWarning,
match="will be removed in release 2.0.0"):
density.notwithin_coordinates_factory(u, 'resname SOL',
'protein', 2,
not_within=False,
use_kdtree=True)()


class TestBfactor2RMSF(object):
@pytest.mark.parametrize('B', [20, 20.0, 20.123,
np.array([45.678, 67.1, 0.0, 100.2])])
def test_Bfactor2RMSF(self, B):
values = density.Bfactor2RMSF(B)
ref = np.sqrt(3. * B / 8.) / np.pi # original equation
assert_almost_equal(values, ref)

@pytest.mark.parametrize('B', [20, 20.0, 20.123,
np.array([45.678, 67.1, 0.0, 100.2])])
def test_Bfactor2RMSF(B):
values = density.Bfactor2RMSF(B)
ref = np.sqrt(3. * B / 8.) / np.pi # original equation
assert_almost_equal(values, ref)
def test_has_DeprecationWarning(self):
with pytest.warns(DeprecationWarning,
match="will be removed in release 2.0.0"):
density.Bfactor2RMSF(100.0)

class Test_density_from_PDB(object):

@pytest.mark.parametrize('sigma,ref_shape,ref_gridsum',
[(None, (21, 26, 29), 108.710149),
(1.5, (21, 26, 29), 160.050167)])
def test_density_from_PDB(self, sigma, ref_shape, ref_gridsum):
# simple test that the code runs: use B-factors from PDB or fixed sigma
D = density.density_from_PDB(PDB_full, delta=2.0, sigma=sigma)

assert isinstance(D, density.Density)
assert_equal(D.grid.shape, ref_shape)
assert_almost_equal(D.grid.sum(), ref_gridsum, decimal=6)

def test_has_DeprecationWarning(self):
with pytest.warns(DeprecationWarning,
match="will be removed in release 2.0.0"):
density.density_from_PDB(PDB_full, delta=5.0, sigma=2)

0 comments on commit 1154532

Please sign in to comment.