Skip to main content

The rcdk and cluster R packages applied to drug candidate selection


The aim of this article is to show how thevpower of statistics and cheminformatics can be combined, in R, using two packages: rcdk and cluster.

We describe the role of clustering methods for identifying similar structures in a group of 23 molecules according to their fingerprints. The most commonly used method is to group the molecules using a “score” obtained by measuring the average distance between them. This score reflects the similarity/non-similarity between compounds and helps us identify active or potentially toxic substances through predictive studies.

Clustering is the process by which the common characteristics of a particular class of compounds are identified. For clustering applications, we are generally measure the molecular fingerprint similarity with the Tanimoto coefficient. Based on the molecular fingerprints, we calculated the molecular distances between the methotrexate molecule and the other 23 molecules in the group, and organized them into a matrix. According to the molecular distances and Ward ’s method, the molecules were grouped into 3 clusters. We can presume structural similarity between the compounds and their locations in the cluster map. Because only 5 molecules were included in the methotrexate cluster, we considered that they might have similar properties and might be further tested as potential drug candidates.


Discovery, synthesis and production of new drugs is still challenging for researchers because of the complex structures of endogenous molecules involved in the pathogenesis of diseases such as AIDS, cancer and autism [16]. Modern drug research is characterized by the growing number of lead molecules and the need to examine and characterize all of these compounds over short periods [14, 39].

Chemical database mining based on the similar compounds search is an in silico method widely used in the drug discovery process [28, 33]. It can be used in the initial stages of drug discovery and speeds up the entire process [10]. The requirement to store, manage and analyse these rapidly growing resources has given rise to a relatively new field known as computer-assisted drug design (CADD) [22, 39, 40].

Computational chemistry is a very effective approach in drug design for the identification of lead compounds. Various virtual screening techniques can be used to reduce the cost and time required to identify a potential drug [2]. As a computational method in drug discovery and virtual screening, clustering of chemical compounds by the similarity of their molecular fingerprints can be used to identify similar structures in a large set of similar data [38, 41]. Their virtual screening performance is comparable to other,more complex, methods. There are many types of fingerprints, each of which represents a different aspect of the molecule [37, 42, 43].

Clustering is an unsupervised machine learning technique that groups data with similar properties. This technique for statistical data analysis is widely used in cheminformatics [19].

A cluster is, in this case, a collection of molecules which are organized in groups, according to their molecular fingerprints [3, 17].

Despite the large number of clustering methods, only a few of them are widely used in practice. In this paper, only two of them were used, which proved to be suitable for chemical structure analysis: hierarchical clustering and K-means method. Regardless of the method, the results were the same.


All the software used for this article can be installed on Windows, Linux or macOS operating systems.

Initially, built as an environment for statistical computing, R, a GNU project, provides a wide variety of packages for cheminformatics that are suitable for calculating molecular fingerprints and clustering [13, 25, 46]. The latest version of R can be downloaded form the CRAN repository. R Studio is considered one of the best IDEs (integrated development environments) for R and was also used for this article. In this paper, R version 3.6.0 and RStudio version 1.2.1335 were used.

Marvin Sketch version 17.3, from ChemAxon, an academic software package, was used to draw, display and characterize the chemical structures [8]. The molecules were imported in Marvin Sketch using their IUPAC names (International Union of Pure and Applied Chemistry) and then saved as SMILES and SDF formats [31]. Then, they were imported and processed in R [12, 28].

R applications for cheminformatics and computational chemistry

Its flexibility and wide application fields have made the R programming environment a popular choice in a large number of areas.

In the field of cheminformatics, R offers several tools that are able to treat a large variety of issues related to the statistical modelling of chemical information. The rcdk package, version:, used in the present work, provides direct access from the R environment to the CDK (Chemistry Development Kit), a powerful Java framework for cheminformatics [6, 47].

CDK is a collection of free Java libraries that supports a wide variety of cheminformatics functionality. This platform allows us to read different molecular formats, calculate molecular descriptors and evaluate molecular fingerprints.

