# Computer-Aided Classification of New Psychoactive Substances

Authors: Alina Bărbulescu; Lucica Barbeș; Cristian-Ştefan Dumitriu
Journal: Journal of Chemistry
 (2021)

---

## Abstract

The appearance on the free market of synthetic cannabinoids raised the researchers’ interest in establishing their molecular similarity by QSAR analysis. A rigorous criterion for classifying drugs is their chemical structure. Therefore, this article presents the structural similarity of two groups of drugs: benzoylindoles and phenylacetylindoles. Statistical analysis and clustering of the molecules are performed based on their numerical characteristics extracted using Cheminformatics methods. Their similarities/dissimilarities are emphasized using the dendrograms and heat map. The highest discrepancies are found in the phenylacetylindoles group.


---

## Body

## 1. Introduction

The consumption of psychoactive substances is among the principal causes of young people’s mortality and health issues. A significant proportion of adults consumed some drugs at a certain period in their lives or presented with drug addiction. Synthetic cannabinoids are found in “herbal smoking mixes” and mixtures with other psychoactive substances such as sedatives/hypnotics, hallucinogens, and stimulants [1]. Manufactured in clandestine laboratories and widely sold, the class of new psychoactive substances known as different spices poses a growing threat to public health, such as Spice, Spice Gold, and Silver. In some countries, acute poisonings of users’ groups [2–4] have been recorded. Synthetic cannabinoids have been missold as THC (Δ9-tetrahydrocannabinol) or CBD (cannabidiol) and are packed in tablets containing powder or utilized in devices as liquids [5].Another recent development is the popularity of “vaping” among young people, synthetic cannabinoids being discovered in cartridges filled with liquid for use in e-cigarettes [6].Psychoactive drugs increase the dopamine available in the brain (called the reward circuit) with participation in pleasure modulation. Synthetic cannabinoids (SCs) are much stronger and more toxic than users can often think. Their increasing prevalence has led to severe health problems among young people and adults. Consumption of these new psychoactive substances poses a health risk and generates severe issues in adapting to reality, education, work, family life, incidents, and road accidents [7, 8]. Research has been performed on the cannabinoids detected in different products to determine how cannabinoids affect the human body [9]. They also aim to analyze the possibilities of using these substances to treat certain neurodegenerative diseases, drug addiction, emotional disorders, or cancer [10]. However, separating the psychoactive effects and the medicinal properties of these drugs is challenging [11].In about 8–10 years, researchers have tried different dynamic simulation models on rapid identification and establishing the relationships of addiction between the chemical structure and effects on consumers to detect possible ways for curing addiction and overdose [12]. Functionally similar to THC (the cannabis active principle), SCs bind to the same receptors in the human body [6].The group of synthetic cannabinoids or synthetic cannabimimetic is the largest monitored by the EU Early Warning System. In December 2016, 169 synthetic cannabinoids were reported to the EMCDDA and increased to 190 substances reported in 2018 [8]. About 280 synthetic cannabinoids were known at the end of 2019 [6].SCs fall into different structural groups [13, 14], such as aminoalkylindoles (which include naphthoylindoles, phenylacetylindoles, naphthylmethylindoles, and benzoylindoles), classical cannabinoids, nonclassical cannabinoids, hybrid cannabinoids, eicosanoids, benzimidazole, indazole-based SCs, and others (naphthoylpyrroles, naphthylmethylindenes, carboxamide-type synthetic cannabinoids, etc.) [15, 16]. The first three groups are presented in Figure 1.Figure 1
Synthetic cannabinoids.In the investigation of the various chemical synthesis compounds or any chemical synthesis preparation, the chemical analyst establishes the relationships between physicochemical parameters and biological (pharmacological) properties, the degree of structural superposability of the chemical products, using filters based on computational structure and quantitative relationship models of structure-activity (QSARs models) [11].The search for similarity using molecular fingerprints are among the most common approaches in identifying active chemical compounds proposed for analysis in this study. Many structural descriptors or 2D and 3D ownership are encoded in different representations of fingerprints. All fingerprints eventually transform molecular likeness analysis into a paired comparison of template compound patterns and the entire database [17].Several fingerprint types are known, each representing a difference in the molecule. Artificial intelligence techniques have been employed in the last period for data analysis in various domains [18–21]. Cheminformatics is a tool used to analyze statistical data based on the digital footprint, which is the dynamic field of modern design for pharmaceuticals, in particular. Cheminformatics plays an essential role in accumulating, grouping, examining chemical data, and finding new entities based on which new chemical structures can build active molecules [22].Swandana et al. [23] showed that the pharmacokinetic parameters could be predicted using the in silico method. The authors emphasize some software (admetSAR, Chemicalize, Molinspiration, QikProp, SwissADME, and pkCSM software) to predict oral-systemic drugs’ pharmacokinetics characteristics. Other scientists [24, 25] emphasized different computational methods utilized for discovering drugs and their limitations.Finding similarities of different molecules has been of interest since the 90s [26] that can be easily done using the new capabilities of the new software. The rcdk, ChemmineR, and rpubchem packages of R are powerful tools in the field of cheminformatics [27–32]. Their utilization helps the cheminformatics-oriented R user to organize chemical information efficiently. Guha et al. [27] showed that one of the main tasks is to include further coverage of the CDK API by developing the CDK itself, such as extending the 3D builder model. The authors noted that the CDK API could be accessed directly via rJava, so the functions provided by the rcdk package make cheminformatics functionality readily available for the end-users.In this study, we used theR software capacity to characterize 14 new cannabinoids from the benzoylindoles and phenylacetylindoles groups and determine the similarities among their molecules. Clustering techniques are applied for grouping the studied molecules based on the determined characteristics to cross-validate similar attributes in the same group of drugs.


