9.1.2. Core objects: Containers — MDAnalysis.core.groups
¶
The Universe
instance contains all
the particles in the system (which MDAnalysis calls
Atom
). Groups of atoms are handled as AtomGroup
instances. The AtomGroup
is probably the most important
object in MDAnalysis because virtually everything can be accessed
through it. AtomGroup instances can be easily created (e.g., from a
AtomGroup.select_atoms()
selection or just by slicing).
For convenience, chemically meaningful groups of atoms such as a
Residue
or a Segment
(typically a whole molecule or
all of the solvent) also exist as containers, as well as groups of
these units ((ResidueGroup
, SegmentGroup
).
9.1.2.1. Classes¶
9.1.2.1.1. Collections¶
-
class
MDAnalysis.core.groups.
AtomGroup
(*args)[source]¶ A group of atoms.
An
AtomGroup
is an ordered collection of atoms. Typically, anAtomGroup
is generated from a selection, or by indexing/slicing theAtomGroup
of all atoms in theUniverse
atMDAnalysis.core.universe.Universe.atoms
.An AtomGroup can be indexed and sliced like a list:
ag[0], ag[-1]
will return the first and the last
Atom
in the group whereas the sliceag[0:6:2]
returns an AtomGroup of every second element, corresponding to indices 0, 2, and 4.
It also supports “advanced slicing” when the argument is a
numpy.ndarray
or alist
:aslice = [0, 3, -1, 10, 3] ag[aslice]
will return a new AtomGroup of atoms with those indices in the old AtomGroup.
Note
AtomGroups originating from a selection are sorted and duplicate elements are removed. This is not true for AtomGroups produced by slicing. Thus slicing can be used when the order of atoms is crucial (for instance, in order to define angles or dihedrals).
AtomGroups can be compared and combined using group operators. For instance, AtomGroups can be concatenated using + or
concatenate()
:ag_concat = ag1 + ag2 # or ag_concat = ag1.concatenate(ag2)
When groups are concatenated, the order of the atoms is conserved. If atoms appear several times in one of the groups, the duplicates are kept in the resulting group. On the contrary to
concatenate()
,union()
treats the AtomGroups as sets, duplicates are removed from the resulting group, and atoms are ordered. The | operator is synomymous tounion()
:ag_union = ag1 | ag2 # or ag_union = ag1.union(ag2)
The opposite operation to
concatenate()
issubtract()
. This method creates a new group with all the atoms of the group that are not in a given other group; the order of the atoms is kept, so as duplicates.difference()
is the set version ofsubtract()
. The resulting group is sorted and deduplicated.All set methods are listed in the table below. These methods treat the groups as sorted and deduplicated sets of atoms.
Operation Equivalent Result s.isdisjoint(t)
True
ifs
andt
do not share elementss.issubset(t)
test if all elements of s
are part oft
s.is_strict_subset(t)
test if all elements of s
are part oft
, ands != t
s.issuperset(t)
test if all elements of t
are part ofs
s.is_strict_superset(t)
test if all elements of t
are part ofs
, ands != t
s.union(t)
s | t
new Group with elements from both s
andt
s.intersection(t)
s & t
new Group with elements common to s
andt
s.difference(t)
s - t
new Group with elements of s
that are not int
s.symmetric_difference(t)
s ^ t
new Group with elements that are part of s
ort
but not bothThe following methods keep the order of the atoms, and keep duplicated atoms.
Operation Equivalent Result len(s)
number of elements (atoms, residues or segment) in the group s == t
test if s
andt
contain the same elements in the same orders.concatenate(t)
s + t
new Group with elements from s
and fromt
s.subtract(t)
new Group with elements from s
that are not int
The in operator allows to test if an
Atom
is in the AtomGroup.AtomGroup instances are always bound to a
MDAnalysis.core.universe.Universe
. They cannot exist in isolation.Deprecated functionality
Instant selectors will be removed in the 1.0 release. See issue #1377 for more details.
Atoms can also be accessed in a Pythonic fashion by using the atom name as an attribute. For instance,
ag.CA
will provide a
AtomGroup
of all CA atoms in the group. These instant selector attributes are auto-generated for each atom name encountered in the group.Notes
The name-attribute instant selector access to atoms is mainly meant for quick interactive work. Thus it either returns a single
Atom
if there is only one matching atom, or a newAtomGroup
for multiple matches. This makes it difficult to use the feature consistently in scripts.See also
Deprecated since version 0.16.2: Instant selectors of AtomGroup will be removed in the 1.0 release. See Instant selectors for details and alternatives.
-
angle
¶ This AtomGroup represented as an Angle object
Returns: Return type: A MDAnalysis.core.topologyobjects.Angle
objectRaises: ValueError if the AtomGroup is not length 3 New in version 0.11.0.
-
atoms
¶ Get another AtomGroup identical to this one.
-
bbox
(**kwargs)¶ Return the bounding box of the selection.
The lengths A,B,C of the orthorhombic enclosing box are
L = AtomGroup.bbox() A,B,C = L[1] - L[0]
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: corners – 2x3 array giving corners of bounding box as [[xmin, ymin, zmin], [xmax, ymax, zmax]]. Return type: array Note
The
MDAnalysis.core.flags
flag use_pbc when set toTrue
allows the pbc flag to be used by default.New in version 0.7.2.
Changed in version 0.8: Added pbc keyword
-
bond
¶ This AtomGroup represented as a Bond object
Returns: Return type: A MDAnalysis.core.topologyobjects.Bond
objectRaises: ValueError if the AtomGroup is not length 2 New in version 0.11.0.
-
bsphere
(**kwargs)¶ Return the bounding sphere of the selection.
The sphere is calculated relative to the centre of geometry.
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: - R (float) – Radius of bounding sphere.
- center (array) – Coordinates of sphere center as
[xcen,ycen,zcen]
.
Note
The
MDAnalysis.core.flags
flag use_pbc when set toTrue
allows the pbc flag to be used by default.New in version 0.7.3.
Changed in version 0.8: Added pbc keyword
-
center
(weights, pbc=None)¶ Calculate center of group given some weights
Parameters: - weights (array_like) – weights to be used
- pbc (boolean, optional) –
True
: Move all atoms within the primary unit cell before calculation [False
]
Returns: center – weighted center of group
Return type: ndarray
Examples
To find the charge weighted center of a given Atomgroup:
>>> sel = u.select_atoms('prop mass > 4.0') >>> sel.center(sel.charges)
Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.
-
center_of_geometry
(pbc=None)¶ Center of geometry (also known as centroid) of the selection.
Parameters: pbc (boolean, optional) – True
: Move all atoms within the primary unit cell before calculation [False
]Returns: center – geometric center of group Return type: ndarray Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.Changed in version 0.8: Added pbc keyword
-
centroid
(pbc=None)¶ Center of geometry (also known as centroid) of the selection.
Parameters: pbc (boolean, optional) – True
: Move all atoms within the primary unit cell before calculation [False
]Returns: center – geometric center of group Return type: ndarray Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.Changed in version 0.8: Added pbc keyword
-
concatenate
(other)¶ Concatenate with another Group or Component of the same level.
Duplicate entries and original order is preserved. It is synomymous to the + operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with elements of self and other concatenated Return type: Group Example
The order of the original contents (including duplicates) are preserved when performing a concatenation.
>>> ag1 = u.select_atoms('name O') >>> ag2 = u.select_atoms('name N') >>> ag3 = ag1 + ag2 # or ag1.concatenate(ag2) >>> ag3[:3].names array(['O', 'O', 'O'], dtype=object) >>> ag3[-3:].names array(['N', 'N', 'N'], dtype=object)
New in version 0.16.0.
-
difference
(other)¶ Elements from this Group that do not appear in another
This method removes duplicate elements and sorts the result. As such, it is different from
subtract()
.difference()
is synomymous to the - operator.Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements of self that are not in other, without duplicate elements Return type: Group See also
-
dihedral
¶ This AtomGroup represented as a Dihedral object
Returns: Return type: A MDAnalysis.core.topologyobjects.Dihedral
objectRaises: ValueError if the AtomGroup is not length 4 New in version 0.11.0.
-
forces
¶ Forces on each atom in the AtomGroup.
The velocities can be changed by assigning an array of the appropriate shape, i.e. either Nx3 to assign individual velocities or 3 to assign the same velocity to all atoms (e.g.
ag.velocity = array([0,0,0])
will give all particles zero velocity).
-
groupby
(topattr)¶ Group together items in this group according to values of topattr
Parameters: topattr (str) – Topology attribute to group components by. Returns: Unique values of the topology attribute as keys, Groups as values. Return type: dict Example
To group atoms with the same mass together:
>>> ag.groupby('masses') {12.010999999999999: <AtomGroup with 462 atoms>, 14.007: <AtomGroup with 116 atoms>, 15.999000000000001: <AtomGroup with 134 atoms>}
New in version 0.16.0.
-
guess_bonds
(vdwradii=None)[source]¶ Guess bonds that exist within this AtomGroup and add to Universe
Parameters: vdwradii (dict, optional) – Dict relating atom type: vdw radii New in version 0.10.0.
-
improper
¶ This AtomGroup represented as an ImproperDihedral object
Returns: Return type: A MDAnalysis.core.topologyobjects.ImproperDihedral
objectRaises: ValueError if the AtomGroup is not length 4 New in version 0.11.0.
-
intersection
(other)¶ Group of elements which are in both this Group and another
This method removes duplicate elements and sorts the result. It is synomymous to the & operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the common elements of self and other, without duplicate elements Return type: Group Example
Intersections can be used when the select atoms string would become too complicated. For example to find the water atoms which are within 4.0A of two segments:
>>> water = u.select_atoms('resname SOL') >>> shell1 = water.select_atoms('around 4.0 segid 1') >>> shell2 = water.select_atoms('around 4.0 segid 2') >>> common = shell1 & shell2 # or shell1.intersection(shell2)
See also
union()
,()
-
is_strict_subset
(other)¶ If this Group is a subset of another Group but not identical
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a strict subset of the other one
- .. versionadded:: 0.16
-
is_strict_superset
(other)¶ If this Group is a superset of another Group but not identical
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a strict superset of the other one
- .. versionadded:: 0.16
-
isdisjoint
(other)¶ If the Group has no elements in common with the other Group
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if the two Groups do not have common elements
- .. versionadded:: 0.16
-
issubset
(other)¶ If all elements of this Group are part of another Group
Note that an empty group is a subset of any group of the same level.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a subset of the other one
- .. versionadded:: 0.16
-
issuperset
(other)¶ If all elements of another Group are part of this Group
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a subset of the other one
- .. versionadded:: 0.16
-
ix
¶ Unique indices of the components in the Group.
- If this Group is an
AtomGroup
, these are the indices of theAtom
instances. - If it is a
ResidueGroup
, these are the indices of theResidue
instances. - If it is a
SegmentGroup
, these are the indices of theSegment
instances.
- If this Group is an
-
ix_array
¶ Unique indices of the components in the Group.
For a Group, ix_array is the same as ix. This method gives a consistent API between components and groups.
See also
-
n_atoms
¶ Number of atoms in AtomGroup.
Equivalent to
len(self)
.
-
n_residues
¶ Number of unique residues represented in the AtomGroup.
Equivalent to
len(self.residues)
.
-
n_segments
¶ Number of unique segments represented in the AtomGroup.
Equivalent to
len(self.segments)
.
-
pack_into_box
(box=None, inplace=True)¶ Shift all atoms in this group to be within the primary unit cell.
Parameters: - box (array_like) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
,[lx, ly, lz, alpha, beta, gamma]
. IfNone
, uses these timestep dimensions. - inplace (bool) –
True
to change coordinates in place.
Returns: coords – Shifted atom coordinates.
Return type: Notes
All atoms will be moved so that they lie between 0 and boxlength \(L_i\) in all dimensions, i.e. the lower left corner of the simulation box is taken to be at (0,0,0):
\[x_i' = x_i - \left\lfloor\frac{x_i}{L_i}\right\rfloor\]The default is to take unit cell information from the underlying
Timestep
instance. The optional argument box can be used to provide alternative unit cell information (in the MDAnalysis standard format[Lx, Ly, Lz, alpha, beta, gamma]
).Works with either orthogonal or triclinic box types.
New in version 0.8.
- box (array_like) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
-
positions
¶ Coordinates of the atoms in the AtomGroup.
The positions can be changed by assigning an array of the appropriate shape, i.e. either Nx3 to assign individual coordinates or 3, to assign the same coordinate to all atoms (e.g.
ag.positions = array([0,0,0])
will move all particles to the origin).Note
Changing the position is not reflected in any files; reading any frame from the trajectory will replace the change with that from the file except if the trajectory is held in memory, e.g., when the
transfer_to_memory
method was used.
-
residues
¶ Get sorted
ResidueGroup
of the (unique) residues represented in the AtomGroup.
-
rotate
(R, point=(0, 0, 0))¶ Apply a rotation matrix R to the selection’s coordinates. \(\mathsf{R}\) is a 3x3 orthogonal matrix that transforms a vector \(\mathbf{x} \rightarrow \mathbf{x}'\):
\[\mathbf{x}' = \mathsf{R}\mathbf{x}\]Parameters: - R (array_like) – 3x3 rotation matrix to use for applying rotation.
- point (array_like, optional) – Center of rotation
Returns: self
Return type: Notes
By default rotates around center of origin
point=(0, 0, 0)
. To rotate around center of geometry of the atomgroup useag.rotate(R, point=ag.centroid)
.See also
rotateby()
- rotate around given axis and angle
MDAnalysis.lib.transformations()
- module of all coordinate transforms
-
rotateby
(angle, axis, point=None)¶ Apply a rotation to the selection’s coordinates.
Parameters: - angle (float) – Rotation angle in degrees.
- axis (array_like) – Rotation axis vector.
- point (array_like, optional) – Center of rotation. If
None
then the center of geometry of this group is used.
Returns: self
Return type: Notes
The transformation from current coordinates \(\mathbf{x}\) to new coordinates \(\mathbf{x}'\) is
\[\mathbf{x}' = \mathsf{R}\,(\mathbf{x}-\mathbf{p}) + \mathbf{p}\]where \(\mathsf{R}\) is the rotation by angle around the axis going through point \(\mathbf{p}\).
See also
MDAnalysis.lib.transformations.rotation_matrix()
- calculate \(\mathsf{R}\)
-
segments
¶ Get sorted
SegmentGroup
of the (unique) segments represented in the AtomGroup.
-
select_atoms
(sel, *othersel, **selgroups)[source]¶ Select atoms using a selection string.
Returns an
AtomGroup
with atoms sorted according to their index in the topology (this is to ensure that there are not any duplicates, which can happen with complicated selections).Existing
AtomGroup
objects can be passed as named arguments, which will then be available to the selection parser.Subselections can be grouped with parentheses.
Selections can be set to update automatically on frame change, by setting the updating named argument to True.
Examples
>>> sel = universe.select_atoms("segid DMPC and not ( name H* or name O* )") >>> sel <AtomGroup with 3420 atoms>
>>> universe.select_atoms("around 10 group notHO", notHO=sel) <AtomGroup with 1250 atoms>
Notes
If exact ordering of atoms is required (for instance, for
angle()
ordihedral()
calculations) then one supplies selections separately in the required order. Also, when multipleAtomGroup
instances are concatenated with the+
operator then the order ofAtom
instances is preserved and duplicates are not removed.See also
Selection syntax
The selection parser understands the following CASE SENSITIVE keywords:
Simple selections
- protein, backbone, nucleic, nucleicbackbone
- selects all atoms that belong to a standard set of residues; a protein is identfied by a hard-coded set of residue names so it may not work for esoteric residues.
- segid seg-name
- select by segid (as given in the topology), e.g.
segid 4AKE
orsegid DMPC
- resid residue-number-range
- resid can take a single residue number or a range of numbers. A
range consists of two numbers separated by a colon (inclusive)
such as
resid 1:5
. A residue number (“resid”) is taken directly from the topology. If icodes are present in the topology, then these will be taken into account. Ie ‘resid 163B’ will only select resid 163 with icode B while ‘resid 163’ will select only residue 163. Range selections will also respect icodes, so ‘resid 162-163B’ will select all residues in 162 and those in 163 up to icode B. - resnum resnum-number-range
- resnum is the canonical residue number; typically it is set to the residue id in the original PDB structure.
- resname residue-name
- select by residue name, e.g.
resname LYS
- name atom-name
- select by atom name (as given in the topology). Often, this is
force field dependent. Example:
name CA
(for Cα atoms) orname OW
(for SPC water oxygen) - type atom-type
- select by atom type; this is either a string or a number and depends on the force field; it is read from the topology file (e.g. the CHARMM PSF file contains numeric atom types). It has non-sensical values when a PDB or GRO file is used as a topology
- atom seg-name residue-number atom-name
- a selector for a single atom consisting of segid resid atomname,
e.g.
DMPC 1 C2
selects the C2 carbon of the first residue of the DMPC segment - altloc alternative-location
- a selection for atoms where alternative locations are available, which is often the case with high-resolution crystal structures e.g. resid 4 and resname ALA and altloc B selects only the atoms of ALA-4 that have an altloc B record.
Boolean
- not
- all atoms not in the selection, e.g.
not protein
selects all atoms that aren’t part of a protein - and, or
- combine two selections according to the rules of boolean
algebra, e.g.
protein and not (resname ALA or resname LYS)
selects all atoms that belong to a protein, but are not in a lysine or alanine residue
Geometric
- around distance selection
- selects all atoms a certain cutoff away from another selection,
e.g.
around 3.5 protein
selects all atoms not belonging to protein that are within 3.5 Angstroms from the protein - point x y z distance
- selects all atoms within a cutoff of a point in space, make sure
coordinate is separated by spaces,
e.g.
point 5.0 5.0 5.0 3.5
selects all atoms within 3.5 Angstroms of the coordinate (5.0, 5.0, 5.0) - prop [abs] property operator value
- selects atoms based on position, using property x, y,
or z coordinate. Supports the abs keyword (for absolute
value) and the following operators: <, >, <=, >=, ==, !=.
For example,
prop z >= 5.0
selects all atoms with z coordinate greater than 5.0;prop abs z <= 5.0
selects all atoms within -5.0 <= z <= 5.0. - sphzone radius selection
- Selects all atoms that are within radius of the center of geometry of selection
- sphlayer inner radius outer radius selection
- Similar to sphzone, but also excludes atoms that are within inner radius of the selection COG
Connectivity
- byres selection
- selects all atoms that are in the same segment and residue as selection, e.g. specify the subselection after the byres keyword
- bonded selection
- selects all atoms that are bonded to selection
eg:
select name H bonded name O
selects only hydrogens bonded to oxygens
Index
- bynum index-range
- selects all atoms within a range of (1-based) inclusive indices,
e.g.
bynum 1
selects the first atom in the universe;bynum 5:10
selects atoms 5 through 10 inclusive. All atoms in theMDAnalysis.Universe
are consecutively numbered, and the index runs from 1 up to the total number of atoms.
Preexisting selections
- group group-name
- selects the atoms in the
AtomGroup
passed to the function as an argument named group-name. Only the atoms common to group-name and the instanceselect_atoms()
was called from will be considered, unlessgroup
is preceded by theglobal
keyword. group-name will be included in the parsing just by comparison of atom indices. This means that it is up to the user to make sure the group-name group was defined in an appropriateUniverse
. - global selection
- by default, when issuing
select_atoms()
from anAtomGroup
, selections and subselections are returned intersected with the atoms of that instance. Prefixing a selection term withglobal
causes its selection to be returned in its entirety. As an example, theglobal
keyword allows forlipids.select_atoms("around 10 global protein")
— wherelipids
is a group that does not contain any proteins. Wereglobal
absent, the result would be an empty selection since theprotein
subselection would itself be empty. When issuingselect_atoms()
from aUniverse
,global
is ignored.
- Dynamic selections
- If
select_atoms()
is invoked with named argument updating set to True, anUpdatingAtomGroup
instance will be returned, instead of a regularAtomGroup
. It behaves just like the latter, with the difference that the selection expressions are re-evaluated every time the trajectory frame changes (this happens lazily, only when theUpdatingAtomGroup
is accessed so that there is no redundant updating going on). Issuing an updating selection from an already updating group will cause later updates to also reflect the updating of the base group. A non-updating selection or a slicing operation made on anUpdatingAtomGroup
will return a staticAtomGroup
, which will no longer update across frames.
Changed in version 0.7.4: Added resnum selection.
Changed in version 0.8.1: Added group and fullgroup selections.
Deprecated since version 0.11: The use of
fullgroup
has been deprecated in favor of the equivalentglobal group
.Changed in version 0.13.0: Added bonded selection
Changed in version 0.16.0: Resid selection now takes icodes into account where present.
New in version 0.16.0: Updating selections now possible by setting the
updating
argument.
-
split
(level)[source]¶ Split AtomGroup into a list of atomgroups by level.
Parameters: level ({'atom', 'residue', 'segment'}) – New in version 0.9.0.
-
subtract
(other)¶ Group with elements from this Group that don’t appear in 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 are removed from the returned Group.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements of self that are not in other, conserves order and duplicates. Return type: Group Example
Unlike
difference()
this method will not sort or remove duplicates.>>> ag1 = u.atoms[[3, 3, 2, 2, 1, 1]] >>> ag2 = u.atoms[2] >>> ag3 = ag1 - ag2 # or ag1.subtract(ag2) >>> ag1.indices array([3, 3, 1, 1])
See also
-
symmetric_difference
(other)¶ Group of elements which are only in one of this Group or another
This method removes duplicate elements and the result is sorted. It is synomym to the ^ operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements that are in self or in other but not in both, without duplicate elements Return type: Group Example
>>> ag1 = u.atoms[[0, 1, 5, 3, 3, 2]] >>> ag2 = u.atoms[[4, 4, 6, 2, 3, 5]] >>> ag3 = ag1 ^ ag2 # or ag1.symmetric_difference(ag2) >>> ag3.indices # 0 and 1 are only in ag1, 4 and 6 are only in ag2 [0, 1, 4, 6]
See also
difference()
,()
-
transform
(M)¶ Apply homogenous transformation matrix M to the coordinates.
Parameters: M (array) – 4x4 matrix, with the rotation in R = M[:3,:3]
and the translation int = M[:3,3]
.Returns: Return type: self See also
MDAnalysis.lib.transformations()
- module of all coordinate transforms
Notes
The rotation \(\mathsf{R}\) is applied before the translation \(\mathbf{t}\):
\[\mathbf{x}' = \mathsf{R}\mathbf{x} + \mathbf{t}\]
-
translate
(t)¶ Apply translation vector t to the selection’s coordinates.
Atom coordinates are translated in-place.
Parameters: t (array_like) – vector to translate coordinates with Returns: Return type: self See also
MDAnalysis.lib.transformations()
- module of all coordinate transforms
Notes
The method applies a translation to the
AtomGroup
from current coordinates \(\mathbf{x}\) to new coordinates \(\mathbf{x}'\):\[\mathbf{x}' = \mathbf{x} + \mathbf{t}\]
-
ts
¶ Temporary Timestep that contains the selection coordinates.
A
Timestep
instance, which can be passed to a trajectory writer.If
ts
is modified then these modifications will be present until the frame number changes (which typically happens when the underlying trajectory frame changes).It is not possible to assign a new
Timestep
to theAtomGroup.ts
attribute; change attributes of the object.
-
union
(other)¶ Group of elements either in this Group or another
On the contrary to concatenation, this method sort the elements and removes duplicate ones. It is synomymous to the | operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the combined elements of self and other, without duplicate elements Return type: Group Example
In contrast to
concatenate()
, any duplicates are dropped and the result is sorted.>>> ag1 = u.select_atoms('name O') >>> ag2 = u.select_atoms('name N') >>> ag3 = ag1 | ag2 # or ag1.union(ag2) >>> ag3[:3].names array(['N', 'O', 'N'], dtype=object)
See also
-
unique
¶ Return an AtomGroup containing sorted and unique atoms only.
Examples
>>> ag = u.atoms[[2, 1, 2, 2, 1, 0]] >>> ag <AtomGroup with 6 atoms> >>> ag.ix array([2, 1, 2, 2, 1, 0]) >>> ag2 = ag.unique >>> ag2 <AtomGroup with 3 atoms> >>> ag2.ix array([0, 1, 2])
New in version 0.16.0.
-
velocities
¶ Velocities of the atoms in the AtomGroup.
The velocities can be changed by assigning an array of the appropriate shape, i.e. either Nx3 to assign individual velocities or 3 to assign the same velocity to all atoms (e.g.
ag.velocity = array([0,0,0])
will give all particles zero velocity).Raises a
NoDataError
if the underlyingTimestep
does not containvelocities
.
-
wrap
(compound=’atoms’, center=’com’, box=None)¶ Shift the contents of this Group back into the unit cell.
This is a more powerful version of
pack_into_box()
, allowing groups of atoms to be kept together through the process.Parameters: - compound ({'atoms', 'group', 'residues', 'segments', 'fragments'}) – The group which will be kept together through the shifting process.
- center ({'com', 'cog'}) – How to define the center of a given group of atoms.
- box (array) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
,[lx, ly, lz, alpha, beta, gamma]
. IfNone
, uses these timestep dimensions.
Notes
When specifying a compound, the translation is calculated based on each compound. The same translation is applied to all atoms within this compound, meaning it will not be broken by the shift. This might however mean that all atoms from the compound are not inside the unit cell, but rather the center of the compound is.
center allows the definition of the center of each group to be specified. This can be either ‘com’ for center of mass, or ‘cog’ for center of geometry.
box allows a unit cell to be given for the transformation. If not specified, an the dimensions information from the current Timestep will be used.
Note
wrap with all default keywords is identical to
pack_into_box()
New in version 0.9.2.
-
write
(filename=None, file_format=’PDB’, filenamefmt=’{trjname}_{frame}’, **kwargs)[source]¶ Write AtomGroup to a file.
The output can either be a coordinate file or a selection, depending on the format. Only single-frame coordinate files are supported. If you need to write out a trajectory, see
MDAnalysis.coordinates
.Parameters: - filename (str, optional) –
None
: create TRJNAME_FRAME.FORMAT from filenamefmt [None
] - file_format (str, optional) – PDB, CRD, GRO, VMD (tcl), PyMol (pml), Gromacs (ndx) CHARMM (str) Jmol (spt); case-insensitive and can also be supplied as the filename extension [PDB]
- filenamefmt (str, optional) – format string for default filename; use substitution tokens ‘trjname’ and ‘frame’ [“%(trjname)s_%(frame)d”]
- bonds (str, optional) – how to handle bond information, especially relevant for PDBs.
"conect"
: write only the CONECT records defined in the original file."all"
: write out all bonds, both the original defined and those guessed by MDAnalysis.None
: do not write out bonds. Default os"conect"
.
Changed in version 0.9.0: Merged with write_selection. This method can now write both selections out.
- filename (str, optional) –
-
-
class
MDAnalysis.core.groups.
ResidueGroup
(*args)[source]¶ ResidueGroup base class.
This class is used by a
Universe
for generating its Topology-specificResidueGroup
class. All theTopologyAttr
components are obtained fromGroupBase
, so this class only includes ad-hoc methods specific to ResidueGroups.ResidueGroups can be compared and combined using group operators. See the list of these operators on
GroupBase
.Deprecated since version 0.16.2: Instant selectors of Segments will be removed in the 1.0 release. See Instant selectors for details and alternatives.
-
atoms
¶ Get an
AtomGroup
of atoms represented in thisResidueGroup
.The atoms are ordered locally by residue in the
ResidueGroup
. No duplicates are removed.
-
bbox
(**kwargs)¶ Return the bounding box of the selection.
The lengths A,B,C of the orthorhombic enclosing box are
L = AtomGroup.bbox() A,B,C = L[1] - L[0]
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: corners – 2x3 array giving corners of bounding box as [[xmin, ymin, zmin], [xmax, ymax, zmax]]. Return type: array Note
The
MDAnalysis.core.flags
flag use_pbc when set toTrue
allows the pbc flag to be used by default.New in version 0.7.2.
Changed in version 0.8: Added pbc keyword
-
bsphere
(**kwargs)¶ Return the bounding sphere of the selection.
The sphere is calculated relative to the centre of geometry.
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: - R (float) – Radius of bounding sphere.
- center (array) – Coordinates of sphere center as
[xcen,ycen,zcen]
.
Note
The
MDAnalysis.core.flags
flag use_pbc when set toTrue
allows the pbc flag to be used by default.New in version 0.7.3.
Changed in version 0.8: Added pbc keyword
-
center
(weights, pbc=None)¶ Calculate center of group given some weights
Parameters: - weights (array_like) – weights to be used
- pbc (boolean, optional) –
True
: Move all atoms within the primary unit cell before calculation [False
]
Returns: center – weighted center of group
Return type: ndarray
Examples
To find the charge weighted center of a given Atomgroup:
>>> sel = u.select_atoms('prop mass > 4.0') >>> sel.center(sel.charges)
Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.
-
center_of_geometry
(pbc=None)¶ Center of geometry (also known as centroid) of the selection.
Parameters: pbc (boolean, optional) – True
: Move all atoms within the primary unit cell before calculation [False
]Returns: center – geometric center of group Return type: ndarray Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.Changed in version 0.8: Added pbc keyword
-
centroid
(pbc=None)¶ Center of geometry (also known as centroid) of the selection.
Parameters: pbc (boolean, optional) – True
: Move all atoms within the primary unit cell before calculation [False
]Returns: center – geometric center of group Return type: ndarray Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.Changed in version 0.8: Added pbc keyword
-
concatenate
(other)¶ Concatenate with another Group or Component of the same level.
Duplicate entries and original order is preserved. It is synomymous to the + operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with elements of self and other concatenated Return type: Group Example
The order of the original contents (including duplicates) are preserved when performing a concatenation.
>>> ag1 = u.select_atoms('name O') >>> ag2 = u.select_atoms('name N') >>> ag3 = ag1 + ag2 # or ag1.concatenate(ag2) >>> ag3[:3].names array(['O', 'O', 'O'], dtype=object) >>> ag3[-3:].names array(['N', 'N', 'N'], dtype=object)
New in version 0.16.0.
-
difference
(other)¶ Elements from this Group that do not appear in another
This method removes duplicate elements and sorts the result. As such, it is different from
subtract()
.difference()
is synomymous to the - operator.Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements of self that are not in other, without duplicate elements Return type: Group See also
-
groupby
(topattr)¶ Group together items in this group according to values of topattr
Parameters: topattr (str) – Topology attribute to group components by. Returns: Unique values of the topology attribute as keys, Groups as values. Return type: dict Example
To group atoms with the same mass together:
>>> ag.groupby('masses') {12.010999999999999: <AtomGroup with 462 atoms>, 14.007: <AtomGroup with 116 atoms>, 15.999000000000001: <AtomGroup with 134 atoms>}
New in version 0.16.0.
-
intersection
(other)¶ Group of elements which are in both this Group and another
This method removes duplicate elements and sorts the result. It is synomymous to the & operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the common elements of self and other, without duplicate elements Return type: Group Example
Intersections can be used when the select atoms string would become too complicated. For example to find the water atoms which are within 4.0A of two segments:
>>> water = u.select_atoms('resname SOL') >>> shell1 = water.select_atoms('around 4.0 segid 1') >>> shell2 = water.select_atoms('around 4.0 segid 2') >>> common = shell1 & shell2 # or shell1.intersection(shell2)
See also
union()
,()
-
is_strict_subset
(other)¶ If this Group is a subset of another Group but not identical
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a strict subset of the other one
- .. versionadded:: 0.16
-
is_strict_superset
(other)¶ If this Group is a superset of another Group but not identical
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a strict superset of the other one
- .. versionadded:: 0.16
-
isdisjoint
(other)¶ If the Group has no elements in common with the other Group
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if the two Groups do not have common elements
- .. versionadded:: 0.16
-
issubset
(other)¶ If all elements of this Group are part of another Group
Note that an empty group is a subset of any group of the same level.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a subset of the other one
- .. versionadded:: 0.16
-
issuperset
(other)¶ If all elements of another Group are part of this Group
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a subset of the other one
- .. versionadded:: 0.16
-
ix
¶ Unique indices of the components in the Group.
- If this Group is an
AtomGroup
, these are the indices of theAtom
instances. - If it is a
ResidueGroup
, these are the indices of theResidue
instances. - If it is a
SegmentGroup
, these are the indices of theSegment
instances.
- If this Group is an
-
ix_array
¶ Unique indices of the components in the Group.
For a Group, ix_array is the same as ix. This method gives a consistent API between components and groups.
See also
-
n_atoms
¶ Number of atoms represented in
ResidueGroup
, including duplicate residues.Equivalent to
len(self.atoms)
.
-
n_residues
¶ Number of residues in ResidueGroup. Equivalent to
len(self)
.
-
n_segments
¶ Number of unique segments represented in the ResidueGroup.
Equivalent to
len(self.segments)
.
-
pack_into_box
(box=None, inplace=True)¶ Shift all atoms in this group to be within the primary unit cell.
Parameters: - box (array_like) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
,[lx, ly, lz, alpha, beta, gamma]
. IfNone
, uses these timestep dimensions. - inplace (bool) –
True
to change coordinates in place.
Returns: coords – Shifted atom coordinates.
Return type: Notes
All atoms will be moved so that they lie between 0 and boxlength \(L_i\) in all dimensions, i.e. the lower left corner of the simulation box is taken to be at (0,0,0):
\[x_i' = x_i - \left\lfloor\frac{x_i}{L_i}\right\rfloor\]The default is to take unit cell information from the underlying
Timestep
instance. The optional argument box can be used to provide alternative unit cell information (in the MDAnalysis standard format[Lx, Ly, Lz, alpha, beta, gamma]
).Works with either orthogonal or triclinic box types.
New in version 0.8.
- box (array_like) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
-
residues
¶ Get another
ResidueGroup
identical to this one.
-
rotate
(R, point=(0, 0, 0))¶ Apply a rotation matrix R to the selection’s coordinates. \(\mathsf{R}\) is a 3x3 orthogonal matrix that transforms a vector \(\mathbf{x} \rightarrow \mathbf{x}'\):
\[\mathbf{x}' = \mathsf{R}\mathbf{x}\]Parameters: - R (array_like) – 3x3 rotation matrix to use for applying rotation.
- point (array_like, optional) – Center of rotation
Returns: self
Return type: Notes
By default rotates around center of origin
point=(0, 0, 0)
. To rotate around center of geometry of the atomgroup useag.rotate(R, point=ag.centroid)
.See also
rotateby()
- rotate around given axis and angle
MDAnalysis.lib.transformations()
- module of all coordinate transforms
-
rotateby
(angle, axis, point=None)¶ Apply a rotation to the selection’s coordinates.
Parameters: - angle (float) – Rotation angle in degrees.
- axis (array_like) – Rotation axis vector.
- point (array_like, optional) – Center of rotation. If
None
then the center of geometry of this group is used.
Returns: self
Return type: Notes
The transformation from current coordinates \(\mathbf{x}\) to new coordinates \(\mathbf{x}'\) is
\[\mathbf{x}' = \mathsf{R}\,(\mathbf{x}-\mathbf{p}) + \mathbf{p}\]where \(\mathsf{R}\) is the rotation by angle around the axis going through point \(\mathbf{p}\).
See also
MDAnalysis.lib.transformations.rotation_matrix()
- calculate \(\mathsf{R}\)
-
segments
¶ Get sorted SegmentGroup of the (unique) segments represented in the ResidueGroup.
-
subtract
(other)¶ Group with elements from this Group that don’t appear in 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 are removed from the returned Group.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements of self that are not in other, conserves order and duplicates. Return type: Group Example
Unlike
difference()
this method will not sort or remove duplicates.>>> ag1 = u.atoms[[3, 3, 2, 2, 1, 1]] >>> ag2 = u.atoms[2] >>> ag3 = ag1 - ag2 # or ag1.subtract(ag2) >>> ag1.indices array([3, 3, 1, 1])
See also
-
symmetric_difference
(other)¶ Group of elements which are only in one of this Group or another
This method removes duplicate elements and the result is sorted. It is synomym to the ^ operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements that are in self or in other but not in both, without duplicate elements Return type: Group Example
>>> ag1 = u.atoms[[0, 1, 5, 3, 3, 2]] >>> ag2 = u.atoms[[4, 4, 6, 2, 3, 5]] >>> ag3 = ag1 ^ ag2 # or ag1.symmetric_difference(ag2) >>> ag3.indices # 0 and 1 are only in ag1, 4 and 6 are only in ag2 [0, 1, 4, 6]
See also
difference()
,()
-
transform
(M)¶ Apply homogenous transformation matrix M to the coordinates.
Parameters: M (array) – 4x4 matrix, with the rotation in R = M[:3,:3]
and the translation int = M[:3,3]
.Returns: Return type: self See also
MDAnalysis.lib.transformations()
- module of all coordinate transforms
Notes
The rotation \(\mathsf{R}\) is applied before the translation \(\mathbf{t}\):
\[\mathbf{x}' = \mathsf{R}\mathbf{x} + \mathbf{t}\]
-
translate
(t)¶ Apply translation vector t to the selection’s coordinates.
Atom coordinates are translated in-place.
Parameters: t (array_like) – vector to translate coordinates with Returns: Return type: self See also
MDAnalysis.lib.transformations()
- module of all coordinate transforms
Notes
The method applies a translation to the
AtomGroup
from current coordinates \(\mathbf{x}\) to new coordinates \(\mathbf{x}'\):\[\mathbf{x}' = \mathbf{x} + \mathbf{t}\]
-
union
(other)¶ Group of elements either in this Group or another
On the contrary to concatenation, this method sort the elements and removes duplicate ones. It is synomymous to the | operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the combined elements of self and other, without duplicate elements Return type: Group Example
In contrast to
concatenate()
, any duplicates are dropped and the result is sorted.>>> ag1 = u.select_atoms('name O') >>> ag2 = u.select_atoms('name N') >>> ag3 = ag1 | ag2 # or ag1.union(ag2) >>> ag3[:3].names array(['N', 'O', 'N'], dtype=object)
See also
-
unique
¶ Return a ResidueGroup containing sorted and unique residues only.
Examples
>>> rg = u.residues[[2, 1, 2, 2, 1, 0]] >>> rg <ResidueGroup with 6 residues> >>> rg.ix array([2, 1, 2, 2, 1, 0]) >>> rg2 = rg.unique >>> rg2 <ResidueGroup with 3 residues> >>> rg2.ix array([0, 1, 2])
New in version 0.16.0.
-
wrap
(compound=’atoms’, center=’com’, box=None)¶ Shift the contents of this Group back into the unit cell.
This is a more powerful version of
pack_into_box()
, allowing groups of atoms to be kept together through the process.Parameters: - compound ({'atoms', 'group', 'residues', 'segments', 'fragments'}) – The group which will be kept together through the shifting process.
- center ({'com', 'cog'}) – How to define the center of a given group of atoms.
- box (array) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
,[lx, ly, lz, alpha, beta, gamma]
. IfNone
, uses these timestep dimensions.
Notes
When specifying a compound, the translation is calculated based on each compound. The same translation is applied to all atoms within this compound, meaning it will not be broken by the shift. This might however mean that all atoms from the compound are not inside the unit cell, but rather the center of the compound is.
center allows the definition of the center of each group to be specified. This can be either ‘com’ for center of mass, or ‘cog’ for center of geometry.
box allows a unit cell to be given for the transformation. If not specified, an the dimensions information from the current Timestep will be used.
Note
wrap with all default keywords is identical to
pack_into_box()
New in version 0.9.2.
-
-
class
MDAnalysis.core.groups.
SegmentGroup
(*args)[source]¶ SegmentGroup base class.
This class is used by a Universe for generating its Topology-specific SegmentGroup class. All the TopologyAttr components are obtained from GroupBase, so this class only includes ad-hoc methods specific to SegmentGroups.
SegmentGroups can be compared and combined using group operators. See the list of these operators on
GroupBase
.Deprecated since version 0.16.2: Instant selectors of Segments will be removed in the 1.0 release. See Instant selectors for details and alternatives.
-
atoms
¶ Get an AtomGroup of atoms represented in this SegmentGroup.
The atoms are ordered locally by residue, which are further ordered by segment in the SegmentGroup. No duplicates are removed.
-
bbox
(**kwargs)¶ Return the bounding box of the selection.
The lengths A,B,C of the orthorhombic enclosing box are
L = AtomGroup.bbox() A,B,C = L[1] - L[0]
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: corners – 2x3 array giving corners of bounding box as [[xmin, ymin, zmin], [xmax, ymax, zmax]]. Return type: array Note
The
MDAnalysis.core.flags
flag use_pbc when set toTrue
allows the pbc flag to be used by default.New in version 0.7.2.
Changed in version 0.8: Added pbc keyword
-
bsphere
(**kwargs)¶ Return the bounding sphere of the selection.
The sphere is calculated relative to the centre of geometry.
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: - R (float) – Radius of bounding sphere.
- center (array) – Coordinates of sphere center as
[xcen,ycen,zcen]
.
Note
The
MDAnalysis.core.flags
flag use_pbc when set toTrue
allows the pbc flag to be used by default.New in version 0.7.3.
Changed in version 0.8: Added pbc keyword
-
center
(weights, pbc=None)¶ Calculate center of group given some weights
Parameters: - weights (array_like) – weights to be used
- pbc (boolean, optional) –
True
: Move all atoms within the primary unit cell before calculation [False
]
Returns: center – weighted center of group
Return type: ndarray
Examples
To find the charge weighted center of a given Atomgroup:
>>> sel = u.select_atoms('prop mass > 4.0') >>> sel.center(sel.charges)
Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.
-
center_of_geometry
(pbc=None)¶ Center of geometry (also known as centroid) of the selection.
Parameters: pbc (boolean, optional) – True
: Move all atoms within the primary unit cell before calculation [False
]Returns: center – geometric center of group Return type: ndarray Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.Changed in version 0.8: Added pbc keyword
-
centroid
(pbc=None)¶ Center of geometry (also known as centroid) of the selection.
Parameters: pbc (boolean, optional) – True
: Move all atoms within the primary unit cell before calculation [False
]Returns: center – geometric center of group Return type: ndarray Notes
If the
MDAnalysis.core.flags
flag use_pbc is set toTrue
then the pbc keyword is used by default.Changed in version 0.8: Added pbc keyword
-
concatenate
(other)¶ Concatenate with another Group or Component of the same level.
Duplicate entries and original order is preserved. It is synomymous to the + operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with elements of self and other concatenated Return type: Group Example
The order of the original contents (including duplicates) are preserved when performing a concatenation.
>>> ag1 = u.select_atoms('name O') >>> ag2 = u.select_atoms('name N') >>> ag3 = ag1 + ag2 # or ag1.concatenate(ag2) >>> ag3[:3].names array(['O', 'O', 'O'], dtype=object) >>> ag3[-3:].names array(['N', 'N', 'N'], dtype=object)
New in version 0.16.0.
-
difference
(other)¶ Elements from this Group that do not appear in another
This method removes duplicate elements and sorts the result. As such, it is different from
subtract()
.difference()
is synomymous to the - operator.Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements of self that are not in other, without duplicate elements Return type: Group See also
-
groupby
(topattr)¶ Group together items in this group according to values of topattr
Parameters: topattr (str) – Topology attribute to group components by. Returns: Unique values of the topology attribute as keys, Groups as values. Return type: dict Example
To group atoms with the same mass together:
>>> ag.groupby('masses') {12.010999999999999: <AtomGroup with 462 atoms>, 14.007: <AtomGroup with 116 atoms>, 15.999000000000001: <AtomGroup with 134 atoms>}
New in version 0.16.0.
-
intersection
(other)¶ Group of elements which are in both this Group and another
This method removes duplicate elements and sorts the result. It is synomymous to the & operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the common elements of self and other, without duplicate elements Return type: Group Example
Intersections can be used when the select atoms string would become too complicated. For example to find the water atoms which are within 4.0A of two segments:
>>> water = u.select_atoms('resname SOL') >>> shell1 = water.select_atoms('around 4.0 segid 1') >>> shell2 = water.select_atoms('around 4.0 segid 2') >>> common = shell1 & shell2 # or shell1.intersection(shell2)
See also
union()
,()
-
is_strict_subset
(other)¶ If this Group is a subset of another Group but not identical
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a strict subset of the other one
- .. versionadded:: 0.16
-
is_strict_superset
(other)¶ If this Group is a superset of another Group but not identical
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a strict superset of the other one
- .. versionadded:: 0.16
-
isdisjoint
(other)¶ If the Group has no elements in common with the other Group
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if the two Groups do not have common elements
- .. versionadded:: 0.16
-
issubset
(other)¶ If all elements of this Group are part of another Group
Note that an empty group is a subset of any group of the same level.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a subset of the other one
- .. versionadded:: 0.16
-
issuperset
(other)¶ If all elements of another Group are part of this Group
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: - Bool – True if this Group is a subset of the other one
- .. versionadded:: 0.16
-
ix
¶ Unique indices of the components in the Group.
- If this Group is an
AtomGroup
, these are the indices of theAtom
instances. - If it is a
ResidueGroup
, these are the indices of theResidue
instances. - If it is a
SegmentGroup
, these are the indices of theSegment
instances.
- If this Group is an
-
ix_array
¶ Unique indices of the components in the Group.
For a Group, ix_array is the same as ix. This method gives a consistent API between components and groups.
See also
-
n_atoms
¶ Number of atoms represented in SegmentGroup, including duplicate segments.
Equivalent to
len(self.atoms)
.
-
n_residues
¶ Number of residues represented in SegmentGroup, including duplicate segments.
Equivalent to
len(self.residues)
.
-
n_segments
¶ Number of segments in SegmentGroup. Equivalent to
len(self)
.
-
pack_into_box
(box=None, inplace=True)¶ Shift all atoms in this group to be within the primary unit cell.
Parameters: - box (array_like) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
,[lx, ly, lz, alpha, beta, gamma]
. IfNone
, uses these timestep dimensions. - inplace (bool) –
True
to change coordinates in place.
Returns: coords – Shifted atom coordinates.
Return type: Notes
All atoms will be moved so that they lie between 0 and boxlength \(L_i\) in all dimensions, i.e. the lower left corner of the simulation box is taken to be at (0,0,0):
\[x_i' = x_i - \left\lfloor\frac{x_i}{L_i}\right\rfloor\]The default is to take unit cell information from the underlying
Timestep
instance. The optional argument box can be used to provide alternative unit cell information (in the MDAnalysis standard format[Lx, Ly, Lz, alpha, beta, gamma]
).Works with either orthogonal or triclinic box types.
New in version 0.8.
- box (array_like) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
-
residues
¶ Get a ResidueGroup of residues represented in this SegmentGroup.
The residues are ordered locally by segment in the SegmentGroup. No duplicates are removed.
-
rotate
(R, point=(0, 0, 0))¶ Apply a rotation matrix R to the selection’s coordinates. \(\mathsf{R}\) is a 3x3 orthogonal matrix that transforms a vector \(\mathbf{x} \rightarrow \mathbf{x}'\):
\[\mathbf{x}' = \mathsf{R}\mathbf{x}\]Parameters: - R (array_like) – 3x3 rotation matrix to use for applying rotation.
- point (array_like, optional) – Center of rotation
Returns: self
Return type: Notes
By default rotates around center of origin
point=(0, 0, 0)
. To rotate around center of geometry of the atomgroup useag.rotate(R, point=ag.centroid)
.See also
rotateby()
- rotate around given axis and angle
MDAnalysis.lib.transformations()
- module of all coordinate transforms
-
rotateby
(angle, axis, point=None)¶ Apply a rotation to the selection’s coordinates.
Parameters: - angle (float) – Rotation angle in degrees.
- axis (array_like) – Rotation axis vector.
- point (array_like, optional) – Center of rotation. If
None
then the center of geometry of this group is used.
Returns: self
Return type: Notes
The transformation from current coordinates \(\mathbf{x}\) to new coordinates \(\mathbf{x}'\) is
\[\mathbf{x}' = \mathsf{R}\,(\mathbf{x}-\mathbf{p}) + \mathbf{p}\]where \(\mathsf{R}\) is the rotation by angle around the axis going through point \(\mathbf{p}\).
See also
MDAnalysis.lib.transformations.rotation_matrix()
- calculate \(\mathsf{R}\)
-
segments
¶ Get another SegmentGroup identical to this one.
-
subtract
(other)¶ Group with elements from this Group that don’t appear in 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 are removed from the returned Group.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements of self that are not in other, conserves order and duplicates. Return type: Group Example
Unlike
difference()
this method will not sort or remove duplicates.>>> ag1 = u.atoms[[3, 3, 2, 2, 1, 1]] >>> ag2 = u.atoms[2] >>> ag3 = ag1 - ag2 # or ag1.subtract(ag2) >>> ag1.indices array([3, 3, 1, 1])
See also
-
symmetric_difference
(other)¶ Group of elements which are only in one of this Group or another
This method removes duplicate elements and the result is sorted. It is synomym to the ^ operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the elements that are in self or in other but not in both, without duplicate elements Return type: Group Example
>>> ag1 = u.atoms[[0, 1, 5, 3, 3, 2]] >>> ag2 = u.atoms[[4, 4, 6, 2, 3, 5]] >>> ag3 = ag1 ^ ag2 # or ag1.symmetric_difference(ag2) >>> ag3.indices # 0 and 1 are only in ag1, 4 and 6 are only in ag2 [0, 1, 4, 6]
See also
difference()
,()
-
transform
(M)¶ Apply homogenous transformation matrix M to the coordinates.
Parameters: M (array) – 4x4 matrix, with the rotation in R = M[:3,:3]
and the translation int = M[:3,3]
.Returns: Return type: self See also
MDAnalysis.lib.transformations()
- module of all coordinate transforms
Notes
The rotation \(\mathsf{R}\) is applied before the translation \(\mathbf{t}\):
\[\mathbf{x}' = \mathsf{R}\mathbf{x} + \mathbf{t}\]
-
translate
(t)¶ Apply translation vector t to the selection’s coordinates.
Atom coordinates are translated in-place.
Parameters: t (array_like) – vector to translate coordinates with Returns: Return type: self See also
MDAnalysis.lib.transformations()
- module of all coordinate transforms
Notes
The method applies a translation to the
AtomGroup
from current coordinates \(\mathbf{x}\) to new coordinates \(\mathbf{x}'\):\[\mathbf{x}' = \mathbf{x} + \mathbf{t}\]
-
union
(other)¶ Group of elements either in this Group or another
On the contrary to concatenation, this method sort the elements and removes duplicate ones. It is synomymous to the | operator.
Parameters: other (Group or Component) – Group or Component with other.level same as self.level Returns: Group with the combined elements of self and other, without duplicate elements Return type: Group Example
In contrast to
concatenate()
, any duplicates are dropped and the result is sorted.>>> ag1 = u.select_atoms('name O') >>> ag2 = u.select_atoms('name N') >>> ag3 = ag1 | ag2 # or ag1.union(ag2) >>> ag3[:3].names array(['N', 'O', 'N'], dtype=object)
See also
-
unique
¶ Return a SegmentGroup containing sorted and unique segments only.
Examples
>>> sg = u.segments[[2, 1, 2, 2, 1, 0]] >>> sg <SegmentGroup with 6 segments> >>> sg.ix array([2, 1, 2, 2, 1, 0]) >>> sg2 = sg.unique >>> sg2 <SegmentGroup with 3 segments> >>> sg2.ix array([0, 1, 2])
New in version 0.16.0.
-
wrap
(compound=’atoms’, center=’com’, box=None)¶ Shift the contents of this Group back into the unit cell.
This is a more powerful version of
pack_into_box()
, allowing groups of atoms to be kept together through the process.Parameters: - compound ({'atoms', 'group', 'residues', 'segments', 'fragments'}) – The group which will be kept together through the shifting process.
- center ({'com', 'cog'}) – How to define the center of a given group of atoms.
- box (array) – Box dimensions, can be either orthogonal or triclinic information.
Cell dimensions must be in an identical to format to those returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
,[lx, ly, lz, alpha, beta, gamma]
. IfNone
, uses these timestep dimensions.
Notes
When specifying a compound, the translation is calculated based on each compound. The same translation is applied to all atoms within this compound, meaning it will not be broken by the shift. This might however mean that all atoms from the compound are not inside the unit cell, but rather the center of the compound is.
center allows the definition of the center of each group to be specified. This can be either ‘com’ for center of mass, or ‘cog’ for center of geometry.
box allows a unit cell to be given for the transformation. If not specified, an the dimensions information from the current Timestep will be used.
Note
wrap with all default keywords is identical to
pack_into_box()
New in version 0.9.2.
-
-
class
MDAnalysis.core.groups.
UpdatingAtomGroup
(base_group, selections, strings)[source]¶ AtomGroup
subclass that dynamically updates its selected atoms.Accessing any attribute/method of an
UpdatingAtomGroup
instance triggers a check for the last frame the group was updated. If the last frame matches the current trajectory frame, the attribute is returned normally; otherwise the group is updated (the stored selections are re-applied), and only then is the attribute returned.New in version 0.16.0.
Parameters: - base_group (
AtomGroup
) – group on which selections are to be applied. - selections (a tuple of
Selection
) – instances selections ready to be applied to base_group.
-
is_uptodate
¶ Checks whether the selection needs updating based on frame number only.
Modifications to the coordinate data that render selections stale are not caught, and in those cases
is_uptodate
may return an erroneous value.Returns: True if the group’s selection is up-to-date, False otherwise. Return type: bool
- base_group (
9.1.2.1.2. Chemical units¶
-
class
MDAnalysis.core.groups.
Atom
(ix, u)[source]¶ Atom base class.
This class is used by a Universe for generating its Topology-specific Atom class. All the TopologyAttr components are obtained from ComponentBase, so this class only includes ad-hoc methods specific to Atoms.
-
force
¶ Force on the atom.
The force can be changed by assigning an array of shape (3,).
Note
changing the force is not reflected in any files; reading any frame from the trajectory will replace the change with that from the file
A
NoDataError
is raised if the trajectory does not contain forces.
-
ix
¶ Unique index of this component.
If this component is an Atom, this is the index of the atom. If it is a Residue, this is the index of the residue. If it is a Segment, this is the index of the segment.
-
ix_array
¶ Unique index of this component as an array.
This method gives a consistent API between components and groups.
See also
-
position
¶ Coordinates of the atom.
The position can be changed by assigning an array of length (3,).
Note
changing the position is not reflected in any files; reading any frame from the trajectory will replace the change with that from the file
-
velocity
¶ Velocity of the atom.
The velocity can be changed by assigning an array of shape (3,).
Note
changing the velocity is not reflected in any files; reading any frame from the trajectory will replace the change with that from the file
A
NoDataError
is raised if the trajectory does not contain velocities.
-
-
class
MDAnalysis.core.groups.
Residue
(ix, u)[source]¶ Residue base class.
This class is used by a Universe for generating its Topology-specific Residue class. All the TopologyAttr components are obtained from ComponentBase, so this class only includes ad-hoc methods specific to Residues.
-
ix
¶ Unique index of this component.
If this component is an Atom, this is the index of the atom. If it is a Residue, this is the index of the residue. If it is a Segment, this is the index of the segment.
-
-
class
MDAnalysis.core.groups.
Segment
(ix, u)[source]¶ Segment base class.
This class is used by a Universe for generating its Topology-specific Segment class. All the TopologyAttr components are obtained from ComponentBase, so this class only includes ad-hoc methods specific to Segments.
Deprecated since version 0.16.2: Instant selectors of Segments will be removed in the 1.0 release. See Instant selectors for details and alternatives.
-
ix
¶ Unique index of this component.
If this component is an Atom, this is the index of the atom. If it is a Residue, this is the index of the residue. If it is a Segment, this is the index of the segment.
-
9.1.2.1.3. Levels¶
Each of the above classes has a level attribute. This can be used to verify that two objects are of the same level, or to access a particular class:
u = mda.Universe()
ag = u.atoms[:10]
at = u.atoms[11]
ag.level == at.level # Returns True
ag.level.singular # Returns Atom class
at.level.plural # Returns AtomGroup class