diff --git a/package/CHANGELOG b/package/CHANGELOG index fa9d2be518d..7b215b398da 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -26,7 +26,8 @@ Enhancements class itself. * Added boolean flag at lib.distances.USED_OPENMP to reveal if OpenMP was used in compilation (Issue #883) - + * Added Principal Component Analysis module for linear dimensionality + reduction. Fixes * GROWriter resids now truncated properly (Issue #886) * reading/writing lambda value in trr files (Issue #859) diff --git a/package/MDAnalysis/analysis/pca.py b/package/MDAnalysis/analysis/pca.py new file mode 100644 index 00000000000..2acce16c4b7 --- /dev/null +++ b/package/MDAnalysis/analysis/pca.py @@ -0,0 +1,353 @@ +# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; coding:utf-8 -*- +# vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4 +# +# MDAnalysis --- http://www.MDAnalysis.org +# Copyright (c) 2006-2015 Naveen Michaud-Agrawal, Elizabeth J. Denning, Oliver Beckstein +# and contributors (see AUTHORS for the full list) +# +# Released under the GNU Public Licence, v2 or any higher version +# +# Please cite your use of MDAnalysis in published work: +# +# N. Michaud-Agrawal, E. J. Denning, T. B. Woolf, and O. Beckstein. +# MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations. +# J. Comput. Chem. 32 (2011), 2319--2327, doi:10.1002/jcc.21787 +# + +""" +Principal Component Analysis (PCA) --- :mod:`MDAnalysis.analysis.pca` +===================================================================== + +:Authors: John Detlefs +:Year: 2016 +:Copyright: GNU Public License v3 + +This module contains the linear dimensions reduction method Principal Component +Analysis (PCA). PCA sorts a simulation into 3N directions of descending +variance, with N being the number of atoms. These directions are called +the principal components. The dimensions to be analyzed are reduced by only +looking at a few projections of the first principal components. To learn how to +run a Principal Component Analysis, please refer to the :ref:`PCA-tutorial`. + +The PCA problem is solved by solving the eigenvalue problem of the covariance +matrix, a 3N x 3N matrix where the element (i, j) is the covariance between +coordinates i and j. The principal components are the eigenvectors of this +matrix. + +For each eigenvector, its eigenvalue is the variance that the eigenvector +explains. Stored in :attribute:`cumulated_variance`, a ratio for each number of +eigenvectors up to index i is provided to quickly find out how many principal +components are needed to explain the amount of variance reflected by those i +eigenvectors. For most data, :attribute:`cumulated_variance` will be +approximately equal to one for some n that is significantly smaller than the +total number of components, these are the components of interest given +by Principal Component Analysis. + +From here, we can project a trajectory onto these principal components and +attempt to retrieve some structure from our high dimensional data. + +For a basic introduction to the module, the :ref:`PCA-tutorial` shows how +to perform Principal Component Analysis. + +.. _PCA-tutorial: +The example uses files provided as part of the MDAnalysis test suite +(in the variables :data:`~MDAnalysis.tests.datafiles.PSF` and +:data:`~MDAnalysis.tests.datafiles.DCD`). This tutorial shows how to use the +PCA class. + +First load all modules and test data + >>> import MDAnalysis as mda + >>> import MDAnalysis.analysis.pca as pca + >>> from MDAnalysis.tests.datafiles import PSF, DCD + +Given a universe containing trajectory data we can perform Principal Component +Analyis by using the class :class:`PCA` and retrieving the principal +components. + >>> u = mda.Universe(PSF,DCD) + >>> PSF_pca = pca.PCA(u, select='backbone') + >>> PSF_pca.run() + +Inspect the components to determine the principal components you would like +to retain. The choice is arbitrary, but I will stop when 95 percent of the +variance is explained by the components. This cumulated variance by the +components is conveniently stored in the one-dimensional array attribute +`cumulated_variance`. The value at the ith index of `cumulated_variance` is the +sum of the variances from 0 to i. + + >>> n_pcs = np.where(PSF_pca.cumulated_var > 0.95)[0][0] + >>> atomgroup = u.select_atoms('backbone') + >>> pca_space = PSF_pca.transform(atomgroup, n_components=n_pcs) + +From here, inspection of the pca_space and conclusions to be drawn from the +data are left to the user. + +Functions +--------- +.. autoclass:: PCA +.. autofunction:: cosine_content +""" +from six.moves import range +import logging +import warnings + +import numpy as np + +from scipy.integrate import simps +from MDAnalysis.core import AtomGroup +from MDAnalysis import Universe +from MDAnalysis.analysis.align import _fit_to +from MDAnalysis.lib.log import ProgressMeter + +from .base import AnalysisBase + + +class PCA(AnalysisBase): + """Principal component analysis on an MD trajectory. + + After initializing and calling method with a universe or an atom group, + principal components ordering the atom coordinate data by decreasing + variance will be available for analysis. As an example: + >>> pca = PCA(atomgroup, select='backbone').run() + >>> pca_space = pca.transform(atomgroup.select_atoms('backbone'), 3) + generates the principal components of the backbone of the atomgroup and + then transforms those atomgroup coordinates by the direction of those + variances. Please refer to the :ref:`PCA-tutorial` for more detailed + instructions. + + Attributes + ---------- + p_components: array, (n_components, n_atoms * 3) + The principal components of the feature space, + representing the directions of maximum variance in the data. + variance : array (n_components, ) + The raw variance explained by each eigenvector of the covariance + matrix. + cumulated_variance : array, (n_components, ) + Percentage of variance explained by the selected components and the sum + of the components preceding it. If a subset of components is not chosen + then all components are stored and the cumulated variance will converge + to 1. + pca_space : array (n_frames, n_components) + After running :method:`pca.transform(atomgroup)` the projection of the + positions onto the principal components will exist here. + mean_atoms: MDAnalyis atomgroup + After running :method:`PCA.run()`, the mean position of all the atoms + used for the creation of the covariance matrix will exist here. + + Methods + ------- + transform(atomgroup, n_components=None) + Take an atomgroup or universe with the same number of atoms as was + used for the calculation in :method:`PCA.run()` and project it onto the + principal components. + """ + + def __init__(self, atomgroup, select='all', align=False, mean=None, + n_components=None, **kwargs): + """ + Parameters + ---------- + atomgroup: MDAnalysis atomgroup + AtomGroup to be used for PCA. + select: string, optional + A valid selection statement for choosing a subset of atoms from + the atomgroup. + align: boolean, optional + If True, the trajectory will be aligned to a reference + structure. + mean: MDAnalysis atomgroup, optional + An optional reference structure to be used as the mean of the + covariance matrix. + n_components : int, optional + The number of principal components to be saved, default saves + all principal components, Default: None + start : int, optional + First frame of trajectory to use for generation + of covariance matrix, Default: None + stop : int, optional + Last frame of trajectory to use for generation + of covariance matrix, Default: None + step : int, optional + Step between frames of trajectory to use for generation + of covariance matrix, Default: None + """ + super(PCA, self).__init__(atomgroup.universe.trajectory, + **kwargs) + if self.n_frames == 1: + raise ValueError('No covariance information can be gathered from a' + 'single trajectory frame.\n') + + self._u = atomgroup.universe + + if self._quiet: + logging.disable(logging.WARN) + # for transform function + self.align = align + # access 0th index + self._u.trajectory[0] + # reference will be 0th index + self._reference = self._u.select_atoms(select) + self._atoms = self._u.select_atoms(select) + self.n_components = n_components + self._n_atoms = self._atoms.n_atoms + self._calculated = False + + if mean is None: + warnings.warn('In order to demean to generate the covariance ' + 'matrix the frames have to be iterated over twice. ' + 'To avoid this slowdown, provide an atomgroup for ' + 'demeaning.') + self.mean = np.zeros(self._n_atoms*3) + self._calc_mean = True + else: + self.mean = mean.positions + self._calc_mean = False + + def _prepare(self): + n_dim = self._n_atoms * 3 + self.cov = np.zeros((n_dim, n_dim)) + self._ref_atom_positions = self._reference.positions + self._ref_cog = self._reference.center_of_geometry() + self._ref_atom_positions -= self._ref_cog + + if self._calc_mean: + interval = int(self.n_frames // 100) + interval = interval if interval > 0 else 1 + format = ("Mean Calculation Step" + "%(step)5d/%(numsteps)d [%(percentage)5.1f%%]\r") + mean_pm = ProgressMeter(self.n_frames if self.n_frames else 1, + interval=interval, quiet=self._quiet, + format=format) + for i, ts in enumerate(self._u.trajectory[self.start:self.stop: + self.step]): + if self.align: + mobile_cog = self._atoms.center_of_geometry() + mobile_atoms, old_rmsd = _fit_to(self._atoms.positions, + self._ref_atom_positions, + self._atoms, + mobile_com=mobile_cog, + ref_com=self._ref_cog) + else: + self.mean += self._atoms.positions.ravel() + mean_pm.echo(i) + self.mean /= self.n_frames + + atom_positions = self.mean.reshape(self._n_atoms, 3) + self.mean_atoms = AtomGroup.AtomGroup(self._atoms) + self.mean_atoms.positions = atom_positions + + def _single_frame(self): + if self.align: + mobile_cog = self._atoms.center_of_geometry() + mobile_atoms, old_rmsd = _fit_to(self._atoms.positions, + self._ref_atom_positions, + self._atoms, + mobile_com=mobile_cog, + ref_com=self._ref_cog) + # now all structures are aligned to reference + x = mobile_atoms.positions.ravel() + else: + x = self._atoms.positions.ravel() + x -= self.mean + self.cov += np.dot(x[:, np.newaxis], x[:, np.newaxis].T) + + def _conclude(self): + self.cov /= self.n_frames - 1 + e_vals, e_vects = np.linalg.eig(self.cov) + sort_idx = np.argsort(e_vals)[::-1] + self.variance = e_vals[sort_idx] + self.variance = self.variance[:self.n_components] + self.p_components = e_vects[:self.n_components, sort_idx] + self.cumulated_variance = (np.cumsum(self.variance) / + np.sum(self.variance)) + self._calculated = True + + def transform(self, atomgroup, n_components=None, start=None, stop=None, + step=None): + """Apply the dimensionality reduction on a trajectory + + Parameters + ---------- + atomgroup: MDAnalysis atomgroup/ Universe + The atomgroup or universe containing atoms to be PCA transformed. + n_components: int, optional + The number of components to be projected onto, Default none: maps + onto all components. + start: int, optional + The frame to start on for the PCA transform. Default: None becomes + 0, the first frame index. + stop: int, optional + Frame index to stop PCA transform. Default: None becomes n_frames. + Iteration stops *before* this frame number, which means that the + trajectory would be read until the end. + step: int, optional + Number of frames to skip over for PCA transform. Default: None + becomes 1. + + Returns + ------- + pca_space : array, shape (number of frames, number of components) + """ + if not self._calculated: + self.run() + + if isinstance(atomgroup, Universe): + atomgroup = atomgroup.atoms + + if(self._n_atoms != atomgroup.n_atoms): + raise ValueError('PCA has been fit for' + '{} atoms. Your atomgroup' + 'has {} atoms'.format(self._n_atoms, + atomgroup.n_atoms)) + if not (self._atoms.types == atomgroup.types).all(): + warnings.warn('Atom types do not match with types used to fit PCA') + + traj = atomgroup.universe.trajectory + start, stop, step = traj.check_slice_indices(start, stop, step) + n_frames = len(range(start, stop, step)) + + dim = (n_components if n_components is not None else + self.p_components.shape[1]) + + dot = np.zeros((n_frames, dim)) + + for i, ts in enumerate(traj[start:stop:step]): + xyz = atomgroup.positions.ravel() - self.mean + dot[i] = np.dot(xyz, self.p_components[:, :n_components]) + + return dot + + +def cosine_content(pca_space, i): + """Measure the cosine content of the PCA projection. + + The cosine content of pca projections can be used as an indicator if a + simulation is converged. Values close to 1 are an indicator that the + simulation isn't converged. For values below 0.7 no statement can be made. + If you use this function please cite [BerkHess1]_. + + + Parameters + ---------- + pca_space: array, shape (number of frames, number of components) + The PCA space to be analyzed. + i: int + The index of the pca_component projectection to be analyzed. + + Returns + ------- + A float reflecting the cosine content of the ith projection in the PCA + space. The output is bounded by 0 and 1, with 1 reflecting an agreement + with cosine while 0 reflects complete disagreement. + + References + ---------- + .. [BerkHess1] + Berk Hess. Convergence of sampling in protein simulations. Phys. Rev. E + 65, 031910 (2002). + """ + t = np.arange(len(pca_space)) + T = len(pca_space) + cos = np.cos(np.pi * t * (i + 1) / T) + return ((2.0 / T) * (simps(cos*pca_space[:, i])) ** 2 / + simps(pca_space[:, i] ** 2)) diff --git a/package/doc/sphinx/source/documentation_pages/analysis/pca.rst b/package/doc/sphinx/source/documentation_pages/analysis/pca.rst new file mode 100644 index 00000000000..dc90112ad4b --- /dev/null +++ b/package/doc/sphinx/source/documentation_pages/analysis/pca.rst @@ -0,0 +1 @@ +.. automodule:: MDAnalysis.analysis.pca diff --git a/package/doc/sphinx/source/documentation_pages/analysis_modules.rst b/package/doc/sphinx/source/documentation_pages/analysis_modules.rst index d4e241cd6a9..6c3e41a581c 100644 --- a/package/doc/sphinx/source/documentation_pages/analysis_modules.rst +++ b/package/doc/sphinx/source/documentation_pages/analysis_modules.rst @@ -113,3 +113,4 @@ Dimension Reduction :maxdepth: 1 analysis/diffusionmap + analysis/pca diff --git a/testsuite/MDAnalysisTests/analysis/test_pca.py b/testsuite/MDAnalysisTests/analysis/test_pca.py new file mode 100644 index 00000000000..29044fca5ac --- /dev/null +++ b/testsuite/MDAnalysisTests/analysis/test_pca.py @@ -0,0 +1,89 @@ +# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; coding:utf-8 -*- +# vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4 fileencoding=utf-8 +# +# MDAnalysis --- http://www.MDAnalysis.org +# Copyright (c) 2006-2015 Naveen Michaud-Agrawal, Elizabeth J. Denning, Oliver +# Beckstein and contributors (see AUTHORS for the full list) +# +# Released under the GNU Public Licence, v2 or any higher version +# +# Please cite your use of MDAnalysis in published work: +# +# N. Michaud-Agrawal, E. J. Denning, T. B. Woolf, and O. Beckstein. +# MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations. +# J. Comput. Chem. 32 (2011), 2319--2327, doi:10.1002/jcc.21787 +# +from __future__ import print_function +import numpy as np +import MDAnalysis +import MDAnalysis.analysis.pca as pca +from MDAnalysis.analysis.align import _fit_to + +from numpy.testing import (assert_almost_equal, assert_equal, + assert_array_almost_equal, raises) + +from MDAnalysisTests.datafiles import (PDB, XTC, RANDOM_WALK, RANDOM_WALK_TOPO, + waterPSF, waterDCD) + + +class TestPCA(object): + """ Test the PCA class """ + def setUp(self): + self.u = MDAnalysis.Universe(PDB, XTC) + self.pca = pca.PCA(self.u.atoms, select='backbone and name CA', + align=False) + self.pca.run() + self.n_atoms = self.u.select_atoms('backbone and name CA').n_atoms + + def test_cov(self): + atoms = self.u.select_atoms('backbone and name CA') + xyz = np.zeros((self.pca.n_frames, self.pca._n_atoms*3)) + for i, ts in enumerate(self.u.trajectory): + xyz[i] = atoms.positions.ravel() + + cov = np.cov(xyz, rowvar=0) + assert_array_almost_equal(self.pca.cov, cov, 4) + + def test_cum_var(self): + assert_almost_equal(self.pca.cumulated_variance[-1], 1) + l = self.pca.cumulated_variance + l = np.sort(l) + assert_almost_equal(self.pca.cumulated_variance, l, 5) + + def test_pcs(self): + assert_equal(self.pca.p_components.shape, + (self.n_atoms*3, self.n_atoms*3)) + + def test_different_steps(self): + dot = self.pca.transform(self.u.select_atoms('backbone and name CA'), + start=5, stop=7, step=1) + assert_equal(dot.shape, (2, self.n_atoms*3)) + + def test_transform(self): + self.ag = self.u.select_atoms('backbone and name CA') + self.pca_space = self.pca.transform(self.ag, n_components=1) + assert_equal(self.pca_space.shape, + (self.u.trajectory.n_frames, 1)) + + # Accepts universe as input, but shapes are not aligned due to n_atoms + @raises(ValueError) + def test_transform_mismatch(self): + self.ag = self.u.select_atoms('backbone and name CA') + self.pca_space = self.pca.transform(self.u, n_components=1) + assert_equal(self.pca_space.shape, + (self.u.trajectory.n_frames, 1)) + + @staticmethod + def test_transform_universe(): + u1 = MDAnalysis.Universe(waterPSF, waterDCD) + u2 = MDAnalysis.Universe(waterPSF, waterDCD) + pca_test = pca.PCA(u1).run() + pca_test.transform(u2) + + @staticmethod + def test_cosine_content(): + rand = MDAnalysis.Universe(RANDOM_WALK_TOPO, RANDOM_WALK) + pca_random = pca.PCA(rand.atoms).run() + dot = pca_random.transform(rand.atoms) + content = pca.cosine_content(dot, 0) + assert_almost_equal(content, .99, 1) diff --git a/testsuite/MDAnalysisTests/data/RANDOM_WALK_TOPO.pdb b/testsuite/MDAnalysisTests/data/RANDOM_WALK_TOPO.pdb new file mode 100644 index 00000000000..ace4d7dd241 --- /dev/null +++ b/testsuite/MDAnalysisTests/data/RANDOM_WALK_TOPO.pdb @@ -0,0 +1,105 @@ +TITLE MDANALYSIS FRAMES FROM 0, SKIP 1: Created by PDBWriter +CRYST1 0.000 0.000 0.000 0.00 0.00 0.00 P 1 1 +MODEL 1 +ATOM 1 X SYST 1 0.192 -0.797 -1.042 1.00 0.00 X +ATOM 2 X SYST 1 -0.150 0.308 0.957 1.00 0.00 X +ATOM 3 X SYST 1 -0.752 -0.125 -0.563 1.00 0.00 X +ATOM 4 X SYST 1 0.483 -1.342 -1.400 1.00 0.00 X +ATOM 5 X SYST 1 0.550 -0.672 -1.462 1.00 0.00 X +ATOM 6 X SYST 1 -0.644 -0.858 0.751 1.00 0.00 X +ATOM 7 X SYST 1 1.092 1.182 -1.172 1.00 0.00 X +ATOM 8 X SYST 1 -1.028 1.086 0.098 1.00 0.00 X +ATOM 9 X SYST 1 -2.109 0.854 -0.326 1.00 0.00 X +ATOM 10 X SYST 1 -0.430 3.140 0.412 1.00 0.00 X +ATOM 11 X SYST 1 -1.051 0.127 -0.654 1.00 0.00 X +ATOM 12 X SYST 1 1.934 1.939 -0.671 1.00 0.00 X +ATOM 13 X SYST 1 -1.622 2.113 -1.297 1.00 0.00 X +ATOM 14 X SYST 1 -0.964 -1.267 -1.626 1.00 0.00 X +ATOM 15 X SYST 1 0.079 -1.222 -1.151 1.00 0.00 X +ATOM 16 X SYST 1 1.328 1.413 -1.109 1.00 0.00 X +ATOM 17 X SYST 1 2.707 -3.728 2.036 1.00 0.00 X +ATOM 18 X SYST 1 -1.879 -1.751 -0.008 1.00 0.00 X +ATOM 19 X SYST 1 0.509 0.497 0.657 1.00 0.00 X +ATOM 20 X SYST 1 -1.859 -1.078 1.447 1.00 0.00 X +ATOM 21 X SYST 1 -1.484 0.366 1.113 1.00 0.00 X +ATOM 22 X SYST 1 1.004 -0.138 -1.316 1.00 0.00 X +ATOM 23 X SYST 1 -1.070 -1.564 -0.288 1.00 0.00 X +ATOM 24 X SYST 1 -1.634 -1.837 1.273 1.00 0.00 X +ATOM 25 X SYST 1 -0.578 -1.060 -1.521 1.00 0.00 X +ATOM 26 X SYST 1 0.855 1.954 -0.705 1.00 0.00 X +ATOM 27 X SYST 1 1.089 -2.214 0.386 1.00 0.00 X +ATOM 28 X SYST 1 0.654 -0.547 0.125 1.00 0.00 X +ATOM 29 X SYST 1 1.196 0.560 -1.884 1.00 0.00 X +ATOM 30 X SYST 1 1.537 -1.519 1.516 1.00 0.00 X +ATOM 31 X SYST 1 -0.333 0.382 -1.087 1.00 0.00 X +ATOM 32 X SYST 1 0.758 -0.508 0.637 1.00 0.00 X +ATOM 33 X SYST 1 0.912 -0.132 -0.114 1.00 0.00 X +ATOM 34 X SYST 1 0.514 -2.753 -1.568 1.00 0.00 X +ATOM 35 X SYST 1 -2.796 0.934 1.568 1.00 0.00 X +ATOM 36 X SYST 1 -0.316 -0.149 -1.460 1.00 0.00 X +ATOM 37 X SYST 1 1.589 -1.918 1.128 1.00 0.00 X +ATOM 38 X SYST 1 0.893 0.243 1.651 1.00 0.00 X +ATOM 39 X SYST 1 -1.718 -1.066 1.226 1.00 0.00 X +ATOM 40 X SYST 1 0.329 0.702 0.037 1.00 0.00 X +ATOM 41 X SYST 1 -0.093 -1.030 -0.573 1.00 0.00 X +ATOM 42 X SYST 1 -0.245 0.664 1.448 1.00 0.00 X +ATOM 43 X SYST 1 1.946 -2.568 0.803 1.00 0.00 X +ATOM 44 X SYST 1 -2.092 -1.243 -0.670 1.00 0.00 X +ATOM 45 X SYST 1 -0.315 0.386 -0.957 1.00 0.00 X +ATOM 46 X SYST 1 -2.788 -1.573 -0.351 1.00 0.00 X +ATOM 47 X SYST 1 -0.680 -1.665 0.872 1.00 0.00 X +ATOM 48 X SYST 1 1.844 1.731 0.719 1.00 0.00 X +ATOM 49 X SYST 1 -0.008 2.494 -1.046 1.00 0.00 X +ATOM 50 X SYST 1 0.594 0.757 -0.859 1.00 0.00 X +ATOM 51 X SYST 1 -1.617 -0.093 2.723 1.00 0.00 X +ATOM 52 X SYST 1 -0.160 1.120 2.035 1.00 0.00 X +ATOM 53 X SYST 1 -0.241 -1.276 1.334 1.00 0.00 X +ATOM 54 X SYST 1 -0.650 -0.517 0.589 1.00 0.00 X +ATOM 55 X SYST 1 0.140 -1.181 -2.026 1.00 0.00 X +ATOM 56 X SYST 1 0.263 0.025 -0.591 1.00 0.00 X +ATOM 57 X SYST 1 -2.602 0.590 -0.103 1.00 0.00 X +ATOM 58 X SYST 1 1.164 -1.947 -0.433 1.00 0.00 X +ATOM 59 X SYST 1 0.335 -2.007 -1.063 1.00 0.00 X +ATOM 60 X SYST 1 3.312 -0.641 0.481 1.00 0.00 X +ATOM 61 X SYST 1 0.318 -0.649 -1.398 1.00 0.00 X +ATOM 62 X SYST 1 0.846 1.393 -1.362 1.00 0.00 X +ATOM 63 X SYST 1 2.293 1.473 0.462 1.00 0.00 X +ATOM 64 X SYST 1 -0.713 -1.557 1.726 1.00 0.00 X +ATOM 65 X SYST 1 2.920 -1.216 -0.189 1.00 0.00 X +ATOM 66 X SYST 1 -1.824 1.101 0.334 1.00 0.00 X +ATOM 67 X SYST 1 0.151 -0.205 -0.040 1.00 0.00 X +ATOM 68 X SYST 1 -1.562 0.608 -0.604 1.00 0.00 X +ATOM 69 X SYST 1 2.523 1.642 0.167 1.00 0.00 X +ATOM 70 X SYST 1 -0.219 0.873 -0.553 1.00 0.00 X +ATOM 71 X SYST 1 1.307 -0.766 0.101 1.00 0.00 X +ATOM 72 X SYST 1 -0.034 -2.097 0.004 1.00 0.00 X +ATOM 73 X SYST 1 0.315 0.406 -0.990 1.00 0.00 X +ATOM 74 X SYST 1 -1.116 -0.050 -0.500 1.00 0.00 X +ATOM 75 X SYST 1 1.444 0.385 -0.008 1.00 0.00 X +ATOM 76 X SYST 1 1.172 1.610 -0.969 1.00 0.00 X +ATOM 77 X SYST 1 0.605 -0.021 -1.358 1.00 0.00 X +ATOM 78 X SYST 1 -2.316 0.837 0.168 1.00 0.00 X +ATOM 79 X SYST 1 0.430 -1.242 0.087 1.00 0.00 X +ATOM 80 X SYST 1 2.290 -0.587 -1.829 1.00 0.00 X +ATOM 81 X SYST 1 0.144 -0.724 -0.735 1.00 0.00 X +ATOM 82 X SYST 1 0.756 -0.873 -0.661 1.00 0.00 X +ATOM 83 X SYST 1 1.778 1.564 1.385 1.00 0.00 X +ATOM 84 X SYST 1 -1.414 1.198 0.396 1.00 0.00 X +ATOM 85 X SYST 1 -0.075 1.270 1.709 1.00 0.00 X +ATOM 86 X SYST 1 -3.111 -1.338 -1.163 1.00 0.00 X +ATOM 87 X SYST 1 0.021 3.037 -0.199 1.00 0.00 X +ATOM 88 X SYST 1 -1.968 -2.604 0.683 1.00 0.00 X +ATOM 89 X SYST 1 -0.082 1.604 -1.400 1.00 0.00 X +ATOM 90 X SYST 1 -0.772 1.781 -1.497 1.00 0.00 X +ATOM 91 X SYST 1 0.159 1.882 0.874 1.00 0.00 X +ATOM 92 X SYST 1 -2.455 -0.851 0.206 1.00 0.00 X +ATOM 93 X SYST 1 -0.461 1.583 3.403 1.00 0.00 X +ATOM 94 X SYST 1 1.318 -1.974 -0.323 1.00 0.00 X +ATOM 95 X SYST 1 0.742 -1.308 -1.363 1.00 0.00 X +ATOM 96 X SYST 1 1.760 0.902 1.184 1.00 0.00 X +ATOM 97 X SYST 1 1.740 -0.107 -1.524 1.00 0.00 X +ATOM 98 X SYST 1 0.084 -1.835 -0.533 1.00 0.00 X +ATOM 99 X SYST 1 2.022 2.084 0.112 1.00 0.00 X +ATOM 100 X SYST 1 -0.454 0.563 2.780 1.00 0.00 X +ENDMDL +END diff --git a/testsuite/MDAnalysisTests/data/xyz_random_walk.xtc b/testsuite/MDAnalysisTests/data/xyz_random_walk.xtc new file mode 100644 index 00000000000..39a68df7b7d Binary files /dev/null and b/testsuite/MDAnalysisTests/data/xyz_random_walk.xtc differ diff --git a/testsuite/MDAnalysisTests/datafiles.py b/testsuite/MDAnalysisTests/datafiles.py index df4432a62bf..ec40b57e00a 100644 --- a/testsuite/MDAnalysisTests/datafiles.py +++ b/testsuite/MDAnalysisTests/datafiles.py @@ -116,7 +116,9 @@ "COORDINATES_TOPOLOGY", "NUCLsel", "GRO_empty_atom", "GRO_missing_atomname", # for testing GROParser exception raise - "ENT" #for testing ENT file extension + "ENT", #for testing ENT file extension + "RANDOM_WALK", + "RANDOM_WALK_TOPO" # garbage topology to go along with XTC positions above ] from pkg_resources import resource_filename @@ -329,6 +331,8 @@ # Contains one of each residue in 'nucleic' selections NUCLsel = resource_filename(__name__, 'data/nucl_res.pdb') +RANDOM_WALK = resource_filename(__name__, 'data/xyz_random_walk.xtc') +RANDOM_WALK_TOPO = resource_filename(__name__, 'data/RANDOM_WALK_TOPO.pdb') # This should be the last line: clean up namespace del resource_filename