5.26. Guessing unknown Topology information — MDAnalysis.topology.guessers

In general guess_atom_X returns the guessed value for a single value, while guess_Xs will work on an array of many atoms.

5.26.1. Example uses of guessers

5.26.1.1. Guessing elements from atom names

Currently, it is possible to guess elements from atom names using guess_atom_element() (or the synonymous guess_atom_type()). This can be done in the following manner:

import MDAnalysis as mda
from MDAnalysis.topology.guessers import guess_atom_element
from MDAnalysisTests.datafiles import PRM7

u = mda.Universe(PRM7)

print(u.atoms.names[1])  # returns the atom name H1

element = guess_atom_element(u.atoms.names[1])

print(element)  # returns element H

In the above example, we take an atom named H1 and use guess_atom_element() to guess the element hydrogen (i.e. H). It is important to note that element guessing is not always accurate. Indeed in cases where the atom type is not recognised, we may end up with the wrong element. For example:

import MDAnalysis as mda
from MDAnalysis.topology.guessers import guess_atom_element
from MDAnalysisTests.datafiles import PRM19SBOPC

u = mda.Universe(PRM19SBOPC)

print(u.atoms.names[-1])  # returns the atom name EPW

element = guess_atom_element(u.atoms.names[-1])

print(element)  # returns element P

Here we find that virtual site atom ‘EPW’ was given the element P, which would not be an expected result. We therefore always recommend that users carefully check the outcomes of any guessers.

In some cases, one may want to guess elements for an entire universe and add this guess as a topology attribute. This can be done using guess_types() in the following manner:

import MDAnalysis as mda
from MDAnalysis.topology.guessers import guess_types
from MDAnalysisTests.datafiles import PRM7

u = mda.Universe(PRM7)

guessed_elements = guess_types(u.atoms.names)

u.add_TopologyAttr('elements', guessed_elements)

print(u.atoms.elements)  # returns an array of guessed elements

More information on adding topology attributes can found in the user guide.

MDAnalysis.topology.guessers.get_atom_mass(element)[source]

Return the atomic mass in u for element.

Masses are looked up in MDAnalysis.topology.tables.masses.

Warning

Unknown masses are set to 0.0

Changed in version 0.20.0: Try uppercase atom type name as well

MDAnalysis.topology.guessers.guess_angles(bonds)[source]

Given a list of Bonds, find all angles that exist between atoms.

Works by assuming that if atoms 1 & 2 are bonded, and 2 & 3 are bonded, then (1,2,3) must be an angle.

Returns:

List of tuples defining the angles. Suitable for use in u._topology

Return type:

list of tuples

See also

guess_bonds()

MDAnalysis.topology.guessers.guess_aromaticities(atomgroup)[source]

Guess aromaticity of atoms using RDKit

Parameters:

atomgroup (mda.core.groups.AtomGroup) – Atoms for which the aromaticity will be guessed

Returns:

aromaticities – Array of boolean values for the aromaticity of each atom

Return type:

numpy.ndarray

New in version 2.0.0.

MDAnalysis.topology.guessers.guess_atom_charge(atomname)[source]

Guess atom charge from the name.

Warning

Not implemented; simply returns 0.

MDAnalysis.topology.guessers.guess_atom_element(atomname)[source]

Guess the element of the atom from the name.

Looks in dict to see if element is found, otherwise it uses the first character in the atomname. The table comes from CHARMM and AMBER atom types, where the first character is not sufficient to determine the atom type. Some GROMOS ions have also been added.

MDAnalysis.topology.guessers.guess_atom_mass(atomname)[source]

Guess a mass based on the atom name.

guess_atom_element() is used to determine the kind of atom.

Warning

Anything not recognized is simply set to 0; if you rely on the masses you might want to double check.

MDAnalysis.topology.guessers.guess_atom_type(atomname)[source]

Guess atom type from the name.

At the moment, this function simply returns the element, as guessed by guess_atom_element().

