diff --git a/package/MDAnalysis/core/groups.py b/package/MDAnalysis/core/groups.py index 95e38e600a3..e3db00c298f 100644 --- a/package/MDAnalysis/core/groups.py +++ b/package/MDAnalysis/core/groups.py @@ -100,11 +100,15 @@ import os import warnings +import Bio.Seq +import Bio.SeqRecord +import Bio.Alphabet + from numpy.lib.utils import deprecate from .. import _ANCHOR_UNIVERSES from ..lib import util -from ..lib.util import cached, warn_if_not_unique, unique_int_1d +from ..lib.util import cached, warn_if_not_unique, unique_int_1d, convert_aa_code from ..lib import distances from ..lib import transformations from ..lib import mdamath @@ -891,6 +895,67 @@ def center_of_geometry(self, pbc=None, compound='group', unwrap=False): centroid = center_of_geometry + @warn_if_not_unique + @check_pbc_and_unwrap + def center_of_mass(self, pbc=None, compound='group', unwrap=False): + """Center of mass of (compounds of) the group. + + Computes the center of mass of :class:`Atoms` in the group. + Centers of mass per :class:`Residue`, :class:`Segment`, molecule, or + fragment can be obtained by setting the `compound` parameter + accordingly. If the masses of a compound sum up to zero, the + center of mass coordinates of that compound will be ``nan`` (not a + number). + + Parameters + ---------- + pbc : bool, optional + If ``True`` and `compound` is ``'group'``, move all atoms to the + primary unit cell before calculation. + If ``True`` and `compound` is ``'segments'`` or ``'residues'``, the + centers of mass of each compound will be calculated without moving + any atoms to keep the compounds intact. Instead, the resulting + center-of-mass position vectors will be moved to the primary unit + cell after calculation. + compound : {'group', 'segments', 'residues', 'molecules', 'fragments'},\ + optional + If ``'group'``, the center of mass of all atoms in the self will + be returned as a single position vector. Otherwise, the centers of + mass of each :class:`Segment`, :class:`Residue`, molecule, or + fragment will be returned as an array of position vectors, i.e. a 2d + array. + Note that, in any case, *only* the positions of :class:`Atoms` + *belonging to the group* will be taken into account. + unwrap : bool, optional + If ``True``, compounds will be unwrapped before computing their centers. + + Returns + ------- + center : numpy.ndarray + Position vector(s) of the center(s) of mass of the self. + If `compound` was set to ``'group'``, the output will be a single + position vector. + If `compound` was set to ``'segments'`` or ``'residues'``, the + output will be a 2d coordinate array of shape ``(n, 3)`` where ``n`` + is the number of compounds. + + Note + ---- + * This method can only be accessed if the underlying topology has + information about atomic masses. + * The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to + ``True`` allows the *pbc* flag to be used by default. + + + .. versionchanged:: 0.8 Added `pbc` parameter + .. versionchanged:: 0.19.0 Added `compound` parameter + .. versionchanged:: 0.20.0 Added ``'molecules'`` and ``'fragments'`` + compounds + .. versionchanged:: 0.20.0 Added `unwrap` parameter + """ + atoms = self.atoms + return atoms.center(weights=atoms.masses, pbc=pbc, compound=compound, unwrap=unwrap) + @warn_if_not_unique def accumulate(self, attribute, function=np.sum, compound='group'): """Accumulates the attribute associated with (compounds of) the group. @@ -899,7 +964,7 @@ def accumulate(self, attribute, function=np.sum, compound='group'): The accumulation per :class:`Residue`, :class:`Segment`, molecule, or fragment can be obtained by setting the `compound` parameter accordingly. By default, the method sums up all attributes per compound, - but any function that takes an array and returns an acuumulation over a + but any function that takes an array and returns an accumulation over a given axis can be used. For multi-dimensional input arrays, the accumulation is performed along the first axis. @@ -927,7 +992,7 @@ def accumulate(self, attribute, function=np.sum, compound='group'): Returns ------- float or numpy.ndarray - Acuumulation of the `attribute`. + Accumulation of the `attribute`. If `compound` is set to ``'group'``, the first dimension of the `attribute` array will be contracted to a single value. If `compound` is set to ``'segments'``, ``'residues'``, @@ -1036,6 +1101,356 @@ def accumulate(self, attribute, function=np.sum, compound='group'): accumulation[compound_mask] = _accumulation return accumulation + @warn_if_not_unique + def total_mass(self, compound='group'): + """Total mass of (compounds of) the group. + + Computes the total mass of :class:`Atoms` in the group. + Total masses per :class:`Residue`, :class:`Segment`, molecule, or + fragment can be obtained by setting the `compound` parameter + accordingly. + + Parameters + ---------- + compound : {'group', 'segments', 'residues', 'molecules', 'fragments'},\ + optional + If ``'group'``, the total mass of all atoms in the self will be + returned as a single value. Otherwise, the total masses per + :class:`Segment`, :class:`Residue`, molecule, or fragment will be + returned as a 1d array. + Note that, in any case, *only* the masses of :class:`Atoms` + *belonging to the self* will be taken into account. + + Returns + ------- + float or numpy.ndarray + Total mass of (compounds of) the self. + If `compound` was set to ``'group'``, the output will be a single + value. Otherwise, the output will be a 1d array of shape ``(n,)`` + where ``n`` is the number of compounds. + + + .. versionchanged:: 0.20.0 Added `compound` parameter + """ + return self.accumulate("masses", compound=compound) + + @warn_if_not_unique + def total_charge(self, compound='group'): + """Total charge of (compounds of) the self. + + Computes the total charge of :class:`Atoms` in the group. + Total charges per :class:`Residue`, :class:`Segment`, molecule, or + fragment can be obtained by setting the `compound` parameter + accordingly. + + Parameters + ---------- + compound : {'group', 'segments', 'residues', 'molecules', 'fragments'},\ + optional + If 'group', the total charge of all atoms in the group will + be returned as a single value. Otherwise, the total charges per + :class:`Segment`, :class:`Residue`, molecule, or fragment + will be returned as a 1d array. + Note that, in any case, *only* the charges of :class:`Atoms` + *belonging to the self* will be taken into account. + + Returns + ------- + float or numpy.ndarray + Total charge of (compounds of) the group. + If `compound` was set to ``'group'``, the output will be a single + value. Otherwise, the output will be a 1d array of shape ``(n,)`` + where ``n`` is the number of compounds. + + + .. versionchanged:: 0.20.0 Added `compound` parameter + """ + return self.accumulate("charges", compound=compound) + + @warn_if_not_unique + def radius_of_gyration(self, **kwargs): + """Radius of gyration. + + Parameters + ---------- + pbc : bool, optional + If ``True``, move all atoms within the primary unit cell before + calculation. [``False``] + + Note + ---- + The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to + ``True`` allows the *pbc* flag to be used by default. + + + .. versionchanged:: 0.8 Added *pbc* keyword + + """ + atomgroup = self.atoms + pbc = kwargs.pop('pbc', flags['use_pbc']) + masses = atomgroup.masses + + com = atomgroup.center_of_mass(pbc=pbc) + if pbc: + recenteredpos = atomgroup.pack_into_box(inplace=False) - com + else: + recenteredpos = atomgroup.positions - com + + rog_sq = np.sum(masses * np.sum(recenteredpos**2, + axis=1)) / atomgroup.total_mass() + + return np.sqrt(rog_sq) + + @warn_if_not_unique + def shape_parameter(self, **kwargs): + """Shape parameter. + + See [Dima2004a]_ for background information. + + Parameters + ---------- + pbc : bool, optional + If ``True``, move all atoms within the primary unit cell before + calculation. [``False``] + + Note + ---- + The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to + ``True`` allows the *pbc* flag to be used by default. + + + References + ---------- + .. [Dima2004a] Dima, R. I., & Thirumalai, D. (2004). Asymmetry + in the shapes of folded and denatured states of + proteins. *J Phys Chem B*, 108(21), + 6564-6570. doi:`10.1021/jp037128y + `_ + + + .. versionadded:: 0.7.7 + .. versionchanged:: 0.8 Added *pbc* keyword + + """ + atomgroup = self.atoms + pbc = kwargs.pop('pbc', flags['use_pbc']) + masses = atomgroup.masses + + com = atomgroup.center_of_mass(pbc=pbc) + if pbc: + recenteredpos = atomgroup.pack_into_box(inplace=False) - com + else: + recenteredpos = atomgroup.positions - com + tensor = np.zeros((3, 3)) + + for x in range(recenteredpos.shape[0]): + tensor += masses[x] * np.outer(recenteredpos[x, :], + recenteredpos[x, :]) + tensor /= atomgroup.total_mass() + eig_vals = np.linalg.eigvalsh(tensor) + shape = 27.0 * np.prod(eig_vals - np.mean(eig_vals)) / np.power(np.sum(eig_vals), 3) + + return shape + + @warn_if_not_unique + @check_pbc_and_unwrap + def asphericity(self, pbc=None, unwrap=None, compound='group'): + """Asphericity. + + See [Dima2004b]_ for background information. + + Parameters + ---------- + pbc : bool, optional + If ``True``, move all atoms within the primary unit cell before + calculation. If ``None`` use value defined in + MDAnalysis.core.flags['use_pbc'] + unwrap : bool, optional + If ``True``, compounds will be unwrapped before computing their centers. + compound : {'group', 'segments', 'residues', 'molecules', 'fragments'}, optional + Which type of component to keep together during unwrapping. + + Note + ---- + The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to + ``True`` allows the *pbc* flag to be used by default. + + + References + ---------- + + .. [Dima2004b] Dima, R. I., & Thirumalai, D. (2004). Asymmetry + in the shapes of folded and denatured states of + proteins. *J Phys Chem B*, 108(21), + 6564-6570. doi:`10.1021/jp037128y + `_ + + + + .. versionadded:: 0.7.7 + .. versionchanged:: 0.8 Added *pbc* keyword + .. versionchanged:: 0.20.0 Added *unwrap* and *compound* parameter + + """ + atomgroup = self.atoms + if pbc is None: + pbc = flags['use_pbc'] + masses = atomgroup.masses + + com = atomgroup.center_of_mass(pbc=pbc, unwrap=unwrap, compound=compound) + if compound is not 'group': + com = (com * self.masses[:, None]).sum(axis=0) / self.masses.sum() + + if pbc: + recenteredpos = (atomgroup.pack_into_box(inplace=False) - com) + elif unwrap: + recenteredpos = (atomgroup.unwrap(inplace=False) - com) + else: + recenteredpos = (atomgroup.positions - com) + + tensor = np.zeros((3, 3)) + for x in range(recenteredpos.shape[0]): + tensor += masses[x] * np.outer(recenteredpos[x], + recenteredpos[x]) + + tensor /= atomgroup.total_mass() + eig_vals = np.linalg.eigvalsh(tensor) + shape = (3.0 / 2.0) * (np.sum((eig_vals - np.mean(eig_vals))**2) / + np.sum(eig_vals)**2) + + return shape + + @warn_if_not_unique + def principal_axes(self, pbc=None): + """Calculate the principal axes from the moment of inertia. + + e1,e2,e3 = AtomGroup.principal_axes() + + The eigenvectors are sorted by eigenvalue, i.e. the first one + corresponds to the highest eigenvalue and is thus the first principal + axes. + + Parameters + ---------- + pbc : bool, optional + If ``True``, move all atoms within the primary unit cell before + calculation. If ``None`` use value defined in setup flags. + + Returns + ------- + axis_vectors : array + 3 x 3 array with ``v[0]`` as first, ``v[1]`` as second, and + ``v[2]`` as third eigenvector. + + Note + ---- + The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to + ``True`` allows the *pbc* flag to be used by default. + + + .. versionchanged:: 0.8 Added *pbc* keyword + + """ + atomgroup = self.atoms + if pbc is None: + pbc = flags['use_pbc'] + e_val, e_vec = np.linalg.eig(atomgroup.moment_of_inertia(pbc=pbc)) + + # Sort + indices = np.argsort(e_val)[::-1] + # Return transposed in more logical form. See Issue 33. + return e_vec[:, indices].T + + def align_principal_axis(self, axis, vector): + """Align principal axis with index `axis` with `vector`. + + Parameters + ---------- + axis : {0, 1, 2} + Index of the principal axis (0, 1, or 2), as produced by + :meth:`~principal_axes`. + vector : array_like + Vector to align principal axis with. + + Notes + ----- + To align the long axis of a channel (the first principal axis, i.e. + *axis* = 0) with the z-axis:: + + u.atoms.align_principal_axis(0, [0,0,1]) + u.atoms.write("aligned.pdb") + """ + p = self.principal_axes()[axis] + angle = np.degrees(mdamath.angle(p, vector)) + ax = transformations.rotaxis(p, vector) + # print "principal[%d] = %r" % (axis, p) + # print "axis = %r, angle = %f deg" % (ax, angle) + return self.rotateby(angle, ax) + + @warn_if_not_unique + @check_pbc_and_unwrap + def moment_of_inertia(self, **kwargs): + """Tensor moment of inertia relative to center of mass as 3x3 numpy + array. + + Parameters + ---------- + pbc : bool, optional + If ``True``, move all atoms within the primary unit cell before + calculation. [``False``] + + Note + ---- + The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to + ``True`` allows the *pbc* flag to be used by default. + + + .. versionchanged:: 0.8 Added *pbc* keyword + .. versionchanged:: 0.20.0 Added `unwrap` parameter + + """ + atomgroup = self.atoms + pbc = kwargs.pop('pbc', flags['use_pbc']) + unwrap = kwargs.pop('unwrap', False) + compound = kwargs.pop('compound', 'group') + + com = atomgroup.center_of_mass(pbc=pbc, unwrap=unwrap, compound=compound) + if compound is not 'group': + com = (com * self.masses[:, None]).sum(axis=0) / self.masses.sum() + + if pbc: + pos = atomgroup.pack_into_box(inplace=False) - com + elif unwrap: + pos = atomgroup.unwrap(compound=compound, inplace=False) - com + else: + pos = atomgroup.positions - com + + masses = atomgroup.masses + # Create the inertia tensor + # m_i = mass of atom i + # (x_i, y_i, z_i) = pos of atom i + # Ixx = sum(m_i*(y_i^2+z_i^2)); + # Iyy = sum(m_i*(x_i^2+z_i^2)); + # Izz = sum(m_i*(x_i^2+y_i^2)) + # Ixy = Iyx = -1*sum(m_i*x_i*y_i) + # Ixz = Izx = -1*sum(m_i*x_i*z_i) + # Iyz = Izy = -1*sum(m_i*y_i*z_i) + tens = np.zeros((3, 3), dtype=np.float64) + # xx + tens[0][0] = (masses * (pos[:, 1] ** 2 + pos[:, 2] ** 2)).sum() + # xy & yx + tens[0][1] = tens[1][0] = - (masses * pos[:, 0] * pos[:, 1]).sum() + # xz & zx + tens[0][2] = tens[2][0] = - (masses * pos[:, 0] * pos[:, 2]).sum() + # yy + tens[1][1] = (masses * (pos[:, 0] ** 2 + pos[:, 2] ** 2)).sum() + # yz + zy + tens[1][2] = tens[2][1] = - (masses * pos[:, 1] * pos[:, 2]).sum() + # zz + tens[2][2] = (masses * (pos[:, 0] ** 2 + pos[:, 1] ** 2)).sum() + + return tens + def bbox(self, **kwargs): """Return the bounding box of the selection. @@ -1447,7 +1862,7 @@ def wrap(self, compound="atoms", center="com", box=None, inplace=True): atoms = _atoms comp = compound.lower() - if comp not in ('atoms', 'group', 'segments', 'residues', 'molecules', \ + if comp not in ('atoms', 'group', 'segments', 'residues', 'molecules', 'fragments'): raise ValueError("Unrecognized compound definition '{}'. " "Please use one of 'atoms', 'group', 'segments', " @@ -1461,7 +1876,7 @@ def wrap(self, compound="atoms", center="com", box=None, inplace=True): positions = distances.apply_PBC(atoms.positions, box) else: ctr = center.lower() - if ctr == 'com': + if ctr == 'com': # Don't use hasattr(self, 'masses') because that's incredibly # slow for ResidueGroups or SegmentGroups if not hasattr(self._u._topology, 'masses'): @@ -1890,7 +2305,7 @@ def subtract(self, other): The original order of this group is kept, as well as any duplicate elements. If an element of this Group is duplicated and appears in - the other Group or Component, then all the occurences of that element + the other Group or Component, then all the occurrences of that element are removed from the returned Group. Parameters @@ -3168,6 +3583,76 @@ def write(self, filename=None, file_format=None, raise ValueError("No writer found for format: {}".format(filename)) + @property + def fragindices(self): + """The + :class:`fragment indices` + of all :class:`Atoms` in this + :class:`~MDAnalysis.core.groups.AtomGroup`. + + A :class:`numpy.ndarray` with + :attr:`~numpy.ndarray.shape`\ ``=(``\ :attr:`~AtomGroup.n_atoms`\ ``,)`` + and :attr:`~numpy.ndarray.dtype`\ ``=numpy.int64``. + + Note + ---- + This property is only accessible if the underlying topology contains + bond information. + + + .. versionadded:: 0.20.0 + """ + fragdict = self.universe._fragdict + return np.array([fragdict[aix].ix for aix in self.ix], dtype=np.int64) + + @property + def fragments(self): + """Read-only :class:`tuple` of + :class:`fragments`. + + Contains all fragments that + any :class:`~MDAnalysis.core.groups.Atom` in this + :class:`~MDAnalysis.core.groups.AtomGroup` is part of. + + A fragment is a + :class:`group of atoms` which are + interconnected by :class:`~MDAnalysis.core.topologyattrs.Bonds`, i.e., + there exists a path along one + or more :class:`~MDAnalysis.core.topologyattrs.Bonds` between any pair + of :class:`Atoms` + within a fragment. Thus, a fragment typically corresponds to a molecule. + + Note + ---- + * This property is only accessible if the underlying topology contains + bond information. + * The contents of the fragments may extend beyond the contents of this + :class:`~MDAnalysis.core.groups.AtomGroup`. + + + .. versionadded:: 0.9.0 + """ + fragdict = self.universe._fragdict + return tuple(sorted(set(fragdict[aix].fragment for aix in self.ix), + key=lambda x: x[0].ix)) + + @property + def n_fragments(self): + """The number of unique + :class:`~MDAnalysis.core.topologyattrs.Bonds.fragments` the + :class:`Atoms` of this + :class:`~MDAnalysis.core.groups.AtomGroup` are part of. + + Note + ---- + This property is only accessible if the underlying topology contains + bond information. + + + .. versionadded:: 0.20.0 + """ + return len(unique_int_1d(self.fragindices)) + class ResidueGroup(GroupBase): """ResidueGroup base class. @@ -3331,6 +3816,99 @@ def unique(self): _unique._cache['unique'] = _unique return _unique + def sequence(self, **kwargs): + """Returns the amino acid sequence. + + The format of the sequence is selected with the keyword *format*: + + ============== ============================================ + *format* description + ============== ============================================ + 'SeqRecord' :class:`Bio.SeqRecord.SeqRecord` (default) + 'Seq' :class:`Bio.Seq.Seq` + 'string' string + ============== ============================================ + + The sequence is returned by default (keyword ``format = 'SeqRecord'``) + as a :class:`Bio.SeqRecord.SeqRecord` instance, which can then be + further processed. In this case, all keyword arguments (such as the + *id* string or the *name* or the *description*) are directly passed to + :class:`Bio.SeqRecord.SeqRecord`. + + If the keyword *format* is set to ``'Seq'``, all *kwargs* are ignored + and a :class:`Bio.Seq.Seq` instance is returned. The difference to the + record is that the record also contains metadata and can be directly + used as an input for other functions in :mod:`Bio`. + + If the keyword *format* is set to ``'string'``, all *kwargs* are + ignored and a Python string is returned. + + .. rubric:: Example: Write FASTA file + + Use :func:`Bio.SeqIO.write`, which takes sequence records:: + + import Bio.SeqIO + + # get the sequence record of a protein component of a Universe + protein = u.select_atoms("protein") + record = protein.sequence(id="myseq1", name="myprotein") + + Bio.SeqIO.write(record, "single.fasta", "fasta") + + A FASTA file with multiple entries can be written with :: + + Bio.SeqIO.write([record1, record2, ...], "multi.fasta", "fasta") + + Parameters + ---------- + format : string, optional + - ``"string"``: return sequence as a string of 1-letter codes + - ``"Seq"``: return a :class:`Bio.Seq.Seq` instance + - ``"SeqRecord"``: return a :class:`Bio.SeqRecord.SeqRecord` + instance + + Default is ``"SeqRecord"`` + id : optional + Sequence ID for SeqRecord (should be different for different + sequences) + name : optional + Name of the protein. + description : optional + Short description of the sequence. + kwargs : optional + Any other keyword arguments that are understood by + class:`Bio.SeqRecord.SeqRecord`. + + Raises + ------ + :exc:`ValueError` if a residue name cannot be converted to a + 1-letter IUPAC protein amino acid code; make sure to only + select protein residues. + + :exc:`TypeError` if an unknown *format* is selected. + + + .. versionadded:: 0.9.0 + """ + formats = ('string', 'Seq', 'SeqRecord') + + format = kwargs.pop("format", "SeqRecord") + if format not in formats: + raise TypeError("Unknown format='{0}': must be one of: {1}".format( + format, ", ".join(formats))) + try: + sequence = "".join([convert_aa_code(r) for r in self.residues.resnames]) + except KeyError as err: + six.raise_from(ValueError("AtomGroup contains a residue name '{0}' that " + "does not have a IUPAC protein 1-letter " + "character".format(err.message)), None) + if format == "string": + return sequence + seq = Bio.Seq.Seq(sequence, alphabet=Bio.Alphabet.IUPAC.protein) + if format == "Seq": + return seq + return Bio.SeqRecord.SeqRecord(seq, **kwargs) + class SegmentGroup(GroupBase): """:class:`SegmentGroup` base class. @@ -3714,6 +4292,53 @@ def force(self, values): except (AttributeError, NoDataError): raise_from(NoDataError("Timestep does not contain forces"), None) + @property + def bonded_atoms(self): + """An :class:`~MDAnalysis.core.groups.AtomGroup` of all + :class:`Atoms` bonded to this + :class:`~MDAnalysis.core.groups.Atom`.""" + idx = [b.partner(self).index for b in self.bonds] + return self.universe.atoms[idx] + + @property + def fragindex(self): + """The index (ID) of the + :class:`~MDAnalysis.core.topologyattrs.Bonds.fragment` this + :class:`~MDAnalysis.core.groups.Atom` is part of. + + Note + ---- + This property is only accessible if the underlying topology contains + bond information. + + + .. versionadded:: 0.20.0 + """ + return self.universe._fragdict[self.ix].ix + + @property + def fragment(self): + """An :class:`~MDAnalysis.core.groups.AtomGroup` representing the + fragment this :class:`~MDAnalysis.core.groups.Atom` is part of. + + A fragment is a + :class:`self of atoms` which are + interconnected by :class:`~MDAnalysis.core.topologyattrs.Bonds`, i.e., + there exists a path along one + or more :class:`~MDAnalysis.core.topologyattrs.Bonds` between any pair + of :class:`Atoms` + within a fragment. Thus, a fragment typically corresponds to a molecule. + + Note + ---- + This property is only accessible if the underlying topology contains + bond information. + + + .. versionadded:: 0.9.0 + """ + return self.universe._fragdict[self.ix].fragment + class Residue(ComponentBase): """:class:`Residue` base class. @@ -3756,6 +4381,93 @@ def segment(self, new): "".format(type(new))) self.universe._topology.tt.move_residue(self.ix, new.segindex) + def phi_selection(self): + """AtomGroup corresponding to the phi protein backbone dihedral + C'-N-CA-C. + + Returns + ------- + AtomGroup + 4-atom selection in the correct order. If no C' found in the + previous residue (by resid) then this method returns ``None``. + """ + # TODO: maybe this can be reformulated into one selection string without + # the additions later + sel_str = "segid {} and resid {} and name C".format( + self.segment.segid, self.resid - 1) + sel = (self.universe.select_atoms(sel_str) + + self.atoms.select_atoms('name N', 'name CA', 'name C')) + + # select_atoms doesnt raise errors if nothing found, so check size + if len(sel) == 4: + return sel + else: + return None + + def psi_selection(self): + """AtomGroup corresponding to the psi protein backbone dihedral + N-CA-C-N'. + + Returns + ------- + AtomGroup + 4-atom selection in the correct order. If no N' found in the + following residue (by resid) then this method returns ``None``. + """ + sel_str = "segid {} and resid {} and name N".format( + self.segment.segid, self.resid + 1) + + sel = (self.atoms.select_atoms('name N', 'name CA', 'name C') + + self.universe.select_atoms(sel_str)) + + if len(sel) == 4: + return sel + else: + return None + + def omega_selection(self): + """AtomGroup corresponding to the omega protein backbone dihedral + CA-C-N'-CA'. + + omega describes the -C-N- peptide bond. Typically, it is trans (180 + degrees) although cis-bonds (0 degrees) are also occasionally observed + (especially near Proline). + + Returns + ------- + AtomGroup + 4-atom selection in the correct order. If no C' found in the + previous residue (by resid) then this method returns ``None``. + + """ + nextres = self.resid + 1 + segid = self.segment.segid + sel = (self.atoms.select_atoms('name CA', 'name C') + + self.universe.select_atoms( + 'segid {} and resid {} and name N'.format(segid, nextres), + 'segid {} and resid {} and name CA'.format(segid, nextres))) + if len(sel) == 4: + return sel + else: + return None + + def chi1_selection(self): + """AtomGroup corresponding to the chi1 sidechain dihedral N-CA-CB-CG. + + Returns + ------- + AtomGroup + 4-atom selection in the correct order. If no CB and/or CG is found + then this method returns ``None``. + + .. versionadded:: 0.7.5 + """ + ag = self.atoms.select_atoms('name N', 'name CA', 'name CB', 'name CG') + if len(ag) == 4: + return ag + else: + return None + class Segment(ComponentBase): """:class:`Segment` base class. diff --git a/package/MDAnalysis/core/topologyattrs.py b/package/MDAnalysis/core/topologyattrs.py index 60b49998df7..6de0287f568 100644 --- a/package/MDAnalysis/core/topologyattrs.py +++ b/package/MDAnalysis/core/topologyattrs.py @@ -35,9 +35,6 @@ import six from six.moves import zip, range -import Bio.Seq -import Bio.SeqRecord -import Bio.Alphabet from collections import defaultdict import copy import functools @@ -46,19 +43,13 @@ import numpy as np import warnings -from numpy.lib.utils import deprecate - -from . import flags -from ..lib.util import (cached, convert_aa_code, iterable, warn_if_not_unique, - unique_int_1d) -from ..lib import transformations, mdamath -from ..exceptions import NoDataError, SelectionError +from ..lib.util import (cached, iterable) +from ..exceptions import NoDataError from .topologyobjects import TopologyGroup from . import selection -from .groups import (ComponentBase, GroupBase, +from .groups import (ComponentBase, Atom, Residue, Segment, - AtomGroup, ResidueGroup, SegmentGroup, - check_pbc_and_unwrap) + AtomGroup, ResidueGroup, SegmentGroup) from .. import _TOPOLOGY_ATTRS @@ -355,7 +346,7 @@ class Resindices(TopologyAttr): Resetting these values changes the residue membership of the atoms. If the group is a ResidueGroup or SegmentGroup, then this gives the - resindices of each residue represented in the group in a 1-D array, in the + resindices of each residue represented in the self in a 1-D array, in the order of the elements in that group. """ @@ -540,102 +531,6 @@ def _get_named_atom(group, name): transplants[Residue].append( ('_get_named_atom', _get_named_atom)) - def phi_selection(residue): - """AtomGroup corresponding to the phi protein backbone dihedral - C'-N-CA-C. - - Returns - ------- - AtomGroup - 4-atom selection in the correct order. If no C' found in the - previous residue (by resid) then this method returns ``None``. - """ - # TODO: maybe this can be reformulated into one selection string without - # the additions later - sel_str = "segid {} and resid {} and name C".format( - residue.segment.segid, residue.resid - 1) - sel = (residue.universe.select_atoms(sel_str) + - residue.atoms.select_atoms('name N', 'name CA', 'name C')) - - # select_atoms doesnt raise errors if nothing found, so check size - if len(sel) == 4: - return sel - else: - return None - - transplants[Residue].append(('phi_selection', phi_selection)) - - def psi_selection(residue): - """AtomGroup corresponding to the psi protein backbone dihedral - N-CA-C-N'. - - Returns - ------- - AtomGroup - 4-atom selection in the correct order. If no N' found in the - following residue (by resid) then this method returns ``None``. - """ - sel_str = "segid {} and resid {} and name N".format( - residue.segment.segid, residue.resid + 1) - - sel = (residue.atoms.select_atoms('name N', 'name CA', 'name C') + - residue.universe.select_atoms(sel_str)) - - if len(sel) == 4: - return sel - else: - return None - - transplants[Residue].append(('psi_selection', psi_selection)) - - def omega_selection(residue): - """AtomGroup corresponding to the omega protein backbone dihedral - CA-C-N'-CA'. - - omega describes the -C-N- peptide bond. Typically, it is trans (180 - degrees) although cis-bonds (0 degrees) are also occasionally observed - (especially near Proline). - - Returns - ------- - AtomGroup - 4-atom selection in the correct order. If no C' found in the - previous residue (by resid) then this method returns ``None``. - - """ - nextres = residue.resid + 1 - segid = residue.segment.segid - sel = (residue.atoms.select_atoms('name CA', 'name C') + - residue.universe.select_atoms( - 'segid {} and resid {} and name N'.format(segid, nextres), - 'segid {} and resid {} and name CA'.format(segid, nextres))) - if len(sel) == 4: - return sel - else: - return None - - transplants[Residue].append(('omega_selection', omega_selection)) - - def chi1_selection(residue): - """AtomGroup corresponding to the chi1 sidechain dihedral N-CA-CB-CG. - - Returns - ------- - AtomGroup - 4-atom selection in the correct order. If no CB and/or CG is found - then this method returns ``None``. - - .. versionadded:: 0.7.5 - """ - ag = residue.atoms.select_atoms('name N', 'name CA', - 'name CB', 'name CG') - if len(ag) == 4: - return ag - else: - return None - - transplants[Residue].append(('chi1_selection', chi1_selection)) - # TODO: update docs to property doc class Atomtypes(AtomAttr): @@ -771,408 +666,6 @@ def get_segments(self, sg): return masses - @warn_if_not_unique - @check_pbc_and_unwrap - def center_of_mass(group, pbc=None, compound='group', unwrap=False): - """Center of mass of (compounds of) the group. - - Computes the center of mass of :class:`Atoms` in the group. - Centers of mass per :class:`Residue`, :class:`Segment`, molecule, or - fragment can be obtained by setting the `compound` parameter - accordingly. If the masses of a compound sum up to zero, the - center of mass coordinates of that compound will be ``nan`` (not a - number). - - Parameters - ---------- - pbc : bool, optional - If ``True`` and `compound` is ``'group'``, move all atoms to the - primary unit cell before calculation. - If ``True`` and `compound` is ``'segments'`` or ``'residues'``, the - centers of mass of each compound will be calculated without moving - any atoms to keep the compounds intact. Instead, the resulting - center-of-mass position vectors will be moved to the primary unit - cell after calculation. - compound : {'group', 'segments', 'residues', 'molecules', 'fragments'},\ - optional - If ``'group'``, the center of mass of all atoms in the group will - be returned as a single position vector. Otherwise, the centers of - mass of each :class:`Segment`, :class:`Residue`, molecule, or - fragment will be returned as an array of position vectors, i.e. a 2d - array. - Note that, in any case, *only* the positions of :class:`Atoms` - *belonging to the group* will be taken into account. - unwrap : bool, optional - If ``True``, compounds will be unwrapped before computing their centers. - - Returns - ------- - center : numpy.ndarray - Position vector(s) of the center(s) of mass of the group. - If `compound` was set to ``'group'``, the output will be a single - position vector. - If `compound` was set to ``'segments'`` or ``'residues'``, the - output will be a 2d coordinate array of shape ``(n, 3)`` where ``n`` - is the number of compounds. - - Note - ---- - * This method can only be accessed if the underlying topology has - information about atomic masses. - * The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to - ``True`` allows the *pbc* flag to be used by default. - - - .. versionchanged:: 0.8 Added `pbc` parameter - .. versionchanged:: 0.19.0 Added `compound` parameter - .. versionchanged:: 0.20.0 Added ``'molecules'`` and ``'fragments'`` - compounds - .. versionchanged:: 0.20.0 Added `unwrap` parameter - """ - atoms = group.atoms - return atoms.center(weights=atoms.masses, pbc=pbc, compound=compound, unwrap=unwrap) - - transplants[GroupBase].append( - ('center_of_mass', center_of_mass)) - - @warn_if_not_unique - def total_mass(group, compound='group'): - """Total mass of (compounds of) the group. - - Computes the total mass of :class:`Atoms` in the group. - Total masses per :class:`Residue`, :class:`Segment`, molecule, or - fragment can be obtained by setting the `compound` parameter - accordingly. - - Parameters - ---------- - compound : {'group', 'segments', 'residues', 'molecules', 'fragments'},\ - optional - If ``'group'``, the total mass of all atoms in the group will be - returned as a single value. Otherwise, the total masses per - :class:`Segment`, :class:`Residue`, molecule, or fragment will be - returned as a 1d array. - Note that, in any case, *only* the masses of :class:`Atoms` - *belonging to the group* will be taken into account. - - Returns - ------- - float or numpy.ndarray - Total mass of (compounds of) the group. - If `compound` was set to ``'group'``, the output will be a single - value. Otherwise, the output will be a 1d array of shape ``(n,)`` - where ``n`` is the number of compounds. - - - .. versionchanged:: 0.20.0 Added `compound` parameter - """ - return group.accumulate("masses", compound=compound) - - transplants[GroupBase].append( - ('total_mass', total_mass)) - - @warn_if_not_unique - @check_pbc_and_unwrap - def moment_of_inertia(group, **kwargs): - """Tensor moment of inertia relative to center of mass as 3x3 numpy - array. - - Parameters - ---------- - pbc : bool, optional - If ``True``, move all atoms within the primary unit cell before - calculation. [``False``] - - Note - ---- - The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to - ``True`` allows the *pbc* flag to be used by default. - - - .. versionchanged:: 0.8 Added *pbc* keyword - .. versionchanged:: 0.20.0 Added `unwrap` parameter - - """ - atomgroup = group.atoms - pbc = kwargs.pop('pbc', flags['use_pbc']) - unwrap = kwargs.pop('unwrap', False) - compound = kwargs.pop('compound', 'group') - - com = atomgroup.center_of_mass(pbc=pbc, unwrap=unwrap, compound=compound) - if compound is not 'group': - com = (com * group.masses[:, None]).sum(axis=0) / group.masses.sum() - - if pbc: - pos = atomgroup.pack_into_box(inplace=False) - com - elif unwrap: - pos = atomgroup.unwrap(compound=compound, inplace=False) - com - else: - pos = atomgroup.positions - com - - masses = atomgroup.masses - # Create the inertia tensor - # m_i = mass of atom i - # (x_i, y_i, z_i) = pos of atom i - # Ixx = sum(m_i*(y_i^2+z_i^2)); - # Iyy = sum(m_i*(x_i^2+z_i^2)); - # Izz = sum(m_i*(x_i^2+y_i^2)) - # Ixy = Iyx = -1*sum(m_i*x_i*y_i) - # Ixz = Izx = -1*sum(m_i*x_i*z_i) - # Iyz = Izy = -1*sum(m_i*y_i*z_i) - tens = np.zeros((3, 3), dtype=np.float64) - # xx - tens[0][0] = (masses * (pos[:, 1] ** 2 + pos[:, 2] ** 2)).sum() - # xy & yx - tens[0][1] = tens[1][0] = - (masses * pos[:, 0] * pos[:, 1]).sum() - # xz & zx - tens[0][2] = tens[2][0] = - (masses * pos[:, 0] * pos[:, 2]).sum() - # yy - tens[1][1] = (masses * (pos[:, 0] ** 2 + pos[:, 2] ** 2)).sum() - # yz + zy - tens[1][2] = tens[2][1] = - (masses * pos[:, 1] * pos[:, 2]).sum() - # zz - tens[2][2] = (masses * (pos[:, 0] ** 2 + pos[:, 1] ** 2)).sum() - - return tens - - transplants[GroupBase].append( - ('moment_of_inertia', moment_of_inertia)) - - @warn_if_not_unique - def radius_of_gyration(group, **kwargs): - """Radius of gyration. - - Parameters - ---------- - pbc : bool, optional - If ``True``, move all atoms within the primary unit cell before - calculation. [``False``] - - Note - ---- - The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to - ``True`` allows the *pbc* flag to be used by default. - - - .. versionchanged:: 0.8 Added *pbc* keyword - - """ - atomgroup = group.atoms - pbc = kwargs.pop('pbc', flags['use_pbc']) - masses = atomgroup.masses - - com = atomgroup.center_of_mass(pbc=pbc) - if pbc: - recenteredpos = atomgroup.pack_into_box(inplace=False) - com - else: - recenteredpos = atomgroup.positions - com - - rog_sq = np.sum(masses * np.sum(recenteredpos**2, - axis=1)) / atomgroup.total_mass() - - return np.sqrt(rog_sq) - - transplants[GroupBase].append( - ('radius_of_gyration', radius_of_gyration)) - - @warn_if_not_unique - def shape_parameter(group, **kwargs): - """Shape parameter. - - See [Dima2004a]_ for background information. - - Parameters - ---------- - pbc : bool, optional - If ``True``, move all atoms within the primary unit cell before - calculation. [``False``] - - Note - ---- - The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to - ``True`` allows the *pbc* flag to be used by default. - - - References - ---------- - .. [Dima2004a] Dima, R. I., & Thirumalai, D. (2004). Asymmetry - in the shapes of folded and denatured states of - proteins. *J Phys Chem B*, 108(21), - 6564-6570. doi:`10.1021/jp037128y - `_ - - - .. versionadded:: 0.7.7 - .. versionchanged:: 0.8 Added *pbc* keyword - - """ - atomgroup = group.atoms - pbc = kwargs.pop('pbc', flags['use_pbc']) - masses = atomgroup.masses - - com = atomgroup.center_of_mass(pbc=pbc) - if pbc: - recenteredpos = atomgroup.pack_into_box(inplace=False) - com - else: - recenteredpos = atomgroup.positions - com - tensor = np.zeros((3, 3)) - - for x in range(recenteredpos.shape[0]): - tensor += masses[x] * np.outer(recenteredpos[x, :], - recenteredpos[x, :]) - tensor /= atomgroup.total_mass() - eig_vals = np.linalg.eigvalsh(tensor) - shape = 27.0 * np.prod(eig_vals - np.mean(eig_vals)) / np.power(np.sum(eig_vals), 3) - - return shape - - transplants[GroupBase].append( - ('shape_parameter', shape_parameter)) - - @warn_if_not_unique - @check_pbc_and_unwrap - def asphericity(group, pbc=None, unwrap=None, compound='group'): - """Asphericity. - - See [Dima2004b]_ for background information. - - Parameters - ---------- - pbc : bool, optional - If ``True``, move all atoms within the primary unit cell before - calculation. If ``None`` use value defined in - MDAnalysis.core.flags['use_pbc'] - unwrap : bool, optional - If ``True``, compounds will be unwrapped before computing their centers. - compound : {'group', 'segments', 'residues', 'molecules', 'fragments'}, optional - Which type of component to keep together during unwrapping. - - Note - ---- - The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to - ``True`` allows the *pbc* flag to be used by default. - - - References - ---------- - - .. [Dima2004b] Dima, R. I., & Thirumalai, D. (2004). Asymmetry - in the shapes of folded and denatured states of - proteins. *J Phys Chem B*, 108(21), - 6564-6570. doi:`10.1021/jp037128y - `_ - - - - .. versionadded:: 0.7.7 - .. versionchanged:: 0.8 Added *pbc* keyword - .. versionchanged:: 0.20.0 Added *unwrap* and *compound* parameter - - """ - atomgroup = group.atoms - if pbc is None: - pbc = flags['use_pbc'] - masses = atomgroup.masses - - com = atomgroup.center_of_mass(pbc=pbc, unwrap=unwrap, compound=compound) - if compound is not 'group': - com = (com * group.masses[:, None]).sum(axis=0) / group.masses.sum() - - if pbc: - recenteredpos = (atomgroup.pack_into_box(inplace=False) - com) - elif unwrap: - recenteredpos = (atomgroup.unwrap(inplace=False) - com) - else: - recenteredpos = (atomgroup.positions - com) - - tensor = np.zeros((3, 3)) - for x in range(recenteredpos.shape[0]): - tensor += masses[x] * np.outer(recenteredpos[x], - recenteredpos[x]) - - tensor /= atomgroup.total_mass() - eig_vals = np.linalg.eigvalsh(tensor) - shape = (3.0 / 2.0) * (np.sum((eig_vals - np.mean(eig_vals))**2) / - np.sum(eig_vals)**2) - - return shape - - transplants[GroupBase].append( - ('asphericity', asphericity)) - - @warn_if_not_unique - def principal_axes(group, pbc=None): - """Calculate the principal axes from the moment of inertia. - - e1,e2,e3 = AtomGroup.principal_axes() - - The eigenvectors are sorted by eigenvalue, i.e. the first one - corresponds to the highest eigenvalue and is thus the first principal - axes. - - Parameters - ---------- - pbc : bool, optional - If ``True``, move all atoms within the primary unit cell before - calculation. If ``None`` use value defined in setup flags. - - Returns - ------- - axis_vectors : array - 3 x 3 array with ``v[0]`` as first, ``v[1]`` as second, and - ``v[2]`` as third eigenvector. - - Note - ---- - The :class:`MDAnalysis.core.flags` flag *use_pbc* when set to - ``True`` allows the *pbc* flag to be used by default. - - - .. versionchanged:: 0.8 Added *pbc* keyword - - """ - atomgroup = group.atoms - if pbc is None: - pbc = flags['use_pbc'] - e_val, e_vec = np.linalg.eig(atomgroup.moment_of_inertia(pbc=pbc)) - - # Sort - indices = np.argsort(e_val)[::-1] - # Return transposed in more logical form. See Issue 33. - return e_vec[:, indices].T - - transplants[GroupBase].append( - ('principal_axes', principal_axes)) - - def align_principal_axis(group, axis, vector): - """Align principal axis with index `axis` with `vector`. - - Parameters - ---------- - axis : {0, 1, 2} - Index of the principal axis (0, 1, or 2), as produced by - :meth:`~principal_axes`. - vector : array_like - Vector to align principal axis with. - - Notes - ----- - To align the long axis of a channel (the first principal axis, i.e. - *axis* = 0) with the z-axis:: - - u.atoms.align_principal_axis(0, [0,0,1]) - u.atoms.write("aligned.pdb") - """ - p = group.principal_axes()[axis] - angle = np.degrees(mdamath.angle(p, vector)) - ax = transformations.rotaxis(p, vector) - # print "principal[%d] = %r" % (axis, p) - # print "axis = %r, angle = %f deg" % (ax, angle) - return group.rotateby(angle, ax) - - transplants[GroupBase].append( - ('align_principal_axis', align_principal_axis)) - # TODO: update docs to property doc class Charges(AtomAttr): @@ -1212,42 +705,6 @@ def get_segments(self, sg): return charges - @warn_if_not_unique - def total_charge(group, compound='group'): - """Total charge of (compounds of) the group. - - Computes the total charge of :class:`Atoms` in the group. - Total charges per :class:`Residue`, :class:`Segment`, molecule, or - fragment can be obtained by setting the `compound` parameter - accordingly. - - Parameters - ---------- - compound : {'group', 'segments', 'residues', 'molecules', 'fragments'},\ - optional - If 'group', the total charge of all atoms in the group will - be returned as a single value. Otherwise, the total charges per - :class:`Segment`, :class:`Residue`, molecule, or fragment - will be returned as a 1d array. - Note that, in any case, *only* the charges of :class:`Atoms` - *belonging to the group* will be taken into account. - - Returns - ------- - float or numpy.ndarray - Total charge of (compounds of) the group. - If `compound` was set to ``'group'``, the output will be a single - value. Otherwise, the output will be a 1d array of shape ``(n,)`` - where ``n`` is the number of compounds. - - - .. versionchanged:: 0.20.0 Added `compound` parameter - """ - return group.accumulate("charges", compound=compound) - - transplants[GroupBase].append( - ('total_charge', total_charge)) - # TODO: update docs to property doc class Bfactors(AtomAttr): @@ -1364,7 +821,7 @@ def _get_named_residue(group, resname): For more than one residue it returns a :class:`MDAnalysis.core.groups.ResidueGroup` instance. A single - :class:`MDAnalysis.core.group.Residue` is returned for a single match. + :class:`MDAnalysis.core.groups.Residue` is returned for a single match. If no residues are found, a :exc:`SelectionError` is raised. .. versionadded:: 0.9.2 @@ -1402,102 +859,6 @@ def _get_named_residue(group, resname): transplants[ResidueGroup].append( ('_get_named_residue', _get_named_residue)) - def sequence(self, **kwargs): - """Returns the amino acid sequence. - - The format of the sequence is selected with the keyword *format*: - - ============== ============================================ - *format* description - ============== ============================================ - 'SeqRecord' :class:`Bio.SeqRecord.SeqRecord` (default) - 'Seq' :class:`Bio.Seq.Seq` - 'string' string - ============== ============================================ - - The sequence is returned by default (keyword ``format = 'SeqRecord'``) - as a :class:`Bio.SeqRecord.SeqRecord` instance, which can then be - further processed. In this case, all keyword arguments (such as the - *id* string or the *name* or the *description*) are directly passed to - :class:`Bio.SeqRecord.SeqRecord`. - - If the keyword *format* is set to ``'Seq'``, all *kwargs* are ignored - and a :class:`Bio.Seq.Seq` instance is returned. The difference to the - record is that the record also contains metadata and can be directly - used as an input for other functions in :mod:`Bio`. - - If the keyword *format* is set to ``'string'``, all *kwargs* are - ignored and a Python string is returned. - - .. rubric:: Example: Write FASTA file - - Use :func:`Bio.SeqIO.write`, which takes sequence records:: - - import Bio.SeqIO - - # get the sequence record of a protein component of a Universe - protein = u.select_atoms("protein") - record = protein.sequence(id="myseq1", name="myprotein") - - Bio.SeqIO.write(record, "single.fasta", "fasta") - - A FASTA file with multiple entries can be written with :: - - Bio.SeqIO.write([record1, record2, ...], "multi.fasta", "fasta") - - Parameters - ---------- - format : string, optional - - ``"string"``: return sequence as a string of 1-letter codes - - ``"Seq"``: return a :class:`Bio.Seq.Seq` instance - - ``"SeqRecord"``: return a :class:`Bio.SeqRecord.SeqRecord` - instance - - Default is ``"SeqRecord"`` - id : optional - Sequence ID for SeqRecord (should be different for different - sequences) - name : optional - Name of the protein. - description : optional - Short description of the sequence. - kwargs : optional - Any other keyword arguments that are understood by - class:`Bio.SeqRecord.SeqRecord`. - - Raises - ------ - :exc:`ValueError` if a residue name cannot be converted to a - 1-letter IUPAC protein amino acid code; make sure to only - select protein residues. - - :exc:`TypeError` if an unknown *format* is selected. - - - .. versionadded:: 0.9.0 - """ - formats = ('string', 'Seq', 'SeqRecord') - - format = kwargs.pop("format", "SeqRecord") - if format not in formats: - raise TypeError("Unknown format='{0}': must be one of: {1}".format( - format, ", ".join(formats))) - try: - sequence = "".join([convert_aa_code(r) for r in self.residues.resnames]) - except KeyError as err: - six.raise_from(ValueError("AtomGroup contains a residue name '{0}' that " - "does not have a IUPAC protein 1-letter " - "character".format(err.message)), None) - if format == "string": - return sequence - seq = Bio.Seq.Seq(sequence, alphabet=Bio.Alphabet.IUPAC.protein) - if format == "Seq": - return seq - return Bio.SeqRecord.SeqRecord(seq, **kwargs) - - transplants[ResidueGroup].append( - ('sequence', sequence)) - # TODO: update docs to property doc class Resnums(ResidueAttr): @@ -1606,7 +967,7 @@ def _get_named_segment(group, segid): For more than one residue it returns a :class:`MDAnalysis.core.groups.SegmentGroup` instance. A single - :class:`MDAnalysis.core.group.Segment` is returned for a single match. + :class:`MDAnalysis.core.groups.Segment` is returned for a single match. If no residues are found, a :exc:`SelectionError` is raised. .. versionadded:: 0.9.2 @@ -1758,141 +1119,6 @@ class Bonds(_Connection): transplants = defaultdict(list) _n_atoms = 2 - def bonded_atoms(self): - """An :class:`~MDAnalysis.core.groups.AtomGroup` of all - :class:`Atoms` bonded to this - :class:`~MDAnalysis.core.groups.Atom`.""" - idx = [b.partner(self).index for b in self.bonds] - return self.universe.atoms[idx] - - transplants[Atom].append( - ('bonded_atoms', property(bonded_atoms, None, None, - bonded_atoms.__doc__))) - - def fragindex(self): - """The index (ID) of the - :class:`~MDAnalysis.core.topologyattrs.Bonds.fragment` this - :class:`~MDAnalysis.core.groups.Atom` is part of. - - Note - ---- - This property is only accessible if the underlying topology contains - bond information. - - - .. versionadded:: 0.20.0 - """ - return self.universe._fragdict[self.ix].ix - - def fragindices(self): - """The - :class:`fragment indices` - of all :class:`Atoms` in this - :class:`~MDAnalysis.core.groups.AtomGroup`. - - A :class:`numpy.ndarray` with - :attr:`~numpy.ndarray.shape`\ ``=(``\ :attr:`~AtomGroup.n_atoms`\ ``,)`` - and :attr:`~numpy.ndarray.dtype`\ ``=numpy.int64``. - - Note - ---- - This property is only accessible if the underlying topology contains - bond information. - - - .. versionadded:: 0.20.0 - """ - fragdict = self.universe._fragdict - return np.array([fragdict[aix].ix for aix in self.ix], dtype=np.int64) - - def fragment(self): - """An :class:`~MDAnalysis.core.groups.AtomGroup` representing the - fragment this :class:`~MDAnalysis.core.groups.Atom` is part of. - - A fragment is a - :class:`group of atoms` which are - interconnected by :class:`~MDAnalysis.core.topologyattrs.Bonds`, i.e., - there exists a path along one - or more :class:`~MDAnalysis.core.topologyattrs.Bonds` between any pair - of :class:`Atoms` - within a fragment. Thus, a fragment typically corresponds to a molecule. - - Note - ---- - This property is only accessible if the underlying topology contains - bond information. - - - .. versionadded:: 0.9.0 - """ - return self.universe._fragdict[self.ix].fragment - - def fragments(self): - """Read-only :class:`tuple` of - :class:`fragments`. - - Contains all fragments that - any :class:`~MDAnalysis.core.groups.Atom` in this - :class:`~MDAnalysis.core.groups.AtomGroup` is part of. - - A fragment is a - :class:`group of atoms` which are - interconnected by :class:`~MDAnalysis.core.topologyattrs.Bonds`, i.e., - there exists a path along one - or more :class:`~MDAnalysis.core.topologyattrs.Bonds` between any pair - of :class:`Atoms` - within a fragment. Thus, a fragment typically corresponds to a molecule. - - Note - ---- - * This property is only accessible if the underlying topology contains - bond information. - * The contents of the fragments may extend beyond the contents of this - :class:`~MDAnalysis.core.groups.AtomGroup`. - - - .. versionadded:: 0.9.0 - """ - fragdict = self.universe._fragdict - return tuple(sorted(set(fragdict[aix].fragment for aix in self.ix), - key=lambda x: x[0].ix)) - - def n_fragments(self): - """The number of unique - :class:`~MDAnalysis.core.topologyattrs.Bonds.fragments` the - :class:`Atoms` of this - :class:`~MDAnalysis.core.groups.AtomGroup` are part of. - - Note - ---- - This property is only accessible if the underlying topology contains - bond information. - - - .. versionadded:: 0.20.0 - """ - return len(unique_int_1d(self.fragindices)) - - transplants[Atom].append( - ('fragment', property(fragment, None, None, - fragment.__doc__))) - - transplants[Atom].append( - ('fragindex', property(fragindex, None, None, - fragindex.__doc__))) - - transplants[AtomGroup].append( - ('fragments', property(fragments, None, None, - fragments.__doc__))) - - transplants[AtomGroup].append( - ('fragindices', property(fragindices, None, None, - fragindices.__doc__))) - - transplants[AtomGroup].append( - ('n_fragments', property(n_fragments, None, None, - n_fragments.__doc__))) - class Angles(_Connection): """Angles between three atoms