Do electrostatic interactions make a difference in physics-based AutoDock4 scoring function?

Physics-based scoring function AutoDock4 is one of the most successfully applied tools in the area of structure-based drug design. However, current scoring functions are still far from being perfect. In a recent work highlighting the strengths and deficiencies of current scoring functions, we discovered that the residual error of ΔGbind predictions made by AutoDock4 is highly correlated to the presence of formally charged fragments in a ligand. In this work, we study how the use of the high-quality atomic charges, applied for contemporary force fields calculation, affects the quality of the experimental ΔGbind prediction by means of AutoDock4. We initially expected that the previously found discrepancy could be attributed to the Gasteiger charges used within AutoDock4. We show that AutoDock4 is, surprisingly, not sensitive to the charges used, and the use of QC-derived atomic charges does not lead to any statistical improvements. We also briefly discuss the role of the explicit empirical hydrogen bond term along with the electrostatic term.

On hidden anisotropy of formally charged fragments.

The proper and precise reproduction of the molecular electrostatic potential (MEP) is crucial to describe correctly electrostatic interactions in molecular modeling. Most of the classical molecular mechanics force fields for biomolecules and drug-like molecules use the atom-centered point charges to describe MEP. However, it has been systematically pointed out in literature that such an approximation is not always enough, and some groups, like amino group or heavy halogens, require the use of anisotropic model for better description of their MEP. At the same time, the formally charged groups have not been as extensively and systematically studied as their neutral counterparts. In this report, we demonstrate that the anisotropic models for formally charged groups do bring improvements in the reference MEP reproduction, that are comparable in magnitude to those for neutral groups.

Investigation of the tight binding mechanism of a new anticoagulant DD217 to factor Xa by means of molecular docking and molecular dynamics.

A new promising drug candidate DD217 has been proposed recently as a potent anticoagulant acting on factor Xa (fXa) target. It exhibits the lowest concentration of doubling the prothrombin time among the known anticoagulants. In order to explain the efficacy of DD217 in terms of molecular interactions with its target we studied the hypothesis of the tight binding mechanism by means of molecular dynamics simulations and statistical analysis of the trajectory. The conducted analysis confirms the significant contributions to the MM/GBSA estimated binding free energy of the S4 pocket residues as well the crucial role of establishing the hydrogen bonds between the ligand and the backbone amides of Gly216 and Gly218 of the target. The simulation results support the hypothesis of the tight binding mechanism of DD217 to fXa.Communicated by Ramaswamy H. Sarma.

Assessing How Residual Errors of Scoring Functions Correlate to Ligand Structural Features.

Scoring functions (SFs) are ubiquitous tools for early stage drug discovery. However, their accuracy currently remains quite moderate. Despite a number of successful target-specific SFs appearing recently, up until now, no ideas on how to systematically improve the general scope of SFs have been formulated. In this work, we hypothesized that the specific features of ligands, corresponding to interactions well appreciated by medicinal chemists (e.g., hydrogen bonds, hydrophobic and aromatic interactions), might be responsible, in part, for the remaining SF errors. The latter provides direction to efforts aimed at the rational and systematic improvement of SF accuracy. In this proof-of-concept work, we took a CASF-2016 coreset of 285 ligands as a basis for comparison and calculated the values of scores for a representative panel of SFs (including AutoDock 4.2, AutoDock Vina, X-Score, NNScore2.0, ΔVina RF20, and DSX). The residual error of linear correlation of each SF value, with the experimental values of affinity and activity, was then analyzed in terms of its correlation with the presence of the fragments responsible for certain medicinal chemistry defined interactions. We showed that, despite the fact that SFs generally perform reasonably, there is room for improvement in terms of better parameterization of interactions involving certain fragments in ligands. Thus, this approach opens a potential way for the systematic improvement of SFs without their significant complication. However, the straightforward application of the proposed approach is limited by the scarcity of reliable available data for ligand-receptor complexes, which is a common problem in the field.

Theoretical Studies of Leu-Pro-Arg-Asp-Ala Pentapeptide (LPRDA) Binding to Sortase A of Staphylococcus aureus.