The cluster package, version 2.1.0, can be used to find groups of molecules that share similar chemical properties [2, 23].

The packages can be installed using the function "install.packages()". The general syntax is listed below:

figure a

To use a package, it must be loaded in the R environment using the function library().

In addition to rcdk, some other packages were also needed:

figure b

Importing and viewing the "drug candidate" molecules in R

In order to manipulate the chemical structures in R, we assigned them a code, starting with CMP1 for methotrexate and ending with CMP24 for the last structure. The Methotrexate molecule (coded as "CMP1") was downloaded from ZINC15, a free database of commercially available compounds, in both SMILE and SDF file format.

In SDF format, the molecule of methotrexate can be imported and visualized in R using the code listed below: [13, 44]

figure c

The result is depicted in Fig. 1:

Fig. 1
figure 1

Methotrexate molecule visualisation in R

All the molecules were imported in SDF format and visualized in R as a grid.

figure d

The result is depicted in Fig. 2:

Fig. 2
figure 2

Molecule set visualization

Computation of the molecular descriptors (physicochemical properties of the molecules)

The rcdk package can also be used to calculate a set of physicochemical properties of the molecules:

The number of atoms:

figure e

The number of chemical bonds:

figure f

The coordinates of the first atom:

figure g

It is also possible to calculate the coordinates for all the atoms present in the molecule:

figure h

R can compute a set of molecular descriptors, grouped into 5 different categories:

figure i

Category 2 (constitutional), important in QSAR, contains 15 descriptors, which are listed below: [14].

figure j

Regarding drug design, the evaluation of AlogP is given a higher importance than that of other descriptors:

figure k

Computation of the molecular fingerprints

Molecular fingerprints can be computed by several methods, but in the case of aromatic compounds the "extended" method is preferred. "Extended" fingerprints have a length (the number of bits ) of 1024, compared to 166 for "maccs "type fingerprints [11, 24, 36].

figure l

Similarly we computed the molecular fingerprints for the entire set of molecules:

figure m

Computation of the intermolecular distances by the Tanimoto index

The Tanimoto coefficient can be expressed as:

$$\begin{aligned} S_{A,B}=c/[a+b-c] \end{aligned}$$

where S is the similarity, a is the number of on bits in molecule A, b is number of on bits in molecule B, and c is the number of on bits in both molecules [49].

Based on molecular fingerprints calculated using the Tanimoto method, the molecular distances between the methotrexate molecule and the other 23 molecules in the group can be evaluated: [24, 30].

figure n

Using this method, a complete set of distances, in matrix form, between each of the 23 molecules of interest was obtained. By analysing these results, it is possible to identify all the molecules located at a certain distance from the target molecule (0.5 in our example): [48].

figure o

From the data presented above we can conclude that only molecules 5 and 7 meet our criteria. This method is the basis for fingerprint-based clustering.

Results and Discussion

In the present study, we used a group of 23 newly synthesized molecules. All of them share the following characteristics: they are pyrazole derivatives, that have never been synthesized, there is no data about them in the literature or in chemical databases, and they have the potential to be drug candidates, such as purine derivatives. Our intention was to check whether the studied chemical compounds can be considered as possible lead molecules [26]. Because the costs of clinical trials are high, even in the preclinical phase, pre-sorting these candidates by computational chemistry and cheminformatics methods would be beneficial [2, 45]. According to the similarity property principle (SPP), which says that drugs with similar molecular structures are likely to have the same properties, a new drug candidate can be identified upon its similarity with another known drug, regardless of how the similarity is evaluated [5]. As a screening criterion, we used the comparison with the traditional methotrexate molecule [32].

Methotrexate is a cytotoxic substance widely used in cancer therapy. It was one of the first purine-inhibiting antimetabolites on the market, and it interferes with the growth of different molecules present in human body, such as like highly reproductive cancer cells. Even though this molecule cannot be considered a "gold standard" for this compound class, the above arguments have contributed to this choice.

Clustering the dataset of molecules

