Subscriber access provided by University of Winnipeg Library
Chemical Information
Three-Dimensional Activity Landscape Models of Different Design and Their Application to Compound Mapping and Potency Prediction Tomoyuki Miyao, Kimito Funatsu, and Jürgen Bajorath J. Chem. Inf. Model., Just Accepted Manuscript • DOI: 10.1021/acs.jcim.8b00661 • Publication Date (Web): 28 Nov 2018 Downloaded from http://pubs.acs.org on December 4, 2018
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Three-Dimensional Activity Landscape Models of Different Design and Their Application to Compound Mapping and Potency Prediction
Tomoyuki Miyao,1 Kimito Funatsu,1,2 and Jürgen Bajorath3*
1Data
Science Center and Graduate School of Science and Technology, Nara Institute of Science and Technology, 8916-5 Takayama-cho, Ikoma, Nara, 630-0192, Japan.
2Department
of Chemical System Engineering, School of Engineering, The University of Tokyo, 73-1 Hongo. Bunkyo-ku, Tokyo 113-8656, Japan.
3Department
of Life Science Informatics, B-IT, LIMES Program Unit Chemical Biology and
Medicinal Chemistry, Endenicher Allee 19c, Rheinische Friedrich-Wilhelms-Universität, D-53115 Bonn, Germany.
*To whom correspondence should be addressed: Tel: +49-228-2699-306, Fax: +49-228-2699-341, E-mail:
[email protected] ACS Paragon Plus Environment
1
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 2 of 49
Abstract
Activity landscapes (ALs) integrate structural and potency data of active compounds and provide graphical access to structure-activity relationships (SARs) contained in compound data sets. Threedimensional (3D) ALs can be conceptualized as a two-dimensional (2D) projection of chemical space with an interpolated activity surface added as a third dimension. Such 3D ALs are particularly intuitive for SAR visualization. In this work, 3D ALs were generated on the basis of different projection methods and fingerprint descriptors and their topologies were compared. Moreover, going beyond qualitative analysis, the use of 3D ALs for semi-quantitative and quantitative potency predictions was investigated. Neuroscale, a neural network variant of multi-dimensional scaling (MDS), combined with Gaussian process regression (GPR) was identified as a preferred approach for generating 3D ALs that accounted with high accuracy for training compounds and their SAR characteristics. On the other hand, GPR-induced overfitting generally limited the accuracy of potency value predictions, regardless of the projection methods that were applied. However, 3D ALs enabled reliable mapping of test compounds with varying potency levels to corresponding AL regions. Most accurate mapping was achieved with neuroscale models. Taken together, the results of our analysis indicate high potential of 3D ALs for graphical SAR exploration and the identification of potent test compounds.
ACS Paragon Plus Environment
2
Page 3 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Introduction
In chemoinformatics and medicinal chemistry, the concept of activity landscapes (ALs) is employed to study the distribution of active compounds in chemical space, elucidate structure-activity relationships (SARs), and identify SAR determinants and key compounds.1-8 An important aspect of the AL concept is that it provides integrated graphical access to compounds and activity data and makes it possible to visualize SARs.3,5,7 For example, large contributions to the topology of ALs and underlying SARs are made by structurally similar compounds with large potency variations that form activity cliffs,9-14 which are prominent features of AL representations.2,3 An AL can principally be envisioned as a section of biologically relevant chemical space to which an activity hypersurface is added that accounts for differences in compound potency. In such extended chemical space representations, structurally similar compounds are close to each other, increasing distance correlates with dissimilarity, and –through the activity hypersurface– compound positions are further differentiated according to varying potencies. In practice, models of chemical space are generated using descriptors to capture compound distance relationships as a measure of similarity/dissimilarity and then transformed into representations that are suitable for AL design. An AL has formally been defined as a graphical representation that integrates structural similarity and potency difference relationships between compounds having the same biological activity.2 Accordingly, a variety of AL representations of different design and complexity have been
ACS Paragon Plus Environment
3
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 4 of 49
introduced.2,4,8 These include simple 2D plots such as structure-activity similarity (SAS) maps15-17 that capture SARs in compound data sets on the basis of pairwise compound comparisons or cellbased SAR maps,18 molecular network and graph representations,19-25 and 2D ALs generated using machine learning approaches for projecting multi-dimensional feature spaces.26-28 In addition, 3D ALs have been constructed.29-31 Some of these 3D AL designs are particularly intuitive for SAR analysis.4,5 In principle, 3D ALs can be rationalized as a 2D projection of chemical feature space with an activity surface added as a third dimension. Such ALs are reminiscent of geographical maps because they consist of regions that can be interpreted as mountains, planes, or valleys. These topological features correlate with SAR characteristics. For example, in smooth and gently sloped regions, gradual changes in compound structure are accompanied by small to moderate changes in potency, which is indicative of SAR continuity.32,33 By contrast, in rugged regions, small chemical changes cause large potency variations, which corresponds to SAR discontinuity.32,33 These regions often contain mountains and steep activity cliffs, which represent the pinnacle of SAR discontinuity. Such 3D ALs can be generated for compound data sets via a two-step procedure.29 First, dimension reduction methods are applied to transform multi-dimensional chemical feature spaces into 2D representations (i.e., projections onto an x-y plane). Second, an activity surface is generated from compound potency values (z-axis). Because potency values are spatially separated, interpolation functions must be applied to generate coherent activity surfaces.29,30 Alternatively, 3D AL views can
ACS Paragon Plus Environment
4
Page 5 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
be directly extracted from high-dimensional feature spaces and annotated with potency information.31 3D ALs of geographical resemblance visualize SARs in an intuitive manner and are readily interpretable.3,29,31 However, such models are mostly descriptive in nature, rather than predictive. Although 2D ALs obtained by generative topographic mapping (GTM)27,28 can be used to predict compounds that fall into regions of activity in a map, the use of 2D and 3D ALs introduced thus far is typically confined to visualization and qualitative analysis, regardless of their design principles and level of complexity.2-4 In this work, we have further investigated 3D ALs of different design and focused on their predictive ability. Hence, compared to earlier studies, we have gone further and evaluated the potential of 3D ALs for semi-quantitative and quantitative predictions of compound activity. 3D ALs were generated using different methods and their topologies were compared. Compounds distance relationships in original and latent space were analyzed and test compounds with varying potency levels were mapped onto 3D ALs. Furthermore, compound potency values were predicted. Our analysis is presented in the following.
ACS Paragon Plus Environment
5
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 49
Materials and Methods
Compound data sets Sets of compounds with high-confidence activity data were taken from a previous study analyzing activity cliffs.34 Ten target sets were selected that originated from ChEMBL (version 23)35 and displayed large intra-set potency variations, giving rise to the frequent formation of activity cliffs.34 Compound potency was numerically specified in the form of equilibrium constants (Ki values; reported as pKi values). The composition of these target sets is reported in Table 1.
Training and test data For each target set, 60% of the compounds were randomly selected and used as a training set and the remaining 40% of the compounds were used as a test set. In test sets, the 10% most potent compounds were classified as highly potent and the 10% least potent compounds as weakly potent and used as subsets for mapping. The composition of training and test sets is summarized in Table 2.
ACS Paragon Plus Environment
6
Page 7 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Table 1. Target sets.
CHEMBL ID
Activity class
237
Kappa opioid receptor ligands
1879
244
Coagulation factor X inhibitors
245 4860 1862
# CPDs Max
Potency (pKi)
SARI Score
Min
Mean
Std
MACCS ECFP4
11.5
4.1
7.4
1.3
1.51
0.95
1608
11.4
3.6
7.5
1.8
1.68
1.03
Muscarinic acetylcholine receptor M3 ligands
652
10.5
4.1
7.8
1.5
1.49
1.00
Apoptosis regulator Bcl-2 inhibitors
630
11.3
3.3
8.1
2.0
1.87
1.15
Tyrosine-protein kinase ABL 553 inhibitors
10.7
4.2
8.1
1.4
1.26
0.90
1983
Serotonin 1d (5-HT1d) receptor ligands
478
10.7
4.4
7.9
1.4
1.58
1.00
2954
Cathepsin S inhibitors
435
9.9
3.5
6.9
1.4
1.47
0.94
3798
Calcitonin gene-related peptide type 1 receptor ligands
416
11.3
4.7
7.9
1.5
1.74
1.06
220
Acetylcholinesterase inhibitors
380
11.8
1.0
6.7
1.9
1.54
1.12
249
Neurokinin 1 receptor ligands
350
12.0
3.8
8.2
1.6
1.75
1.18
For each target set, the CHEMBL ID, number of compounds (# CPDs), and the largest and smallest potency (pKi) value are reported. Mean and standard deviations (Std) are also given. In addition, SARI discontinuity scores54 were calculated applying Tc similarity thresholds of 0.65 for MACCS and 0.43 for ECFP4.55
ACS Paragon Plus Environment
7
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 8 of 49
Table 2. Training and test sets.
Potency threshold CHEMB L ID
# Test CPDs
# Training CPDs
Similarity between test CPDs and NN training CPDs
Highly potent
Weakly potent
Highly potent
Weakly potent
Mean
Std
237
1127
752
76
76
9.31
5.77
0.78
0.13
244
964
644
65
65
9.61
5.13
0.79
0.11
245
391
261
27
27
9.74
5.59
0.73
0.19
4860
378
252
26
26
10.30
4.94
0.83
0.12
1862
331
222
24
23
9.70
6.13
0.76
0.11
1983
286
192
21
20
9.75
5.96
0.72
0.17
2954
261
174
18
18
8.67
4.71
0.78
0.10
3798
249
167
17
17
10.15
6.07
0.81
0.08
220
228
152
16
16
8.76
3.93
0.70
0.18
249
210
140
14
14
9.94
5.83
0.76
0.15
For each target set, the composition of training and test sets is reported. The 10% most potent and 10% least potent test compounds are selected as subsets for mapping and the corresponding potency (pKi) threshold values are given. Tanimoto similarity statistics for each test compound and the corresponding nearest neighbor in the training data set are reported. ECFP4 was used as a molecular representation.
ACS Paragon Plus Environment
8
Page 9 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Molecular representation As molecular representations, two standard fingerprints of different design and complexity were used. The extended connectivity fingerprints of bond diameter 4 (ECFP4),36 a topological feature set fingerprint, was folded into a 1024-bit vector by modulo operation. MACCS structural keys37 consisted of 166-bit vector representing structural fragments and patterns. Both fingerprints were generated using in-house Python scripts with the aid of the OEChem toolkit.38
Construction of 3D activity landscape models 3D ALs were designed to consist of a 2D manifold (x-y plane) capturing structural relationships between compounds and a third dimension (z-axis) accounting for observed or predicted potency values. The 2D manifold represented a 2D projection of chemical space. For our analysis, a coordinate-free chemical reference space was generated on the basis of pairwise fingerprint Tanimoto distances between training compounds. In 2D projections, compounds were assigned coordinates. From compound coordinates and sparsely distributed potency values, coherent activity surfaces were interpolated. For each target set, 3D ALs were generated for training compounds. Test compounds were used for mapping and potency predictions.
ACS Paragon Plus Environment
9
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 10 of 49
Two-dimensional projections Different dimension reduction and machine learning methods were applied to generate 2D projections of fingerprint spaces. Multi-dimensional scaling Multi-dimensional scaling (MDS)39 aims to preserve distance relationships in original and latent feature spaces. When first generating 3D ALs for compound data sets29 MDS was applied for dimension reduction of feature spaces and differences between Euclidian distances in original space and 2D projections were minimized. Here, MDS was applied by minimizing an objective function accounting for differences between Tanimoto distances in fingerprint space and Euclidian distances in the x-y projection. Principal component analysis As a conceptually different dimension reduction approach, principal component analysis (PCA)40 was also used, which generates linear combinations of original descriptors as orthogonal axes of lower dimensional space. The proportion of the variance of original data that was accounted for by chosen principal components was determined. Generative topographic mapping Generative topographic mapping (GTM)41 enables non-linear mapping of a set of molecules from high-dimensional space into low-dimensional latent space. In 2D latent space, evenly distributed lattice points are placed and the points are re-projected into the high-dimensional space via radial
ACS Paragon Plus Environment
10
Page 11 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
basis function (RBF) networks. Projected points form centers of Gaussian distributions and their coordinates in high-dimensional space are derived by maximizing a logarithmic likelihood function. Since GTM makes use of probability density functions such as Gaussians, a fingerprint vector x of a molecule in the high-dimensional space can be connected to coordinates (t) in latent space via probabilistic operations. Projected coordinates of a compound represent expected coordinates in latent space weighted by the posterior density of t given x:
, where M is the
number of lattice points. Although the expectation-maximization (EM) algorithm can efficiently identify local optima of the logarithmic likelihood function (i.e., a proper network weight matrix and the Gaussian variance) the number of lattice points (termed map size), the number of RBFs, the width of RBFs, and the regularization term (determining the smoothness of the mapping function in the RBF network) should be determined in advance. Here, map size was set to 225 (15×15) and the remaining hyper-parameters were optimized by two-fold cross validation based on the average distance between original and re-projected coordinates of compounds in high-dimensional space. Latent trait model The latent trait model (LTM)42 represents a specific embodiment of GTM, which facilitates the generation of different noise models for representing the probability distribution of a compound set in high-dimensional space. For representing the probability density of binary data (such as fingerprints), LTM utilizes an independent Bernoulli noise model.42 Accordingly, density of fingerprint vectors is represented by a mixture of independent Bernoulli distributions. Each lattice
ACS Paragon Plus Environment
11
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 12 of 49
point in low- dimensional latent space corresponds to the input of an inverse link function for Bernoulli distributions. For Gaussians, this is identity, which means that the center of a Gaussian in high-dimensional space corresponds to a lattice point in latent space. Although LTM can directly assess distributions in fingerprint space (different from GTM), each iteration of the EM algorithm requires non-linear optimization, which increases training costs. For LTM, map size was also set to 225 and the number of RBFs and their width were determined by two-fold cross validation based on the average distance between the original and re-projected coordinates of compounds in highdimensional space, corresponding to GTM. Neuroscale Neuroscale43 represents an RBF neural network (NN) that starts from coordinates in highdimensional space to generate coordinates in latent space. Neuroscale can be understood as
an NN
variant of MDS since it uses a stress function as the loss function during training. An efficient algorithm44 enables fast convergence of stress function minimization. Here, coordinates of RBF centers were initially determined by optimizing a mixture of Gaussians in fingerprint space. The number of RBFs is the only hyper-parameter of neuroscale, which was determined by seven-fold cross validation on the basis of the stress function.
ACS Paragon Plus Environment
12
Page 13 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
T-distributed stochastic neighbor embedding So-called t-distributed stochastic embedding (t-SNE) is a non-linear manifold mapping method that emphasizes local structure.45 It is based upon the Kullback-Leibler (KL) divergence formalism.46 By minimizing the KL divergence, similar compounds should be projected to similar coordinates. We note that lower-dimensional latent spaces have recently also been derived using autoencoder deep neural networks, albeit not for visualization purposes.47
Surface interpolation From 2D projections, coherent activity surfaces were interpolated using the so-called Kriging or Gaussian process regression (GPR) approach48,49 originating from geostatistics, which interpolates intermediate values by a Gaussian process on the basis of prior covariances of known data points. An exponential kernel function was used assuming a mean of zero such that the covariance kernel established relationships between observed potency values. GPR can be carried out in the absence or presence of Gaussian noise factors. Without introducing a noise factor, an interpolated surface exactly passes through known potency values when deriving the activity surface. However, applying noise factors of increasing noise level tolerance, varying z-values are generated for a point on the x-y plane and the best fit of surfaces to training data can be determined.
ACS Paragon Plus Environment
13
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 14 of 49
Compound mapping Test compounds can be mapped onto a given 3D AL if 2D projection methods possess mapping functions: t = f (x), where x denotes a vector in the original fingerprint space and t the corresponding coordinates in the x-y plane. Mapping functions can be derived using PCA, GTM, LTM, and neuroscale. By contrast, conventional MDS and t-SNE are not capable of deriving such mapping functions.
Potency prediction Since GPR uses coordinates in the x-y plane to predict surface (potency) values, such predictions can be carried out on the basis of projection methods that include or generate mapping functions. For potency predictions, Gaussian noise models were introduced and the noise level tolerance was controlled to accurately reproduce potency values of training compounds. Therefore, as a parameter for model optimization, a range of noise values was used for the diagonal kernel. Optimization was based on the logarithmic likelihood function for training data and constrained by given noise values. As reference calculations for potency predictions, k-nearest neighbor (k-NN) regression models were also generated. Potency values were predicted as the weighted average of the k-NN training compounds on the basis of 1-Tc distances in fingerprint space.
Evaluation metrics
ACS Paragon Plus Environment
14
Page 15 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Pearson correlation coefficient (r) between pairwise Tanimoto distances in fingerprint space and Euclidian distances in 2D space was calculated as a metric for evaluating projection methods. In addition, the explanatory and predictive capacity of GPR models was assessed by calculating the coefficient of determination (R2) and mean absolute error (MAE) for both training and test data.
Libraries and implementation For PCA, MDS, and t-SNE, Scikit-learn version 0.19.1 was used 50 and GPR was carried out with GPy.51 As t-SNE parameters, default settings were used perplexity: 30 and the number of steps: 1000. As a stress function for MDS, the sum (original distance – mapped distance)2 was calculated and minimized. Furthermore, GTM, LTM, and neuroscales were implemented in Python based on available MATLAB libraries.52 Compound data sets and scripts are available upon request.
ACS Paragon Plus Environment
15
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 16 of 49
Results and Discussion
Study design Given their intuitive nature, “geographical” 3D ALs are attractive for visualizing SARs. So far, 3D ALs have mostly been used for qualitative graphical analysis and comparison of data sets. For qualitative purposes, approximate 2D projections and 3D surface interpolations are usually sufficient. In this work, we have compared 3D ALs of different design and investigated the use of these models for semi-quantitative applications (mapping of test compounds with varying potency) and quantitative applications (potency value predictions). Such predictions require accurate projections and interpolation. By design, 3D ALs represent non-linear quantitative SAR (QSAR) models because they combine SAR continuity and SAR discontinuity including activity cliffs within the same reference frame. ALs are set apart from other QSAR models because they provide graphical SAR views that also enable predictions. It is emphasized that SAR continuity is based upon the similarity-property principle and chemical neighborhood behavior and represents the paradigm upon which classical QSAR analysis is based.2,4,53 By contrast, SAR discontinuity falls outside the applicability domain of classical QSAR.9,10 Accordingly, we have explored to what extent 3D ALs, which combine these different SAR features, might be useful for QSAR-like predictions. For AL analysis, compound data sets were selected that had been used for SAR discontinuity and activity cliff analyses. As a measure of SAR discontinuity in fingerprint spaces, SAR Index (SARI)
ACS Paragon Plus Environment
16
Page 17 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
discontinuity score components53 were calculated, as reported in Table 1, confirming the presence of substantial SAR discontinuity in the target sets. Such data sets generally represent challenging test cases for potency predictions. A variety of methods for generating 2D projections of chemical space have been evaluated for AL design. It would also be possible to apply PCA regression or partial least squares techniques. The resulting models were tested by quantifying distance relationships in original and latent space and by 3D compound mapping. Furthermore, potency values of test compounds were predicted and Gaussian process models with different noise level tolerance were generated to analyze the influence of different tolerance levels on potency predictions.
Landscape topology Figure 1 shows exemplary 3D ALs for the largest target set generated using different projection methods. In each case, the same target set, molecular representation, and interpolation protocol was used. All models showed coexistence of continuous (flat green/yellow) and discontinuous (rugged green/yellow and red) regions containing many activity cliffs (red peaks mark positions of potent cliff compounds). However, the comparison also revealed topological variations among 3D ALs, dependent on the methods that were used. The conceptually related GTM/LTM approaches on the one hand and MDS/neuroscale on the other produced pairs of similar landscapes that, however, differed from each other. Notably, GTM and LTM produced 3D ALs that were overall smoother than
ACS Paragon Plus Environment
17
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 18 of 49
the MDS and neuroscale models, which exhibited steeper cliffs. For MDS, the root mean square deviation (RMSD) of distances between two spaces ranged from 0.20 (Bcl-2 inhibitors) to 0.31 (coagulation factor X inhibitors) for ECFP4 and from 0.07 (Bcl-2 inhibitors) to 0.17 (acetylcholinesterase inhibitors) for MACCS. Furthermore, the PCA-based 3D AL also showed pronounced ruggedness in some regions, but rather low resolution (fuzziness) of landscape features in others. Variance proportions captured by the first two PCs range from 0.15 (coagulation factor X inhibitors) to 0.34 (tyrosine-protein kinase ABL inhibitors) for ECFP4, and from 0.26 (kappa opioid receptor ligands) to 0.43 (tyrosine-protein kinase ABL inhibitors) for MACCS. Moreover, the t-SNE 3D AL was characterized by a clearly evident loss in information compared to the GTM/LTM and MDS/neuroscale models, with a predominance of rather unstructured regions and extrapolated intermediate potency. Accordingly, the t-SNE landscape view resulted from close co-projection of compounds with varying potency in the y-x plane, giving rise to an averaging effect in interpolation. Molecular representations also affected the topologies of ALs. In Figure 2, 3D ALs derived on the basis of MACCS and ECFP4 are compared for kappa opioid receptor ligands and neurokinin 1 receptor ligands. The ALs displayed notable topology variations, especially in regions containing highly potent compounds. The resulting 3D AL did not account for the potency distribution within the target set. Hence, graphical analysis of alternative 3D ALs revealed interesting commonalities and differences between
ACS Paragon Plus Environment
18
Page 19 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
models generated with related or unrelated projection methods. These differences were further explored in quantitative terms.
ACS Paragon Plus Environment
19
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 20 of 49
Figure 1. Alternative activity landscapes.
For kappa opioid receptor ligands (CHEMBL ID 237, Table 1), ALs were generated using different methods for generating 2D projections on the basis of MACCS fingerprint distances. (a) GTM, (b) LTM, (c) MDS, (d) Neuroscale, (e) PCA, (f) t-SNE.
ACS Paragon Plus Environment
20
Page 21 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Figure 2. Activity landscapes using different molecular representations.
For two compound data sets, a side-by-side comparison of MACCS- and ECFP4-based 3D ALs is shown that were generated using neuroscale.
ACS Paragon Plus Environment
21
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 22 of 49
Distance correlation for training compounds First, the correlation of distances between training compounds in original descriptor space and 2D projections was determined. The results are reported in Table 3 for ECFP4 and MACCS. For both representations, highest correlation across target sets was achieved by MDS, followed by neuroscale. For MACCS distances, MDS yielded highest correlation for all the targets sets and for ECFP4, for eight out of 10 sets. For MACCS, consistently high correlation coefficients were obtained, ranging from 0.78 to 0.96; for ECFP4, some lower values were observed, which was likely due to the higher dimensionality of this fingerprint. Overall strongest correlation of distances between training compounds obtained for MDS and neuroscale was consistent with the fact that both methods utilized a stress function as an objective function to preserve distances in original and latent space. Lower correlation for neuroscale compared to MDS might be attributable to an intrinsic restriction imposed by the RBF network and the need to optimize a hyper-parameter. As pointed out by a reviewer, we also note that PCA on binary fingerprints is equivalent to MDS using Hamming distance, rather than 1-Tc. Hamming distance assigns equal weights to 0 and 1 bit settings, different from fingerprint similarity calculations.
ACS Paragon Plus Environment
22
Page 23 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Table 3. Distance correlation for training compounds. CHEMBL ID
GTM
LTM
MDS
Neuroscale
PCA
t-SNE
237
0.51 (0.48)
0.45 (0.50)
0.59 (0.81)
0.55 (0.71)
0.40 (0.61)
0.53 (0.68)
244
0.50 (0.59)
0.51 (0.57)
0.57 (0.78)
0.57 (0.71)
0.48 (0.66)
0.54 (0.47)
245
0.37 (0.52)
0.34 (0.50)
0.59 (0.80)
0.55 (0.73)
0.30 (0.55)
0.49 (0.56)
4860
0.53 (0.52)
0.65 (0.71)
0.83 (0.96)
0.82 (0.93)
0.75 (0.89)
0.53 (0.62)
1862
0.60 (0.69)
0.75 (0.79)
0.84 (0.90)
0.82 (0.88)
0.83 (0.86)
0.79 (0.78)
1983
0.50 (0.54)
0.47 (0.63)
0.64 (0.85)
0.65 (0.82)
0.38 (0.77)
0.58 (0.61)
2954
0.62 (0.62)
0.53 (0.57)
0.66 (0.78)
0.63 (0.76)
0.49 (0.63)
0.61 (0.68)
3798
0.76 (0.56)
0.72 (0.59)
0.73 (0.82)
0.70 (0.73)
0.74 (0.61)
0.64 (0.61)
220
0.35 (0.6)
0.42 (0.55)
0.61 (0.82)
0.59 (0.82)
0.40 (0.64)
0.54 (0.69)
249
0.44 (0.55)
0.33 (0.44)
0.62 (0.83)
0.57 (0.76)
0.23 (0.68)
0.56 (0.56)
For training compounds, Pearson correlation coefficients are reported for Tanimoto distances in ECFP4 (MACCS) fingerprint space and Euclidean distances in 2D projections using different methods. For each target set and representation, the largest value is given in bold.
ACS Paragon Plus Environment
23
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 24 of 49
Overall lowest and in part poor correlation, especially for ECFP4, was observed for the PCA and tSNE models. In the case of PCA, 2D projections only utilized the first two principal components, which typically leads to a significant loss in distance information, as indicated by the presence of region of low resolution in the PCA-based 3D AL. On the other hand, t-SNE is known for extensive clustering of similar compounds in high-dimensional space, resulting in short inter-compound distances, which translate into projections. This makes it difficult to distinguish between regions of varying potency through interpolation with the characteristics of the t-SNE model. The relative smoothness of the GTM- and LTM models compared to MDS and neuroscale observed in Figure 1 was further investigated. In GTM and LTM landscapes, multiple training compounds were found to be projected onto the same coordinates (grid points) in the x-y plane, which clearly rationalized the smoothing effect relative to the MDS and neuroscale models.
Mapping of test compounds Test compounds not involved in landscape generation could only be projected onto 3D ALs using methods for which mapping functions were available. This constituted an important difference between the related MDS and neuroscale approaches since mapping was only possible using neuroscale, but not MDS. Correlation of distances between test and training compounds in original space and in 2D projections generated with GTM, LTM, neuroscale, and PCA was analyzed. The results are reported in Table 4 for ECFP4 and MACCS. In this case, neuroscale produced strongest
ACS Paragon Plus Environment
24
Page 25 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
correlation across target sets, consistently on the basis of MACCS and in eight of 10 cases on the basis of ECFP4. For test compounds, the RMSD of distances between two spaces using neuroscale range from 0.24 (Bcl-2 inhibitors) to 0.32 (kappa opioid receptor ligands) for ECFP4 and from 0.09 (calcitonin gene-related peptide type 1 receptor ligands) to 0.18 (acetylcholinesterase inhibitors) for MACCS. Figure 3 depicts the mapping of exemplary highly potent acetylcholinesterase inhibitors (test compounds) on 2D projections for a MACCS- and ECFP4-based landscape generated using neuroscale. In Figure 3a (MACCS), all compounds were correctly mapped to local regions of high potency and in Figure 3b (ECFP4), only one mismatch was observed (compound D).
ACS Paragon Plus Environment
25
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 26 of 49
Table 4. Distance correlation for test versus training compounds. CHEMBL ID
GTM
LTM
Neuroscale
PCA
237
0.49 (0.47)
0.43 (0.47)
0.53 (0.71)
0.38 (0.58)
244
0.48 (0.58)
0.5 (0.57)
0.55 (0.71)
0.47 (0.67)
245
0.34 (0.51)
0.31 (0.48)
0.52 (0.71)
0.29 (0.53)
4860
0.52 (0.50)
0.62 (0.71)
0.79 (0.92)
0.73 (0.87)
1862
0.59 (0.67)
0.74 (0.77)
0.81 (0.85)
0.82 (0.84)
1983
0.50 (0.53)
0.46 (0.60)
0.65 (0.79)
0.38 (0.76)
2954
0.61 (0.63)
0.53 (0.57)
0.62 (0.73)
0.47 (0.64)
3798
0.76 (0.55)
0.73 (0.56)
0.69 (0.71)
0.74 (0.62)
220
0.39 (0.62)
0.41 (0.58)
0.58 (0.81)
0.40 (0.66)
249
0.46 (0.55)
0.35 (0.42)
0.56 (0.74)
0.24 (0.68)
For each test compound, ECFP4 (MACCS) Tanimoto distances and 2D Euclidean distances to all training compounds were determined and the Pearson correlation coefficient was calculated for all distances per target set. For each target set and representation, the largest value is given in bold.
ACS Paragon Plus Environment
26
Page 27 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Figure 3. Activity landscapes and representative test compounds.
For acetylcholinesterase inhibitors (ID 220), two ALs were generated with neuroscale on the basis of (a) MACCS and (b) ECFP4. Corresponding 3D and 2D views are shown. Five exemplary test compounds
were
mapped
onto
ACS Paragon Plus Environment
the
landscapes.
27
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 28 of 49
In Figure 4, subsets of highly potent and weakly potent kappa opioid receptor ligands (test compounds) are mapped to 2D projections of all models. Comparison of the projections generated with GTM/LTM, neuroscale, and PCA illustrate the high local resolution of the neuroscale projection (displaying a fine grained structure) compared to GTM/LTM as well as the relative loss of information in the PCA projection compared to the others. However, despite these differences, for both highly and weakly potent test compounds, mapping to corresponding regions was generally accurate across the projections, with only few misplaced compounds. This is well illustrated by comparing the distributions of highly and weakly potent compounds in both Figure 4a (MACCS) and Figure 4b (ECFP4). Mapping of compounds with different potency levels to 3D ALs is a semiquantitative exercise that goes well beyond qualitative analysis and comparison of landscape topology. Thus, high mapping accuracy was considered an encouraging finding. Tables S1-S4 of the Supporting Information report RMSD and coefficient of variation (CV) values of distances between fingerprint and projected spaces for all training and test sets.
ACS Paragon Plus Environment
28
Page 29 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Figure 4a. Mapping of highly and weakly potent test compounds.
ACS Paragon Plus Environment
29
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 30 of 49
Figure 4b Mapping of highly and weakly potent test compounds.
For kappa opioid receptor ligands (ID 237), 2D AL views are shown that were generated using different methods on the basis of (a) MACCS and (b) ECFP4. Subsets of highly potent (left) and
ACS Paragon Plus Environment
30
Page 31 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
weakly potent (right) test compounds (according to Table 2) were mapped onto the landscapes (black dots).
Prediction of potency values In light of the findings discussed above, quantitative prediction of potency values was attempted. Quantitative potency prediction was worth exploring because 3D ALs can be rationalized as nonlinear QSAR models, as discussed above. Table 5 reports the results of potency predictions on the basis of MACCS and ECFP4, respectively. Different from distance correlations, which were dominated by MDS and neuroscale, overall best potency predictions were achieved using GTM and LTM models. However, prediction accuracy was generally limited, with best R2 values ranging from about 0.50 to 0.80 across target sets. For neuroscale and PCA, only one R2 value greater than 0.50 was obtained on the basis of MACCS and two (neuroscale) or three (PCA) such values on the basis of ECFP4. Nonetheless, limited prediction accuracy still resulted in relatively low MAE values. For example, for GTM and LTM, MAE values reported in Table 5 were maximally 0.80, i.e., well within one order of magnitude. These results were comparable with k-NN based predictions in fingerprint space. Hence, given the limited potency prediction accuracy, there was no advantage of AL over k-NN based predictions. In addition, for neuroscale and PCA, MAE values close to or slightly exceeding one order of magnitude were obtained. For the subset of highly potent compounds, LTM produced overall best predictions, with MAE values within one order of
ACS Paragon Plus Environment
31
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 32 of 49
magnitude for eight of 10 target sets. Mean errors in potency values within or around an order of magnitude make it possible to distinguish between highly and weakly potent compounds.
Table 5. Potency predictions.
CHEMB L ID 237 244 245 4860 1862 1983 2954 3798 220 249
GTM
LTM
Neuroscale
5-NN (original space)
PCA
R2
MAE
R2
MAE
R2
MAE
R2
MAE
R2
MAE
0.55
0.68
0.61
0.64
0.31
0.85
0.24
0.88
0.69
0.55
(0.50)
(0.72)
(0.48)
(0.73)
(0.04)
(0.98)
(0.09)
(0.96)
(0.60)
(0.64)
0.73
0.69
0.74
0.67
0.41
1.00
0.40
1.03
0.82
0.56
(0.69)
(0.73)
(0.66)
(0.76)
(0.41)
(1.01)
(0.32)
(1.04)
(0.77)
(0.64)
0.66
0.68
0.59
0.72
0.42
0.85
0.18
1.01
0.66
0.67
(0.56)
(0.77)
(0.57)
(0.74)
(0.34)
(0.94)
(0.30)
(0.95)
(0.64)
(0.67)
0.79
0.70
0.78
0.70
0.72
0.81
0.69
0.86
0.83
0.60
(0.77)
(0.71)
(0.80)
(0.66)
(0.69)
(0.84)
(0.73)
(0.78)
(0.81)
(0.66)
0.69
0.62
0.68
0.60
0.53
0.73
0.63
0.66
0.63
0.61
(0.61)
(0.65)
(0.62)
(0.66)
(0.29)
(0.86)
(0.39)
(0.82)
(0.53)
(0.69)
0.58
0.64
0.66
0.60
0.38
0.77
0.20
0.91
0.65
0.61
(0.48)
(0.70)
(0.41)
(0.77)
(0.28)
(0.89)
(0.12)
(0.96)
(0.52)
(0.69)
0.64
0.68
0.73
0.59
0.49
0.75
0.51
0.80
0.74
0.59
(0.58)
(0.70)
(0.63)
(0.68)
(0.12)
(1.01)
(0.22)
(0.99)
(0.64)
(0.68)
0.67
0.66
0.63
0.67
0.61
0.70
0.54
0.78
0.74
0.58
(0.62)
(0.72)
(0.51)
(0.80)
(0.29)
(0.95)
(0.19)
(1.03)
(0.59)
(0.74)
0.72
0.70
0.64
0.81
0.29
1.20
0.05
1.38
0.66
0.76
(0.68)
(0.75)
(0.70)
(0.72)
(0.43)
(1.01)
(0.36)
(1.10)
(0.71)
(0.73)
0.57
0.70
0.67
0.65
0.48
0.84
0.15
1.12
0.64
0.65
(0.59)
(0.69)
(0.59)
(0.71)
(0.13)
(1.14)
(0.20)
(1.04)
(0.60)
(0.68)
For test compounds, potency predictions based on ECFP4 (MACCS) using ALs generated with different projection methods are compared based on the coefficient of determination (R2) and mean absolute error (MAE). For each target set, the largest R2 value is given in bold. As a reference for
ACS Paragon Plus Environment
32
Page 33 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
potency
Journal of Chemical Information and Modeling
predictions,
5-NN
calculations
ACS Paragon Plus Environment
are
also
reported.
33
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 34 of 49
The limited accuracy of quantitative predictions was further investigated. It was found to be mainly due to a strong tendency of GPR to overfit interpolated models, especially in the absence of noise level tolerance. Overfitting was strongest for neuroscale, which provided a rationale for better potency value predictions achieved with GTM and LTM models having smoother activity surfaces than neuroscale models. Despite achieving overall higher R2 and lower MAE values, it must also be considered that relative smoothing or “flattening” of potency surfaces corresponds to a loss in SAR information such that activity cliffs can no longer be identified or predicted. Table 6 reports potency predictions using neuroscale-based GPR models with varying noise level tolerance. The results clearly illustrate the negative influence of overfitting on prediction accuracy. Under the condition of small noise level tolerance, the models achieved a perfect or nearly perfect fit to training data, with only minute prediction errors. However, for nine (MACCS) and seven (ECFP4) of 10 target sets, these models yielded R2 values smaller than 0.50. When 3D ALs were derived under the condition of large noise level tolerance, this situation changed. The quality of the fit to training data was clearly reduced and the prediction accuracy of so derived models on test data consistently increased, as reported in Table 6. For ECFP4, models for seven of 10 target sets yielded R2 values greater than 0.50 and nine of 10 models produced MAE values smaller than 1.00.
ACS Paragon Plus Environment
34
Page 35 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Table 6. Predictions with Gaussian process models of varying noise level tolerance.
CHEMB L ID
Small noise level tolerance Training R2
MAE
Test R2
Large noise level tolerance Training
MAE
R2
MAE
Test R2
MAE
237
0.99 (0.97) 0.06 (0.14) 0.31 (0.04) 0.85 (0.98) 0.62 (0.48) 0.67 (0.79) 0.48 (0.23) 0.75 (0.93)
244
1.00 (0.98) 0.05 (0.16) 0.41 (0.41) 1.00 (1.01) 0.76 (0.61) 0.67 (0.86) 0.53 (0.48) 0.94 (0.98)
245
0.99 (0.98) 0.06 (0.11) 0.42 (0.34) 0.85 (0.94) 0.69 (0.80) 0.65 (0.52) 0.53 (0.43) 0.84 (0.93)
4860
1.00 (0.97) 0.05 (0.17) 0.72 (0.69) 0.81 (0.84) 0.83 (0.76) 0.61 (0.72) 0.76 (0.70) 0.77 (0.88)
1862
1.00 (0.98) 0.03 (0.11) 0.53 (0.29) 0.73 (0.86) 0.80 (0.81) 0.50 (0.48) 0.60 (0.49) 0.66 (0.74)
1983
1.00 (0.97) 0.03 (0.11) 0.38 (0.28) 0.77 (0.89) 0.87 (0.61) 0.40 (0.69) 0.49 (0.43) 0.73 (0.83)
2954
1.00 (0.99) 0.03 (0.06) 0.49 (0.12) 0.75 (1.01) 0.86 (0.84) 0.36 (0.41) 0.53 (0.21) 0.76 (1.00)
3798
1.00 (0.95) 0.04 (0.19) 0.61 (0.29) 0.70 (0.95) 0.84 (0.62) 0.45 (0.72) 0.63 (0.37) 0.69 (0.94)
220
1.00 (0.99) 0.04 (0.07) 0.29 (0.43) 1.2 (1.01)
249
1.00 (0.99) 0.02 (0.06) 0.48 (0.13) 0.84 (1.14) 0.88 (0.69) 0.43 (0.69) 0.52 (0.25) 0.83 (1.10)
0.61 (0.81) 0.94 (0.63) 0.35 (0.45) 1.19 (1.03)
Gaussian process optimization was carried out with a small or large range of allowed noise values. 2D projections of ECFP4 (MACCS) fingerprint space were generated using neuroscale. For training and test compounds, potency predictions are compared based on the coefficient of determination (R2) and mean absolute error (MAE). For each target set and representation, the largest R2 value obtained for the test set is given in bold.
ACS Paragon Plus Environment
35
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 36 of 49
Figure 5 shows a comparison of two neuroscale-GPR models with different noise level tolerance. Both 3D ALs displayed a clear separation of subsets of compounds with high or low potency. While training compounds were similarity distributed over the two landscapes, the model generated under the condition of large noise level tolerance (Figure 5, left) had a smoother surface than the model generated with small tolerance (right). Accordingly, the latter model exhibited more activity cliffs on its surface and achieved a highly accurate fit to the training data, as exemplified by correct representation of the steepest activity cliff in the data set. This cliff was formed by two close structural analogs with a potency difference of more than three orders of magnitude (compounds A and B in Figure 5). However, the activity cliff “disappeared” on the smoother surface of the model generated with higher noise level tolerance, which corresponded to a loss in SAR information. However, the smoother landscape favored potency predictions for many test compounds falling into densely populated areas on the activity surface. Importantly, accurate 3D ALs can be built for given compound data sets via low tolerance GPR based on neuroscale projections. These models are not suitable for potency value predictions but highly desirable for other applications, as discussed.
ACS Paragon Plus Environment
36
Page 37 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Figure 5. Landscapes for Gaussian process models with different noise level tolerance.
For Bcl-2 inhibitors (ID 4860), two neuroscale-based 3D ALs were generated via GPR with different noise level tolerance (left, large noise level tolerance; right, small tolerance). Training compounds were mapped onto the landscapes (black dots). Two compounds (A and B) forming an activity cliff are highlighted.
ACS Paragon Plus Environment
37
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 38 of 49
Conclusions In this work, we have further investigated the design of 3D ALs reminiscent of geographical maps. The study was motivated by the intuitive nature of these landscape models for SAR visualization. Going beyond qualitative analysis of landscape views, we were interested in exploring as to whether such 3D ALs could also be used for semi-quantitative prediction of compounds with varying potency and quantitative prediction of actual potency values. Therefore, 3D ALs were generated starting from conventional fingerprint spaces using different projection methods in combination with interpolative GPR modeling. It was shown that different projection approaches yielded 3D ALs with varying topological features. AL topology was also influenced by fingerprint representations. For generation of detailed landscape views, preference may be given to high-resolution fingerprints such as ECFPs. We found that AL views were informative and interpretable for all target sets. A major difference between GTM/LTM and MDS/neuroscale 3D ALs was the smoother activity surface of GTM- and LTM-based models. We found that 3D ALs with a highly accurate fit to training data could be derived on the basis of MDS and neuroscale, due to the strong tendency of overfitting. The resulting MDS or neuroscale ALs are excellent tools for SAR visualization and analysis of given data sets, but not suitable for potency value prediction. An important distinction between MDS and neuroscale is that neuroscale is capable of deriving a mapping function, whereas MDS is not. This made it possible to map test compounds onto neuroscale-based 3D ALs and globally distinguish between highly and weakly potent compounds. This represents a promising application of ALs. By contrast, for
ACS Paragon Plus Environment
38
Page 39 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
prediction of potency values, GTM and LTM model with smoother activity surfaces were preferred. However, the accuracy of potency value predictions remained limited in all instances. The smoothness of GTM- and LTM-based activity surfaces led to a reduction or loss of activity cliff information. On the other hand, absolute potency prediction errors for many compounds were within an order of magnitude, which would not preclude practical applications. A major reason for the overall limited accuracy of potency predictions was the strong tendency of GPR with no or low noise level tolerance to overfit models. However, GPR overfitting also enabled the derivation of highly accurate 3D ALs for given data sets, representing a trade-off. These models provided the basis for mapping of highly and weakly potent compounds to corresponding landscape regions. Taken together, these findings provided strong support for the use of 3D ALs for qualitative and semiquantitative SAR exploration and exploitation.
Supporting Information Tables S1-S4 of the Supporting Information report RMSD and CV values of distances between fingerprint and projected spaces for all training and test sets.
Acknowledgement We are grateful to Dagmar Stumpfe, University of Bonn, for providing target sets. We also thank OpenEye Scientific Software, Inc., for providing a free academic license of the OpenEye toolkit.
ACS Paragon Plus Environment
39
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 40 of 49
Conflict of Interest The authors declare no conflict of interest.
ACS Paragon Plus Environment
40
Page 41 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
References
1. Bajorath, J.; Peltason, L.; Wawer, M.; Guha, R.; Lajiness, M. S.; Van Drie, J. H. Navigating Structure-Activity Landscapes. Drug Discovery Today 2009, 14, 698-705. 2. Wassermann, A. M.; Wawer, M.; Bajorath, J. Activity Landscape Representations for StructureActivity Relationship Analysis. J. Med. Chem. 2010, 53, 8209−8223. 3. Stumpfe, D.; Bajorath J. Methods for SAR Visualization. RSC Adv. 2012, 2, 369-378. 4. Bajorath, J. Modeling of Activity Landscapes for Drug Discovery. Expert Opin. Drug Discovery 2012, 7, 463-473. 5. Hu, Y.; Stumpfe, D.; Bajorath, J. Visualization of Activity Landscapes and Chemogenomics Data. Mol. Inf. 2013, 32, 954-963. 6. Medina-Franco, J.; Navarrete-Vázquez, G.; Méndez-Lucio, O. Activity and Property Landscape Modeling is at the Interface of Chemoinformatics and Medicinal Chemistry. Future Med. Chem. 2015, 7, 1197-1211. 7. Stumpfe, D.; Bajorath, J. Recent Developments in SAR Visualization. Med. Chem. Commun. 2016, 7, 1045-1055. 8. Vogt, M. Recent Progress with Modeling Activity Landscapes in Drug Design. Expert Opin. Drug Discovery 2018, 13, 605-615. 9. Maggiora, G. M. On Outliers and Activity Cliffs − Why QSAR Often Disappoints. J. Chem. Inf. Model. 2006, 46, 1535-1535.
ACS Paragon Plus Environment
41
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 42 of 49
10. Stumpfe, D.; Bajorath, J. Exploring Activity Cliffs in Medicinal Chemistry. J. Med. Chem. 2012, 55, 2932-2942. 11. Stumpfe, D.; Hu, Y.; Dimova, D.; Bajorath, J. Recent Progress in Understanding Activity Cliffs and their Utility in Medicinal Chemistry. J. Med. Chem. 2014, 57, 18-28. 12. Cruz-Monteagudo, M.; Medina-Franco, J. L.; Perez-Castillo, Y.; Nicolotti, O.; Cordeiro, M. N. D.; Borges, F.. Activity Cliffs in Drug Discovery: Dr. Jekyll or Mr. Hyde? Drug Discovery Today 2014, 19, 1069-1080. 13. Dimova, D.; Bajorath, J. Advances in Activity Cliff Research. Mol. Inf. 2016, 35, 181-191. 14. Bajorath, J. Representation and identification of activity cliffs. Expert Opin. Drug Discovery 2017, 12, 897-883. 15. Shanmugasundaram, V.; Maggiora, G. M. Characterizing Property and Activity Landscapes Using an Information-Theoretic Approach. Proceedings of 222nd American Chemical Society National Meeting, Division of Chemical Information, Chicago, IL, August 26-30, 2001; American Chemical Society: Washington, D.C., 2001; Abstract no. 77. 16. Perez-Villanueva, J.; Santos, R.; Hernandez-Campos, A.; Giulianotti, M. A.; Castillo, R.; Medina-Franco, J. L. Structure-activity Relationships of Benzimidazole Derivatives as Antiparasitic Agents: Dual-activity Difference (DAD) Maps. Med. Chem. Commun. 2011, 2, 44-49.
17. Yongye, A. B.; Byler, K.; Santos, R.; Martínez-Mayorga, K.; Maggiora, G. M.; Medina-Franco,
ACS Paragon Plus Environment
42
Page 43 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
J. L. Consensus Models of Activity Landscapes with Multiple Chemical, Conformer, and Property Representations. J. Chem. Inf. Model. 2011, 51, 2427-2439. 18. Agrafiotis, D.; Shemanarev, M.; Connolly, P.; Farnum, M.; Lobanov, V. SAR Maps: A New SAR Visualization Technique for Medicinal Chemists. J. Med. Chem. 2007, 50, 5926-5937. 19. Guha, R.; Van Drie, J. H. Structure-Activity Landscape Index: Identifying and Quantifying Activity Cliffs. J. Chem. Inf. Model 2008, 48, 646-658. 20. Peltason, L.; Weskamp, N.; Teckentrup, A.; Bajorath, J. Exploration of Structure-Activity Relationship Determinants in Analogue Series. J. Med. Chem. 2009, 52, 3212-3224. 21. Wawer, M.; Bajorath, J. Similarity-Potency Trees: A Method to Search for SAR Information in Compound Data Sets and Derive SAR Rules. J. Chem. Inf. Model. 2010, 50, 1395-1409. 22. Dimova, D.; Wawer, M.; Wassermann, A. M.; Bajorath J. Design of Multi-Target Activity Landscapes that Capture Hierarchical Activity Cliff Distributions. J. Chem. Inf. Model. 2011, 51, 256-288. 23. Iyer, P.; Dimova, D.; Vogt, M.; Bajorath, J. Navigating High-Dimensional Activity Landscapes: Design and Application of the Ligand-Target Differentiation Map. J. Chem. Inf. Model. 2012, 52, 1962-1969. 24. Maggiora, G. M.; Bajorath, J. Chemical Space Networks – A Powerful New Paradigm for the Description of Chemical Space. J. Comput.-Aided Mol. Des. 2014, 28, 795-802.
ACS Paragon Plus Environment
43
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 44 of 49
25. Zwierzyna, M.; Vogt, M.; Maggiora, G. M.; Bajorath, J. Design and Characterization of Chemical Space Networks for Different Compound Data Sets. J. Comput.-Aided Mol. Des. 2015, 29, 113-125. 26. Mishima, K.; Kaneko, H.; Funatsu, K. Development of a New De Novo Design Algorithm for Exploring Chemical Space. Mol. Inf. 2014, 33, 779-789. 27. Gaspar, H. A.; Baskin, I. I.; Marcou, G.; Horvath, D.; Varnek, A. Stargate GTM: Bridging Descriptor and Activity Spaces. J. Chem. Inf. Model. 2015, 55, 2403-2410. 28. Klimenko, K.; Marcou, G.; Horvath, D.; Varnek, A. Chemical Space Mapping and StructureActivity Analysis of the ChEMBL Antiviral Compound Set. J. Chem. Inf. Model. 2016, 56, 1438-1454. 29. Peltason, L.; Iyer, P.; Bajorath, J. Rationalizing Three-Dimensional Activity Landscapes and the Influence of Molecular Representations on Landscape Topology and the Formation of Activity Cliffs. J. Chem. Inf. Model. 2010, 50, 1021-1033. 30. Iyer, P.; Wawer, M.; Bajorath, J. Comparison of Two- and Three-Dimensional Activity Landscape Representations for Different Compound Data sets. Med. Chem. Commun. 2011, 2, 113-118. 31. Reutlinger, M.; Guba, W.; Martin, R. E.; Alanine, A. I.; Hoffmann, T.; Klenner, A.; Hiss, J. A.; Schneider, P.; Schneider, G. Neighborhood‐Preserving Visualization of Adaptive Structure–
ACS Paragon Plus Environment
44
Page 45 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
Activity Landscapes: Application to Drug Discovery. Angew. Chem. Intl. Ed. 2011, 123, 1183711840. 32. Peltason L.; Bajorath J. SAR Index: Quantifying the Nature of Structure-Activity Relationships. J. Med. Chem. 2007, 50, 5571-5578. 33. Peltason, L.; Bajorath, J. Systematic Computational Analysis of Structure-Activity Relationships: Concepts, Challenges and Recent Advances. Future Med. Chem. 2009, 1, 451-466. 34. Hu, H.; Stumpfe, D.; Bajorath, J. Rationalizing the Formation of Activity Cliffs in Different Compound Data Sets. ACS Omega 2018, 3, 7736–7744. 35. Bento, A. P.; Gaulton, A.; Hersey, A.; Bellis, L. J.; Chambers, J.; Davies, M.; Krüger, F. A.; Light, Y.; Mak, L.; McGlinchey, S.; Nowotka, M.; Papadatos, G.; Santos, R.; Overington J. P. The ChEMBL Bioactivity Database: An Update. Nucleic Acids Res. 2014, 42, D1083-D1090. 36. Rogers, D.; Hahn, M. Extended-Connectivity Fingerprints. J. Chem. Inf. Model. 2010, 50, 742– 754. 37. MACCS Structural Keys; Accelrys: San Diego, CA, 2011. 38. OEChem TK, version 2.1.4; OpenEye Scientific Software, Santa Fe, NM. 39. Borg, I.; Groenen, P. J. F. Modern Multidimensional Scaling : Theory and Applications; Springer: New York, 2005. 40. Abdi, H.; Williams, L. J. Principal Component Analysis. WIRES: Comput. Stat. 2010, 2, 433459.
ACS Paragon Plus Environment
45
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 46 of 49
41. Bishop, C. M.; Svensén, M.; Williams, C. K. I. GTM: The Generative Topographic Mapping. Neural Comput. 1998, 10, 215–234. 42. Kabán, A.; Girolami, M. A Combined Latent Class and Trait Model for the Analysis and Visualization of Discrete Data. IEEE Trans. Pattern Anal. Mach. Intell. 2001, 23, 859–872. 43. Lowe, D.; Tipping, M. E. NeuroScale: Novel Topographic Feature Extraction Using RBF Networks. In Advances in Neural Information Processing Systems 9; Mozer, M. C., Jordan, M. I., Petsche, T., Eds.; MIT Press, 1997; pp 543–549. 44. Tipping, M. E.; Lowe, D. Shadow Targets: A Novel Algorithm for Topographic Projections by Radial Basis Functions. Neurocomputing 1998, 19, 211–222. 45. Maaten, L. van der; Hinton, G. Visualizing Data Using T-SNE. J. Mach. Learn. Res. 2008, 9, 2579–2605. 46. Joyce, J. M. Kullback-Leibler Divergence. In International Encyclopedia of Statistical Science; Springer: Berlin, Heidelberg, 2011, pp. 720-722. 47. Gómez-Bombarelli, R.; Wei, J. N.; Duvenaud, D.; Hernández-Lobato, J. M.; Sánchez-Lengeling, B.; Sheberla, D.; Aguilera-Iparraguirre, J.; Hirzel, T. D.; Adams, R. P.; Aspuru-Guzik, A., 2018. Automatic Chemical Design Using a Data-driven Continuous Representation of Molecules. ACS Centr, Sci. 2018, 4, 268-276. 48. Cressie, N. Statistics for Spatial Data; Wiley: New York, 1993.
ACS Paragon Plus Environment
46
Page 47 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
49. Rasmussen, C. E. Gaussian Processes in Machine Learning. In Advanced Lectures on Machine Learning; Springer: Berlin, Heidelberg, 2005, pp. 63-71. 50. Pedregosa, F.; Varoquaux, G.; Gramfort, A.; Michel, V.; Thirion, B.; Grisel, O.; Blondel, M.; Prettenhofer, P.; Weiss, R.; Dubourg, V.; Vanderplas, J.; Passos, A.; Cournapeau, D.; Brucher, M.; Perrot, M.; Duchesnay, É. Scikit-Learn: Machine Learning in Python. J. Mach. Learn. Res. 2011, 12, 2825–2830. 51. GPy. {GPy}: A Gaussian Process Framework in Python. http://github.com/SheffieldML/GPy 52. Owen, J. R.; Nabney, I. T.; Medina-Franco, J. L.; López-Vallejo, F. Visualization of Molecular Fingerprints. J. Chem. Inf. Model. 2011, 51, 1552–1563. 53. Cherkasov, A.; Muratov, E. N.; Fourches, D.; Varnek, A.; Baskin, I. I.; Cronin, M.; Dearden, J.; Gramatica, P.; Martin, Y. C.; Todeschini, R.; Consonni, V.; Kuz’min, V. E.; Cramer, R.; Benigni, R.; Yang, C.; Rathman, J.; Terfloth, L.; Gasteiger, J.; Richard, A.; Tropsha, A. QSAR Modeling: Where Have You Been? Where Are You Going To? J. Med. Chem. 2014, 57, 49775010. 54. Wawer, M.; Peltason, L.; Weskamp, N.; Teckentrup, A.; Bajorath, J. Structure−Activity Relationship Anatomy by Network-like Similarity Graphs and Local Structure−Activity Relationship Indices. J. Med. Chem. 2008, 51, 6075–6084.
ACS Paragon Plus Environment
47
Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 48 of 49
55. de la Vega de León, A.; Bajorath, J. Formation of Activity Cliffs Is Accompanied by Systematic Increases in Ligand Efficiency from Lowly to Highly Potent Compounds. AAPS J. 2014, 16, 335–341.
ACS Paragon Plus Environment
48
Page 49 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Journal of Chemical Information and Modeling
TOC
ACS Paragon Plus Environment
49