Surely you are joking, Mr Docking!

In the wake of recent COVID-19 pandemics scientists around the world rushed to deliver numerous CADD (Computer-Aided Drug Discovery) methods and tools that could be reliably used to discover novel drug candidates against the SARS-CoV-2 virus. With that, there emerged a trend of a significant democratization of CADD that contributed to the rapid development of various COVID-19 drug candidates currently undergoing different stages of validation. On the other hand, this democratization also inadvertently led to the surge rapidly performed molecular docking studies to nominate multiple scores of novel drug candidates supported by computational arguments only. Albeit driven by best intentions, most of such studies also did not follow best practices in the field that require experience and expertise learned through multiple rigorously designed benchmarking studies and rigorous experimental validation. In this Viewpoint we reflect on recent disbalance between small number of rigorous and comprehensive studies and the proliferation of purely computational studies enabled by the ease of docking software availability. We further elaborate on the hyped oversale of CADD methods' ability to rapidly yield viable drug candidates and reiterate the critical importance of rigor and adherence to the best practices of CADD in view of recent emergence of AI and Big Data in the field.

Beware of R(2): Simple, Unambiguous Assessment of the Prediction Accuracy of QSAR and QSPR Models.

The statistical metrics used to characterize the external predictivity of a model, i.e., how well it predicts the properties of an independent test set, have proliferated over the past decade. This paper clarifies some apparent confusion over the use of the coefficient of determination, R(2), as a measure of model fit and predictive power in QSAR and QSPR modeling. R(2) (or r(2)) has been used in various contexts in the literature in conjunction with training and test data for both ordinary linear regression and regression through the origin as well as with linear and nonlinear regression models. We analyze the widely adopted model fit criteria suggested by Golbraikh and Tropsha ( J. Mol. Graphics Modell. 2002 , 20 , 269 - 276 ) in a strict statistical manner. Shortcomings in these criteria are identified, and a clearer and simpler alternative method to characterize model predictivity is provided. The intent is not to repeat the well-documented arguments for model validation using test data but rather to guide the application of R(2) as a model fit statistic. Examples are used to illustrate both correct and incorrect uses of R(2). Reporting the root-mean-square error or equivalent measures of dispersion, which are typically of more practical importance than R(2), is also encouraged, and important challenges in addressing the needs of different categories of users such as computational chemists, experimental scientists, and regulatory decision support specialists are outlined.

QSAR analysis of the toxicity of nitroaromatics in Tetrahymena pyriformis: structural factors and possible modes of action.

The Hierarchical Technology for Quantitative Structure-Activity Relationships (HiT QSAR) was applied to 95 diverse nitroaromatic compounds (including some widely known explosives) tested for their toxicity (50% inhibition growth concentration, IGC50) against the ciliate Tetrahymena pyriformis. The dataset was divided into subsets according to putative mechanisms of toxicity. The Classification and Regression Trees (CART) approach implemented within HiT QSAR has been used for prediction of mechanism of toxicity for new compounds. The resulting models were shown to have ~80% accuracy for external datasets indicating that the mechanistic dataset division was sensible. The Partial Least Squares (PLS) statistical approach was then used to develop 2D QSAR models. Validated PLS models were explored to: (1) elucidate the effects of different substituents in nitroaromatic compounds on toxicity; (2) differentiate compounds by probable mechanisms of toxicity based on their structural descriptors; and (3) analyse the role of various physical-chemical factors responsible for compounds' toxicity. Models were interpreted in terms of molecular fragments promoting or interfering with toxicity. It was also shown that mutual influence of substituents in benzene ring plays the determining role in toxicity variation. Although chemical mechanism based models were statistically significant and externally predictive (r²(ext) = 0.64 for the external set of 63 nitroaromatics identified after all calculations have been completed), they were also shown to have limited coverage (57% for modelling and 76% for external set).

QSAR modeling of GPCR ligands: methodologies and examples of applications.

GPCR ligands represent not only one of the major classes of current drugs but the major continuing source of novel potent pharmaceutical agents. Because 3D structures of GPCRs as determined by experimental techniques are still unavailable, ligand-based drug discovery methods remain the major computational molecular modeling approaches to the analysis of growing data sets of tested GPCR ligands. This paper presents an overview of modern Quantitative Structure Activity Relationship (QSAR) modeling. We discuss the critical issue of model validation and the strategy for applying the successfully validated QSAR models to virtual screening of available chemical databases. We present several examples of applications of validated QSAR modeling approaches to GPCR ligands. We conclude with the comments on exciting developments in the QSAR modeling of GPCR ligands that focus on the study of emerging data sets of compounds with dual or even multiple activities against two or more of GPCRs.

QSAR modeling of datasets with enantioselective compounds using chirality sensitive molecular descriptors.

Shape descriptors used in 3D QSAR studies naturally take into account chirality; however, for flexible and structurally diverse molecules such studies require extensive conformational searching and alignment. QSAR modeling studies of two datasets of fragrance compounds with complex stereochemistry using simple alignment-free chirality sensitive descriptors developed in our laboratories are presented. In the first investigation, 44 alpha-campholenic derivatives with sandalwood odor were represented as derivatives of several common structural templates with substituents numbered according to their relative spatial positions in the molecules. Both molecular and substituent descriptors were used as independent variables in MLR calculations, and the best model was characterized by the training set q2 of 0.79 and external test set r2 of 0.95. In the second study, several types of chirality descriptors were employed in combinatorial QSAR modeling of 98 ambergris fragrance compounds. Among 28 possible combinations of seven types of descriptors and four statistical modeling techniques, k nearest neighbor classification with CoMFA descriptors was initially found to generate the best models with the internal and external accuracies of 76 and 89%, respectively. The same dataset was then studied using novel atom pair chirality descriptors (cAP). The cAP are based on a modified definition of the atomic chirality, in which the seniority of the substituents is defined by their relative partial charge values: higher values correspond to higher seniorities. The resulting models were found to have higher predictive power than those developed with CoMFA descriptors; the best model was characterized by the internal and external accuracies of 82 and 94%, respectively. The success of modeling studies using simple alignment free chirality descriptors discussed in this paper suggests that they should be applied broadly to QSAR studies of many datasets when compound stereochemistry plays an important role in defining their activity.

Accurate classification of protein structural families using coherent subgraph analysis.

Protein structural annotation and classification is an important problem in bioinformatics. We report on the development of an efficient subgraph mining technique and its application to finding characteristic substructural patterns within protein structural families. In our method, protein structures are represented by graphs where the nodes are residues and the edges connect residues found within certain distance from each other. Application of subgraph mining to proteins is challenging for a number reasons: (1) protein graphs are large and complex, (2) current protein databases are large and continue to grow rapidly, and (3) only a small fraction of the frequent subgraphs among the huge pool of all possible subgraphs could be significant in the context of protein classification. To address these challenges, we have developed an information theoretic model called coherent subgraph mining. From information theory, the entropy of a random variable X measures the information content carried by X and the Mutual Information (MI) between two random variables X and Y measures the correlation between X and Y. We define a subgraph X as coherent if it is strongly correlated with every sufficiently large sub-subgraph Y embedded in it. Based on the MI metric, we have designed a search scheme that only reports coherent subgraphs. To determine the significance of coherent protein subgraphs, we have conducted an experimental study in which all coherent subgraphs were identified in several protein structural families annotated in the SCOP database (Murzin et al, 1995). The Support Vector Machine algorithm was used to classify proteins from different families under the binary classification scheme. We find that this approach identifies spatial motifs unique to individual SCOP families and affords excellent discrimination between families.

Diversity and coverage of structural sublibraries selected using the SAGE and SCA algorithms.

It is often impractical to synthesize and test all compounds in a large exhaustive chemical library. Herein, we discuss rational approaches to selecting representative subsets of virtual libraries that help direct experimental synthetic efforts for diverse library design. We compare the performance of two stochastic sampling algorithms, Simulating Annealing Guided Evaluation (SAGE; Zheng, W.; Cho, S. J.; Waller, C. L.; Tropsha, A. J. Chem. Inf. Comput. Sci. 1999, 39, 738-746.) and Stochastic Cluster Analysis (SCA; Reynolds, C. H.; Druker, R.; Pfahler, L. B. Lead Discovery Using Stochastic Cluster Analysis (SCA): A New Method for Clustering Structurally Similar Compounds J. Chem. Inf. Comput. Sci. 1998, 38, 305-312.) for their ability to select both diverse and representative subsets of the entire chemical library space. The SAGE and SCA algorithms were compared using u- and s-optimal metrics as an independent assessment of diversity and coverage. This comparison showed that both algorithms were capable of generating sublibraries in descriptor space that are diverse and give reasonable coverage (i.e. are representative) of the original full library. Tests were carried out using simulated two-dimensional data sets and a 27 000 compound proprietary structural library as represented by computed Molconn-Z descriptors. One of the key observations from this work is that the algorithmically simple SCA method is capable of selecting subsets that are comparable to the more computationally intensive SAGE method.

Four-body potentials reveal protein-specific correlations to stability changes caused by hydrophobic core mutations.

Mutational experiments show how changes in the hydrophobic cores of proteins affect their stabilities. Here, we estimate these effects computationally, using four-body likelihood potentials obtained by simplicial neighborhood analysis of protein packing (SNAPP). In this procedure, the volume of a known protein structure is tiled with tetrahedra having the center of mass of one amino acid side-chain at each vertex. Log-likelihoods are computed for the 8855 possible tetrahedra with equivalent compositions from structural databases and amino acid frequencies. The sum of these four-body potentials for tetrahedra present in a given protein yields the SNAPP score. Mutations change this sum by changing the compositions of tetrahedra containing the mutated residue and their related potentials. Linear correlation coefficients between experimental mutational stability changes, Delta(DeltaG(unfold)), and those based on SNAPP scoring range from 0.70 to 0.94 for hydrophobic core mutations in five different proteins. Accurate predictions for the effects of hydrophobic core mutations can therefore be obtained by virtual mutagenesis, based on changes to the total SNAPP likelihood potential. Significantly, slopes of the relation between Delta(DeltaG(unfold)) and DeltaSNAPP for different proteins are statistically distinct, and we show that these protein-specific effects can be estimated using the average SNAPP score per residue, which is readily derived from the analysis itself. This result enhances the predictive value of statistical potentials and supports previous suggestions that "comparable" mutations in different proteins may lead to different Delta(DeltaG(unfold)) values because of differences in their flexibility and/or conformational entropy.

Accurate prediction of the bound conformation of galanthamine in the active site of Torpedo californica acetylcholinesterase using molecular docking.

The alkaloid (-)-galanthamine is known to produce significant improvement of cognitive performances in patients with the Alzheimer's disease. Its mechanism of action involves competitive and reversible inhibition of acetylcholinesterase (AChE). Herein, we correctly predict the orientation and conformation of the galanthamine molecule in the active site of AChE from Torpedo californica (TcAChE) using a combination of rigid docking and flexible geometry optimization with a molecular mechanics force field. The quality of the predicted model is remarkable, as indicated by the value of the RMS deviation of approximately 0.5A when compared with the crystal structure of the TcAChE-galanthamine complex. A molecular model of the complex between TcAChE and a galanthamine derivative, SPH1107, with a long chain substituent on the nitrogen has been generated as well. The side chain of this ligand is predicted to extend along the enzyme active site gorge from the anionic subsite, at the bottom, to the peripheral anionic site, at the top. The docking procedure described in this paper can be applied to produce models of ligand-receptor complexes for AChE and other macromolecular targets of drug design.

Identification of the descriptor pharmacophores using variable selection QSAR: applications to database mining.

The pharmacophore concept is central to the rational drug design and discovery process. Traditionally, a pharmacophore is defined as a specific three-dimensional (3D) arrangement of chemical functional groups found in active molecules, which are characteristic of a certain pharmacological class of compounds. Herein, by analogy with 3D pharmacophores, a more general concept of descriptor pharmacophore is introduced. The descriptor pharmacophores are defined by the means of variable selection QSAR as a subset of molecular descriptors that afford the most statistically significant structure-activity correlation. The two variable selection QSAR methods developed in this laboratory are discussed; these include Genetic Algorithms--Partial Least Squares (GA-PLS) and K-Nearest Neighbors (KNN). Both methods employ multiple topological descriptors of chemical structures such as molecular connectivity indices or atom pairs (AP), and stochastic optimization algorithms to achieve a robust QSAR model, which is characterized by the highest value of cross-validated R2 (q2). By default, the descriptor pharmacophore represents an invariant selection of descriptor types however, descriptor values are generally different for different molecules. We demonstrate that chemical similarity searches using descriptor pharmacophores as opposed to using all descriptors afford more efficient mining of chemical databases or virtual libraries to discover compounds with a desired biological activity.

Lattice protein folding with two and four-body statistical potentials.

The cooperative folding of proteins implies a description by multibody potentials. Such multibody potentials can be generalized from common two-body statistical potentials through a relation to probability distributions of residue clusters via the Boltzmann condition. In this exploratory study, we compare a four-body statistical potential, defined by the Delaunay tessellation of protein structures, to the Miyazawa-Jernigan (MJ) potential for protein structure prediction, using a lattice chain growth algorithm. We use the four-body potential as a discriminatory function for conformational ensembles generated with the MJ potential and examine performance on a set of 22 proteins of 30-76 residues in length. We find that the four-body potential yields comparable results to the two-body MJ potential, namely, an average coordinate root-mean-square deviation (cRMSD) value of 8 A for the lowest energy configurations of all-alpha proteins, and somewhat poorer cRMSD values for other protein classes. For both two and four-body potentials, superpositions of some predicted and native structures show a rough overall agreement. Formulating the four-body potential using larger data sets and direct, but costly, generation of conformational ensembles with multibody potentials may offer further improvements. Proteins 2001;43:161-174.

