Chemical space

Embedding, diversity, clustering, nearest-neighbour and scaffold analysis of compound collections.

These are the questions you ask before modelling: does this library cover one region or several, does the test set sit inside the training set’s cloud, is this screening hit an outlier.

Embedding

>>> from qsarkit.chemspace import ChemicalSpaceAnalyzer
>>> analyzer = ChemicalSpaceAnalyzer(random_state=0).fit(demo_mols)
>>> analyzer.embedding_.shape
(24, 2)
>>> type(analyzer.plot()).__name__
'Figure'

Warning

The three methods are not interchangeable. PCA preserves global variance and its axes are interpretable. t-SNE and UMAP preserve local neighbourhoods and give the familiar island plots — but between-cluster distances in those plots are not meaningful. Reading them as chemical distance is the commonest misuse of the technique.

Diversity

>>> from qsarkit.chemspace import DiversityAnalyzer
>>> report = DiversityAnalyzer().analyze(demo_mols)
>>> round(report["internal_diversity"], 3)
0.763
>>> report["n_scaffolds"], report["acyclic_fraction"]
(3.0, 0.0)

Comparing libraries needs no extra code:

>>> frame = DiversityAnalyzer().compare({
...     "benzoic acids": demo_mols[:6],
...     "benzimidazoles": demo_mols[18:],
... })
>>> len(frame)
2

Scaffolds

>>> from qsarkit.chemspace import ScaffoldAnalyzer
>>> scaffolds = ScaffoldAnalyzer().fit(demo_mols)
>>> scaffolds.n_scaffolds
3
>>> scaffolds.most_common(2)[0][1]
12

Acyclic molecules have no Bemis-Murcko framework at all. They are reported separately rather than counted as a scaffold of their own — otherwise a library of straight chains would look scaffold-diverse:

>>> scaffolds.summary()["acyclic_fraction"]
0.0

Novelty and coverage

>>> from qsarkit.chemspace import NearestNeighborAnalyzer
>>> nn = NearestNeighborAnalyzer().fit(demo_mols[:18])
>>> novelty = nn.novelty(demo_mols[18:])
>>> novelty.shape
(6,)
>>> bool((novelty > 0).all())
True

redundancy compares a library against itself, excluding each molecule’s own row — so an exact duplicate is detected rather than masked:

>>> library = demo_mols + [demo_mols[0]]
>>> analyzer = NearestNeighborAnalyzer().fit(library)
>>> round(analyzer.redundancy(library, threshold=0.99), 3)
0.08

>>> from qsarkit.chemspace import ChemicalSpaceCoverage
>>> coverage = ChemicalSpaceCoverage(threshold=0.5).compare(
...     demo_mols[18:], demo_mols[:18])
>>> round(coverage["novel_fraction"], 3)
1.0

The benzimidazoles are entirely novel relative to the other three series — which is exactly why a scaffold split holding them out is a hard test.

API

Chemical space analysis: fingerprints, similarity and scaffolds.

References

  • Bajusz, D., Racz, A. & Heberger, K. (2015). “Why is Tanimoto Index an Appropriate Choice for Fingerprint-Based Similarity Calculations?” J. Cheminform., 7, 20. https://doi.org/10.1186/s13321-015-0069-3

  • Bemis, G. W. & Murcko, M. A. (1996). “The Properties of Known Drugs. 1. Molecular Frameworks.” J. Med. Chem., 39(15), 2887-2893. https://doi.org/10.1021/jm9602928

class qsarkit.chemspace.ChemicalSpaceAnalyzer(method='pca', n_components=2, metric='jaccard', random_state=None, **kwargs)[source]

Bases: object

Project a compound collection into two dimensions for inspection.

Chemical space is high-dimensional and sparse, so any 2D picture of it is a lossy projection — but the right projection answers real questions: whether a library covers one region or several, whether the test set sits inside the training set’s cloud, whether a screening hit is an outlier.

The three methods answer different questions and are not interchangeable. PCA preserves global variance and its axes are interpretable, but it flattens the non-linear structure fingerprints actually have. t-SNE and UMAP preserve local neighbourhoods and give the familiar island plots, but between-cluster distances in those plots are not meaningful — reading them as chemical distance is the commonest misuse of the technique.