Sortase A (SrtA) of Staphylococcus aureus is a well-defined molecular target to combat the virulence of these clinically important bacteria. However up to now no efficient drugs or even clinical candidates are known, hence the search for such drugs is still relevant and necessary. SrtA is a complex target, so many straight-forward techniques for modeling using the structure-based drug design (SBDD) fail to produce the results they used to bring for other, simpler, targets. In this work we conduct theoretical studies of the binding/activity of Leu-Pro-Arg-Asp-Ala (LPRDA) polypeptide, which was recently shown to possess antivirulence activity against S. aureus. Our investigation was aimed at establishing a framework for the estimation of the key interactions and subsequent modification of LPRDA, targeted at non-peptide molecules, with better drug-like properties than the original polypeptide. Firstly, the available PDB structures are critically analyzed and the criteria to evaluate the quality of the ligand-SrtA complex geometry are proposed. Secondly, the docking protocol was investigated to establish its applicability to the LPRDA-SrtA complex prediction. Thirdly, the molecular dynamics studies were carried out to refine the geometries and estimate the stability of the complexes, predicted by docking. The main finding is that the previously reported partially chaotic movement of the ß6/ß7 and ß7/ß8 loops of SrtA (being the intrinsically disordered parts related to the SrtA binding site) is exaggerated when SrtA is complexed with LPRDA, which in turn reveals all the signs of the flexible and structurally disordered molecule. As a result, a wealth of plausible LPRDA-SrtA complex conformations are hard to distinguish using simple modeling means, such as docking. The use of more elaborate modeling approaches may help to model the system reliably but at the cost of computational efficiency.

Is an Inductive Effect Explicit Account Required for Atomic Charges Aimed at Use within the Force Fields?

Polarization and inductive effects are the concepts that have been widely used in qualitative and even quantitative descriptions of experimentally observed properties in chemistry. The polarization effect has proven to be important in cases of biomolecular modeling though still the vast majority of molecular simulations use the classical non-polarizable force fields. In the last few decades, a lot of effort has been put into promoting the polarization effect and incorporating it into modern force fields and charge calculation methods. In contrast, the inductive effect has not attracted such attention and is effectively absent in both classic and modern force fields. Thus, a question is whether this difference corresponds to the difference in the physical significance of the effects and their explicit account, or is an artifact that should be corrected in the next generation of force fields. The significance of the electronic effects is studied in this paper through the prism of performance of specific models for atomic charge calculation that take into explicit account a nested set of effects: the formal charge, the nearest neighbors, the inductive effect, and finally the model, which takes into account all effects, which are possible to account for using atomic charges. The specific choice for the methods is the following: formal charges, MMFF94 bond charge increments, Dynamic Electronegativity Relaxation (DENR), and RESP. We propose a special scheme for the separate estimation of each particular effect contribution. By pairwise comparing the residual molecular electrostatic potential (MEP) errors of those charge models (aimed at best reproducing the quantum chemical reference MEP), we sequentially revealed how the account of each effect contributes to the better-quality MEP reproduction. The following relative importance of effects was estimated; thus, the natural hierarchy of the effects was established. First, the account of formal charges is of primordial importance. Second, the nearest neighbors account is the next in significance. Third, the explicit account of inductive effect in empirical charge calculation schemes was shown to significantly─both qualitatively and quantitatively─improve the quality of MEP reproduction. Fourth, the contribution of polarization is indirectly assessed. Surprisingly, it is of the order of magnitude of the inductive effect even for the molecular systems, for which it is anticipated to be more significant. Finally, the relative importance of anisotropic effects in neutral molecules was additionally reviewed.

In Silico Structure-Based Approach for Group Efficiency Estimation in Fragment-Based Drug Design Using Evaluation of Fragment Contributions.

The notion of a contribution of a specific group in an organic molecule's property and/or activity is both common in our thinking and is still not strictly correct due to the inherent non-additivity of free energy with respect to molecular fragments composing a molecule. The fragment- based drug discovery (FBDD) approach has proven to be fruitful in addressing the above notions. The main difficulty of the FBDD, however, is in its reliance on the low throughput and expensive experimental means of determining the fragment-sized molecules binding. In this article we propose a way to enhance the throughput and availability of the FBDD methods by judiciously using an in silico means of assessing the contribution to ligand-receptor binding energy of fragments of a molecule under question using a previously developed in silico Reverse Fragment Based Drug Discovery (R-FBDD) approach. It has been shown that the proposed structure-based drug discovery (SBDD) type of approach fills in the vacant niche among the existing in silico approaches, which mainly stem from the ligand-based drug discovery (LBDD) counterparts. In order to illustrate the applicability of the approach, our work retrospectively repeats the findings of the use case of an FBDD hit-to-lead project devoted to the experimentally based determination of additive group efficiency (GE)-an analog of ligand efficiency (LE) for a group in the molecule-using the Free-Wilson (FW) decomposition. It is shown that in using our in silico approach to evaluate fragment contributions of a ligand and to estimate GE one can arrive at similar decisions as those made using the experimentally determined activity-based FW decomposition. It is also shown that the approach is rather robust to the choice of the scoring function, provided the latter demonstrates a decent scoring power. We argue that the proposed approach of in silico assessment of GE has a wider applicability domain and expect that it will be widely applicable to enhance the net throughput of drug discovery based on the FBDD paradigm.

