ligand and submodule contents

ligand

Ligand support functions.

Todo

Some of the definitions in this module belong in a configuration file other than here.

Code author: Jens Erik Nielsen

pdb2pqr.ligand.ELEMENT_BY_GROUP = {1: ['H', 'Li', 'Na', 'K', 'Rb', 'Cs', 'Fr'], 2: ['Be', 'Mg', 'Ca', 'Sr', 'Ba', 'Ra'], 13: ['B', 'Al', 'Ga', 'In', 'Tl', 'Nh'], 14: ['C', 'Si', 'Ge', 'Sn', 'Pb', 'Fl'], 15: ['N', 'P', 'As', 'Sb', 'Bi', 'Mc'], 16: ['O', 'S', 'Se', 'Te', 'Po', 'Lv'], 17: ['F', 'Cl', 'Br', 'I', 'At', 'Ts'], 18: ['He', 'Ne', 'Ar', 'Kr', 'Xe', 'Rn', 'Og']}

Groups of the periodic table

pdb2pqr.ligand.NONBONDED_BY_TYPE = {'Al': 0, 'Br': 6, 'C.1': 0, 'C.2': 0, 'C.3': 0, 'C.ar': 0, 'Ca': 0, 'Cl': 6, 'F': 6, 'H': 0, 'I': 6, 'K': 0, 'Li': 0, 'N.1': 2, 'N.2': 2, 'N.3': 2, 'N.4': 0, 'N.am': 0, 'N.ar': 2, 'N.pl3': 0, 'Na': 0, 'O.2': 4, 'O.3': 4, 'O.co2': 4.5, 'P.3': 0, 'S.2': 4, 'S.3': 4, 'S.o': 2, 'S.o2': 0, 'Si': 0}

Numbers of non-bonded electrons for Sybyl-type atoms. Adapted from https://onlinelibrary.wiley.com/doi/abs/10.1002/jcc.540100804 (Table I).

pdb2pqr.ligand.RADII = {'bondi': {'Ar': 1.88, 'As': 1.85, 'Br': 1.85, 'C': 1.7, 'Cl': 1.75, 'F': 1.47, 'H': 1.2, 'He': 1.4, 'I': 1.98, 'Kr': 2.02, 'N': 1.55, 'Ne': 1.54, 'O': 1.52, 'P': 1.8, 'S': 1.8, 'Se': 1.9, 'Si': 2.1, 'Te': 2.06, 'Xe': 2.16}, 'bondi-zap': {'C': 1.7, 'Cl': 1.75, 'F': 1.47, 'H': 1.2, 'I': 1.98, 'N': 1.55, 'O.co2': 1.52, 'S': 1.8}, 'not parse - do not use': {'Br': 2.5, 'C': 1.7, 'Cl': 1.75, 'F': 1.2, 'H': 1.0, 'N': 1.5, 'O': 1.4, 'P': 1.9, 'S': 1.85}, 'parse': {'C': 1.7, 'C.1': 2.0, 'C.2': 2.0, 'C.3': 2.0, 'H': 1.0, 'N': 1.5, 'O': 1.4, 'S': 1.85}, 'zap9': {'C': 1.87, 'Cl': 1.82, 'F': 2.4, 'H': 1.1, 'I': 2.65, 'N': 1.4, 'O.co2': 1.76, 'S': 2.15}}

Radii for different atom types. When using these tables, the most specific Sybyl atom type should be used first and then the generic element should be used

pdb2pqr.ligand.VALENCE_BY_ELEMENT = {'Al': 3, 'Ar': 8, 'As': 5, 'At': 7, 'B': 3, 'Ba': 2, 'Be': 2, 'Bi': 5, 'Br': 7, 'C': 4, 'Ca': 2, 'Cl': 7, 'Cs': 1, 'F': 7, 'Fl': 4, 'Fr': 1, 'Ga': 3, 'Ge': 4, 'H': 1, 'He': 8, 'I': 7, 'In': 3, 'K': 1, 'Kr': 8, 'Li': 1, 'Lv': 6, 'Mc': 5, 'Mg': 2, 'N': 5, 'Na': 1, 'Ne': 8, 'Nh': 3, 'O': 6, 'Og': 8, 'P': 5, 'Pb': 4, 'Po': 6, 'Ra': 2, 'Rb': 1, 'Rn': 8, 'S': 6, 'Sb': 5, 'Se': 6, 'Si': 4, 'Sn': 4, 'Sr': 2, 'Te': 6, 'Tl': 3, 'Ts': 7, 'Xe': 8}

Valence electrons by element

pdb2pqr.ligand.VALENCE_BY_GROUP = {1: 1, 2: 2, 13: 3, 14: 4, 15: 5, 16: 6, 17: 7, 18: 8}

Numbers of valence electrons for the groups of the periodic table

hydrogens.ligand.mol2

Support molecules in Tripos MOL2 format.

For further information look at (web page exists: 25 August 2005): http://www.tripos.com/index.php?family=modules,SimplePage,,,&page=sup_mol2&s=0

class pdb2pqr.ligand.mol2.Mol2Atom[source]

MOL2 molecule atoms.

__init__()[source]

Initialize self. See help(type(self)) for accurate signature.

assign_radius(primary_dict, secondary_dict)[source]

Assign radius to atom.

Todo

It seems inconsistent that this function pulls radii from a dictionary and the biomolecule routines use force field files.

Parameters:
  • primary_dict (dict) – primary dictionary of radii indexed by atom type or element
  • secondary_dict (dict) – backup dictionary for radii not found in primary dictionary
bond_order

Total number of electrons in bonds with other atoms.

Returns:total number of electrons in bonds with other atoms
Return type:int
bonded_atom_names

Bonded atom names.

Returns:bonded atom names
Return type:list
coords

Coordinates.

Returns:coordinates
Return type:numpy.ndarray
distance(other)[source]

Get distance between two atoms.

Parameters:other (Mol2Atom) – other atom for distance measurement
Returns:distance
Return type:float
element

Element for this atom (uppercase).

Returns:element for this atom
Return type:str
formal_charge

Formal charge for this atom

Returns:formal charge for this atom
Return type:int
num_bonded_heavy

Number of heavy atoms bonded to this atom.

Returns:number of heavy atoms
Return type:int
num_bonded_hydrogen

Number of hydrogen atoms bonded to this atom.

Returns:number of hydrogen atoms
Return type:int
class pdb2pqr.ligand.mol2.Mol2Bond(atom1, atom2, bond_type, bond_id=0)[source]

MOL2 molecule bonds.

__init__(atom1, atom2, bond_type, bond_id=0)[source]

Initialize bond.

Parameters:
  • atom1 (str) – name of first atom in bond
  • atom2 (str) – name of second atom in bond
  • bond_type (int) – type of bond: 1 (single), 2 (double), or ar (aromatic)
  • bond_id (int) – integer ID of bond
atom_names

Get atom names in bond.

Returns:tuple with names of atoms in bond.
Return type:(str, str)
length

Get bond length.

Returns:bond length
Return type:float
class pdb2pqr.ligand.mol2.Mol2Molecule[source]

Tripos MOL2 molecule.

__init__()[source]

Initialize self. See help(type(self)) for accurate signature.

assign_charges()[source]

Assign charges to atoms in molecule.

assign_parameters(primary_dict={'C': 1.87, 'Cl': 1.82, 'F': 2.4, 'H': 1.1, 'I': 2.65, 'N': 1.4, 'O.co2': 1.76, 'S': 2.15}, secondary_dict={'Ar': 1.88, 'As': 1.85, 'Br': 1.85, 'C': 1.7, 'Cl': 1.75, 'F': 1.47, 'H': 1.2, 'He': 1.4, 'I': 1.98, 'Kr': 2.02, 'N': 1.55, 'Ne': 1.54, 'O': 1.52, 'P': 1.8, 'S': 1.8, 'Se': 1.9, 'Si': 2.1, 'Te': 2.06, 'Xe': 2.16})[source]

Assign charges and radii to atoms in molecule.

Parameters:
  • primary_dict – primary dictionary of radii indexed by atom type or element
  • secondary_dict – backup dictionary for radii not found in primary dictionary
assign_radii(primary_dict, secondary_dict)[source]

Assign radii to atoms in molecule.

Parameters:
  • primary_dict (dict) – primary dictionary of radii indexed by atom type or element
  • secondary_dict (dict) – backup dictionary for radii not found in primary dictionary
find_atom_torsions(start_atom)[source]

Set the torsion angles that start with this atom (name).

Parameters:start_atom (str) – starting atom name
Returns:list of 4-tuples containing atom names comprising torsions
find_new_rings(path, rings, level=0)[source]

Find new rings in molecule.

This was borrowed from StackOverflow: https://j.mp/2AHaukj

Parameters:
  • path (list of str) – list of atom names
  • rings (list of str) – current list of rings
  • level (int) – recursion level
Returns:

new list of rings

Return type:

int

parse_atoms(mol2_file)[source]

Parse @<TRIPOS>ATOM section of file.

Parameters:

mol2_file – file-like object with MOL2 data

