Source code for schrodinger.application.matsci.etarotamers
"""
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'
# this is a factor that multiplies the number of atoms in the
# ring to determine the number of rotamers sampled
N_ROTAMERS_F = 2
ROTATE_BOTH_LIGANDS = True
[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):
        """
        Create all the rotated ligand structures by rotating the ring.  The
        number of rotations is equal to N_ROTAMERS_F * num_ring_atoms and
        each rotation is 360 / (N_ROTAMERS_F * num_ring_atoms) degrees.
        :type base_title: str
        :param base_title: The base title for the structures
        :rtype: list
        :return: Each item of the list is a rotated structure
        """
        structs = []
        num_ringat = len(self.ring_atoms)
        stepsize = 360.0 / (N_ROTAMERS_F * num_ringat)
        for step in range(N_ROTAMERS_F * num_ringat):
            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])
            # Clean up the new geometry
            buildcomplex.minimize_complex(struct)
            struct.title = base_title + '_%d' % rotation
            structs.append(struct)
        return structs 
[docs]    def getRotamers(self, base_title):
        """
        Return the rotamers.
        :type base_title: str
        :param base_title: The base title for the structures
        :rtype: list
        :return: Each item of the list is a rotated structure
        """
        self.addRotationAtoms()
        self.defineTorsionIndexes()
        return self.createRotatedStructures(base_title)  
[docs]def create_rotated_complexes(struct,
                             lignums,
                             base_title,
                             rotate_both_ligands=ROTATE_BOTH_LIGANDS):
    """
    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
    :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)
    smaller_idx = int(ligands[0].numat >= ligands[1].numat)
    rotamers = ligands[smaller_idx].getRotamers(base_title)
    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))
    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):
    """
    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.
    :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)
    # 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