Parameters:

  • method (Literal['pca', 'tsne', 'mds', 'umap']) – Projection method. "umap" needs the optional umap-learn package.

  • n_components (int) – Output dimensionality.

  • metric (Literal['jaccard', 'euclidean']) – Distance used by the neighbourhood methods. Jaccard/Tanimoto is correct for fingerprints.

  • random_state (Optional[int]) – Seed.

  • **kwargs (Any) – Forwarded to the underlying estimator (perplexity, n_neighbors, …).

Variables:

embedding (ndarray of shape (n_molecules, n_components)) – The projected coordinates.

Examples

>>> from rdkit import Chem
>>> mols = [Chem.MolFromSmiles(s) for s in ("CCO", "CCN", "c1ccccc1", "CCC")]
>>> analyzer = ChemicalSpaceAnalyzer(random_state=0).fit(mols)
>>> analyzer.embedding_.shape
(4, 2)

References

embedding_: ndarray[tuple[Any, ...], dtype[float64]]
fit(mols, y=None)[source]

Project the molecules.

Parameters:
Return type:

ChemicalSpaceAnalyzer

fit_transform(mols, y=None)[source]

Project and return the coordinates.

Return type:

ndarray[tuple[Any, ...], dtype[double]]

plot(color=None, labels=None, title='Chemical space')[source]

Scatter the embedding, optionally coloured by a property.

Parameters:
Return type:

Figure

class qsarkit.chemspace.DiversityAnalyzer(radius=2, n_bits=2048)[source]

Bases: object

Quantify how diverse a compound collection is.

“Diverse” needs a definition before it can be measured, and the available ones disagree. Mean pairwise distance rewards a few far-out outliers; scaffold count rewards structural variety regardless of distance; internal diversity is bounded and comparable between sets of different sizes. All three are reported, because a library can score well on one and poorly on another and the difference is informative.

Parameters:

  • radius (int) – Morgan radius.

  • n_bits (int) – Fingerprint length.

Examples

>>> from rdkit import Chem
>>> mols = [Chem.MolFromSmiles(s) for s in ("CCO", "c1ccccc1", "CCCCCC")]
>>> report = DiversityAnalyzer().analyze(mols)
>>> 0.0 <= report["mean_pairwise_distance"] <= 1.0
True

References

  • Waldman, M., Li, H. & Hassan, M. (2000). “Novel Algorithms for the Optimization of Molecular Diversity of Combinatorial Libraries.” J. Mol. Graph. Model., 18(4-5), 412-426. https://doi.org/10.1016/S1093-3263(00)00071-2

  • Bemis, G. W. & Murcko, M. A. (1996). “The Properties of Known Drugs. 1. Molecular Frameworks.” J. Med. Chem., 39(15), 2887-2893. https://doi.org/10.1021/jm9602928

  • Benhenda, M. (2017). “ChemGAN Challenge for Drug Discovery: Can AI Reproduce Natural Chemical Diversity?” arXiv:1708.08227. https://arxiv.org/abs/1708.08227

  • Shannon, C. E. (1948). “A Mathematical Theory of Communication.” Bell Syst. Tech. J., 27(3), 379-423. https://doi.org/10.1002/j.1538-7305.1948.tb01338.x

analyze(mols)[source]

Compute the diversity measures.

Parameters:

mols (Sequence[Any])

Returns:

n_molecules, mean_pairwise_distance, median_pairwise_distance, internal_diversity, n_scaffolds (acyclic molecules excluded – they have no framework), acyclic_fraction, scaffold_diversity (scaffolds per molecule), scaffold_entropy and bit_entropy.

Return type:

Dict[str, float]

compare(libraries)[source]

Compare several libraries on every diversity measure.

Parameters:

libraries (Dict[str, Sequence[Any]]) – name -> molecules.

Returns:

One row per library.

Return type:

DataFrame

class qsarkit.chemspace.ClusterAnalyzer(method='butina', cutoff=0.35, n_clusters=None, radius=2, n_bits=2048, random_state=None)[source]

Bases: object

Cluster a compound collection and describe the result.

Delegates to qsarkit.cluster for the cheminformatics methods and to scikit-learn for the general ones, then summarizes what came out: how many clusters, how big, and how many singletons. The singleton fraction is the number worth watching — a library that clusters into mostly singletons is either genuinely diverse or being clustered at too tight a cutoff.

