oddt package¶
Subpackages¶
Submodules¶
oddt.datasets module¶
oddt.interactions module¶
Module calculates interactions between two molecules (proein-protein, protein-ligand, small-small). Currently following interacions are implemented:
- hydrogen bonds
- halogen bonds
- pi stacking (parallel and perpendicular)
- salt bridges
- hydrophobic contacts
- pi-cation
- metal coordination
- pi-metal
-
oddt.interactions.
close_contacts
(x, y, cutoff, x_column='coords', y_column='coords')[source]¶ Returns pairs of atoms which are within close contac distance cutoff.
Parameters: x, y : atom_dict-type numpy array
Atom dictionaries generated by oddt.toolkit.Molecule objects.
- cutoff : float
Cutoff distance for close contacts
- x_column, ycolumn : string, (default=’coords’)
Column containing coordinates of atoms (or pseudo-atoms, i.e. ring centroids)
Returns: x_, y_ : atom_dict-type numpy array
Aligned pairs of atoms in close contact for further processing.
-
oddt.interactions.
hbond_acceptor_donor
(mol1, mol2, cutoff=3.5, base_angle=120, tolerance=30)[source]¶ Returns pairs of acceptor-donor atoms, which meet H-bond criteria
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute H-bond acceptor and H-bond donor pairs
- cutoff : float, (default=3.5)
Distance cutoff for A-D pairs
- base_angle : int, (default=120)
Base angle determining allowed direction of hydrogen bond formation, which is devided by the number of neighbors of acceptor atom to establish final directional angle
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) in which H-bonds are considered as strict.
Returns: a, d : atom_dict-type numpy array
Aligned arrays of atoms forming H-bond, firstly acceptors, secondly donors.
- strict : numpy array, dtype=bool
Boolean array align with atom pairs, informing whether atoms form ‘strict’ H-bond (pass all angular cutoffs). If false, only distance cutoff is met, therefore the bond is ‘crude’.
-
oddt.interactions.
hbond
(mol1, mol2, *args, **kwargs)[source]¶ Calculates H-bonds between molecules
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute H-bond acceptor and H-bond donor pairs
- cutoff : float, (default=3.5)
Distance cutoff for A-D pairs
- base_angle : int, (default=120)
Base angle determining allowed direction of hydrogen bond formation, which is devided by the number of neighbors of acceptor atom to establish final directional angle
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) in which H-bonds are considered as strict.
Returns: mol1_atoms, mol2_atoms : atom_dict-type numpy array
Aligned arrays of atoms forming H-bond
- strict : numpy array, dtype=bool
Boolean array align with atom pairs, informing whether atoms form ‘strict’ H-bond (pass all angular cutoffs). If false, only distance cutoff is met, therefore the bond is ‘crude’.
-
oddt.interactions.
halogenbond_acceptor_halogen
(mol1, mol2, base_angle_acceptor=120, base_angle_halogen=180, tolerance=30, cutoff=4)[source]¶ Returns pairs of acceptor-halogen atoms, which meet halogen bond criteria
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute halogen bond acceptor and halogen pairs
- cutoff : float, (default=4)
Distance cutoff for A-H pairs
- base_angle_acceptor : int, (default=120)
Base angle determining allowed direction of halogen bond formation, which is devided by the number of neighbors of acceptor atom to establish final directional angle
- base_angle_halogen : int (default=180)
Ideal base angle between halogen bond and halogen-neighbor bond
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) in which halogen bonds are considered as strict.
Returns: a, h : atom_dict-type numpy array
Aligned arrays of atoms forming halogen bond, firstly acceptors, secondly halogens
- strict : numpy array, dtype=bool
Boolean array align with atom pairs, informing whether atoms form ‘strict’ halogen bond (pass all angular cutoffs). If false, only distance cutoff is met, therefore the bond is ‘crude’.
-
oddt.interactions.
halogenbond
(mol1, mol2, **kwargs)[source]¶ Calculates halogen bonds between molecules
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute halogen bond acceptor and halogen pairs
- cutoff : float, (default=4)
Distance cutoff for A-H pairs
- base_angle_acceptor : int, (default=120)
Base angle determining allowed direction of halogen bond formation, which is devided by the number of neighbors of acceptor atom to establish final directional angle
- base_angle_halogen : int (default=180)
Ideal base angle between halogen bond and halogen-neighbor bond
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) in which halogen bonds are considered as strict.
Returns: mol1_atoms, mol2_atoms : atom_dict-type numpy array
Aligned arrays of atoms forming halogen bond
- strict : numpy array, dtype=bool
Boolean array align with atom pairs, informing whether atoms form ‘strict’ halogen bond (pass all angular cutoffs). If false, only distance cutoff is met, therefore the bond is ‘crude’.
-
oddt.interactions.
pi_stacking
(mol1, mol2, cutoff=5, tolerance=30)[source]¶ Returns pairs of rings, which meet pi stacking criteria
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute ring pairs
- cutoff : float, (default=5)
Distance cutoff for Pi-stacking pairs
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (parallel or perpendicular) in which pi-stackings are considered as strict.
Returns: r1, r2 : ring_dict-type numpy array
Aligned arrays of rings forming pi-stacking
- strict_parallel : numpy array, dtype=bool
Boolean array align with ring pairs, informing whether rings form ‘strict’ parallel pi-stacking. If false, only distance cutoff is met, therefore the stacking is ‘crude’.
- strict_perpendicular : numpy array, dtype=bool
Boolean array align with ring pairs, informing whether rings form ‘strict’ perpendicular pi-stacking (T-shaped, T-face, etc.). If false, only distance cutoff is met, therefore the stacking is ‘crude’.
-
oddt.interactions.
salt_bridge_plus_minus
(mol1, mol2, cutoff=4)[source]¶ Returns pairs of plus-mins atoms, which meet salt bridge criteria
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute plus and minus pairs
- cutoff : float, (default=4)
Distance cutoff for A-H pairs
Returns: plus, minus : atom_dict-type numpy array
Aligned arrays of atoms forming salt bridge, firstly plus, secondly minus
-
oddt.interactions.
salt_bridges
(mol1, mol2, *args, **kwargs)[source]¶ Calculates salt bridges between molecules
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute plus and minus pairs
- cutoff : float, (default=4)
Distance cutoff for plus-minus pairs
Returns: mol1_atoms, mol2_atoms : atom_dict-type numpy array
Aligned arrays of atoms forming salt bridges
-
oddt.interactions.
hydrophobic_contacts
(mol1, mol2, cutoff=4)[source]¶ Calculates hydrophobic contacts between molecules
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute hydrophobe pairs
- cutoff : float, (default=4)
Distance cutoff for hydrophobe pairs
Returns: mol1_atoms, mol2_atoms : atom_dict-type numpy array
Aligned arrays of atoms forming hydrophobic contacts
-
oddt.interactions.
pi_cation
(mol1, mol2, cutoff=5, tolerance=30)[source]¶ Returns pairs of ring-cation atoms, which meet pi-cation criteria
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute ring-cation pairs
- cutoff : float, (default=5)
Distance cutoff for Pi-cation pairs
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (perpendicular) in which pi-cation are considered as strict.
Returns: r1 : ring_dict-type numpy array
Aligned rings forming pi-stacking
- plus2 : atom_dict-type numpy array
Aligned cations forming pi-cation
- strict_parallel : numpy array, dtype=bool
Boolean array align with ring-cation pairs, informing whether they form ‘strict’ pi-cation. If false, only distance cutoff is met, therefore the interaction is ‘crude’.
-
oddt.interactions.
acceptor_metal
(mol1, mol2, base_angle=120, tolerance=30, cutoff=4)[source]¶ Returns pairs of acceptor-metal atoms, which meet metal coordination criteria Note: This function is directional (mol1 holds acceptors, mol2 holds metals)
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute acceptor and metal pairs
- cutoff : float, (default=4)
Distance cutoff for A-M pairs
- base_angle : int, (default=120)
Base angle determining allowed direction of metal coordination, which is devided by the number of neighbors of acceptor atom to establish final directional angle
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) in metal coordination are considered as strict.
Returns: a, d : atom_dict-type numpy array
Aligned arrays of atoms forming metal coordination, firstly acceptors, secondly metals.
- strict : numpy array, dtype=bool
Boolean array align with atom pairs, informing whether atoms form ‘strict’ metal coordination (pass all angular cutoffs). If false, only distance cutoff is met, therefore the interaction is ‘crude’.
-
oddt.interactions.
pi_metal
(mol1, mol2, cutoff=5, tolerance=30)[source]¶ Returns pairs of ring-metal atoms, which meet pi-metal criteria
Parameters: mol1, mol2 : oddt.toolkit.Molecule object
Molecules to compute ring-metal pairs
- cutoff : float, (default=5)
Distance cutoff for Pi-metal pairs
- tolerance : int, (default=30)
Range (+/- tolerance) from perfect direction (perpendicular) in which pi-metal are considered as strict.
Returns: r1 : ring_dict-type numpy array
Aligned rings forming pi-metal
- m : atom_dict-type numpy array
Aligned metals forming pi-metal
- strict_parallel : numpy array, dtype=bool
Boolean array align with ring-metal pairs, informing whether they form ‘strict’ pi-metal. If false, only distance cutoff is met, therefore the interaction is ‘crude’.
oddt.metrics module¶
oddt.spatial module¶
Spatial functions included in ODDT Mainly used by other modules, but can be accessed directly.
-
oddt.spatial.
angle
(p1, p2, p3)[source]¶ Returns an angle from a series of 3 points (point #2 is centroid).Angle is returned in degrees.
Parameters: p1,p2,p3 : numpy arrays, shape = [n_points, n_dimensions]
Triplets of points in n-dimensional space, aligned in rows.
Returns: angles : numpy array, shape = [n_points]
Series of angles in degrees
-
oddt.spatial.
angle_2v
(v1, v2)[source]¶ Returns an angle between two vecors.Angle is returned in degrees.
Parameters: v1,v2 : numpy arrays, shape = [n_vectors, n_dimensions]
Pairs of vectors in n-dimensional space, aligned in rows.
Returns: angles : numpy array, shape = [n_vectors]
Series of angles in degrees
-
oddt.spatial.
dihedral
(p1, p2, p3, p4)[source]¶ Returns an dihedral angle from a series of 4 points. Dihedral is returned in degrees. Function distingishes clockwise and antyclockwise dihedrals.
Parameters: p1,p2,p3,p4 : numpy arrays, shape = [n_points, n_dimensions]
Quadruplets of points in n-dimensional space, aligned in rows.
Returns: angles : numpy array, shape = [n_points]
Series of angles in degrees
-
oddt.spatial.
distance
(XA, XB, metric='euclidean', p=2, V=None, VI=None, w=None)¶ Computes distance between each pair of the two collections of inputs.
The following are common calling conventions:
Y = cdist(XA, XB, 'euclidean')
Computes the distance between \(m\) points using Euclidean distance (2-norm) as the distance metric between the points. The points are arranged as \(m\) \(n\)-dimensional row vectors in the matrix X.
Y = cdist(XA, XB, 'minkowski', p)
Computes the distances using the Minkowski distance \(||u-v||_p\) (\(p\)-norm) where \(p \geq 1\).
Y = cdist(XA, XB, 'cityblock')
Computes the city block or Manhattan distance between the points.
Y = cdist(XA, XB, 'seuclidean', V=None)
Computes the standardized Euclidean distance. The standardized Euclidean distance between two n-vectors
u
andv
is\[\sqrt{\sum {(u_i-v_i)^2 / V[x_i]}}.\]V is the variance vector; V[i] is the variance computed over all the i’th components of the points. If not passed, it is automatically computed.
Y = cdist(XA, XB, 'sqeuclidean')
Computes the squared Euclidean distance \(||u-v||_2^2\) between the vectors.
Y = cdist(XA, XB, 'cosine')
Computes the cosine distance between vectors u and v,
\[1 - \frac{u \cdot v} {{||u||}_2 {||v||}_2}\]where \(||*||_2\) is the 2-norm of its argument
*
, and \(u \cdot v\) is the dot product of \(u\) and \(v\).Y = cdist(XA, XB, 'correlation')
Computes the correlation distance between vectors u and v. This is
\[1 - \frac{(u - \bar{u}) \cdot (v - \bar{v})} {{||(u - \bar{u})||}_2 {||(v - \bar{v})||}_2}\]where \(\bar{v}\) is the mean of the elements of vector v, and \(x \cdot y\) is the dot product of \(x\) and \(y\).
Y = cdist(XA, XB, 'hamming')
Computes the normalized Hamming distance, or the proportion of those vector elements between two n-vectors
u
andv
which disagree. To save memory, the matrixX
can be of type boolean.Y = cdist(XA, XB, 'jaccard')
Computes the Jaccard distance between the points. Given two vectors,
u
andv
, the Jaccard distance is the proportion of those elementsu[i]
andv[i]
that disagree where at least one of them is non-zero.Y = cdist(XA, XB, 'chebyshev')
Computes the Chebyshev distance between the points. The Chebyshev distance between two n-vectors
u
andv
is the maximum norm-1 distance between their respective elements. More precisely, the distance is given by\[d(u,v) = \max_i {|u_i-v_i|}.\]Y = cdist(XA, XB, 'canberra')
Computes the Canberra distance between the points. The Canberra distance between two points
u
andv
is\[d(u,v) = \sum_i \frac{|u_i-v_i|} {|u_i|+|v_i|}.\]Y = cdist(XA, XB, 'braycurtis')
Computes the Bray-Curtis distance between the points. The Bray-Curtis distance between two points
u
andv
is\[d(u,v) = \frac{\sum_i (u_i-v_i)} {\sum_i (u_i+v_i)}\]Y = cdist(XA, XB, 'mahalanobis', VI=None)
Computes the Mahalanobis distance between the points. The Mahalanobis distance between two pointsu
andv
is \((u-v)(1/V)(u-v)^T\) where \((1/V)\) (theVI
variable) is the inverse covariance. IfVI
is not None,VI
will be used as the inverse covariance matrix.Y = cdist(XA, XB, 'yule')
Computes the Yule distance between the boolean vectors. (see yule function documentation)Y = cdist(XA, XB, 'matching')
Computes the matching distance between the boolean vectors. (see matching function documentation)Y = cdist(XA, XB, 'dice')
Computes the Dice distance between the boolean vectors. (see dice function documentation)Y = cdist(XA, XB, 'kulsinski')
Computes the Kulsinski distance between the boolean vectors. (see kulsinski function documentation)Y = cdist(XA, XB, 'rogerstanimoto')
Computes the Rogers-Tanimoto distance between the boolean vectors. (see rogerstanimoto function documentation)Y = cdist(XA, XB, 'russellrao')
Computes the Russell-Rao distance between the boolean vectors. (see russellrao function documentation)Y = cdist(XA, XB, 'sokalmichener')
Computes the Sokal-Michener distance between the boolean vectors. (see sokalmichener function documentation)Y = cdist(XA, XB, 'sokalsneath')
Computes the Sokal-Sneath distance between the vectors. (see sokalsneath function documentation)Y = cdist(XA, XB, 'wminkowski')
Computes the weighted Minkowski distance between the vectors. (see sokalsneath function documentation)Y = cdist(XA, XB, f)
Computes the distance between all pairs of vectors in X using the user supplied 2-arity function f. For example, Euclidean distance between the vectors could be computed as follows:
dm = cdist(XA, XB, lambda u, v: np.sqrt(((u-v)**2).sum()))
Note that you should avoid passing a reference to one of the distance functions defined in this library. For example,:
dm = cdist(XA, XB, sokalsneath)
would calculate the pair-wise distances between the vectors in X using the Python function sokalsneath. This would result in sokalsneath being called \({n \choose 2}\) times, which is inefficient. Instead, the optimized C version is more efficient, and we call it using the following syntax.:
dm = cdist(XA, XB, 'sokalsneath')
Parameters: XA : ndarray
An \(m_A\) by \(n\) array of \(m_A\) original observations in an \(n\)-dimensional space.
XB : ndarray
An \(m_B\) by \(n\) array of \(m_B\) original observations in an \(n\)-dimensional space.
metric : string or function
The distance metric to use. The distance function can be ‘braycurtis’, ‘canberra’, ‘chebyshev’, ‘cityblock’, ‘correlation’, ‘cosine’, ‘dice’, ‘euclidean’, ‘hamming’, ‘jaccard’, ‘kulsinski’, ‘mahalanobis’, ‘matching’, ‘minkowski’, ‘rogerstanimoto’, ‘russellrao’, ‘seuclidean’, ‘sokalmichener’, ‘sokalsneath’, ‘sqeuclidean’, ‘wminkowski’, ‘yule’.
w : ndarray
The weight vector (for weighted Minkowski).
p : double
The p-norm to apply (for Minkowski, weighted and unweighted)
V : ndarray
The variance vector (for standardized Euclidean).
VI : ndarray
The inverse of the covariance matrix (for Mahalanobis).
Returns: Y : ndarray
A \(m_A\) by \(m_B\) distance matrix is returned. For each \(i\) and \(j\), the metric
dist(u=XA[i], v=XB[j])
is computed and stored in the \(ij\) th entry.Raises: An exception is thrown if ``XA`` and ``XB`` do not have
the same number of columns.
oddt.virtualscreening module¶
Module contents¶
Open Drug Discovery Toolkit¶
Universal and easy to use resource for various drug discovery tasks, ie docking, virutal screening, rescoring.
- toolkit : module,
- Toolkits backend module, currenlty OpenBabel [ob] and RDKit [rdk]. This setting is toolkit-wide, and sets given toolkit as default