Selection of Promising Novel Fragment Sized S. aureus SrtA Noncovalent Inhibitors Based on QSAR and Docking Modeling Studies.

Sortase A (SrtA) of Staphylococcus aureus has been identified as a promising target to a new type of antivirulent drugs, and therefore, the design of lead molecules with a low nanomolar range of activity and suitable drug-like properties is important. In this work, we aimed at identifying new fragment-sized starting points to design new noncovalent S. aureus SrtA inhibitors by making use of the dedicated molecular motif, 5-arylpyrrolidine-2-carboxylate, which has been previously shown to be significant for covalent binding SrtA inhibitors. To this end, an in silico approach combining QSAR and molecular docking studies was used. The known SrtA inhibitors from the ChEMBL database with diverse scaffolds were first employed to derive descriptors and interpret their significance and correlation to activity. Then, the classification and regression QSAR models were built, which were used for rough ranking of the virtual library of the synthetically feasible compounds containing the dedicated motif. Additionally, the virtual library compounds were docked into the "activated" model of SrtA (PDB:2KID). The consensus ranking of the virtual library resulted in the most promising structures, which will be subject to further synthesis and experimental testing in order to establish new fragment-like molecules for further development into antivirulent drugs.

Computational characterization of the glutamate receptor antagonist perampanel and its close analogs: density functional exploration of conformational space and molecular docking study.

Perampanel approved by FDA in 2012 is a first-in-class antiepileptic drug which inhibits α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid (AMPA) receptor currents. It is markedly more active than many of its close analogs, and the reasons for this activity difference are not quite clear. Recent crystallographic studies allowed the authors to identify the location of its binding site. Unfortunately, the resolution is low, and the detailed description of perampanel binding mode is still in part speculative. Here we provide a detailed DFT-level conformational analysis of perampanel in a vacuum and in the solvents, mimicking the protein environment, followed by quantum theory of atoms in molecules (QTAIM), non-covalent interactions (NCI), and natural bond orbital (NBO) analyses. The findings indicate the electrostatic nature of the intramolecular interactions which contribute to energy differences of the conformations in a vacuum whereas the increase of dielectric constant leads to the energy equalization of conformations. Based on these results, the docking study was performed to investigate possible binding modes of perampanel and its close analogs in AMPA receptors. The influence of the pyridine nitrogen and cyano group position was explained based on the results of conformational analysis and molecular docking. These findings may contribute to the design of novel antiepileptic drugs and the development of novel approaches to treat neurodegenerative diseases and major depressive disorder.

Iterative Solvers for Empirical Partial Atomic Charges: Breaking the Curse of Cubic Numerical Complexity.

Rational drug design involves a vast amount of computations to get thermodynamically reliable results and often relies on atomic charges as a means to model electrostatic interactions within the system. Computational inefficiency often hampers the development of new and wider dissemination of the known methods; thus, any source to speed up the calculations without a sacrifice in quality is warranted. At the heart of many empirical methods of calculating atomic charges is the solution of a system of linear algebraic equations (SLAE). The classical method of solving SLAE-the Gauss method-has in general case a cubic computational complexity. It is shown that the use of iterative methods for solving SLAE, characteristic to typical empirical atomic charge calculation methods, makes it possible to significantly reduce the amount of calculations and to obtain a computational complexity approaching a quadratic one. Despite the fact that this phenomenon is well-known in numerical methods, iterative solvers surprisingly do not seem to have been systematically applied to calculation of atomic charges via empirical schemes. Another finding is the relative values of the matrix elements, determined by the physical grounds of the interactions within the empirical system, generally lead to SLAE's with well-defined matrices, suited to use with iterative solvers to fasten computation compared to using the noniterative solvers. This finding broadens the applicability range of atomic charges obtained with empirical methods for such cases as, e.g., account of polarizability via "on-the-fly" recalculation of charges in changing surroundings within the force fields in molecular dynamics settings.