Parameters:

  • method (Literal['butina', 'sphere_exclusion', 'hierarchical', 'kmeans', 'dbscan']) – Clustering algorithm.

  • cutoff (float) – Distance cutoff for the cheminformatics methods.

  • n_clusters (Optional[int]) – Cluster count for "hierarchical" and "kmeans".

  • radius (int) – Morgan radius.

  • n_bits (int) – Fingerprint length.

  • random_state (Optional[int]) – Seed for k-means.

Variables:

labels (ndarray of shape (n_molecules,)) – Cluster assignment.

Examples

>>> from rdkit import Chem
>>> mols = [Chem.MolFromSmiles(s) for s in
...         ("CCO", "CCN", "c1ccccc1", "c1ccccc1C")]
>>> analyzer = ClusterAnalyzer(cutoff=0.5).fit(mols)
>>> analyzer.labels_.shape
(4,)

References

  • Butina, D. (1999). “Unsupervised Data Base Clustering Based on Daylight’s Fingerprint and Tanimoto Similarity.” J. Chem. Inf. Comput. Sci., 39(4), 747-750. https://doi.org/10.1021/ci9803381

  • Downs, G. M. & Barnard, J. M. (2002). “Clustering Methods and Their Uses in Computational Chemistry.” Rev. Comput. Chem., 18, 1-40. https://doi.org/10.1002/0471433519.ch1

labels_: ndarray[tuple[Any, ...], dtype[int64]]
fit(mols, y=None)[source]

Cluster the molecules.

Parameters:
Return type:

ClusterAnalyzer

fit_predict(mols, y=None)[source]

Cluster and return the labels.

Return type:

ndarray[tuple[Any, ...], dtype[int_]]

summary()[source]

Describe the clustering.

Returns:

n_clusters, n_singletons, singleton_fraction, largest_cluster, mean_cluster_size. DBSCAN’s noise label (-1) is counted as singletons.

Return type:

Dict[str, float]

cluster_members()[source]

Molecule indices grouped by cluster.

Returns:

cluster label -> list of molecule indices.

Return type:

Dict[int, List[int]]

class qsarkit.chemspace.NearestNeighborAnalyzer(n_neighbors=1, radius=2, n_bits=2048)[source]

Bases: object

Nearest-neighbour statistics for a compound collection.

Underpins two everyday questions: how close is this new compound to anything we already have (novelty), and how self-similar is this library (redundancy). A screening set whose members are all each other’s near neighbours is smaller than its compound count suggests.

Parameters:

  • n_neighbors (int) – Neighbours considered.

  • radius (int) – Morgan radius.

  • n_bits (int) – Fingerprint length.

Examples

>>> from rdkit import Chem
>>> library = [Chem.MolFromSmiles(s) for s in ("CCO", "CCN", "c1ccccc1")]
>>> analyzer = NearestNeighborAnalyzer().fit(library)
>>> float(analyzer.nearest_similarity([Chem.MolFromSmiles("CCO")])[0])
1.0

References

  • Sheridan, R. P. et al. (2004). “Similarity to Molecules in the Training Set Is a Good Discriminator for Prediction Accuracy in QSAR.” J. Chem. Inf. Comput. Sci., 44(6), 1912-1928. https://doi.org/10.1021/ci049782w

  • Willett, P. (2006). “Similarity-Based Virtual Screening Using 2D Fingerprints.” Drug Discov. Today, 11(23-24), 1046-1053. https://doi.org/10.1016/j.drudis.2006.10.005

fit(mols, y=None)[source]

Store the reference library.

Parameters:
Return type:

NearestNeighborAnalyzer

nearest_similarity(mols, exclude_self=False)[source]

Similarity of each query to its nearest library member.

Parameters:

  • mols (Sequence[Any]) – Query molecules.

  • exclude_self (bool) – Ignore each query’s own row, which is what you want when the queries are the library and you are measuring internal redundancy. Only the matching position is masked, not every perfect match – masking those would hide exact duplicates, which is precisely what such a measurement is looking for. Requires the query set to be the fitted library.

Return type:

ndarray[tuple[Any, ...], dtype[double]]

Raises:

ValueError – If exclude_self is set but the query set is not the same size as the library, so there is no “self” to exclude.

novelty(mols)[source]

