The outline of our analysis is presented in Figure 1. Each node in the diagram represents a data source and the method followed to compute the similarity between different compounds. The psychoactive drugs investigated in this study are shown in Figure 2 (Ray, 2010) and Figure 3 (Rickli et al., 2015, 2016), with the phenethylamine and tryptamine functional groups highlighted in red and blue, respectively. The upper node in Figure 1 represents data on the reported subjective effects of the compounds, which was provided by Erowid Center’s Experience Vaults (Erowid et al., 2017). Latent Semantic Analysis (LSA; Landauer, 2006) was used to compute the semantic similarity between sets of reports of subjective effects corresponding to the different psychedelic compounds, following our analysis in a previously published article (Sanz et al., 2018). The left node in Figure 1 corresponds to data on the molecular structure of compounds included in this study (identified using simplified molecular-input line-entry system [SMILES] strings). The 2D molecular structure similarity was computed using ChemMine Tools² (Backman et al., 2011). Finally, the right node in Figure 1 corresponds to data on the binding affinity profiles of the molecules (Ray, 2010; Rickli et al., 2015, 2016), obtained from assays performed by the NIMH Psychoactive Drug Screening Program (PDSP)³ and by researchers from the University of Basel in Switzerland, as described in Rickli et al. (2015, 2016). The similarity between the binding affinity profiles of each pair of compounds was obtained as the Pearson’s linear correlation coefficient of their sign-inverted and log-scaled equilibrium dissociation constants (pK_i).

Data on the binding affinity of 18 psychedelic molecules was obtained from tables in the supplementary material of Ray (2010). This data was complemented with binding affinity measurements of 19 compounds published in Rickli et al. (2015, 2016), which were used to provide independent validation of the main results of our analysis.

Two primary sources of data comprise the binding affinities published in Ray (2010). Of the drugs investigated in this work, data for LSD and ibogaine were obtained from the literature. For LSD this information came from NIMH-PDSP assays (Nichols et al., 2002), while the data for ibogaine was compiled from different sources (Ray, 2010). New PDSP binding assays were conducted for the remaining molecules using the methodology described in Ray (2010) and Glennon et al. (2000). Of all sites tested only 42 presented at least one “hit” (i.e., equilibrium dissociation constant K_i < 10,000 nm). These sites included 5-HT receptors (5-HT_2A, 5-HT_2C, 5-HT_2B, 5-HT_1A, 5-HT_1B, 5-HT_1D, 5-HT_1E, 5-HT_5A, 5-HT₆, 5-HT₇), dopamine receptors (D₁, D₂, D₃, D₄, D₅), adrenergic receptors (α_1A, α_1B, α_2A, α_2B, α_2C, β₁, β₂), serotonin transporter (SERT), DA transporter (DAT), norepinephrine transporter (NET), imidazoline1 receptor (I₁), sigma receptors (σ₁, σ₂), delta opioid receptor (DOR), kappa opioid receptor (KOR), mu opioid receptor (MOR), muscarinic receptors (M₁, M₂, M₃, M₄, M₅), histamine receptors (H₁, H₂), cannabinoid receptors (CB₁, CB₂), calcium+ ion channel (Ca+) and N-methyl D-aspartate (NMDA) glutamate receptor. The procedure described in Ray (2010) was followed for the normalization of the K_i values⁴. A logarithmic transformation with a sign inversion was applied (since higher affinities result in lower K_i values). This transformation was given by pK_i = −log₁₀(K_i; Neubig et al., 2003). Since the K_i values generated by the NIMH-PDSP assays were capped at 10,000, the lowest possible pK_i value was −4 (note that pK_i = −4 was also assigned in other circumstances, as described in Ray (2010); e.g., when the primary assay did not produce >50% inhibition). All the K_i values are available as online Supplementary Material in Ray (2010); as well as in Supplementary Table S1 (raw K_i values) and Supplementary Table S2 (normalized pK_i values) of the Supplementary Material.