Different types of methods can be used for clustering, including partitioning methods (K-means), hierarchical clustering, fuzzy clustering, density-based clustering and model-based clustering. The K-means and hierarchical clustering were chosen because they are suitable for our goal [29].

The number of clusters

The optimal number of clusters can be estimated using the NbClust package. The function fviz.nbclust is used for visualizing the result [1]. The R code for the elbow method is presented below:

figure p

The result is depicted in Fig. 3:

Fig. 3
figure 3

Estimation of the optimal number of clusters

The optimal number of clusters is 3.

All the considered molecules were then grouped into clusters by taking into account the calculated intermolecular distances.

Hierarchical clustering with the hclust package

The hierarchical clustering algorithm creates clusters with sets of data that are similar internally but different from each other externally [30]. The most common and useful graphical representation of molecular clusters is hierarchical clustering (dendrogram). We performed hierarchical clusterization using Ward’s method [30, 37].

Ward’s method is based on an ANOVA approach and its goal is to maximize the \(r^{2}\) value.

To obtain this dendogram, we used the following R code:

figure q

The graphical representation of the dendrogram is depicted in Fig. 4.

Fig. 4
figure 4

Dendrogram−hierarchical clustering using Ward’s method

K-means Clustering

K-means clustering is one of the most commonly used clustering algorithms because it is easy to code and implement. Each cluster has a centre, which is called a centroid. The algorithm combines the distances between points and centroids [27]. The R code for K-means clustering is shown below.

figure r

The result is presented in Fig. 5:

Fig. 5
figure 5

Polygonal clusters

The molecules included in cluster 1, containing methotrexate (CMP1), were visualized using the following R code:

figure s

The result are depicted in Fig. 6:

Fig. 6
figure 6

Molecules similar to methotrexate

Clustering validation

Statistics for K-means clustering

The R code for the cluster statistics is listed below:

figure t

The most important information for the cluster analysis provided by this function can be considered the silhouette index and Dunn index: [7, 34].

figure u

The Dunn index is equal to the ratio of the smallest inter-cluster distance divided by the largest intra-cluster distance [20].

It takes a value between zero and infinity, and a higher DI means that clusters are compact and well separated. A larger distance between clusters means a better separation, and smaller cluster sizes lead to a higher Dunn Index [1, 21].

The Silhouette coefficient is a method of cluster validation that combines both cohesion and separation [35]. It measures, for each point \(M_{i}\), the mean distance to each cluster, and the mean distance to the other points in its cluster. Silhouette values range between − 1 and 1 [4]. A Silhouette coefficient with a value near +1 indicates that the point is far from its neighbouring cluster and very close to the cluster to which it is assigned. These values are preferred.

The R code for visualizing a Silhouette plot for K-means clustering:

figure v

The plot is visualized in Fig. 7:

Fig. 7
figure 7

Silhouette plot for K-means clustering

The Silhouette plot for K-means clustering reveals a coefficient of 0.4 for the first cluster and a mean value of 0.49 for all the other clusters. These values can be considered acceptable.


Cheminformatics is a dynamic and powerful field that is considered the heart of modern drug design [9, 15]. It plays an important role in collecting, storing and analysing chemical data [18]. It is also an emerging interdisciplinary field that aims to discover new chemical entities that ultimately result in the design of a new active molecule (chemical data) [14, 22].

This work focused on cheminformatics and its application in the discovery and testing of new active molecules. In addition we focused on modern data mining techniques that help chemists and medical researchers to discover, produce and test new active molecules for the treatment of certain diseases.

Our goal was to test in vitro a set of 23 newly synthesized molecules, about which we do not have enough experimental data. The lack of information about the physicochemical properties of the respective molecules, especially those related to QSAR, was supplemented by the computational methods offered by the rcdk software package [12, 25].

Because we studied 23 compounds from the pyrazole class, we expected that at least some of them would behave similar to the like cytotoxic antimetabolite class (purine inhibitors). The clusters were obtained using hierarchical and K-means clustering methods. The results of clustering were confirmed using the Dunn index and Silhouette coefficient.