Novelty of each query: 1 - nearest similarity.

Parameters:

mols (Sequence[Any])

Returns:

0 means an exact match in the library; 1 means nothing alike.

Return type:

ndarray[tuple[Any, ...], dtype[double]]

redundancy(mols, threshold=0.9)[source]

Fraction of a set having a near-duplicate elsewhere in it.

Parameters:

  • mols (Sequence[Any])

  • threshold (float) – Similarity above which two molecules count as duplicates.

Returns:

In [0, 1].

Return type:

float

class qsarkit.chemspace.ScaffoldAnalyzer(generic=False)[source]

Bases: object

Bemis-Murcko scaffold analysis of a compound collection.

Scaffolds are how medicinal chemists actually partition a library, and the scaffold distribution says something a compound count cannot: a 10,000-compound set built on twelve scaffolds is a different asset from one built on three thousand.

Parameters:

  • generic (bool) – Reduce scaffolds to their carbon skeleton, which merges

  • system. (heteroatom-substituted variants of the same ring)

Examples

>>> from rdkit import Chem
>>> mols = [Chem.MolFromSmiles(s) for s in
...         ("c1ccccc1C", "c1ccccc1CC", "c1ccncc1C", "CCO")]
>>> analyzer = ScaffoldAnalyzer().fit(mols)
>>> analyzer.n_scaffolds
2

Toluene and ethylbenzene share the benzene framework; picoline contributes a second. Ethanol is acyclic, so it has no Bemis-Murcko framework at all and is not counted as a scaffold of its own – summary() reports it separately instead:

>>> analyzer.summary()["acyclic_fraction"]
0.25

References

scaffolds_: List[str]
fit(mols, y=None)[source]

Extract the scaffold of every molecule.

Parameters:
Return type:

ScaffoldAnalyzer

property n_scaffolds: int

Number of distinct scaffolds, excluding the acyclic empty one.

most_common(n=10)[source]

The most frequent scaffolds.

Parameters:

n (int)

Return type:

List[Tuple[str, int]]

groups()[source]

Molecule indices grouped by scaffold.

Returns:

scaffold SMILES -> list of molecule indices.

Return type:

Dict[str, List[int]]

summary()[source]

Describe the scaffold distribution.

Returns:

n_molecules, n_scaffolds, scaffold_diversity, n_singleton_scaffolds, largest_scaffold_group, acyclic_fraction and top_scaffold_share (the fraction of the library sitting on its single most common scaffold — the quickest way to spot a library that is really one series).

Return type:

Dict[str, float]

to_dataframe()[source]

Scaffold frequency table.

Returns:

Columns scaffold, count, fraction, descending.

Return type:

DataFrame

class qsarkit.chemspace.ChemicalSpaceCoverage(threshold=0.7, radius=2, n_bits=2048)[source]

Bases: object

Compare the chemical space covered by two collections.

The practical question behind a library purchase or a virtual screen: how much of the target space does this set actually reach, and how much of it is already covered by what we own?

Parameters:

  • threshold (float) – Similarity at which a reference compound counts as covered. 0.7 on ECFP4 is the conventional “similar enough” cutoff.

  • radius (int) – Morgan radius.

  • n_bits (int) – Fingerprint length.

Examples

>>> from rdkit import Chem
>>> a = [Chem.MolFromSmiles(s) for s in ("CCO", "CCN")]
>>> b = [Chem.MolFromSmiles(s) for s in ("CCO", "c1ccccc1")]
>>> report = ChemicalSpaceCoverage().compare(a, b)
>>> 0.0 <= report["coverage_of_reference"] <= 1.0
True

References

compare(query, reference)[source]

Measure how well query covers reference.

Parameters:

  • query (Sequence[Any]) – The collection being assessed.

  • reference (Sequence[Any]) – The space to be covered.

Returns:

coverage_of_reference (fraction of reference compounds with a similar query compound), mean_nearest_similarity, n_novel_in_query (query compounds unlike anything in the reference) and novel_fraction.

Return type:

Dict[str, float]

qsarkit.chemspace.morgan_generator(radius=2, n_bits=2048)[source]

Return a configured RDKit Morgan fingerprint generator.

Parameters:

  • radius (int) – Morgan radius. radius=2 corresponds to ECFP4.

  • n_bits (int) – Folded bit-vector length.