Returns:

file-like object advanced to bonds section

Raises:
parse_bonds(mol2_file)[source]

Parse @<TRIPOS>BOND section of file.

Atoms must already have been parsed. Also sets up torsions and rings.

Parameters:mol2_file – file-like object with MOL2 data
Returns:file-like object advanced to SUBSTRUCTURE section
read(mol2_file)[source]

Routines for reading MOL2 file.

Parameters:mol2_file – file-like object with MOL2 data
static rotate_to_smallest(path)[source]

Rotate cycle path so that it begins with the smallest node.

This was borrowed from StackOverflow: https://j.mp/2AHaukj

Parameters:path (list of str) – list of atom names
Returns:rotated path
Return type:list of str
set_rings()[source]

Set all rings in molecule.

This was borrowed from StackOverflow: https://j.mp/2AHaukj

set_torsions()[source]

Set all torsions in molecule.

hydrogens.ligand.peoe

Implements the PEOE method.

The PEOE method is described in: Paul Czodrowski, Ingo Dramburg, Christoph A. Sotriffer Gerhard Klebe. Development, validation, and application of adapted PEOE charges to estimate pKa values of functional groups in protein–ligand complexes. Proteins, 65, 424-437, 2006. https://doi.org/10.1002/prot.21110

pdb2pqr.ligand.peoe.assign_terms(atoms, term_dict)[source]

Assign polynomial terms to each atom.

Parameters:
  • atoms (list) – list of Mol2Atom atoms
  • term_dict (dict) – dictionary of polynomial terms
Returns:

modified list of atoms

Return type:

list

pdb2pqr.ligand.peoe.electronegativity(charge, poly_terms, atom_type)[source]

Calculate the electronegativity.

Calculation is based on a third-order polynomial in the atomic charge as described in Equation 2 of https://doi.org/10.1002/prot.21110.

Parameters:
  • charge (float) – charge of atom
  • poly_terms – polynomial terms ordered from 0th- to 3rd-order
  • atom_type (str) – string with atom type
Returns:

electronegativity value

Return type:

float

Raises:

IndexError – if incorrect number of poly_terms given

pdb2pqr.ligand.peoe.equilibrate(atoms, damp=0.778, scale=1.56, num_cycles=6, term_dict={'BR': (10.88, 8.47, 1.16, 19.71), 'C.1': (10.39, 9.45, 0.73, 20.57), 'C.2': (9.29, 9.32, 1.51, 19.62), 'C.3': (7.98, 9.18, 1.88, 19.04), 'C.AR': (8.530000000000001, 9.18, 1.88, 19.04), 'C.CAT': (7.98, 9.18, 1.88, 19.04), 'CL': (10.38, 9.69, 1.35, 22.04), 'F': (12.36, 13.85, 2.31, 30.82), 'H': (7.17, 6.24, -0.56, 12.85), 'I': (10.9, 7.96, 0.96, 18.82), 'N.1': (15.68, 11.7, -0.27, 27.11), 'N.2': (12.87, 11.15, 0.85, 24.87), 'N.3': (17.54, 10.28, 1.36, 28.0), 'N.4': (17.54, 10.28, 1.36, 28.0), 'N.AM': (16.369999999999997, 11.15, 0.85, 24.87), 'N.AR': (11.579999999999998, 11.15, 0.85, 24.87), 'N.PL3': (13.37, 11.15, 0.85, 24.87), 'O.2': (14.18, 12.92, 1.39, 28.49), 'O.3': (11.08, 12.92, 1.39, 28.49), 'O.CO2': (15.25, 13.79, 0.47, 31.33), 'O.OH': (14.98, 12.92, 1.39, 28.49), 'P.3': (10.63, 9.13, 1.38, 20.65), 'S.2': (10.63, 9.13, 1.38, 20.65), 'S.3': (10.63, 9.13, 1.38, 20.65), 'S.O2': (10.63, 9.13, 1.38, 20.65)})[source]

Equilibrate the atomic charges.

Parameters:
  • atoms (list) – list of Mol2Atom atoms to equilibrate
  • damp (float) – damping factor for equilibration process
  • scale (float) – scaling factor for equilibration process
  • num_cycles (int) – number of PEOE cycles
  • term_dict (dict) – dictionary of polynomial terms
Returns:

revised list of atoms

Return type:

list

hydrogens.ligand.topology

Ligand topology classes.

class pdb2pqr.ligand.topology.Topology(molecule)[source]

Ligand topology class.

__init__(molecule)[source]

Initialize with molecule.

Parameters:molecule (Mol2Molecule) – Mol2Molecule object