Novel chirality descriptors derived from molecular topology.

Several series of novel chirality descriptors of chemical organic molecules have been introduced. The descriptors have been developed on the basis of conventional topological descriptors of molecular graphs. They include modified molecular connectivity indices, Zagreb group indices, extended connectivity, overall connectivity, and topological charge indices. These modified descriptors make use of an additional term called chirality correction, which is added to the vertex degrees of asymmetric atoms in a molecular graph. Chirality descriptors can be real or complex numbers. Advantages and drawbacks of different series of chirality descriptors are discussed. These descriptors circumvent the inability of conventional topological indices to distinguish chiral or enantiomeric isomers, which so far has been the major drawback of 2D descriptors as compared to true 3D descriptors (e.g., shape, molecular fields) of molecular structure. These novel chirality descriptors have been implemented in a quantitative structure-activity releationship (QSAR) study of a set of ecdysteroids with a high content of chiral and enantiomeric compounds using the k nearest neighbor QSAR method (kNN) recently developed in this laboratory. We show that the results of this study compare favorably with those obtained with the comparative molecular field analysis (CoMFA) applied to the same dataset. The novel chirality descriptors of molecular structure should find their applications in QSAR studies and related investigations of molecular sdatasets.

Novel variable selection quantitative structure--property relationship approach based on the k-nearest-neighbor principle

A novel automated variable selection quantitative structure-activity relationship (QSAR) method, based on the kappa-nearest neighbor principle (kNN-QSAR) has been developed. The kNN-QSAR method explores formally the active analogue approach, which implies that similar compounds display similar profiles of pharmacological activities. The activity of each compound is predicted as the average activity of K most chemically similar compounds from the data set. The robustness of a QSAR model is characterized by the value of cross-validated R2 (q2) using the leave-one-out cross-validation method. The chemical structures are characterized by multiple topological descriptors such as molecular connectivity indices or atom pairs. The chemical similarity is evaluated by Euclidean distances between compounds in multidimensional descriptor space, and the optimal subset of descriptors is selected using simulated annealing as a stochastic optimization algorithm. The application of the kNN-QSAR method to 58 estrogen receptor ligands as well as to several other groups of pharmacologically active compounds yielded QSAR models with q2 values of 0.6 or higher. Due to its relative simplicity, high degree of automation, nonlinear nature, and computational efficiency, this method could be applied routinely to a large variety of experimental data sets.

An efficient projection protocol for chemical databases: singular value decomposition combined with truncated-newton minimization.

A rapid algorithm for visualizing large chemical databases in a low-dimensional space (2D or 3D) is presented as a first step in database analysis and design applications. The projection mapping of the compound database (described as vectors in the high-dimensional space of chemical descriptors) is based on the singular value decomposition (SVD) combined with a minimization procedure implemented with the efficient truncated-Newton program package (TNPACK). Numerical experiments on four chemical datasets with real-valued descriptors (ranging from 58 to 27 255 compounds) show that the SVD/TNPACK projection duo achieves a reasonable accuracy in 2D, varying from 30% to about 100% of pairwise distance segments that lie within 10% of the original distances. The lowest percentages, corresponding to scaled datasets, can be made close to 100% with projections onto a 10-dimensional space. We also show that the SVD/TNPACK duo is efficient for minimizing the distance error objective function (especially for scaled datasets), and that TNPACK is much more efficient than a current popular approach of steepest descent minimization in this application context. Applications of our projection technique to similarity and diversity sampling in drug design can be envisioned.

Antitumor agents. 199. Three-dimensional quantitative structure-activity relationship study of the colchicine binding site ligands using comparative molecular field analysis.