Returns:

A generator object exposing GetFingerprint.

Return type:

Any

Examples

>>> from rdkit import Chem
>>> gen = morgan_generator()
>>> fp = gen.GetFingerprint(Chem.MolFromSmiles("c1ccccc1"))
>>> fp.GetNumBits()
2048

References

qsarkit.chemspace.compute_fingerprints(mols, radius=2, n_bits=2048)[source]

Compute ECFP bit vectors for a sequence of molecules.

Parameters:

  • mols (Sequence[Any]) – Input molecules. None entries are not allowed.

  • radius (int) – Morgan radius.

  • n_bits (int) – Folded bit-vector length.

Returns:

One fingerprint per input molecule, positionally aligned.

Return type:

List[Any]

Examples

>>> from rdkit import Chem
>>> fps = compute_fingerprints([Chem.MolFromSmiles("CCO")])
>>> len(fps)
1

References

qsarkit.chemspace.fingerprints_to_array(fps)[source]

Convert a list of RDKit bit vectors into a dense (n, n_bits) array.

Parameters:

fps (Sequence[Any]) – Fingerprints to densify.

Returns:

The dense binary matrix.

Return type:

ndarray

Examples

>>> from rdkit import Chem
>>> arr = fingerprints_to_array(compute_fingerprints([Chem.MolFromSmiles("CCO")]))
>>> arr.shape
(1, 2048)

References

qsarkit.chemspace.tanimoto_matrix(fps, other=None)[source]

Full pairwise Tanimoto similarity matrix.

Parameters:

  • fps (Sequence[Any]) – Query fingerprints (rows of the output).

  • other (Optional[Sequence[Any]]) – Reference fingerprints (columns). Defaults to fps itself, in which case the returned matrix is symmetric with a unit diagonal.

Returns:

Tanimoto (Jaccard) similarities in [0, 1].

Return type:

ndarray

Examples

>>> from rdkit import Chem
>>> fps = compute_fingerprints([Chem.MolFromSmiles(s) for s in ("CCO", "CCO")])
>>> float(tanimoto_matrix(fps)[0, 1])
1.0

References

  • Tanimoto, T. T. (1958). IBM Internal Report.

  • Bajusz, D. et al. (2015). “Why is Tanimoto index an appropriate choice for fingerprint-based similarity calculations?” J. Cheminform., 7, 20. https://doi.org/10.1186/s13321-015-0069-3

qsarkit.chemspace.bemis_murcko_smiles(mol, generic=False)[source]

Canonical SMILES of a molecule’s Bemis-Murcko scaffold.

Parameters:

  • mol (Any) – Input molecule.

  • generic (bool) – If True, strip element and bond-order information to obtain the cyclic skeleton (“graph framework”).

Returns:

Canonical scaffold SMILES (empty string for acyclic molecules).

Return type:

str

Examples

>>> from rdkit import Chem
>>> bemis_murcko_smiles(Chem.MolFromSmiles("CC(=O)Oc1ccccc1C(=O)O"))
'c1ccccc1'

References

  • Bemis, G. W. & Murcko, M. A. (1996). “The Properties of Known Drugs. 1. Molecular Frameworks.” J. Med. Chem., 39(15), 2887-2893. https://doi.org/10.1021/jm9602928

References

  • van der Maaten, L. & Hinton, G. (2008). “Visualizing Data Using t-SNE.” J. Mach. Learn. Res., 9, 2579-2605. https://jmlr.org/papers/v9/vandermaaten08a.html

  • McInnes, L., Healy, J. & Melville, J. (2018). “UMAP.” arXiv:1802.03426

  • Wattenberg, M., Viegas, F. & Johnson, I. (2016). “How to Use t-SNE Effectively.” Distill. doi:10.23915/distill.00002

  • Bemis, G. W. & Murcko, M. A. (1996). “The Properties of Known Drugs. 1. Molecular Frameworks.” J. Med. Chem., 39(15), 2887-2893. doi:10.1021/jm9602928

  • Sheridan, R. P. et al. (2004). “Similarity to Molecules in the Training Set Is a Good Discriminator for Prediction Accuracy in QSAR.” J. Chem. Inf. Comput. Sci., 44(6), 1912-1928. doi:10.1021/ci049782w