"""
Utilities for creating rotamers of eta-complexes.
Copyright Schrodinger, LLC. All rights reserved.
"""
from schrodinger import structure
from schrodinger.application.matsci import buildcomplex
from schrodinger.application.matsci import etatoggle as tesu
from schrodinger.structutils import analyze
from schrodinger.structutils import rings
from schrodinger.structutils import transform
METAL_INDEX_PROP = 'i_matsci_metal_index'
BOND_PROP = 'i_matsci_metal_neighbor'
ZERO = structure.BondType.Zero
ONE = structure.BondType.Single
LIGNUM_PROP = 'i_matsci_ligand_number'
INTERVENING_ATOM_PROP = 'b_matsci_intervening_atom'
[docs]class EtaRotamersException(Exception):
pass
[docs]class HapticLigand(object):
"""
Manages manipulation of a haptic ligand
"""
[docs] def __init__(self, struct, num):
"""
Create a HapticLigand object
:type struct: `schrodinger.structure.Structure`
:param struct: The structure containing the ligand
:type num: int
:param num: The value of the LIGNUM_PROP property for atoms in this
ligand
"""
# All atoms in this ligand
self.atoms = []
# The haptic ring atoms in this ligand
self.ring_atoms = []
for atom in struct.atom:
if atom.property.get(LIGNUM_PROP) == num:
self.atoms.append(atom)
if atom.property.get(BOND_PROP) is not None:
self.ring_atoms.append(atom)
self.numat = len(self.atoms)
self.struct = struct
[docs] def addRotationAtoms(self):
"""
We add two atoms that help define the torsion for rotation.
1) The centroid of the ring. We'll rotate the ring about the
centroid-metal axis
2) A fake atom as the same location as one of the ring atoms. We'll
rotate the ring a number of degrees relative to this atom
"""
ring_indexes = [x.index for x in self.ring_atoms]
centroidxyz = transform.get_centroid(self.struct,
atom_list=ring_indexes)[:3]
self.cent_index = self.struct.addAtom('C', *centroidxyz).index
self.marker_index = self.struct.addAtom('C',
*self.ring_atoms[0].xyz).index
# Must add a bond to a ring atom so that the ring rotates when the
# metal-marker axis rotates via the struct.adjust method
self.struct.addBond(self.cent_index, self.ring_atoms[0], ZERO)
[docs] def defineTorsionIndexes(self):
"""
Define the torsion as the four atoms in the following order:
1) An arbitrary haptic atom in this ligand
2) The centroid of this ligand
3) The metal atom
4) An arbitrary haptic atom in the other ligand
"""
self.torsion_indexes = [
self.ring_atoms[0].index, 1, self.cent_index, self.marker_index
]
[docs] def createRotatedStructures(self, base_title, n_rotamers=None):
"""
Create all the rotated ligand structures by rotating the ring.
:type base_title: str
:param base_title: The base title for the structures
:type n_rotamers: int or None
:param n_rotamers: the number of rotamers to return, if None
then it is determined as twice the number of ring atoms
:rtype: list
:return: Each item of the list is a rotated structure
"""
if not n_rotamers:
n_rotamers = 2 * len(self.ring_atoms)
stepsize = 360.0 / n_rotamers
structs = []
for step in range(n_rotamers):
rotation = step * stepsize
struct = self.struct.copy()
# We have to remove the metal-ligand bonds so that the moving atoms
# aren't considered part of a non-moving ring
bond_orders = []
metal_index = self.struct.property[METAL_INDEX_PROP]
for atom in self.ring_atoms:
struct.deleteBond(atom, metal_index)
# Perform the rotation of the now freely-rotating ligand
struct.adjust(rotation, *self.torsion_indexes)
# Restore the bonds
for atom in self.ring_atoms:
struct.addBond(atom, metal_index, atom.property[BOND_PROP])
# Delete the rotation marker atoms
struct.deleteAtoms([self.cent_index, self.marker_index])
struct.title = base_title + '_%d' % rotation
structs.append(struct)
# Clean up the new geometries
buildcomplex.minimize_complexes(structs)
return structs
[docs] def getRotamers(self, base_title, n_rotamers=None):
"""
Return the rotamers.
:type base_title: str
:param base_title: The base title for the structures
:type n_rotamers: int or None
:param n_rotamers: the number of rotamers to return, if None
then it is determined as twice the number of ring atoms
:rtype: list
:return: Each item of the list is a rotated structure
"""
self.addRotationAtoms()
self.defineTorsionIndexes()
return self.createRotatedStructures(base_title, n_rotamers=n_rotamers)
[docs]def create_rotated_complexes(struct,
lignums,
base_title,
rotate_both_ligands=True,
n_rotamers=None):
"""
Create a series of complexes by rotating one haptic ligand around its
centroid-metal axis
:type struct: `schrodinger.structure.Structure`
:param struct: The complex
:type lignums: list
:param lignums: The LIGNUM_PROP vals for atoms in haptic ligands
:type base_title: str
:param base_title: The base title for the structures
:type rotate_both_ligands: bool
:param rotate_both_ligands: whether to rotate both ligands
:type n_rotamers: int or None
:param n_rotamers: the number of rotamers per ligand, if None
then it is determined as twice the number of ring atoms
:rtype: list
:return: Each item of the list is a rotated structure
"""
ligands = [HapticLigand(struct, x) for x in lignums]
if len(ligands) == 1:
return ligands[0].getRotamers(base_title, n_rotamers=n_rotamers)
smaller_idx = int(ligands[0].numat >= ligands[1].numat)
rotamers = ligands[smaller_idx].getRotamers(base_title,
n_rotamers=n_rotamers)
if not rotate_both_ligands:
return rotamers
larger_idx = int(not smaller_idx)
all_rotamers = []
for rotamer in rotamers:
ligand = HapticLigand(rotamer, lignums[larger_idx])
all_rotamers.extend(
ligand.getRotamers(rotamer.title, n_rotamers=n_rotamers))
return all_rotamers
[docs]def get_rotatable_haptic_ligands(st, only_rings=True):
"""
Return rotatable haptic ligand molecules in the given structure.
:type st: schrodinger.structure.Structure
:param st: the structure
:type only_rings: bool
:param only_rings: if True then only allow rotation of eta-bound rings,
if False then also allow rotation of ligands where the eta-bound
motif is acyclic, for example ethene, etc.
:rtype: list
:return: contains schrodinger.structure._Molecule
"""
haptic_ligs = []
for mol in st.molecule:
coordinating_idxs = [
atom.index
for atom in mol.atom
if atom.property.get(BOND_PROP) is not None
]
# skip metal atom
if not coordinating_idxs:
continue
# skip multi-dentate ligands, for example tethered ligands
coordinating_groups = analyze.group_by_connectivity(
st, coordinating_idxs)
if len(coordinating_groups) > 1:
continue
# optionally skip acyclic ligands, for example ethene
sssr_rings = rings.find_rings(st)
try:
ring_idxs = [
idxs for idxs in sssr_rings
if set(coordinating_idxs).issubset(idxs)
][0]
except IndexError:
ring_idxs = []
if only_rings and not ring_idxs:
continue
# MATSCI-8384 - skip rotation of certain eta-bound rings
if ring_idxs and len(coordinating_idxs) < len(ring_idxs) - 1:
continue
n_zob = 0
for atom in mol.atom:
atom.property[LIGNUM_PROP] = mol.number
if atom.property.get(BOND_PROP) == 0:
n_zob += 1
if n_zob >= 2:
haptic_ligs.append(mol)
return haptic_ligs
[docs]def get_rotamers(struct,
only_rings=True,
rotate_both_ligands=True,
n_rotamers=None):
"""
Return the rotamers for the given eta-complex.
:type struct: schrodinger.structure.Structure
:param struct: the structure
:type only_rings: bool
:param only_rings: if True then only allow rotation of eta-bound rings,
if False then also allow rotation of ligands where the eta-bound
motif is acyclic, for example ethene, etc.
:type rotate_both_ligands: bool
:param rotate_both_ligands: whether to rotate both ligands
:type n_rotamers: int or None
:param n_rotamers: the number of rotamers per ligand, if None
then it is determined as twice the number of ring atoms
:raise EtaRotamersException: if there is an issue
:rtype: list
:return: contains schrodinger.structure.Structure of rotamers
"""
st = struct.copy()
# get metal
metal_index = find_metal(st)
metal_atom = st.atom[metal_index]
# toggle representation
dummy_atom_repr = False
for neighbor in metal_atom.bonded_atoms:
if neighbor.atomic_number in tesu.DUMMY_ATOMIC_NUMBERS:
dummy_atom_repr = True
break
if dummy_atom_repr:
try:
tesu.toggle_structure(st)
except ValueError as msg:
raise EtaRotamersException(str(msg))
# mark the atoms bound to the metal
st.property[METAL_INDEX_PROP] = metal_index
neighbors = []
for neighbor in metal_atom.bonded_atoms:
index = neighbor.index
bond = st.getBond(metal_index, index)
neighbor.property[BOND_PROP] = bond.order
neighbors.append(neighbor)
# deleting all the bonds to the metal leaves each ligand as a separate
# molecule
for neighbor in neighbors:
metal_atom.deleteBond(neighbor)
# define haptic ligands
haptic_ligs = [
mol.number
for mol in get_rotatable_haptic_ligands(st, only_rings=only_rings)
]
# ensure hapticity
if len(haptic_ligs) > 2:
raise EtaRotamersException('More than 2 separate haptic ligands found')
elif not haptic_ligs:
raise EtaRotamersException('No mono-dentate haptic ligands found')
# rebond the metal atom
for neighbor in neighbors:
metal_atom.addBond(neighbor, neighbor.property[BOND_PROP])
# get rotamers
base_title = st.title
rotamers = create_rotated_complexes(st,
haptic_ligs,
st.title,
rotate_both_ligands=rotate_both_ligands,
n_rotamers=n_rotamers)
# toggle representation
if dummy_atom_repr:
for rotamer in rotamers:
try:
tesu.toggle_structure(rotamer)
except ValueError as msg:
raise EtaRotamersException(str(msg))
return rotamers