The radioligand binding assays conducted in Rickli et al. (2015, 2016) are described with detail elsewhere (Hysek et al., 2012; Simmler et al., 2013). Briefly, measurements were performed at membrane preparations of human embryonic kidney cells. Affinity values (Ki) were derived from IC₅₀ values (Cheng-Prusoff equation) obtained from 3–5 independently replicated concentration-response curves; this procedure was performed at the following binding sites: 5-HT receptors (5-HT_1A, 5-HT_2A, 5-HT_2C), dopamine receptors (D₁, D₂, D₃), adrenergic receptors (α_1A, α_2A), the histamine H₁ receptor, trace amine-associated receptors (TAAR₁, from rat and mouse genes), SERT, DAT and NET. A logarithmic transformation was applied as described above to obtain the pK_i values. The supplementary information includes the analysis (Figure A1) of activation potency measurements (EC₅₀) at 5-HT_2A receptors and of monoamine transporter inhibition at NET, DAT and SERT. Details of these experiments are provided in Rickli et al. (2015, 2016).

Finally, to compare the binding affinities of two psychedelic compounds, we computed the Pearson’s linear correlation coefficient between the corresponding sets of pK_i values. When applied to all pairs of compounds, this resulted in matrices containing all correlation coefficients. Whenever the compounds under comparison presented missing data at a subset of sites, the correlation coefficient was computed using only the pK_i values that were simultaneously available for both compounds.

A subset of the molecules investigated in Ray (2010) was included in the present study (Figure 2). These partially overlap with those investigated in Rickli et al. (2015, 2016); those that were assayed in the latter but not in the former are presented in Figure 3. We followed two criteria for inclusion of compounds in our study. First, we discarded molecules lacking binding affinity data at most of the sites. Second, only compounds with a sufficiently high (>10) number of subjective reports in the Erowid database were included. The application of these criteria to the data presented in Ray (2010) resulted in the inclusion of 18 compounds, 10 of which were classified as substituted phenethylamines and seven as substituted tryptamines. The same criteria applied to the data provided in Rickli et al. (2015, 2016) resulted in the inclusion of 19 compounds (11 phenethylamines and eight tryptamines). LSD belongs to the class of ergolines, which can be considered as a subtype of tryptamines (“rigidified tryptamines,” Nichols, 2016). Therefore, for reasons of simplicity, we classified LSD as a tryptamine (note however that the phenethylamine functional group can be found in the tetracyclic ring structure of ergolines). For the purposes of this work, this classification is analogous to that into phenylalkylamines (which includes phenethylamines) and indolealkylamines (which include tryptamines, but also ergolines and β-carbolines; Glennon, 1994).

The 18 compounds from Ray (2010) that were included in this study are the following:

Ibogaine, a naturally occurring tryptamine which can be found in plants of the Apocynaceae family (Schultes and Hofmann, 1979). The synthesis of ibogaine was achieved in 1956 by M. M. Janot and R. Goutarel (Goutarel and Janot, 1956). Missing data for 5-HT_2B, 5-HT_1E, 5-HT_5A, 5-HT₆, 5-HT₇, D₄, D₅, α_1A, α_1B, α_2A, α_2B, α_2C, β₂, NET, I₁, M₄, M₅, H₂, CB₁, CB₂ and Ca+. While ibogaine is used primarily in ritual and therapeutic contexts, the range of receptors involved (including, but not limited to, the 5-HT_2A receptor; Sweetnam et al., 1995) and its equally heterogenous subjective effects allow for variability that is of advantage when investigating the correlation between both dimensions.

The additional compounds assayed in Rickli et al. (2015, 2016) and included in this study are the following:

While many reports of experiences with mushrooms of the Psilocybe genus (containing psilocybin as their main active compound) were available, we decided not to use them as a proxy for the isolated active principle. Mushrooms can contain other psychoactive alkaloids (e.g., baeocystin and norbaeocystin) capable of influencing the subjective effects (Leung and Paul, 1968). They may also contain phenethylamine, which could underlie some of the undesired effects associated with consumption of Psilocybe mushrooms, such as tachycardia, nausea, anxiety (Beck et al., 1998). A survey reported that these undesired effects were frequently observed in users who consumed mushrooms, but not in those who consumed isolated psilocybin (van Amsterdam et al., 2011). Finally, interoception (understood as the sensing of peripheral bodily signals) is known to play a role in anticipatory events that can modulate behavior and subjective effects even before a drug has time to reach the brain and elicit its pharmacological action (Wise et al., 2008). Hofmann and Wasson speculated that the chewing of mushrooms could elicit a faster onset of the subjective effects relative to ingesting a pill containing psilocybin, because in the former case psilocybin can be absorbed by the mucous membrane, while in the latter case the pills must first reach the stomach and dissolve before psilocybin can enter the system (an account of this speculation in the context of a traditional mushroom ceremony carried out with pills containing psiloybin instead of mushrooms can be found in Hofmann, 1980).

Reports of the subjective effects of the psychedelic compounds included in this study were acquired from Erowid’s Experience Vaults⁵. Erowid Center is a nonprofit, member-supported organization providing access to reliable, non-judgmental information about psychoactive plants, chemicals and related issues, containing a large number (>20,000) of published reports of the subjective effects of different psychoactive substances. The full Erowid corpus and the classification of its reports into different categories is extensively described in a previous work (Sanz et al., 2018). A total of 16 reports were obtained for mescaline, 143 for 2C-B, 206 for 2C-E, 101 for 2C-T-2, 36 for DOB, 32 for DOI, 23 for DOM, 19 for TMA-2, 63 for MDA, 770 for MDMA, 236 for DMT, 247 for 5-MeO-DMT, 69 for 5-MeO-MiPT, 45 for DiPT, 182 for 5-MeO-DiPT, 137 for DPT, 718 for LSD, 32 for ibogaine, 64 for 2C-C, 383 for 2C-I, 57 for 2C-P, 16 for 2C-T-4, 171 for 2C-T-7, 48 for 2C-D, 144 for 25I-NBOMe, 51 for 4-OH-MET, 109 for 5-MeO-AMT, 208 for 4-OH-DiPT and 8 for psilocin/psilocybin. All the Erowid users authorize the use of their narratives with research purposes at the time of submission. The data obtained from Erowid’s Experience Vaults was completely anonymous, consisting of lists of terms derived from the tokenization of all the reports associated with each substance, without links to individual reports nor the users who submitted them. The measures followed to secure data confidentiality and privacy at the Erowid Center servers are described on its website. This data is both accessible online and available upon request at the Erowid Center for its use with research purposes.

The preprocessing of the reports was performed using the Natural Language Toolkit (NLTK⁶) in Python 3.4.6 (Bird, 2006). The reports were first separated into individual words after discarding all punctuation marks. Each word was lemmatized using NLTK (i.e., converted to the root from which it is inflected). All words were converted to lowercase and lemmatized words containing less than three characters were discarded. As described in a previous work (Sanz et al., 2018), words including substance names, slang variations and terms related to possible routes of administration were removed from the reports.

To quantify the semantic similarity between the reports of the psychedelic compounds included in this study we applied LSA (Landauer, 2006), a natural language processing technique based on the hypothesis that words with similar meaning appear with similar frequency in texts (Sahlgren, 2008). The rationale behind LSA is comparing the profile of word frequencies after reducing the rank of the word-by-document matrix, thus attenuating the problems of sparseness and lack of context-sensitive term occurrence. Before applying LSA we computed the frequency of the different words using the term frequency-inverse document frequency (tf-idf) transform, as implemented in scikit-learn⁷ (Pedregosa et al., 2011; Leskovec et al., 2014). The tf-idf transform computes a matrix in which rows are unique words in the corpus and columns represent “documents” (for this analysis, a “document” contains all the reports associated with a psychoactive compound contained in the published corpus from the Erowid Experience Vaults, see Sanz et al., 2018). The product of the term frequency and the inverse document frequency determines the entries of this matrix. The term frequency is defined as the number of times each term appears in each document. The inverse document frequency is defined as the logarithmically scaled inverse fraction of the documents containing the term. To eliminate very frequent/rare terms from the corpus, only those terms appearing in more/less than 5%/95% of the documents were retained, resulting in a vocabulary of 1,466 words.

The word-by-document matrix obtained using the tf-idf transform was decomposed into the product of three matrices using Singular Value Decomposition (SDV; Klema and Laub, 1980). Of the three resulting matrices (U, S, V), S contains in its diagonal the matrix of singular values ordered by size, and U and V are real unitary matrices (their size is determined by the number of words and documents, respectively). To reduce the number of linearly independent rows (terms) while preserving the similarity structure among columns (documents), only the first D largest singular values were retained and all others were set to zero, resulting in the reduced matrix of singular values S*. Following previous work (Sanz et al., 2018), we retained the D = 20 largest singular values. Computing the product of the matrices U, S*, V yields a low-rank approximation of the word-by-document frequency matrix.

Finally, the correlation matrix containing the Pearson’s linear correlation coefficient between the columns of the rank-reduced word-by-document matrix (corresponding to each of the drugs) was computed as a metric of the semantic similarity of the reports of each pair of compounds.

To measure the 2D structural similarity between the molecular structures of the compounds we computed atom pair similarities using the Tanimoto coefficient, as implemented in ChemMine Tools⁸ (Backman et al., 2011). Each molecule was represented as a set of atom pairs, comprising pairs of non-hydrogen atom types together with the distances between them (computed as the lengths of the shortest bond-by-bond paths between them, Sheridan et al., 1996). The Tanimoto coefficient provides an index ranging between 0 and 1, obtained as the size of the intersection of the atom pair sets divided by the size of the union of the atom pair sets. This method has been shown to perform more efficiently than other metrics and is a standard tool to quantify 2D molecular similarity (Chen and Reynolds, 2002; Willett, 2006).

We applied principal component analysis (PCA) to reduce the word-by-document matrix into a smaller number of components that capture idiosyncratic features of the reported experiences in the Erowid database. PCA was performed using an algorithm based on SVD as implemented in MATLAB 2014 (Wall et al., 2003). We retained the first five components explaining the highest variance in the tf-idf rank-reduced word frequency matrices. For each component we obtained 1,466 coefficients, which were used to represent in word cloud format,⁹ the words that were most representative of each component.

The scores of each document (i.e., the combined set of subjective reports of a given psychoactive substance in the Erowid corpus) were used to produce radar plots quantifying the presence of the five principal components in the narratives. By correlating the binding affinity of each receptor with the corresponding principal component scores across the 18 compounds, we produced radar plots quantifying the relevance of each receptor/transporter/Ca+ channel for each of the five principal components.

Following the methodology introduced by Vidal and Mestres (2010) suggesting that an unknown interaction of a ligand for a particular target can be approximated by its molecular resemblance to a selection of assayed active compounds, we implemented a similar “in silico” approach for the virtual screening of binding affinity profiles (Vidal and Mestres, 2010). The main objective of this analysis was to obtain an independent set of binding affinity values to compare with those reported in Ray (2010) and evaluate the plausibility given the relatively low number of replicates.

Using the open-source data mining and integration platform KNIME¹⁰ (Berthold et al., 2009), we extracted the molecular structures with their respective assayed binding affinity values from two public access chemical databases: ChEMBL¹¹ (Gaulton et al., 2012) curated by the European Bioinformatics Institute (EBI) of the European Molecular Biology Laboratory (EMBL), and the KiDatabase¹² (Evans et al., 2001) from the PDSP of the National Institute of Mental Health (NIMH). The molecular descriptor summarizing the structure of each compound as a binary (“bit”) vector was calculated using the Fingerprint node of the open-source cheminformatics extension RDKit for KNIME¹³. The similarity between the fingerprints of a pair of molecules was obtained using the Tanimoto distance coefficient (described in a previous section). For this purpose, the Fingerprint Similarity node from the KNIME-CDK Integration was used (Beisken et al., 2013). With the objective of minimizing the noise due to the contribution of compounds without a sufficiently high structural correlation with the drugs included in this study, we implemented the filtration criteria proposed by Paulke et al. (2013) consisting of a cut-off similarity value equivalent to the 95% quantile of the overall similarity distribution. Finally, the K_i values were predicted by inverse distance weighting, assigning different contribution weights to each molecule in the databases that are inversely proportional to the structural difference between the molecules:

where Kij(m) represents the binding affinity value of the molecule m (i.e., one of those presented in Figure 2) at the j-th receptor, u_i represents the binding affinity of the i-th ligand at the databases included in the analysis, and w_i(m) is a weighting term for the affinity value based on an exponent of the inverse distance between the molecule m and the x_i ligand. The distance was computed as one minus the Tanimoto coefficient based on different molecular descriptor algorithms. The best power parameter value “p” and the most suitable molecular descriptor algorithm were obtained after running an iterative optimization procedure benchmarked against the known binding affinity profile of 13 non-selective psychiatric drugs assayed by Roth et al. (2004). Using the optimal value p = 4.0 and the “Layered” algorithm gave a 39% of hits against 34 receptors, an acceptable result compared to the reference value of 59.3%–66.7% achieved by Vidal and Mestres (2010) using the same procedure applied to a larger database of compounds.

We first computed the similarities between the compounds in terms of the reported subjective effects, binding affinity profiles and molecular structures based on the data provided in Ray (2010). These are presented in matrix and graph forms in Figure 4A (graphs were visualized using Gephi¹⁴ with the “force atlas” layout; Bastian et al., 2009). It is clear from the matrix and graph representations of the molecular similarity that the Tanimoto coefficient differentiates phenethylamines from tryptamines. 2C-x and DOx compounds presented high within-group structural similarity, as well as MDA and MDMA. As expected, tryptamines and their methoxylated analogs (e.g., DMT and 5-MeO-DMT) presented a high degree of structural similarity. Ibogaine and LSD presented a moderate similarity to other tryptamines (likely a consequence of their complex multi-cyclic structures) and a lower similarity to phenethylamines. Such clear differentiation between tryptamines and phenethylamines was not as marked for the similarity of binding affinity profiles and the reported subjective effects. Compounds of the 2C-x family, tryptamines and their methoxylated analogs presented the highest similarities in terms of binding affinity profiles, while ibogaine presented a generally low similarity to the other compounds (however, this could be influenced by the relatively high amount of missing data for ibogaine—see the “Materials and Methods” section).

As shown in Figure 4B, we observed positive and significant correlations between the binding affinity profile similarity and the reported subjective effects similarity (R = 0.60, p = 2.2e–16), between the binding affinity profile similarity and the molecular structure similarity (R = 0.46, p = 2.5e–9), and between the reported subjective effects similarity and molecular structure similarity (R = 0.3, p = 1.9e–4). To test whether these correlations were driven by two clusters presenting high/low within/between group similarities, we employed different symbols for phenethylamine-phenethylamine, tryptamine-tryptamine and phenethylamine-tryptamine pairs in the scatter plots of Figure 4B. These clusters are not apparent in the visualization of the data. We observed that correlations restricted to phenethylamine-tryptamine pairs were in all cases significant (p < 1e–4).

We tested how the correlation values and their significance depended on the number of components retained in the SVD step of LSA (see “Materials and Methods” section). As shown in Figure 5A, both sets of correlation coefficients (reported subjective effects vs. molecular structure/binding affinity profiles) remained positive, but declined as a function of the number of components retained. We also employed a previously validated computational framework for the “in silico” estimation of binding affinity profiles (“virtual receptorome screening”). This approach consisted of querying large databases for ligands that have been experimentally assayed at certain sites and that are also structurally close to the molecules in Figure 2. The binding affinities were then estimated by the average of the experimental data of all those ligands weighted by the structural similarity with each of the psychedelic molecules included in the present study (see the “Virtual Screening of Binding Affinity Profiles” section of the “Materials and Methods” for more information). This method was implemented to provide an independent evaluation of the feasibility of the PDSP data provided in Ray (2010); as well as of the data compiled from different sources. The significant positive correlation (R = 0.42, p = 7.3e–8) between the empirical and estimated binding affinity profiles supports such feasibility (see Figure 5B).

