11.3.3. Topology attribute objects — MDAnalysis.core.topologyattrs
¶
Common TopologyAttr
instances that are used by most topology
parsers.
TopologyAttrs are used to contain attributes such as atom names or resids. These are usually read by the TopologyParser.
-
class
MDAnalysis.core.topologyattrs.
AltLocs
(vals, guessed=False)[source]¶ AltLocs for each atom
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Angles
(values, types=None, guessed=False, order=None)[source]¶ Angles between three atoms
Initialise with a list of 3 long tuples E.g., [(0, 1, 2), (1, 2, 3), (2, 3, 4)]
These indices refer to the atom indices.
-
class
MDAnalysis.core.topologyattrs.
Aromaticities
(values, guessed=False)[source]¶ Aromaticity (RDKit)
-
dtype
¶ alias of
builtins.bool
-
-
class
MDAnalysis.core.topologyattrs.
AtomAttr
(values, guessed=False)[source]¶ Base class for atom attributes.
-
get_residues
(rg)[source]¶ By default, the values for each atom present in the set of residues are returned in a single array. This behavior can be overriden in child attributes.
-
-
class
MDAnalysis.core.topologyattrs.
Atomids
(values, guessed=False)[source]¶ ID for each atom.
-
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
Atomindices
[source]¶ Globally unique indices for each atom in the group.
If the group is an AtomGroup, then this gives the index for each atom in the group. This is the unambiguous identifier for each atom in the topology, and it is not alterable.
If the group is a ResidueGroup or SegmentGroup, then this gives the indices of each atom represented in the group in a 1-D array, in the order of the elements in that group.
-
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
Atomnames
(vals, guessed=False)[source]¶ Name for each atom.
-
chi1_selection
(n_name='N', ca_name='CA', cb_name='CB', cg_name='CG CG1 OG OG1 SG')[source]¶ Select AtomGroup corresponding to the chi1 sidechain dihedral
N-CA-CB-*G.
The gamma atom is taken to be the heavy atom in the gamma position. If more than one heavy atom is present (e.g. CG1 and CG2), the one with the lower number is used (CG1).Warning
This numbering of chi1 atoms here in following with the IUPAC 1970 rules. However, it should be noted that analyses which use dihedral angles may have different definitions. For example, the
MDAnalysis.analysis.dihedrals.Janin
class does not incorporate amino acids where the gamma atom is not carbon, into its chi1 selections.Parameters: Returns: - AtomGroup – 4-atom selection in the correct order. If no CB and/or CG is found
then this method returns
None
. - .. versionchanged:: 1.0.0 – Added arguments for flexible atom names and refactored code for faster atom matching with boolean arrays.
- .. versionadded:: 0.7.5
- AtomGroup – 4-atom selection in the correct order. If no CB and/or CG is found
then this method returns
-
chi1_selections
(n_name='N', ca_name='CA', cb_name='CB', cg_name='CG')[source]¶ Select list of AtomGroups corresponding to the chi1 sidechain dihedral N-CA-CB-CG.
Parameters: Returns: - List of AtomGroups – 4-atom selections in the correct order. If no CB and/or CG is found
then the corresponding item in the list is
None
. - .. versionadded:: 1.0.0
- List of AtomGroups – 4-atom selections in the correct order. If no CB and/or CG is found
then the corresponding item in the list is
-
dtype
¶ alias of
builtins.object
-
omega_selection
(c_name='C', n_name='N', ca_name='CA')[source]¶ Select 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).
Parameters: Returns: - AtomGroup – 4-atom selection in the correct order. If no C’ found in the
previous residue (by resid) then this method returns
None
. - .. versionchanged:: 1.0.0 – Added arguments for flexible atom names and refactored code for faster atom matching with boolean arrays.
- AtomGroup – 4-atom selection in the correct order. If no C’ found in the
previous residue (by resid) then this method returns
-
omega_selections
(c_name='C', n_name='N', ca_name='CA')[source]¶ Select list of AtomGroups 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).
Parameters: Returns: - List of AtomGroups – 4-atom selections in the correct order. If no C’ found in the
previous residue (by resid) then the corresponding item in the
list is
None
. - .. versionadded:: 1.0.0
- List of AtomGroups – 4-atom selections in the correct order. If no C’ found in the
previous residue (by resid) then the corresponding item in the
list is
-
phi_selection
(c_name='C', n_name='N', ca_name='CA')[source]¶ Select AtomGroup corresponding to the phi protein backbone dihedral C’-N-CA-C.
Parameters: Returns: - AtomGroup – 4-atom selection in the correct order. If no C’ found in the
previous residue (by resid) then this method returns
None
. - .. versionchanged:: 1.0.0 – Added arguments for flexible atom names and refactored code for faster atom matching with boolean arrays.
- AtomGroup – 4-atom selection in the correct order. If no C’ found in the
previous residue (by resid) then this method returns
-
phi_selections
(c_name='C', n_name='N', ca_name='CA')[source]¶ Select list of AtomGroups corresponding to the phi protein backbone dihedral C’-N-CA-C.
Parameters: Returns: - list of AtomGroups – 4-atom selections in the correct order. If no C’ found in the
previous residue (by resid) then corresponding item in the list
is
None
. - .. versionadded:: 1.0.0
- list of AtomGroups – 4-atom selections in the correct order. If no C’ found in the
previous residue (by resid) then corresponding item in the list
is
-
psi_selection
(c_name='C', n_name='N', ca_name='CA')[source]¶ Select AtomGroup corresponding to the psi protein backbone dihedral N-CA-C-N’.
Parameters: Returns: - AtomGroup – 4-atom selection in the correct order. If no N’ found in the
following residue (by resid) then this method returns
None
. - .. versionchanged:: 1.0.0 – Added arguments for flexible atom names and refactored code for faster atom matching with boolean arrays.
- AtomGroup – 4-atom selection in the correct order. If no N’ found in the
following residue (by resid) then this method returns
-
psi_selections
(c_name='C', n_name='N', ca_name='CA')[source]¶ Select list of AtomGroups corresponding to the psi protein backbone dihedral N-CA-C-N’.
Parameters: Returns: - List of AtomGroups – 4-atom selections in the correct order. If no N’ found in the
following residue (by resid) then the corresponding item in the
list is
None
. - .. versionadded:: 1.0.0
- List of AtomGroups – 4-atom selections in the correct order. If no N’ found in the
following residue (by resid) then the corresponding item in the
list is
-
-
class
MDAnalysis.core.topologyattrs.
Atomtypes
(vals, guessed=False)[source]¶ Type for each atom
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Bonds
(values, types=None, guessed=False, order=None)[source]¶ Bonds between two atoms
Must be initialised by a list of zero based tuples. These indices refer to the atom indices. E.g., ` [(0, 1), (1, 2), (2, 3)]`
Also adds the bonded_atoms, fragment and fragments attributes.
-
fragindices
()[source]¶ The
fragment indices
of allAtoms
in thisAtomGroup
.A
numpy.ndarray
withshape
=(
n_atoms
,)
anddtype
=numpy.int64
.New in version 0.20.0.
-
fragment
()[source]¶ An
AtomGroup
representing the fragment thisAtom
is part of.A fragment is a
group of atoms
which are interconnected byBonds
, i.e., there exists a path along one or moreBonds
between any pair ofAtoms
within a fragment. Thus, a fragment typically corresponds to a molecule.New in version 0.9.0.
-
fragments
()[source]¶ -
Contains all fragments that any
Atom
in thisAtomGroup
is part of.A fragment is a
group of atoms
which are interconnected byBonds
, i.e., there exists a path along one or moreBonds
between any pair ofAtoms
within a fragment. Thus, a fragment typically corresponds to a molecule.Note
- The contents of the fragments may extend beyond the contents of this
AtomGroup
.
New in version 0.9.0.
- The contents of the fragments may extend beyond the contents of this
-
-
class
MDAnalysis.core.topologyattrs.
CMaps
(values, types=None, guessed=False, order=None)[source]¶ A connection between five atoms .. versionadded:: 1.0.0
-
class
MDAnalysis.core.topologyattrs.
ChainIDs
(vals, guessed=False)[source]¶ ChainID per atom
Note
This is an attribute of the Atom, not Residue or Segment
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Charges
(values, guessed=False)[source]¶ -
dtype
¶ alias of
builtins.float
-
get_residues
(rg)[source]¶ By default, the values for each atom present in the set of residues are returned in a single array. This behavior can be overriden in child attributes.
-
get_segments
(sg)[source]¶ By default, the values for each atom present in the set of residues are returned in a single array. This behavior can be overriden in child attributes.
-
total_charge
(compound='group')[source]¶ Total charge of (compounds of) the group.
Computes the total charge of
Atoms
in the group. Total charges perResidue
,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 Segment
,Residue
, molecule, or fragment will be returned as a 1d array. Note that, in any case, only the charges ofAtoms
belonging to the group will be taken into account.Returns: 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,)
wheren
is the number of compounds.Return type: float or numpy.ndarray Changed in version 0.20.0: Added compound parameter
-
-
class
MDAnalysis.core.topologyattrs.
Dihedrals
(values, types=None, guessed=False, order=None)[source]¶ A connection between four sequential atoms
-
class
MDAnalysis.core.topologyattrs.
Elements
(vals, guessed=False)[source]¶ Element for each atom
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Epsilon14s
(values, guessed=False)[source]¶ The epsilon LJ parameter for 1-4 interactions
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
Epsilons
(values, guessed=False)[source]¶ The epsilon LJ parameter
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
GBScreens
(values, guessed=False)[source]¶ Generalized Born screening factor
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
ICodes
(vals, guessed=False)[source]¶ Insertion code for Atoms
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Impropers
(values, types=None, guessed=False, order=None)[source]¶ An imaginary dihedral between four atoms
-
class
MDAnalysis.core.topologyattrs.
Masses
(values, guessed=False)[source]¶ -
align_principal_axis
(axis, vector)[source]¶ 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
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")
- axis ({0, 1, 2}) – Index of the principal axis (0, 1, or 2), as produced by
-
asphericity
(pbc=False, unwrap=None, compound='group')[source]¶ Asphericity.
See [Dima2004b] for background information.
Parameters: - pbc (bool, optional) – If
True
, move all atoms within the primary unit cell before calculation. [False
] - 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.
New in version 0.7.7.
Changed in version 0.8: Added pbc keyword
Changed in version 0.20.0: Added unwrap and compound parameter
- pbc (bool, optional) – If
-
center_of_mass
(pbc=False, compound='group', unwrap=False)[source]¶ Center of mass of (compounds of) the group.
Computes the center of mass of
Atoms
in the group. Centers of mass perResidue
,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 benan
(not a number).Parameters: - pbc (bool, optional) – If
True
and compound is'group'
, move all atoms to the primary unit cell before calculation. IfTrue
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 eachSegment
,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 ofAtoms
belonging to the group will be taken into account. - unwrap (bool, optional) – If
True
, compounds will be unwrapped before computing their centers.
Returns: center – 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)
wheren
is the number of compounds.Return type: Note
- This method can only be accessed if the underlying topology has information about atomic masses.
Changed in version 0.8: Added pbc parameter
Changed in version 0.19.0: Added compound parameter
Changed in version 0.20.0: Added
'molecules'
and'fragments'
compoundsChanged in version 0.20.0: Added unwrap parameter
- pbc (bool, optional) – If
-
dtype
¶ alias of
numpy.float64
-
get_residues
(rg)[source]¶ By default, the values for each atom present in the set of residues are returned in a single array. This behavior can be overriden in child attributes.
-
get_segments
(sg)[source]¶ By default, the values for each atom present in the set of residues are returned in a single array. This behavior can be overriden in child attributes.
-
moment_of_inertia
(pbc=False, **kwargs)[source]¶ 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
]Changed in version 0.8: Added pbc keyword
Changed in version 0.20.0: Added unwrap parameter
-
principal_axes
(pbc=False)[source]¶ 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.
The eigenvectors form a right-handed coordinate system.
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Returns: axis_vectors – 3 x 3 array with v[0]
as first,v[1]
as second, andv[2]
as third eigenvector.Return type: array Changed in version 0.8: Added pbc keyword
Changed in version 1.0.0: Always return principal axes in right-hand convention.
-
radius_of_gyration
(pbc=False, **kwargs)[source]¶ Radius of gyration.
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]Changed in version 0.8: Added pbc keyword
-
shape_parameter
(pbc=False, **kwargs)[source]¶ Shape parameter.
See [Dima2004a] for background information.
Parameters: pbc (bool, optional) – If True
, move all atoms within the primary unit cell before calculation. [False
]New in version 0.7.7.
Changed in version 0.8: Added pbc keyword
-
total_mass
(compound='group')[source]¶ Total mass of (compounds of) the group.
Computes the total mass of
Atoms
in the group. Total masses perResidue
,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 perSegment
,Residue
, molecule, or fragment will be returned as a 1d array. Note that, in any case, only the masses ofAtoms
belonging to the group will be taken into account.Returns: 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,)
wheren
is the number of compounds.Return type: float or numpy.ndarray Changed in version 0.20.0: Added compound parameter
-
-
class
MDAnalysis.core.topologyattrs.
Molnums
(values, guessed=False)[source]¶ Index of molecule from 0
-
dtype
¶ alias of
numpy.int64
-
-
class
MDAnalysis.core.topologyattrs.
Moltypes
(vals, guessed=False)[source]¶ Name of the molecule type
Two molecules that share a molecule type share a common template topology.
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
NonbondedIndices
(values, guessed=False)[source]¶ Nonbonded index (AMBER)
-
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
Occupancies
(values, guessed=False)[source]¶ -
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
RMin14s
(values, guessed=False)[source]¶ The Rmin/2 LJ parameter for 1-4 interactions
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
RMins
(values, guessed=False)[source]¶ The Rmin/2 LJ parameter
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
Radii
(values, guessed=False)[source]¶ Radii for each atom
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
RecordTypes
(vals, guessed=False)[source]¶ For PDB-like formats, indicates if ATOM or HETATM
Defaults to ‘ATOM’
Changed in version 0.20.0: Now stores array of dtype object rather than boolean mapping
-
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Resids
(values, guessed=False)[source]¶ Residue ID
-
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
ResidueAttr
(values, guessed=False)[source]¶
-
class
MDAnalysis.core.topologyattrs.
Resindices
[source]¶ Globally unique resindices for each residue in the group.
If the group is an AtomGroup, then this gives the resindex for each atom in the group. This unambiguously determines each atom’s residue membership. 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 order of the elements in that group.
-
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
Resnames
(vals, guessed=False)[source]¶ -
dtype
¶ alias of
builtins.object
-
sequence
(**kwargs)[source]¶ Returns the amino acid sequence.
The format of the sequence is selected with the keyword format:
format description ‘SeqRecord’ Bio.SeqRecord.SeqRecord
(default)‘Seq’ Bio.Seq.Seq
‘string’ string The sequence is returned by default (keyword
format = 'SeqRecord'
) as aBio.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 toBio.SeqRecord.SeqRecord
.If the keyword format is set to
'Seq'
, all kwargs are ignored and aBio.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 inBio
.If the keyword format is set to
'string'
, all kwargs are ignored and a Python string is returned.Example: Write FASTA file
Use
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 aBio.Seq.Seq
instance"SeqRecord"
: return aBio.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: ValueError
if a residue name cannot be converted to a- 1-letter IUPAC protein amino acid code; make sure to only
- select protein residues.
TypeError
if an unknown format is selected.
New in version 0.9.0.
- format (string, optional) –
-
-
class
MDAnalysis.core.topologyattrs.
Resnums
(values, guessed=False)[source]¶ -
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
Segids
(vals, guessed=False)[source]¶ -
dtype
¶ alias of
builtins.object
-
-
class
MDAnalysis.core.topologyattrs.
Segindices
[source]¶ Globally unique segindices for each segment in the group.
If the group is an AtomGroup, then this gives the segindex for each atom in the group. This unambiguously determines each atom’s segment membership. It is not possible to set these, since membership in a segment is an attribute of each atom’s residue.
If the group is a ResidueGroup or SegmentGroup, then this gives the segindices of each segment represented in the group in a 1-D array, in the order of the elements in that group.
-
dtype
¶ alias of
builtins.int
-
-
class
MDAnalysis.core.topologyattrs.
SegmentAttr
(values, guessed=False)[source]¶ Base class for segment attributes.
-
class
MDAnalysis.core.topologyattrs.
SolventRadii
(values, guessed=False)[source]¶ Intrinsic solvation radius
-
dtype
¶ alias of
builtins.float
-
-
class
MDAnalysis.core.topologyattrs.
Tempfactors
(values, guessed=False)[source]¶ Tempfactor for atoms
-
dtype
¶ alias of
builtins.float
-
group
¶ alias of
MDAnalysis.core.groups.SegmentGroup
-
-
class
MDAnalysis.core.topologyattrs.
TopologyAttr
(values, guessed=False)[source]¶ Base class for Topology attributes.
Note
This class is intended to be subclassed, and mostly amounts to a skeleton. The methods here should be present in all
TopologyAttr
child classes, but by default they raise appropriate exceptions.-
dtype
¶ Type to coerce this attribute to be. For string use ‘object’
Type: int/float/object
-
classmethod
from_blank
(n_atoms=None, n_residues=None, n_segments=None, values=None)[source]¶ Create a blank version of this TopologyAttribute
Parameters:
-
is_guessed
¶ Bool of if the source of this information is a guess
-