## 2. Materials and Methods

The structures subject to the study have been downloaded as .sdf files from PubChem [33]. They belong to the class aminoalkylindoles, subclasses benzoylindoles, and phenylacetylindoles.TheR software has been used for performing the study. For this aim, the following libraries have been installed: ChemmineR, chemometrics, cluster, factoextra, fingerprint, fmcsR, ggplot2, gridExtra, iqspr, NbClust, rcdk, rJava, rgl, and vegan. They are necessary for drawing the structures, finding the characteristics of the chemical structures, or determining the appurtenance of the chemical structures to some clusters.The stages of this study are presented in Figure2.Figure 2
The flowchart.They are(i)
Import the molecule from ChemPub, as SDF files, and visualize them(ii)
Visualize the molecule structure and molecular formula (MF), and determine the molar weight filter (MW), the types and numbers of atoms, and the functional groups together with the corresponding frequencies(iii)
Compute different descriptors utilizing the library ChemmineOB:(i)
HBA1: number of hydrogen bond acceptors 1(ii)
HBA2: number of hydrogen bond acceptors 2(iii)
HBD: number of hydrogen bond donors(iv)
logP: logarithmic of partition coefficients (log P)(v)
MR: molar refractivity(vi)
TPSA: topological polar surface areaThe hydrogen bonds are among the most significant interactions between the solute and solvent. In binding the ligands, the hydrogen bonds have important contributions to the substrates, antagonists, and agonist recognition [34]. The octanol/water partition coefficient is used to characterize the lipophilicity of the compounds and is defined by the ratio [solute]octanol/[solute]water. Generally, the logarithm of this ratio, logP, is reported [35].The molar refractivity (MR) measures the steric factor, which is the volume covered by an atom or a group of atoms. It is correlated to the London dispersive forces acting during the drug-receptor interaction [36].The hydrophilic interaction and the hydrogen bond formation are reflected by the polar surface (PSA) [37], computed by eliminating the area of nonpolar hydrogen, carbon, and halogen atoms from the molecular surface. Determined by the method of Ertl et al. [38], it is referred to as topological PSA (TPSA) [39]. According to the available literature, Turner and Agatonovic-Kustrin [40] showed that most drugs have a molar weight lower than 500 and PSA less than 120 Å.(i)
Search the atom pair similarity by cmp.similarity command from ChemmineR. For this aim, the Tanimoto index was utilized [41] due to its performance [26]. For dichotomous variables, it is computed by(1)SSA,B=ca+b−c.The following relationship exists between the distance and similarity metrics:(2)similarity=11+distance.(ii)
Determine the intermolecular distances, using the fmcs function from the fmscR library and the compound similarities.(iii)
Clustering the molecules’ sets: for this aim, the binning clustering (using cmp.function that implements the Tanimoto and Tversky methods) [41, 42], the Jarvis–Patrick algorithm [43], the k-means, and the hierarchical clustering (with the NbClust package) have been utilized [44]. Comparisons are provided.