To provide independent support for the results in Figure 4, we investigated the correlation between the semantic similarity of the reported subjective effects and the similarity of the binding affinity profiles using the K_i values published in Rickli et al. (2015, 2016). The results of this analysis are presented in Figure 6. Figure 6A shows the semantic (retaining 20 dimensions in the LSA algorithm) and binding affinity similarity matrices. We observed high binding affinity profile similarities for 2C-x compounds, consistent with the results obtained from the PDSP data (Figure 4). 25I-NBOMe presented a relatively lower similarity to other phenethylamines. Also consistent with Figure 4, tryptamines presented lower within-group similarity in terms of binding affinities and of the reported subjective effects. Conversely, phenethylamines presented higher within-group similarity in both senses. Figure 6B shows the correlation coefficient between the semantic and binding affinity profile similarities as a function of the number of dimensions. We observed R > 0.4 except for the lowest number of dimensions. The scatter plot in the inset provides an example of the correlation between both variables (R = 0.48, p = 1e–11).

The correlations reported in Figure 4B correspond to those obtained using the binding affinities measured at all 42 sites (Ray, 2010). It could be the case that higher correlations are achievable when considering only a subset of those sites. This would imply that the binding affinities at those sites are of a higher relative importance to explain the subjective effects of the compounds. We first evaluated whether the binding affinity profiles restricted to 5-HT receptor subtypes yielded a higher correlation than the reported in Figure 4B for all 42 sites, which did not occur (R = 0.31). This suggested the relevance of other neurotransmitters and neuromodulators for the prediction of the similarity of the reported subjective effects.

We explored the combinations of receptor types to find the one that maximized the correlations in Figure 4B. The search space was restricted to combinations including at least one of the following groups of receptors: 5-HT, dopamine (D), adrenergic (α/β) and muscarinic (M). This restriction ensured that sufficient data points were available for the computation of the correlation coefficients. Figure 7A shows the sorted linear correlation coefficients between the binding affinity profile similarity and reported subjective effects similarity (black), and between the binding affinity profile similarity and the molecular structure similarity (gray). The panels below display the combinations of receptors/transporters/Ca+ channel corresponding to each correlation value shown above. The optimal combination of sites (R = 0.63, p < 1e–17) to predict the reported subjective effects similarity from the binding affinity profile similarity was based on data from 5-HT, DA, opioid, muscarinic, NMDA receptors and the Ca+ channel (Figure 7B, left). The optimal prediction (R = 0.52, p < 1e–11) of the molecular structure similarity was based on the binding affinity profile restricted to 5-HT, imidazoline, σ and muscarinic receptors (Figure 7B, right). In Figure 7C we show the matrices containing the pairwise binding affinity profile similarity, restricted to the set of receptors obtained by the optimization procedures described above. The resemblance with the reported subjective effects and molecular structure similarity matrices shown in Figure 4A is apparent upon direct visual inspection.

To assess the robustness of the results presented in Figure 7, we conducted a bootstrap analysis with 1,000 iterations, applying in each iteration the LSA algorithm to 50% of the words in the Erowid reports (randomly chosen, with replacement). Each iteration yielded a correlation coefficient between the binding affinity profile similarity and the reported subjective effects similarity, for each combination of receptors/transporters/Ca+ channel shown in Figure 7A. The average correlation coefficients across all 1,000 bootstrap iterations are shown in Figure 8A, using the same ordering of the permutations as in Figure 7A. Even though small fluctuations occurred, a similar functional dependence was observed, implying that the combinations explored in Figure 7A produced very similar correlation coefficients across all bootstrap iterations. Figure 8B shows the optimal selection of receptors/transporters/Ca+ channel for all bootstrap iterations in matrix form. It is clear that certain receptors were frequently included (e.g., 5-HT receptors) while others were systematically excluded. Figure 8C shows the percentage of times (over all 1,000 iterations) that a receptor type/transporter protein/Ca+ channel appeared in the set that maximized the correlation between the binding affinity profile similarity and the reported subjective effects similarity. 5-HT, DA, opioid, muscarinic, cannabinoid and NMDA receptors and the Ca+ channel were selected for all or almost all bootstrap iterations.