Quadrupole Correction: From Molecular Electrostatic Potential to Free Energies of Halogen Bonding.

Recently several molecular mechanics models of halogen bonding have been published. They describe the electrostatic potential anisotropy near the heavy halogen atom (known as a σ-hole) in different ways, ranging from an all-atom multipole expansion to a single positive extra-point charge. However, the question of a reasonable balance between the accuracy and the simplicity of the model remains open. In this work, we introduce the simplistic RESPQ electrostatics model built on the RESP charges complemented with fixed atomic quadrupoles. We show that it: (1) correctly describes the MEP anisotropy of aromatic halogen atoms, (2) improves the description of the halogen-water interaction energies in both halogen and hydrogen bonding cases, (3) provides an excellent estimation of solvation free energy differences of aromatic halogens, and (4) is compatible (with the help of multipole charge cluster approximation) with contemporary molecular modeling packages.

Perspectives of Halogen Bonding Description in Scoring Functions and QSAR/QSPR: Substituent Effects in Aromatic Core.

Halogen bonding (XB) is a new promising interaction pattern in medicinal chemistry. It has predominantly electrostatic nature - high electrostatic potential anisotropy. However to fully unleash the potential of XB in rational drug design fast and robust empirical methods of XB description should be developed. Current approaches rely heavily on ab initio calculation for each molecule studied. Thus fast prediction of electrostatic parameters for description of XB for arbitrary organic molecules is of paramount importance to promptly establish QSAR/QSPR, virtual screening and molecular docking pipelines suitable for today's agile development requirements. The two most promising approaches to describe anisotropic electrostatic models - the extra point (EP) charge model and the multipole expansion (ME) model - were studied on their ability (1) to describe ab initio molecular electrostatic potential (MEP) and (2) to produce parameters that can be predicted for each molecule empirically rather than estimated via ab initio calculations. The reference ab initio MEP was calculated for a set of 730 substituted halobenzenes. Parameters for anisotropic electrostatics of both empirical models (EP and ME) studied were extracted from ab initio MEP. The FreeWilson and Hansch type QSPR models relating XB parameters with aromatic substituents were built and analyzed, providing the guidelines for further development.

Simulation of intramolecular hydrogen bond dynamics in manzamine A as a sensitive test for charge distribution quality.

Subtle balance of inter- and intramolecular hydrogen bond strength in aqueous solutions often governs the structure and dynamics of molecular species used as potential drugs and in supramolecular applications. In silico molecular dynamics study of water solution of manzamine A has been performed with different atomic charges in order to investigate the influence of charge distribution choice on predicting qualitative and quantitative features of the simulated systems. Various well known charge schemes (MK-ESP, RESP, Mulliken, AMI-BCC, Gasteiger-Hückel, Gasteiger-Marsili, MMFF94, and Dynamic Electronegativity Relaxation - DENR) led to qualitatively different pictures of dynamic behavior of the intramolecular hydrogen bond. The reported calculation framework represents a relatively rare case where differences in charge distributions lead to noticeable differences in simulated properties, thus providing a useful test case for force field and charge distribution development, provided high quality experiments are conducted to use as references.

General Purpose Electronegativity Relaxation Charge Models Applied to CoMFA and CoMSIA Study of GSK-3 Inhibitors.

Two fast empirical charge models, Kirchhoff Charge Model (KCM) and Dynamic Electronegativity Relaxation (DENR), had been developed in our laboratory previously for widespread use in drug design research. Both models are based on the electronegativity relaxation principle (Adv. Quantum Chem. 2006, 51, 139-156) and parameterized against ab initio dipole/quadrupole moments and molecular electrostatic potentials, respectively. As 3D QSAR studies comprise one of the most important fields of applied molecular modeling, they naturally have become the first topic to test our charges and thus, indirectly, the assumptions laid down to the charge model theories in a case study. Here these charge models are used in CoMFA and CoMSIA methods and tested on five glycogen synthase kinase 3 (GSK-3) inhibitor datasets, relevant to our current studies, and one steroid dataset. For comparison, eight other different charge models, ab initio through semiempirical and empirical, were tested on the same datasets. The complex analysis including correlation and cross-validation, charges robustness and predictability, as well as visual interpretability of 3D contour maps generated was carried out. As a result, our new electronegativity relaxation-based models both have shown stable results, which in conjunction with other benefits discussed render them suitable for building reliable 3D QSAR models.