## 3. Results and Discussion

The structures of the molecules from the benzoylindoles group (Group 1) are represented in Figure3, and those from the phenylacetylindoles (Group 2) are represented in Figure 4, together with the identification number from PubChem.Figure 3
The structures of the molecules from the benzoylindoles group.Figure 4
The structures of the molecules from the phenylacetylindoles group.A similar representation, but showing the number of atoms, is presented in Figure5, for the molecule 9889172. It contains an automatic numbering added by R software when extracting the information from PubChem: the atoms are numbered consecutively (from 1 to the total number of atoms), based on the following order: elements in the seventh group, from the highest to the lowest atomic number, elements in the sixth group, from the highest to the lowest atomic number, and so on, from left to right in the Mendeleev table.Figure 5
The structure of the molecule 9889172.Table1 contains the molecules formula, the types and numbers of atoms in each molecule from both groups, the molecular weight filters (MW), the number of rings, and the number of aromatic rings retrieved using ChemmineR. The molar weight filters vary in a larger range in the first group by comparison to the second one.Table 1
The cod identification (CID), molecular formula (MF), molar weight filter (MW), types, and the numbers of atoms, rings, and aromatic rings in the molecules.
CIDFormulaMWCFHINOClRingsAromaticGroup 19889172C20H19FINO435.273820119111033117587582C20H19ClINO451.72842001911113310226340C22H23IN2O458.33532202312104356841530C21H23NO2321.41282102301203357507911C21H23NO2321.41282102301203357507905C20H21NO2307.38622002101203356463C23H26N2O3378.464123026023043Group 244397641C22H25NO2335.43942202501203344397500C21H22ClNO339.85852102201113344397540C22H25NO2335.43942202501203353494930C24H28N2O2376.49132402802204311616723C22H25NO319.44002202501103311493740C22H25NO2335.43942202501203353394099C25H29NO2375.503325029012043All but the first molecule do not contain F; all the molecules in Group 2 do not contain I. Only one molecule in each group contains Cl. The number of N and O is one or two. The number of C and H atoms is lower in the first group than the second one.The analysis shows that no molecule contains groups RNH2, R2NH, ROPO3, ROH, RCHO, RCOOH, RCOOR, RCCH, and RCN, but all contain rings, most of them aromatic.Table2 contains the most important descriptors that characterize the molecules. The values of HBA1 are lower for the first group, HBD is absent for both groups, logP is in the range of 3.4594–5.8860 for Group 1 and 4.4975–6.0457 for Group 2 (showing higher hydrophilicity for the second group compared to the first one). TPSA varies in larger limits for Group 1 than for Group 2. The higher the TPSA is, the higher the drug transport is. In the studied case, the molecule with CID = 56463 has the highest TPSA.Table 2
Molecules’ descriptors.
CIDHBA1HBA2HBDlogPMRTPSAGroup 1988917220205.6167105.070522.0011758758220205.8860109.815522.001022634025304.8991119.330525.245684153025305.071198.794531.235750791125305.071198.794531.235750790523304.681093.987531.235646330403.4594114.349543.70Group 24439764127305.2655103.601531.234439750023205.9103102.119522.004439754027305.2655103.601531.235349493031404.4975117.912534.471161672326205.5653102.075522.001149374027305.2655103.601531.235339409931306.0457115.908531.23The similarity among the pairs of atoms in different molecules is given in Table3.Table 3
Tanimoto coefficients for the atoms similarity.
Group 1CID9889172117587582102263405750791156841530575079055646398891721.00000.84620.67880.61400.53760.49010.44691175875821.00000.67880.53760.61400.49010.4469102263401.00000.51390.57740.48590.5054575079111.00000.72500.84970.5425568415301.00000.66350.5069575079051.00000.4953564631.0000Group 2CID44397641443975004439754053494930116167231149374053394099443976411.00000.66090.72490.57770.70000.76540.5455443975001.00000.78400.57310.84670.68020.5545443975401.00000.68730.85260.76540.6545534949301.00000.61580.58880.6882116167231.00000.73570.5769114937401.00000.5490533940991.0000Most similarity coefficients are above 0.5, the best ones corresponding to the pair of CIDs (57507911, 57507905), (988917, 0.8462) in Group 1 and (44397540, 11616723), (44397500, 11616723) in Group 2.The compounds’ similarity has been studied using the fmcs function. The output shows the size of the query molecule and that of the target molecule (Query_Size and Target_Size, respectively), the Tanimoto coefficient, and the overlap coefficient, measuring the structures’ superposition degree.For example, in Table4, the query molecule has CID 9889172, and the target molecule in row three has CID 117587582. For both, the query and target sizes are equal (24 atoms, the number of hydrogen atoms is not counted), the Tanimoto coefficient is 0.92, and the overlap coefficient is 0.9583. The higher the Tanimoto and the overlap coefficients are, the higher the similarity of the structures is. The highest dissimilarity is that between molecule CID 56463 and the other molecules in Group 1.Table 4
Similarity measures for the molecule CID 9889172 and the molecules in Group 1.
CIDQuery_SizeTarget_SizeTanimoto_CoefficientOverlap_Coefficient988917224241.00001.000011758758224240.92000.95831022634024260.85190.95835750791124240.84690.91675684153024240.84690.91675750790524230.80770.91305646324280.57580.7917The similarity coefficients in Group 2 are higher than in the first group, as shown in Table5. Hence, the compounds in the second group have more similar structural characteristics than those in the first group.Table 5
Similarity measures for the molecule CID 44397641 and the molecules in Group 2.
CIDQuery_SizeTarget_SizeTanimoto_CoefficientOverlap_Coefficient4439764125251.00001.00004439750025240.88460.95834439754025250.85190.92005349493025280.76670.92001161672325240.88460.95831149374025250.85190.92005339409925280.76670.9200The computed similarity coefficients between the pairs of atoms in different molecules from Group 1 and those from Group 2 are shown in Table6.Table 6
Pair-atoms similarities for molecules from different groups.
CID4439764144397500443975405349493011616723114937405339409998891720.48830.52480.49220.45660.52480.48840.45981175875820.48830.54620.49220.45660.52480.48840.4598102263400.47060.49130.48100.60870.50620.47410.4894575079110.63640.58170.68420.56460.66770.66000.5000568415300.66490.52490.61340.51390.57260.71940.4696575079050.60290.49010.59370.48470.53330.63610.4244564630.48030.39150.44560.45950.421730.49340.4766The similarities of atom pairs from 9889172 (Group 1) with those from 44397500 and 1161723 (Group 2) are higher than those from 57507905 and 56463 (Group 1). There are high similarities of atoms pairs from 57507911 (Group 1), on the one hand, and 44397540 and 11616723 (Group 2), on the other hand (higher than those with other molecules from Group 1). Similarities above 0.60, which are significant, are noticed for other pairs of molecules belonging to different groups.A similar result has been obtained using the Tanimoto coefficient for molecules’ similarities. Table7 displays the Tanimoto coefficients and their ranked similarities (from 1 to 7) for pairs of molecules: one from the first group and another from the second group. Blank means a similarity rank over 8. There is only one molecule not listed in Table 7(ID 53394099, Group 2), whose similarity with all the molecules in Group 1 is very low.Table 7
Tanimoto coefficients for the atoms similarity in different groups of molecules.
CID44397641443975004439754053494930116167231149374098891720.5262 (6)0.5262 (7)1175875820.5262 (7)102263400.6096 (4)0.5075 (7)575079110.6374 (7)0.6851 (3)0.6609 (6)568415300.6657 (5)0.6145 (6)0.7202 (4)575079050.6040 (5)0.5948 (6)0.5347 (7)0.6371 (4)564630.4815 (7)0.4945 (6)Based on the structural characteristics, the molecules’ clustering has been performed. The results for Group 1 are presented in the following. The binning clustering, with cutoffs (lower similarity bounds) of 0.3, 0.6, and 0.9 resulted in 1, 2, or 7 clusters, respectively (for both groups). The Jarvis–Patrick algorithm provided 2 clusters when considering four neighbors, and a single cluster when considering 5 and 6 neighbors. The optimal number of clusters (3) in the k-means algorithm has been selected by the elbow method.Running the k-means algorithm resulted in three clusters with the sizes 1, 3, and 3 and between the sum of squares/total sum of squares = 77.0%. The first clusters contain the first three molecules; the second one, the next three molecules; and the last cluster contains molecule CID 56463.For the second group, the algorithms binning and Jarvis–Patrick algorithms provided a single cluster.The heat map and the dendrogram resulted from the hierarchical clustering are presented in Figure6 for Group 1 and in Figure 7 for Group 2.Figure 6
The heat map and dendrogram for Group 1.Figure 7
The heat map and dendrogram for Group 2.In the heat maps, the darker the color is, the lower the distance between the compounds is, so the higher the similarity is. Hence, the pair of compounds whose corresponding squares are dark green are the most similar. A very low similarity is noticed for the pairs of compounds situated in the bottom right-hand side of both heat maps.The dendrograms (right-hand side of Figures6 and 7) confirm the similarity strength. The higher the branch between two compounds is, the lower the similarity is.In Figure6, the highest similarity is found between the groups of molecules (117587582 and 9889172) and (79507911 and 56841530). They both have a Tanimoto coefficient of 0.900 and an overlap coefficient of 0.9583 (listed in Table 4 for the first pair, not listed for the second one for the lack of space). From the viewpoint of similarity intensity, the second place is occupied by the pairs of molecules (10226340, 117587582) and (10226340, 9889172), with a Tanimoto coefficient of 0.8519 and the same overlap coefficient.In Figure7, the most similar molecules are in the pairs (44397500, 11616723), (44397641, 44397500), and (44397641, 11616723), which is in concordance with the results from Table 6, because the second and third pairs have the same Tanimoto coefficient (0.8846) and the overlap coefficient 0.8846. The Tanimoto coefficient of the pair (11616723, 44397500) is 0.9152, and the overlap coefficient is 0.9011. They are not listed in the tables to avoid repeating similar information.Analyzing Figures6 and 7, remark that the molecules in Group 2 have the distance between them, measured on the dendrograms, lower than those from Group 1, confirming the previous findings related to the molecules’ similarity.The same procedures have been utilized, putting together all the 14 molecules. The heat map and the dendrogram are shown in Figure8.Figure 8
The heat map and dendrogram for Groups 1 and 2.The highest similarities are found for the triple (44397500, 44397540, 11616723) from Group 2 and the pairs (957507905, 56481530) and (117587582, 9889172) from the first group.The clusters determined by the k-means algorithm on the 14 molecules are shown in Figure9.Figure 9
The clusters.The first cluster contains all the molecules in Group 2, a result consistent with the previous findings. The first group is divided into two clusters: the first one containing the first three molecules and the second one, the rest. The result is consistent with the dendrogram in Figure6. Cluster 2 (Figure 9) contains the same elements as the top subtree in Figure 7.


## 4. Conclusions

Cheminformatics is a study field that offers many facilities in collecting, organizing, and studying chemical information [45]. Its contribution is essential to drug design [46, 47]. In this view, our article focused on some applications of Cheminformatics to analyze the benzoylindoles and phenylacetylindoles groups of drugs. The structures and molecules functional groups have been emphasized using the facilities provided by rcdk and ChemmineR packages. Structural similarities of different molecules have been determined using the Tanimoto index and clustering techniques. A higher similarity in the second group compared to the first one is observed. Further practical studies should confirm the hypothesis of the similar actions and effects of these drugs and the possible cure using the same inhibitors.



---

*Source: 4816970-2021-12-31.xml*