Inhibitors of tubulin polymerization interacting at the colchicine binding site are potential anticancer agents. We have been involved in the synthesis of a number of colchicine site agents, such as thiocolchicinoids and allocolchicinoids, which are colchicine analogues, and 2-phenyl-quinolones and 2-aryl-naphthyridinones, which are the amino analogues of cytotoxic antimitotic flavonoids. The most cytotoxic of the latter compounds strongly inhibit binding of radiolabeled colchicine to tubulin, and these agents therefore probably bind in the colchicine site of tubulin. We have applied conventional CoMFA and q(2)-GRS CoMFA to identify the essential structural requirements for increasing the ability of these compounds to form tubulin complexes. The CoMFA model for the training set of 51 compounds yielded cross-validated R(2) (q(2)) values of 0.637 for conventional CoMFA and 0.692 for q(2)-GRS CoMFA. The predictive power of this model was confirmed by successful activity prediction for a test set of 53 compounds with known potencies as inhibitors of tubulin polymerization. The activities of 88% of the compounds were predicted with absolute value of residuals of less than 0.5. The predictive q(2) values were 0.546 for conventional CoMFA and 0.426 for q(2)-GRS CoMFA. The conventional CoMFA model with the highest predictive q(2) (0.546) was analyzed in detail in terms of underlying structure-activity relationships.

Recent trends in computer-aided drug discovery.

The field of computer-aided drug design and discovery is experiencing one of the most exciting periods in its history. The changes have been brought about by an extraordinarily rapid expansion in molecular and biomolecular database technology, which have dramatically altered the way we store, process and exploit molecular information to establish a relationship between chemical structure and biological action. Modern approaches to drug design and discovery are characterized by computational tool integration, an increasing role for bioinformatics and a paradigm shift from the design of high-affinity ligands to the discovery of drug-like molecules with optimal ADME (absorption, distribution, metabolism, excretion) and toxicity properties. This review briefly summarizes recent and developing trends in the field that concur with advances in chemo- and bioinformatics.

Automated pharmacophore identification for large chemical data sets.

The identification of three-dimensional pharmacophores from large, heterogeneous data sets is still an unsolved problem. We developed a novel program, SCAMPI (statistical classification of activities of molecules for pharmacophore identification), for this purpose by combining a fast conformation search with recursive partitioning, a data-mining technique, which can easily handle large data sets. The pharmacophore identification process is designed to run recursively, and the conformation spaces are resampled under the constraints of the evolving pharmacophore model. This program is capable of deriving pharmacophores from a data set of 1000-2000 compounds, with thousands of conformations generated for each compound and in less than 1 day of computational time. For two test data sets, the identified pharmacophores are consistent with the known results from the literature.

The "random-coil" state of proteins: comparison of database statistics and molecular simulations.

This study presents a comparison of two models of the random-coil state, one based on statistical distributions from the structural database and the other based on molecular dynamics simulations. The database model relies on the assumption that the random- or statistical-coil state of a particular residue can be described by its conformational distribution in a sufficiently diverse subset of protein structures. The molecular dynamics model is based on distributions from molecular simulations carried out on "dipeptide" models (single residues with N-terminal acetyl and C-terminal N'-methyl amide blocking groups). A comparison of the two models for the residues Ala, Asn, Asp, Gly, and Val indicates that the database distributions are greatly influenced by long-range interactions and dominated by specific recognizable elements of protein structure. In contrast, the limited structural scope of the dipeptide models presents the extreme case of a peptide under the influence of only short-range interactions. The models were evaluated by a comparison of scalar coupling constants calculated from the conformational distributions and compared with experimentally values determined for unstructured peptides. Although the models gave different distributions, there was similar agreement with experiment. This comparison emphasizes the differences and limitations in each model and highlights the difficulty in presenting an accurate picture of the random-coil state. Proteins 1999;36:407- 418.

Quantitative structure-activity relationship modeling of dopamine D(1) antagonists using comparative molecular field analysis, genetic algorithms-partial least-squares, and K nearest neighbor methods.

Several quantitative structure-activity relationship (QSAR) methods were applied to 29 chemically diverse D(1) dopamine antagonists. In addition to conventional 3D comparative molecular field analysis (CoMFA), cross-validated R(2) guided region selection (q(2)-GRS) CoMFA (see ref 1) was employed, as were two novel variable selection QSAR methods recently developed in one of our laboratories. These latter methods included genetic algorithm-partial least squares (GA-PLS) and K nearest neighbor (KNN) procedures (see refs 2-4), which utilize 2D topological descriptors of chemical structures. Each QSAR approach resulted in a highly predictive model, with cross-validated R(2) (q(2)) values of 0.57 for CoMFA, 0.54 for q(2)-GRS, 0.73 for GA-PLS, and 0.79 for KNN. The success of all of the QSAR methods indicates the presence of an intrinsic structure-activity relationship in this group of compounds and affords more robust design and prediction of biological activities of novel D(1) ligands.