MDAnalysis.topology.guessers.guess_bonds(atoms, coords, box=None, **kwargs)[source]

Guess if bonds exist between two atoms based on their distance.

Bond between two atoms is created, if the two atoms are within

\[d < f \cdot (R_1 + R_2)\]

of each other, where \(R_1\) and \(R_2\) are the VdW radii of the atoms and \(f\) is an ad-hoc fudge_factor. This is the same algorithm that VMD uses.

Parameters:
  • atoms (AtomGroup) – atoms for which bonds should be guessed

  • coords (array) – coordinates of the atoms (i.e., AtomGroup.positions))

  • fudge_factor (float, optional) – The factor by which atoms must overlap eachother to be considered a bond. Larger values will increase the number of bonds found. [0.55]

  • vdwradii (dict, optional) – To supply custom vdwradii for atoms in the algorithm. Must be a dict of format {type:radii}. The default table of van der Waals radii is hard-coded as MDAnalysis.topology.tables.vdwradii. Any user defined vdwradii passed as an argument will supercede the table values. [None]

  • lower_bound (float, optional) – The minimum bond length. All bonds found shorter than this length will be ignored. This is useful for parsing PDB with altloc records where atoms with altloc A and B maybe very close together and there should be no chemical bond between them. [0.1]

  • box (array_like, optional) – Bonds are found using a distance search, if unit cell information is given, periodic boundary conditions will be considered in the distance search. [None]

Returns:

List of tuples suitable for use in Universe topology building.

Return type:

list

Warning

No check is done after the bonds are guessed to see if Lewis structure is correct. This is wrong and will burn somebody.

Raises:

ValueError

New in version 0.7.7.

Changed in version 0.9.0: Updated method internally to use more numpy, should work faster. Should also use less memory, previously scaled as \(O(n^2)\). vdwradii argument now augments table list rather than replacing entirely.

MDAnalysis.topology.guessers.guess_dihedrals(angles)[source]

Given a list of Angles, find all dihedrals that exist between atoms.

Works by assuming that if (1,2,3) is an angle, and 3 & 4 are bonded, then (1,2,3,4) must be a dihedral.

Returns:

  • list of tuples – List of tuples defining the dihedrals. Suitable for use in u._topology

  • .. versionadded 0.9.0

MDAnalysis.topology.guessers.guess_gasteiger_charges(atomgroup)[source]

Guess Gasteiger partial charges using RDKit

Parameters:

atomgroup (mda.core.groups.AtomGroup) – Atoms for which the charges will be guessed

Returns:

charges – Array of float values representing the charge of each atom

Return type:

numpy.ndarray

New in version 2.0.0.

MDAnalysis.topology.guessers.guess_improper_dihedrals(angles)[source]

Given a list of Angles, find all improper dihedrals that exist between atoms.

Works by assuming that if (1,2,3) is an angle, and 2 & 4 are bonded, then (2, 1, 3, 4) must be an improper dihedral. ie the improper dihedral is the angle between the planes formed by (1, 2, 3) and (1, 3, 4)

Returns:

  • List of tuples defining the improper dihedrals.

  • Suitable for use in u._topology

MDAnalysis.topology.guessers.guess_masses(atom_types)[source]

Guess the mass of many atoms based upon their type

Parameters:

atom_types – Type of each atom

Returns:

atom_masses

Return type:

np.ndarray dtype float64

MDAnalysis.topology.guessers.guess_types(atom_names)[source]

Guess the atom type of many atoms based on atom name

Parameters:

atom_names – Name of each atom

Returns:

atom_types

Return type:

np.ndarray dtype object

MDAnalysis.topology.guessers.validate_atom_types(atom_types)[source]

Vaildates the atom types based on whether they are available in our tables

Parameters:

atom_types – Type of each atom

Returns:

  • None

  • .. versionchanged:: 0.20.0 – Try uppercase atom type name as well