Package rdkit :: Package Chem :: Module AllChem
[hide private]
[frames] | no frames]

Module AllChem

source code

Import all RDKit chemistry modules

Functions [hide private]
 
TransformMol(mol, tform, confId=-1, keepConfs=False)
Applies the transformation (usually a 4x4 double matrix) to a molecule...
source code
 
ComputeMolShape(mol, confId=-1, boxDim=(20, 20, 20), spacing=0.5, **kwargs)
returns a grid representation of the molecule's shape...
source code
 
ComputeMolVolume(mol, confId=-1, gridSpacing=0.2, boxMargin=2.0)
Calculates the volume of a particular conformer of a molecule based on a grid-encoding of the molecular shape.
source code
 
GenerateDepictionMatching2DStructure(mol, reference, confId=-1, referencePattern=None, acceptFailure=False, **kwargs)
Generates a depiction for a molecule where a piece of the molecule is constrained to have the same coordinates as a reference.
source code
 
GenerateDepictionMatching3DStructure(mol, reference, confId=-1, **kwargs)
Generates a depiction for a molecule where a piece of the molecule is constrained to have coordinates similar to those of a 3D reference structure.
source code
 
GetBestRMS(ref, probe, refConfId=-1, probeConfId=-1, maps=None)
Returns the optimal RMS for aligning two molecules, taking symmetry into account.
source code
 
EnumerateLibraryFromReaction(reaction, sidechainSets)
Returns a generator for the virtual library defined by a reaction and a sequence of sidechain sets
source code
 
ConstrainedEmbed(mol, core, useTethers=True, coreConfId=-1, randomseed=2342, getForceField=<Boost.Python.function object at 0x418bab0>, **kwargs)
generates an embedding of a molecule where part of the molecule is constrained to have particular coordinates
source code
 
AssignBondOrdersFromTemplate(refmol, mol)
assigns bond orders to a molecule based on the bond orders in a template molecule
source code
 
_test() source code
Variables [hide private]
  logger = logger()
  LayeredFingerprint_substructLayers = 7
  __package__ = 'rdkit.Chem'

Imports: SupplierFromFilename, InchiToInchiKey, rdGeometry, DataStructs, rdchem, INCHI_AVAILABLE, QuickSmartsMatch, rdBase, pyPeriodicTable, CanonSmiles, MolFromInchi, MolToInchiAndAuxInfo, RDConfig, FindMolChiralCenters, InchiReadWriteError, MolToInchi, MCFF_GetFeaturesForMol, ForceField, numpy, os, AddHs, AddRecursiveQuery, AlignMol, AlignMolConformers, AssignAtomChiralTagsFromStructure, AssignRadicals, AssignStereochemistry, Atom, AtomMonomerInfo, AtomMonomerType, AtomNumEqualsQueryAtom, AtomNumGreaterQueryAtom, AtomNumLessQueryAtom, AtomPDBResidueInfo, AtomPairsParameters, BRICS, Bond, BondDir, BondStereo, BondType, BuildFeatureFactory, BuildFeatureFactoryFromString, BuildFragmentCatalog, CalcChi0n, CalcChi0v, CalcChi1n, CalcChi1v, CalcChi2n, CalcChi2v, CalcChi3n, CalcChi3v, CalcChi4n, CalcChi4v, CalcChiNn, CalcChiNv, CalcCrippenDescriptors, CalcExactMolWt, CalcFractionCSP3, CalcHallKierAlpha, CalcKappa1, CalcKappa2, CalcKappa3, CalcLabuteASA, CalcMolFormula, CalcNumAliphaticCarbocycles, CalcNumAliphaticHeterocycles, CalcNumAliphaticRings, CalcNumAmideBonds, CalcNumAromaticCarbocycles, CalcNumAromaticHeterocycles, CalcNumAromaticRings, CalcNumHBA, CalcNumHBD, CalcNumHeteroatoms, CalcNumLipinskiHBA, CalcNumLipinskiHBD, CalcNumRings, CalcNumRotatableBonds, CalcNumSaturatedCarbocycles, CalcNumSaturatedHeterocycles, CalcNumSaturatedRings, CalcTPSA, CanonicalizeConformer, CanonicalizeMol, ChemicalFeatures, ChemicalReaction, ChiralType, Cleanup, CombineMols, CompositeQueryType, Compute2DCoords, Compute2DCoordsForReaction, Compute2DCoordsMimicDistmat, ComputeCanonicalTransform, ComputeCentroid, ComputeConfBox, ComputeConfDimsAndOffset, ComputeGasteigerCharges, ComputeUnionBox, Conformer, Crippen, DeleteSubstructs, Descriptors, Draw, EState, EditableMol, EmbedMolecule, EmbedMultipleConfs, EncodeShape, ExplicitDegreeEqualsQueryAtom, ExplicitDegreeGreaterQueryAtom, ExplicitDegreeLessQueryAtom, ExplicitValenceEqualsQueryAtom, ExplicitValenceGreaterQueryAtom, ExplicitValenceLessQueryAtom, FastFindRings, FeatFinderCLI, FindAllPathsOfLengthN, FindAllSubgraphsOfLengthMToN, FindAllSubgraphsOfLengthN, FindAtomEnvironmentOfRadiusN, FindUniqueSubgraphsOfLengthN, FormalChargeEqualsQueryAtom, FormalChargeGreaterQueryAtom, FormalChargeLessQueryAtom, ForwardSDMolSupplier, FragmentCatalog, FragmentMatcher, FragmentOnBRICSBonds, FragmentOnBonds, Fragments, FreeChemicalFeature, FunctionalGroups, Get3DDistanceMatrix, GetAdjacencyMatrix, GetAlignmentTransform, GetAngleDeg, GetAngleRad, GetAtomMatch, GetAtomPairAtomCode, GetAtomPairFingerprint, GetBondLength, GetConnectivityInvariants, GetCrippenO3A, GetDihedralDeg, GetDihedralRad, GetDistanceMatrix, GetFeatureInvariants, GetFormalCharge, GetHashedAtomPairFingerprint, GetHashedAtomPairFingerprintAsBitVect, GetHashedMorganFingerprint, GetHashedTopologicalTorsionFingerprint, GetHashedTopologicalTorsionFingerprintAsBitVect, GetMACCSKeysFingerprint, GetMolFrags, GetMoleculeBoundsMatrix, GetMorganFingerprint, GetMorganFingerprintAsBitVect, GetO3A, GetPeriodicTable, GetSSSR, GetSymmSSSR, GetTopologicalTorsionFingerprint, GraphDescriptors, Graphs, HCountEqualsQueryAtom, HCountGreaterQueryAtom, HCountLessQueryAtom, HybridizationEqualsQueryAtom, HybridizationGreaterQueryAtom, HybridizationLessQueryAtom, HybridizationType, InNRingsEqualsQueryAtom, InNRingsGreaterQueryAtom, InNRingsLessQueryAtom, IsAliphaticQueryAtom, IsAromaticQueryAtom, IsInRingQueryAtom, IsUnsaturatedQueryAtom, IsotopeEqualsQueryAtom, IsotopeGreaterQueryAtom, IsotopeLessQueryAtom, Kekulize, LayeredFingerprint, Lipinski, MCS, MMFFGetMoleculeForceField, MMFFGetMoleculeProperties, MMFFHasAllMoleculeParams, MMFFOptimizeMolecule, MMFFSanitizeMolecule, MQNs_, MassEqualsQueryAtom, MassGreaterQueryAtom, MassLessQueryAtom, MergeQueryHs, MinRingSizeEqualsQueryAtom, MinRingSizeGreaterQueryAtom, MinRingSizeLessQueryAtom, Mol, MolAddRecursiveQueries, MolCatalog, MolChemicalFeature, MolChemicalFeatureFactory, MolFragmentToSmiles, MolFromMol2Block, MolFromMol2File, MolFromMolBlock, MolFromMolFile, MolFromPDBBlock, MolFromPDBFile, MolFromQuerySLN, MolFromSLN, MolFromSmarts, MolFromSmiles, MolFromTPLBlock, MolFromTPLFile, MolSurf, MolToMolBlock, MolToMolFile, MolToPDBBlock, MolToPDBFile, MolToSmarts, MolToSmiles, MolToTPLBlock, MolToTPLFile, MurckoDecompose, O3A, PDBWriter, PEOE_VSA_, ParseMolQueryDefFile, PathToSubmol, PatternFingerprint, PeriodicTable, PropertyMol, PyMol, QueryAtom, RDKFingerprint, RandomTransform, ReactionFromRxnBlock, ReactionFromRxnFile, ReactionFromSmarts, ReactionToRxnBlock, ReactionToSmarts, RemoveHs, RemoveStereochemistry, RenumberAtoms, ReplaceCore, ReplaceSidechains, ReplaceSubstructs, RingBondCountEqualsQueryAtom, RingBondCountGreaterQueryAtom, RingBondCountLessQueryAtom, RingInfo, SATIS, SDMolSupplier, SDWriter, SMR_VSA_, SaltRemover, SanitizeFlags, SanitizeMol, Scaffolds, SetAngleDeg, SetAngleRad, SetAromaticity, SetBondLength, SetConjugation, SetDihedralDeg, SetDihedralRad, SetHybridization, ShapeProtrudeDist, ShapeTanimotoDist, ShowMols, SlogP_VSA_, SmilesMolSupplier, SmilesMolSupplierFromText, SmilesWriter, SplitMolByPDBChainId, SplitMolByPDBResidues, TDTMolSupplier, TDTWriter, TotalDegreeEqualsQueryAtom, TotalDegreeGreaterQueryAtom, TotalDegreeLessQueryAtom, TotalValenceEqualsQueryAtom, TotalValenceGreaterQueryAtom, TotalValenceLessQueryAtom, TransformConformer, UFFGetMoleculeForceField, UFFHasAllMoleculeParams, UFFOptimizeMolecule, WedgeMolBonds, fmcs, inchi, rdChemReactions, rdChemicalFeatures, rdDistGeom, rdForceFieldHelpers, rdMolAlign, rdMolCatalog, rdMolChemicalFeatures, rdMolDescriptors, rdMolTransforms, rdPartialCharges, rdSLNParse, rdShapeHelpers, rdfragcatalog, rdinchi, rdmolfiles, rdmolops, tossit


Function Details [hide private]

TransformMol(mol, tform, confId=-1, keepConfs=False)

source code 
Applies the transformation (usually a 4x4 double matrix) to a molecule
if keepConfs is False then all but that conformer are removed

ComputeMolShape(mol, confId=-1, boxDim=(20, 20, 20), spacing=0.5, **kwargs)

source code 
returns a grid representation of the molecule's shape
  

GenerateDepictionMatching2DStructure(mol, reference, confId=-1, referencePattern=None, acceptFailure=False, **kwargs)

source code 
Generates a depiction for a molecule where a piece of the molecule
   is constrained to have the same coordinates as a reference.

   This is useful for, for example, generating depictions of SAR data
   sets so that the cores of the molecules are all oriented the same
   way.

Arguments:
  - mol:          the molecule to be aligned, this will come back
                  with a single conformer.
  - reference:    a molecule with the reference atoms to align to;
                  this should have a depiction.
  - confId:       (optional) the id of the reference conformation to use
  - referencePattern:  (optional) an optional molecule to be used to
                       generate the atom mapping between the molecule
                       and the reference.
  - acceptFailure: (optional) if True, standard depictions will be generated
                   for molecules that don't have a substructure match to the
                   reference; if False, a ValueError will be raised

GenerateDepictionMatching3DStructure(mol, reference, confId=-1, **kwargs)

source code 
Generates a depiction for a molecule where a piece of the molecule
   is constrained to have coordinates similar to those of a 3D reference
   structure.

Arguments:
  - mol:          the molecule to be aligned, this will come back
                  with a single conformer.
  - reference:    a molecule with the reference atoms to align to;
                  this should have a depiction.
  - confId:       (optional) the id of the reference conformation to use

GetBestRMS(ref, probe, refConfId=-1, probeConfId=-1, maps=None)

source code 
Returns the optimal RMS for aligning two molecules, taking
symmetry into account. As a side-effect, the probe molecule is
left in the aligned state.

Arguments:
  - ref: the reference molecule
  - probe: the molecule to be aligned to the reference
  - refConfId: (optional) reference conformation to use
  - probeConfId: (optional) probe conformation to use
  - maps: (optional) a list of lists of (probeAtomId,refAtomId)
    tuples with the atom-atom mappings of the two molecules.
    If not provided, these will be generated using a substructure
    search.

EnumerateLibraryFromReaction(reaction, sidechainSets)

source code 
Returns a generator for the virtual library defined by
 a reaction and a sequence of sidechain sets

>>> from rdkit import Chem
>>> from rdkit.Chem import AllChem
>>> s1=[Chem.MolFromSmiles(x) for x in ('NC','NCC')]
>>> s2=[Chem.MolFromSmiles(x) for x in ('OC=O','OC(=O)C')]
>>> rxn = AllChem.ReactionFromSmarts('[O:2]=[C:1][OH].[N:3]>>[O:2]=[C:1][N:3]')
>>> r = AllChem.EnumerateLibraryFromReaction(rxn,[s2,s1])
>>> [Chem.MolToSmiles(x[0]) for x in list(r)]
['CNC=O', 'CCNC=O', 'CNC(C)=O', 'CCNC(C)=O']

Note that this is all done in a lazy manner, so "infinitely" large libraries can
be done without worrying about running out of memory. Your patience will run out first:

Define a set of 10000 amines:
>>> amines = (Chem.MolFromSmiles('N'+'C'*x) for x in range(10000))

... a set of 10000 acids
>>> acids = (Chem.MolFromSmiles('OC(=O)'+'C'*x) for x in range(10000))

... now the virtual library (1e8 compounds in principle):
>>> r = AllChem.EnumerateLibraryFromReaction(rxn,[acids,amines])

... look at the first 4 compounds:
>>> [Chem.MolToSmiles(r.next()[0]) for x in range(4)]
['NC=O', 'CNC=O', 'CCNC=O', 'CCCNC=O']

ConstrainedEmbed(mol, core, useTethers=True, coreConfId=-1, randomseed=2342, getForceField=<Boost.Python.function object at 0x418bab0>, **kwargs)

source code 
generates an embedding of a molecule where part of the molecule
is constrained to have particular coordinates

Arguments
  - mol: the molecule to embed
  - core: the molecule to use as a source of constraints
  - useTethers: (optional) if True, the final conformation will be
      optimized subject to a series of extra forces that pull the
      matching atoms to the positions of the core atoms. Otherwise
      simple distance constraints based on the core atoms will be
      used in the optimization.
  - coreConfId: (optional) id of the core conformation to use
  - randomSeed: (optional) seed for the random number generator


  An example, start by generating a template with a 3D structure:
  >>> from rdkit.Chem import AllChem
  >>> template = AllChem.MolFromSmiles("c1nn(Cc2ccccc2)cc1")
  >>> AllChem.EmbedMolecule(template)
  0
  >>> AllChem.UFFOptimizeMolecule(template)
  0

  Here's a molecule:
  >>> mol = AllChem.MolFromSmiles("c1nn(Cc2ccccc2)cc1-c3ccccc3")

  Now do the constrained embedding
  >>> newmol=AllChem.ConstrainedEmbed(mol, template)

  Demonstrate that the positions are the same:
  >>> newp=newmol.GetConformer().GetAtomPosition(0)
  >>> molp=mol.GetConformer().GetAtomPosition(0)
  >>> list(newp-molp)==[0.0,0.0,0.0]
  True
  >>> newp=newmol.GetConformer().GetAtomPosition(1)
  >>> molp=mol.GetConformer().GetAtomPosition(1)
  >>> list(newp-molp)==[0.0,0.0,0.0]
  True

AssignBondOrdersFromTemplate(refmol, mol)

source code 
assigns bond orders to a molecule based on the
    bond orders in a template molecule

Arguments
  - refmol: the template molecule
  - mol: the molecule to assign bond orders to

  An example, start by generating a template from a SMILES
  and read in the PDB structure of the molecule
  >>> from rdkit.Chem import AllChem
  >>> template = AllChem.MolFromSmiles("CN1C(=NC(C1=O)(c2ccccc2)c3ccccc3)N")
  >>> mol = AllChem.MolFromPDBFile(os.path.join(RDConfig.RDCodeDir, 'Chem', 'test_data', '4DJU_lig.pdb'))
  >>> print len([1 for b in template.GetBonds() if b.GetBondTypeAsDouble() == 1.0])
  8
  >>> print len([1 for b in mol.GetBonds() if b.GetBondTypeAsDouble() == 1.0])
  22

  Now assign the bond orders based on the template molecule
  >>> newMol = AllChem.AssignBondOrdersFromTemplate(template, mol)
  >>> print len([1 for b in newMol.GetBonds() if b.GetBondTypeAsDouble() == 1.0])
  8

  Note that the template molecule should have no explicit hydrogens
  else the algorithm will fail.

  It also works if there are different formal charges (this was github issue 235):
  >>> template=AllChem.MolFromSmiles('CN(C)C(=O)Cc1ccc2c(c1)NC(=O)c3ccc(cc3N2)c4ccc(c(c4)OC)[N+](=O)[O-]')
  >>> mol = AllChem.MolFromMolFile(os.path.join(RDConfig.RDCodeDir, 'Chem', 'test_data', '4FTR_lig.mol'))
  >>> AllChem.MolToSmiles(mol)
  'COC1CC(C2CCC3C(C2)NC2CCC(CC(O)N(C)C)CC2NC3O)CCC1N(O)O'
  >>> newMol = AllChem.AssignBondOrdersFromTemplate(template, mol)
  >>> AllChem.MolToSmiles(newMol)
  'COc1cc(-c2ccc3c(c2)Nc2ccc(CC(=O)N(C)C)cc2NC3=O)ccc1[N+](=O)[O-]'