The next step of our analysis was to identify the facets of the reported subjective experiences linked to the degree of binding affinity at the 42 sites published in Ray (2010). For this purpose, we applied PCA to the word-by-document rank-reduced matrix obtained after applying LSA to the tf-idf matrix (see “Materials and Methods” section). This analysis included all the reports of psychoactive substances in the Erowid corpus as described in a previous publication (Sanz et al., 2018). We retained the first five principal components, which explained ≈58% of the variance of the data. The coefficients associated to each word indicated its relevance for each component and were used as the weights to construct the word clouds presented in Figure 9A. The five principal components are summarized below:

The names of the principal components were chosen to indicate the nature of their most relevant terms. The “perception” component presented high weights for terms associated with conscious perception of the self (e.g., “ego,” “fear”) and the environment (e.g., “visual,” “color”). While the “body load” component also included such terms, others associated with stress and negative bodily sensations also presented high weights (e.g., “nausea,” “tension,” “headache,” “jaw,” “stimulate”). The “preparation” component was clearly associated with the processing and preparation of psychoactive materials. The remaining two components (“dependence” and “therapeutic”) were labeled following similar considerations. Words relevant for the “dependence” component included terms such as “addict,” “addiction,” “withdrawal,” and “tolerance,” as well as others that could be related to the adverse reactions observed during the withdrawal syndrome of some addictive drugs (“nausea,” “vomit,” “sick,” “puke”). With respect to the last component (“therapeutic”), many terms of high relevance were directly associated with medical disorders (“depression,” “insomnia,” “disorder,” “anxiety,” “treatment,” “prescribe,” “med,” “symptom”). Words in common with the “dependence” component could represent adverse secondary reactions of medicaments. The inclusion of terms associated with psychedelic or visionary experiences (“universe,” “visuals,” “vision,” “ego”) might indicate the therapeutic use of some of the compounds included in our analysis.

Figure 9B shows the score of each principal component for the 18 psychedelic compounds included from Ray (2010). As expected, the “perception” component presented the highest score for all compounds (except ibogaine), with the highest scores corresponding to 2C-E, DPT and DMT. In terms of the “body load” component, the highest scores were observed for 5-MeO-MIPT, TMA-2 and 2C-T-2; in terms of the “preparation” component were DMT, DiPT and 5-MeO-DMT; in terms of the “dependence” component were MDMA, ibogaine and 5-MeO-DMT and in terms of the “therapeutic” component were the same three as in the “dependence” component. Figure 9C presents radar plots showing the principal component scores of three compounds included in this study (LSD, TMA-2 and ibogaine), while Figure 9D presents an analogous plot for three psychoactive compounds included in the Erowid database but not analyzed in this study (cocaine, venlafaxine and plants of the Datura genus; Sanz et al., 2018). These substances presented scores agreeing with their known subjective effects, “set” and “setting”; e.g., cocaine scored high in the “dependence” component, venlafaxine in the “therapeutic” component, and plants of the Datura genus in the “perception” and “preparation” components.

Knowing the principal component scores for each psychedelic compound, as well as their receptor binding affinities, it was possible to obtain their correlation and thus to determine whether a higher binding affinity for a certain receptor predicted higher principal component scores and vice-versa. The results of this analysis are presented in Figure 10A. The most salient results were the positive correlations between binding affinities at some 5-HT receptor subtypes and the score of the “perception” component. We observed generally negative correlations for the “body load” component, and positive correlations between affinity at transporters (SERT, DAT), opioid (MOR, KOR) and NMDA receptors and the “therapeutic” component scores, while Ca+ affinity correlated positively with the “dependence” component scores. Figure 10B shows a radar plot for the correlation at sites presenting different profiles (5-HT_2A, MOR, Ca+).