RDKit
Open-source cheminformatics and machine learning.
|
Namespaces | |
namespace | details |
Classes | |
class | LAP |
class | MolAlignException |
class | MolHistogram |
class | O3A |
class | O3AConstraint |
class | O3AConstraintVect |
struct | O3AFuncData |
class | SDM |
Enumerations | |
enum | { O3_USE_MMFF_WEIGHTS = (1 << 0) , O3_ACCURACY_MASK = (1 << 0 | 1 << 1) , O3_LOCAL_ONLY = (1 << 2) } |
Variables | |
const int | O3_DUMMY_COST = 100000 |
const int | O3_LARGE_NEGATIVE_WEIGHT = -1e9 |
const int | O3_DEFAULT_CONSTRAINT_WEIGHT = 100.0 |
const unsigned int | O3_MAX_H_BINS = 20 |
const unsigned int | O3_MAX_SDM_ITERATIONS = 100 |
const unsigned int | O3_MAX_SDM_THRESHOLD_ITER = 3 |
const double | O3_RANDOM_TRANS_COEFF = 5.0 |
const double | O3_THRESHOLD_DIFF_DISTANCE = 0.1 |
const double | O3_SDM_THRESHOLD_START = 0.7 |
const double | O3_SDM_THRESHOLD_STEP = 0.3 |
const double | O3_CHARGE_WEIGHT = 10.0 |
const double | O3_CRIPPEN_WEIGHT = 10.0 |
const double | O3_RMSD_THRESHOLD = 1.0e-04 |
const double | O3_SCORE_THRESHOLD = 0.01 |
const double | O3_SCORING_FUNCTION_ALPHA = 5.0 |
const double | O3_SCORING_FUNCTION_BETA = 0.5 |
const double | O3_CHARGE_COEFF = 5.0 |
const double | O3_CRIPPEN_COEFF = 1.0 |
const int | O3_MAX_WEIGHT_COEFF = 5 |
Enumerator | |
---|---|
O3_USE_MMFF_WEIGHTS | |
O3_ACCURACY_MASK | |
O3_LOCAL_ONLY |
Definition at line 125 of file O3AAlignMolecules.h.
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::alignMol | ( | ROMol & | prbMol, |
const ROMol & | refMol, | ||
int | prbCid = -1 , |
||
int | refCid = -1 , |
||
const MatchVectType * | atomMap = nullptr , |
||
const RDNumeric::DoubleVector * | weights = nullptr , |
||
bool | reflect = false , |
||
unsigned int | maxIters = 50 |
||
) |
Optimally (minimum RMSD) align a molecule to another molecule.
The 3D transformation required to align the specified conformation in the probe molecule to a specified conformation in the reference molecule is computed so that the root mean squared distance between a specified set of atoms is minimized. This transform is then applied to the specified conformation in the probe molecule
prbMol | molecule that is to be aligned |
refMol | molecule used as the reference for the alignment |
prbCid | ID of the conformation in the probe to be used for the alignment (defaults to first conformation) |
refCid | ID of the conformation in the ref molecule to which the alignment is computed (defaults to first conformation) |
atomMap | a vector of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If this mapping is not specified an attempt is made to generate on by substructure matching |
weights | Optionally specify weights for each of the atom pairs |
reflect | if true reflect the conformation of the probe molecule |
maxIters | maximum number of iterations used in minimizing the RMSD |
Returns RMSD value
RDKIT_MOLALIGN_EXPORT void RDKit::MolAlign::alignMolConformers | ( | ROMol & | mol, |
const std::vector< unsigned int > * | atomIds = nullptr , |
||
const std::vector< unsigned int > * | confIds = nullptr , |
||
const RDNumeric::DoubleVector * | weights = nullptr , |
||
bool | reflect = false , |
||
unsigned int | maxIters = 50 , |
||
std::vector< double > * | RMSlist = nullptr |
||
) |
Align the conformations of a molecule using a common set of atoms. If the molecules contains queries, then the queries must also match exactly.
mol | The molecule of interest. |
atomIds | vector of atoms to be used to generate the alignment. All atoms will be used is not specified |
confIds | vector of conformations to align - defaults to all |
weights | (optional) weights for each pair of atoms. |
reflect | toggles reflecting (about the origin) the alignment |
maxIters | the maximum number of iterations to attempt |
RMSlist | if nonzero, this will be used to return the RMS values between the reference conformation and the other aligned conformations |
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::CalcRMS | ( | ROMol & | prbMol, |
const ROMol & | refMol, | ||
int | prbCid, | ||
int | refCid, | ||
const std::vector< MatchVectType > & | map, | ||
int | maxMatches, | ||
const RDNumeric::DoubleVector * | weights | ||
) |
Returns the RMS between two molecules, taking symmetry into account. In contrast to getBestRMS, the RMS is computed "in place", i.e. probe molecules are not aligned to the reference ahead of the RMS calculation. This is useful, for example, to compute the RMSD between docking poses and the co-crystallized ligand.
This function will attempt to match all permutations of matching atom orders in both molecules, for some molecules it will lead to 'combinatorial explosion' especially if hydrogens are present.
prbMol | the molecule to be aligned to the reference |
refMol | the reference molecule |
prbCid | (optional) probe conformation to use |
refCid | (optional) reference conformation to use |
map | (optional) a vector of vectors of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If not provided, these will be generated using a substructure search. |
maxMatches | (optional) if map is empty, this will be the max number of matches found in a SubstructMatch(). |
weights | (optional) weights for each pair of atoms. |
Returns Best RMSD value found
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::CalcRMS | ( | ROMol & | prbMol, |
const ROMol & | refMol, | ||
int | prbCid = -1 , |
||
int | refCid = -1 , |
||
const std::vector< MatchVectType > & | map = std::vector< MatchVectType >() , |
||
int | maxMatches = 1e6 , |
||
bool | symmetrizeConjugatedTerminalGroups = true , |
||
const RDNumeric::DoubleVector * | weights = nullptr |
||
) |
Returns the RMS between two molecules, taking symmetry into account. In contrast to getBestRMS, the RMS is computed "in place", i.e. probe molecules are not aligned to the reference ahead of the RMS calculation. This is useful, for example, to compute the RMSD between docking poses and the co-crystallized ligand.
This function will attempt to match all permutations of matching atom orders in both molecules, for some molecules it will lead to 'combinatorial explosion' especially if hydrogens are present.
prbMol | the molecule to be aligned to the reference |
refMol | the reference molecule |
prbCid | (optional) probe conformation to use |
refCid | (optional) reference conformation to use |
map | (optional) a vector of vectors of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If not provided, these will be generated using a substructure search. |
maxMatches | (optional) if map is empty, this will be the max number of matches found in a SubstructMatch(). |
symmetrizeConjugatedTerminalGroups | (optional) if set, conjugated terminal functional groups (like nitro or carboxylate) will be considered symmetrically |
weights | (optional) weights for each pair of atoms. |
Returns Best RMSD value found
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::getAlignmentTransform | ( | const ROMol & | prbMol, |
const ROMol & | refMol, | ||
RDGeom::Transform3D & | trans, | ||
int | prbCid = -1 , |
||
int | refCid = -1 , |
||
const MatchVectType * | atomMap = nullptr , |
||
const RDNumeric::DoubleVector * | weights = nullptr , |
||
bool | reflect = false , |
||
unsigned int | maxIters = 50 |
||
) |
Alignment functions.
Compute the transformation required to align a molecule
The 3D transformation required to align the specified conformation in the probe molecule to a specified conformation in the reference molecule is computed so that the root mean squared distance between a specified set of atoms is minimized
prbMol | molecule that is to be aligned |
refMol | molecule used as the reference for the alignment |
trans | storage for the computed transform |
prbCid | ID of the conformation in the probe to be used for the alignment (defaults to first conformation) |
refCid | ID of the conformation in the ref molecule to which the alignment is computed (defaults to first conformation) |
atomMap | a vector of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If this mapping is not specified an attempt is made to generate on by substructure matching |
weights | Optionally specify weights for each of the atom pairs |
reflect | if true reflect the conformation of the probe molecule |
maxIters | maximum number of iterations used in minimizing the RMSD |
Returns RMSD value
RDKIT_MOLALIGN_EXPORT std::vector< double > RDKit::MolAlign::getAllConformerBestRMS | ( | const ROMol & | mol, |
int | numThreads = 1 , |
||
const std::vector< MatchVectType > & | map = std::vector< MatchVectType >() , |
||
int | maxMatches = 1e6 , |
||
bool | symmetrizeConjugatedTerminalGroups = true , |
||
const RDNumeric::DoubleVector * | weights = nullptr |
||
) |
Returns the symmetric distance matrix between the conformers of a molecule. getBestRMS() is used to calculate the inter-conformer distances
This function will attempt to align all permutations of matching atom orders in both molecules, for some molecules it will lead to 'combinatorial explosion' especially if hydrogens are present.
mol | the molecule to be considered |
numThreads | (optional) number of threads to use during the calculation |
map | (optional) a vector of vectors of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If not provided, these will be generated using a substructure search. |
maxMatches | (optional) if map is empty, this will be the max number of matches found in a SubstructMatch(). |
symmetrizeConjugatedTerminalGroups | (optional) if set, conjugated terminal functional groups (like nitro or carboxylate) will be considered symmetrically |
weights | (optional) weights for each pair of atoms. |
Returns a vector with the RMSD values stored in the order: [(1,0), (2,0), (2,1), (3,0), (3, 2), (3,1), ...]
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::getBestAlignmentTransform | ( | const ROMol & | prbMol, |
const ROMol & | refMol, | ||
RDGeom::Transform3D & | bestTrans, | ||
MatchVectType & | bestMatch, | ||
int | prbCid = -1 , |
||
int | refCid = -1 , |
||
const std::vector< MatchVectType > & | map = std::vector< MatchVectType >() , |
||
int | maxMatches = 1e6 , |
||
bool | symmetrizeConjugatedTerminalGroups = true , |
||
const RDNumeric::DoubleVector * | weights = nullptr , |
||
bool | reflect = false , |
||
unsigned int | maxIters = 50 , |
||
int | numThreads = 1 |
||
) |
Compute the optimal RMS, transformation and atom map for aligning two molecules, taking symmetry into account. Molecule coordinates are left unaltered.
This function will attempt to align all permutations of matching atom orders in both molecules, for some molecules it will lead to 'combinatorial explosion' especially if hydrogens are present. Use 'RDKit::MolAlign::getAlignmentTransform' to align molecules without changing the atom order.
prbMol | the molecule to be aligned to the reference |
refMol | the reference molecule |
bestTrans | storage for the best computed transform |
bestMatch | storage for the MatchVectType corresponding to the best match found. |
prbCid | (optional) probe conformation to use |
refCid | (optional) reference conformation to use |
map | (optional) a vector of vectors of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If not provided, these will be generated using a substructure search. |
maxMatches | (optional) if map is empty, this will be the max number of matches found in a SubstructMatch(). |
symmetrizeConjugatedTerminalGroups | (optional) if set, conjugated terminal functional groups (like nitro or carboxylate) will be considered symmetrically |
weights | (optional) weights for each pair of atoms. |
reflect | if true reflect the conformation of the probe molecule |
maxIters | maximum number of iterations used in minimizing the RMSD |
numThreads | (optional) number of threads to use during the calculation |
Returns Best RMSD value found
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::getBestRMS | ( | ROMol & | prbMol, |
const ROMol & | refMol, | ||
int | prbCid = -1 , |
||
int | refCid = -1 , |
||
const std::vector< MatchVectType > & | map = std::vector< MatchVectType >() , |
||
int | maxMatches = 1e6 , |
||
bool | symmetrizeConjugatedTerminalGroups = true , |
||
const RDNumeric::DoubleVector * | weights = nullptr , |
||
int | numThreads = 1 |
||
) |
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.
This function will attempt to align all permutations of matching atom orders in both molecules, for some molecules it will lead to 'combinatorial explosion' especially if hydrogens are present. Use 'RDKit::MolAlign::alignMol' to align molecules without changing the atom order.
prbMol | the molecule to be aligned to the reference |
refMol | the reference molecule |
trans | storage for the computed transform |
prbCid | (optional) probe conformation to use |
refCid | (optional) reference conformation to use |
map | (optional) a vector of vectors of pairs of atom IDs (probe AtomId, ref AtomId) used to compute the alignments. If not provided, these will be generated using a substructure search. |
maxMatches | (optional) if map is empty, this will be the max number of matches found in a SubstructMatch(). |
symmetrizeConjugatedTerminalGroups | (optional) if set, conjugated terminal functional groups (like nitro or carboxylate) will be considered symmetrically |
weights | (optional) weights for each pair of atoms. |
numThreads | (optional) number of threads to use during the calculation |
Returns Best RMSD value found
RDKIT_MOLALIGN_EXPORT void RDKit::MolAlign::getO3AForProbeConfs | ( | ROMol & | prbMol, |
const ROMol & | refMol, | ||
void * | prbProp, | ||
void * | refProp, | ||
std::vector< boost::shared_ptr< O3A > > & | res, | ||
int | numThreads = 1 , |
||
O3A::AtomTypeScheme | atomTypes = O3A::MMFF94 , |
||
const int | refCid = -1 , |
||
const bool | reflect = false , |
||
const unsigned int | maxIters = 50 , |
||
unsigned int | options = 0 , |
||
const MatchVectType * | constraintMap = nullptr , |
||
const RDNumeric::DoubleVector * | constraintWeights = nullptr |
||
) |
Definition at line 39 of file O3AAlignMolecules.h.
RDKIT_MOLALIGN_EXPORT int RDKit::MolAlign::o3aCrippenCostFunc | ( | const unsigned int | prbIdx, |
const unsigned int | refIdx, | ||
double | hSum, | ||
void * | data | ||
) |
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::o3aCrippenScoringFunc | ( | const unsigned int | prbIdx, |
const unsigned int | refIdx, | ||
void * | data | ||
) |
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::o3aCrippenWeightFunc | ( | const unsigned int | prbIdx, |
const unsigned int | refIdx, | ||
void * | data | ||
) |
RDKIT_MOLALIGN_EXPORT int RDKit::MolAlign::o3aMMFFCostFunc | ( | const unsigned int | prbIdx, |
const unsigned int | refIdx, | ||
double | hSum, | ||
void * | data | ||
) |
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::o3aMMFFScoringFunc | ( | const unsigned int | prbIdx, |
const unsigned int | refIdx, | ||
void * | data | ||
) |
RDKIT_MOLALIGN_EXPORT double RDKit::MolAlign::o3aMMFFWeightFunc | ( | const unsigned int | prbIdx, |
const unsigned int | refIdx, | ||
void * | data | ||
) |
RDKIT_MOLALIGN_EXPORT void RDKit::MolAlign::randomTransform | ( | ROMol & | mol, |
const int | cid = -1 , |
||
const int | seed = -1 |
||
) |
RDKIT_MOLALIGN_EXPORT const RDGeom::POINT3D_VECT * RDKit::MolAlign::reflect | ( | const Conformer & | conf | ) |
Definition at line 122 of file O3AAlignMolecules.h.
Definition at line 116 of file O3AAlignMolecules.h.
Definition at line 123 of file O3AAlignMolecules.h.
Definition at line 117 of file O3AAlignMolecules.h.
Definition at line 108 of file O3AAlignMolecules.h.
Definition at line 106 of file O3AAlignMolecules.h.
Definition at line 107 of file O3AAlignMolecules.h.
Definition at line 109 of file O3AAlignMolecules.h.
Definition at line 110 of file O3AAlignMolecules.h.
Definition at line 111 of file O3AAlignMolecules.h.
Definition at line 124 of file O3AAlignMolecules.h.
Definition at line 112 of file O3AAlignMolecules.h.
Definition at line 118 of file O3AAlignMolecules.h.
Definition at line 119 of file O3AAlignMolecules.h.
Definition at line 120 of file O3AAlignMolecules.h.
Definition at line 121 of file O3AAlignMolecules.h.
Definition at line 114 of file O3AAlignMolecules.h.
Definition at line 115 of file O3AAlignMolecules.h.
Definition at line 113 of file O3AAlignMolecules.h.