To avoid the additional costs that pre-clinical and clinical trials for all these compounds would have involved, we tried to reduce the number of "candidate drugs" by computational methods. This reduction was accomplished by calculating the molecular fingerprints of all the studied molecules and then comparing them with the molecular marker methotrexate, which still has a wide use. As a result of this comparison and after the clusterization of the molecules according to the Tanimoto distances, an optimal number of 3 clusters was obtained. In the cluster containing the methotrexate molecule of, marked with a 1, we can also find the molecules marked with a 2, 5, 6 and 7. The remaining 17 molecules are part of the other two clusters [30]. Therefore, starting from the assumption that "similar chemical structures have similar biological properties and actions", the number of compounds worth considering for further studies has been significantly reduced, from 23 to 4, which will lead to a significant decrease in all future costs [9].

Data and material availability

Data and materials are available on GitHub:

The R code used for this paper:


The cemical structures in sdf file format:



  1. Arbelaitz O, Gurrutxaga I, Muguerza J, PéRez JM, Perona I (2013) An extensive comparative study of cluster validity indices. Pattern Recognit 46(1):243–256

    Article  Google Scholar 

  2. Backman Tyler WH, Yiqun C, Thomas G (2011) Chemmine tools: an online service for analyzing and clustering small molecules. Nucleic Acids Res 39(suppl–2):W486–W491

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Bajusz D, Rácz A, Héberger K (2015) Why is tanimoto index an appropriate choice for fingerprint-based similarity calculations? J Cheminform 7(1):20

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  4. Baridam BB (2012) More work on k-means clustering algorithm: the dimensionality problem. Int J Comput Appl 44(2):23–30

    Google Scholar 

  5. Begam BF, Kumar JS (2012) A study on cheminformatics and its applications on modern drug discovery. Procedia Eng 38:1264–1275

    Article  CAS  Google Scholar 

  6. Beisken S, Meinl T, Wiswedel B, de Figueiredo LF, Berthold M, Steinbeck C (2013) Knime-cdk: workflow-driven cheminformatics. BMC Bioinform 14(1):257

    Article  Google Scholar 

  7. Brock G, Pihur V, Datta S, Datta S et al. (2008) clValid, an R package for cluster validation. J Stat Softw 25(4):1–22

    Article  Google Scholar 

  8. ChemAxon L (2013) Marvinsketch.

  9. DiMasi JA, Hansen RW, Grabowski HG (2003) The price of innovation: new estimates of drug development costs. J Health Econ 22(2):151–185

    Article  PubMed  Google Scholar 

  10. Georgiou KR, Scherer MA, Fan CM, Cool JC, King TJ, Foster BK, Xian CJ (2012) Methotrexate chemotherapy reduces osteogenesis but increases adipogenic potential in the bone marrow. J Cell Physiol 227(3):909–918

    Article  CAS  PubMed  Google Scholar 

  11. Godden JW, Stahura FL, Bajorath J (2005) Anatomy of fingerprint search calculations on structurally diverse sets of active compounds. J Chem Inform Model 45(6):1812–1819

    Article  CAS  Google Scholar 

  12. Guha R, Cherto MR (2017) Integrating the CDK with R. Chemical informatics functionality in R, pp 1–17

  13. Guha R et al (2007) Chemical informatics functionality in r. J Stat Softw 18(5):1–16

    Article  Google Scholar 

  14. Guha R, Gilbert K, Fox G, Pierce M, Wild D, Yuan H (2010) Advances in cheminformatics methodologies and infrastructure to support the data mining of large, heterogeneous chemical datasets. Curr Comput Aided Drug Design 6(1):50–67

    Article  CAS  Google Scholar 

  15. Hassan Baig M, Ahmad K, Roy S, Mohammad Ashraf J, Adil M, Haris Siddiqui M, Khan S, Amjad Kamal M, Provazník I, Choi I (2016) Computer aided drug design: success and limitations. Curr Pharma Design 22(5):572–581

    Article  CAS  Google Scholar 

  16. Hughes JP, Rees S, Kalindjian SB, Philpott KL (2011) Principles of early drug discovery. Br J Pharmacol 162(6):1239–1249

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Jacques Julien, Preda Cristian (2014) Functional data clustering: a survey. Adv Data Anal Classif 8(3):231–255

    Article  Google Scholar 

  18. Karthikeyan M, Vyas R (2014) Machine learning methods in chemoinformatics for drug discovery. In: Karthikeyan M, Vyas R (eds) Practical chemoinformatics. Springer, New Delhi, pp 133–194

    Chapter  Google Scholar 

  19. Kovács F, Legány C, Babos A (2005) Cluster validity measurement techniques. In: 6th International symposium of hungarian researchers on computational intelligence, p 35. Citeseer

  20. Kryszczuk K, Hurley P (2010) Estimation of the number of clusters using multiple clustering validity indices. In: International workshop on multiple classifier systems. Springer, pp 114–123

  21. Liu Y, Li Z, Xiong H, Gao X, Wu J (2010) Understanding of internal clustering validation measures. In: 2010 IEEE international conference on data mining. IEEE, pp 911–916

  22. Macalino SJY, Gosu V, Hong S, Choi S (2015) Role of computer-aided drug design in modern drug discovery. Arch Pharm Res 38(9):1686–1701

    Article  CAS  PubMed  Google Scholar 

  23. MacCuish JD, MacCuish NE (2014) Chemoinformatics applications of cluster analysis. Wiley Interdiscip Rev Comput Mol Sci 4(1):34–48

    Article  CAS  Google Scholar 

  24. Martin E, Cao E (2015) Euclidean chemical spaces from molecular fingerprints: hamming distance. J Comput Aided Mol Design 29(5):387–395

    Article  CAS  Google Scholar 

  25. Mente S, Kuhn M (2012) The use of the r language for medicinal chemistry applications. Curr Topics Med Chem 12(18):1957–1964

    Article  CAS  Google Scholar 

  26. Mioc M, Avram S, Tomescu AB, Chiriac DV, Heghes A, Voicu M, Voicu A, Citu C, Kurunczi L (2017) Docking study of 3-mercapto-1, 2, 4-triazole derivatives as inhibitors for vegfr and egfr. Rev Chim 68(3):500–503

    CAS  Google Scholar 

  27. Morissette L, Chartier S (2013) The k-means clustering technique: general considerations and implementation in mathematica. Tutor Quant Methods Psychol 9(1):15–24

    Article  Google Scholar 

  28. Muchmore SW, Edmunds JJ, Stewart KD, Hajduk PJ (2010) Cheminformatic tools for medicinal chemists. J Med Chem 53(13):4830–4841

    Article  CAS  PubMed  Google Scholar 

  29. Murtagh F, Contreras P (2012) Algorithms for hierarchical clustering: an overview. Wiley Interdiscip Rev Data Min Knowl Discov 2(1):86–97

    Article  Google Scholar 

  30. Murtagh F, Legendre P (2014) Ward’s hierarchical agglomerative clustering method: which algorithms implement ward’s criterion? J Classif 31(3):274–295

    Article  Google Scholar 

  31. O’Boyle NM, Banck M, James CA, Morley C, Vandermeersch T, Hutchison GR (2011) Open babel: an open chemical toolbox. J Cheminform 3(1):33

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  32. OBoyle NM (2012) Towards a universal smiles representation-a standard method to generate canonical smiles based on the inchi. J Cheminform 4(1):22

    Article  CAS  Google Scholar 

  33. Prakash N, Gareja DA (2010) Cheminformatics. J Proteomics Bioinform 3:249–252

    Article  CAS  Google Scholar 

  34. Rendón E, Abundez I, Arizmendi A, Quiroz EM (2011) Internal versus external cluster validation indexes. Int J Comput Commun 5(1):27–34

    Google Scholar 

  35. Rendón E, Abundez IM, Gutierrez C, Zagal SD, Arizmendi A, Quiroz EM, Arzate HE (2011) A comparison of internal and external cluster validation indexes. In: Proceedings of the 5th WSEAS international conference on computer engineering and applications, pp 158–163

  36. Rogers D, Hahn M (2010) Extended-connectivity fingerprints. J Chem Inform Model 50(5):742–754

    Article  CAS  Google Scholar 

  37. Saeed F, Salim N, Abdo A (2012) Voting-based consensus clustering for combining multiple clusterings of chemical structures. J Cheminform 4(1):37

    Article  PubMed  PubMed Central  Google Scholar 

  38. Sliwoski G, Kothiwale S, Meiler J, Lowe EW (2014) Computational methods in drug discovery. Pharmacol Rev 66(1):334–395

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  39. Szymański P, Markowicz M, Mikiciuk-Olasik E (2012) Adaptation of high-throughput screening in drug discovery–toxicological screening tests. Int J Mol Sci 13(1):427–452

    Article  PubMed  CAS  Google Scholar 

  40. Taft CA, Da Silva VB et al (2008) Current topics in computer-aided drug design. J Pharm Sci 97(3):1089–1098

    Article  CAS  PubMed  Google Scholar 

  41. Taguchi Y-H (2017) Identification of candidate drugs using tensor-decomposition-based unsupervised feature extraction in integrated analysis of gene expression between diseases and drugmatrix datasets. Sci Rep 7(1):13733

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  42. Vogt M, Stumpfe D, Geppert H, Bajorath J (2010) Scaffold hopping using two-dimensional fingerprints: true potential, black magic, or a hopeless endeavor? Guidelines for virtual screening. J Med Chem 53(15):5707–5715

    Article  CAS  PubMed  Google Scholar 

  43. Wagener M, van Geerestein VJ (2000) Potential drugs and nondrugs: prediction and identification of important structural features. J Chem Inf Comput Sci 40(2):280–292

    Article  CAS  PubMed  Google Scholar 

  44. Warr WA (2011) Representation of chemical structures. Wiley Interdiscip Rev Comput Mol Sci 1(4):557–579

    Article  CAS  Google Scholar 

  45. Willett P (2009) Similarity methods in chemoinformatics. Annu Rev Inform Sci Technol 43:3–71

    Article  Google Scholar 

  46. Willett Peter (2010) Similarity searching using 2d structural fingerprints. In: Chemoinformatics and computational chemical biology. Springer, pp 133–158

  47. Willighagen EL, Mayfield JW, Alvarsson J, Berg A, Carlsson L, Jeliazkova N, Kuhn S, Pluskal T, Rojas-Chertó M, Spjuth O et al (2017) The chemistry development kit (cdk) v2. 0: atom typing, depiction, molecular formulas, and substructure searching. J Cheminform 9(1):33

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  48. Zhang B, Vogt M, Maggiora GM, Bajorath J (2015) Design of chemical space networks using a tanimoto similarity variant based upon maximum common substructures. J Comput Aided Mol design 29(10):937–950

    Article  CAS  Google Scholar 

  49. Zhang C, Idelbayev Y, Roberts N, Tao Y, Nannapaneni Y, Duggan BM, Min J, Lin EC, Gerwick EC, Cottrell GW et al (2017) Small molecule accurate recognition technology (smart) to enhance natural products research. Sci Rep 7(1):14243

    Article  PubMed  PubMed Central  CAS  Google Scholar 

Download references


We would like to thank Dr. Sabina Nitu who kindly provide us the structures and additional information regarding the studied compounds.

Author information

Authors and Affiliations



AV—clustering and computing in R and Chem Axon—MS, MV, DV and VD—documenting and selection of possible active molecules AV, MV, ND—wrote the first draft of the manuscript. All authors contributed to the interpretation of experimental data, discussion and the preparation of the final manuscript. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Narcis Duteanu.

Ethics declarations

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Voicu, A., Duteanu, N., Voicu, M. et al. The rcdk and cluster R packages applied to drug candidate selection. J Cheminform 12, 3 (2020).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: