Molecular mechanics Poisson–Boltzmann surface area (MM/PBSA) and molecular mechanics generalized Born surface area (MM/GBSA) are arguably very popular methods for binding free energy prediction since ...they are more accurate than most scoring functions of molecular docking and less computationally demanding than alchemical free energy methods. MM/PBSA and MM/GBSA have been widely used in biomolecular studies such as protein folding, protein–ligand binding, protein–protein interaction, etc. In this review, methods to adjust the polar solvation energy and to improve the performance of MM/PBSA and MM/GBSA calculations are reviewed and discussed. The latest applications of MM/GBSA and MM/PBSA in drug design are also presented. This review intends to provide readers with guidance for practically applying MM/PBSA and MM/GBSA in drug design and related research fields.
Abstract
Combustion is a complex chemical system which involves thousands of chemical reactions and generates hundreds of molecular species and radicals during the process. In this work, a neural ...network-based molecular dynamics (MD) simulation is carried out to simulate the benchmark combustion of methane. During MD simulation, detailed reaction processes leading to the creation of specific molecular species including various intermediate radicals and the products are intimately revealed and characterized. Overall, a total of 798 different chemical reactions were recorded and some new chemical reaction pathways were discovered. We believe that the present work heralds the dawn of a new era in which neural network-based reactive MD simulation can be practically applied to simulating important complex reaction systems at ab initio level, which provides atomic-level understanding of chemical reaction processes as well as discovery of new reaction pathways at an unprecedented level of detail beyond what laboratory experiments could accomplish.
The human ether-à-go-go-related gene (hERG) K+ channel plays an important role in cardiac action potentials. The inhibition of the hERG channel may lead to long QT syndrome (LQTS) and even sudden ...cardiac death. Due to severe hERG-related cardiotoxicity, many drugs have been withdrawn from the market. Therefore, it is necessary to estimate the chemical blockade of hERG in the early stage of drug discovery. In this study, we collected 12,850 compounds with hERG inhibition data from the literature and trained a series of hERG blocking classification models based on the MACCS and Morgan fingerprints. A consensus model named HergSPred was generated based on the individual models using voting principles. The accuracy of HergSPred is higher than previous models using identical training and test sets. Moreover, we analyzed the contribution of each input fingerprint to the prediction output to obtain intuitive chemical insights into the hERG inhibition, which allows visualization of warning substructures that may cause cardiotoxicity in the input compound. The model is available at http://www.icdrug.com/ICDrug/T.
COVID-19 has emerged as the most serious international pandemic in early 2020 and the lack of comprehensive knowledge in the recognition and transmission mechanisms of this virus hinders the ...development of suitable therapeutic strategies. The specific recognition during the binding of the spike glycoprotein (S protein) of coronavirus to the angiotensin-converting enzyme 2 (ACE2) in the host cell is widely considered the first step of infection. However, detailed insights on the underlying mechanism of dynamic recognition and binding of these two proteins remain unknown. In this work, molecular dynamics simulation and binding free energy calculation were carried out to systematically compare and analyze the receptor-binding domain (RBD) of six coronavirus’ S proteins. We found that affinity and stability of the RBD from SARS-CoV-2 under the binding state with ACE2 are stronger than those of other coronaviruses. The solvent-accessible surface area (SASA) and binding free energy of different RBD subunits indicate an “anchor-locker” recognition mechanism involved in the binding of the S protein to ACE2. Loop 2 (Y473-F490) acts as an anchor for ACE2 recognition, and Loop 3 (G496-V503) locks ACE2 at the other nonanchoring end. Then, the charged or long-chain residues in the β-sheet 1 (N450-F456) region reinforce this binding. The proposed binding mechanism was supported by umbrella sampling simulation of the dissociation process. The current computational study provides important theoretical insights for the development of new vaccines against SARS-CoV-2.
The molecular mechanics/Poisson-Boltzmann surface area (MM/PBSA) method is constantly used to calculate the binding free energy of protein-ligand complexes, and has been shown to effectively balance ...computational cost against accuracy. The relative binding affinities obtained by the MM/PBSA approach are acceptable, while it usually overestimates the absolute binding free energy. This paper proposes four free energy estimators based on the MM/PBSA for enthalpy change combined with interaction entropy (IE) for entropy change using different weights for individual energy terms. The Δ
G
PBSA_IE
method is determined to be an optimal estimator based on its performance in terms of the correlation between experimental and theoretical values and error estimations. This approach is optimized using high-quality experimental values from a training set containing 84 protein-ligand systems, and the coefficients for the sum of electrostatic energy and polar solvation free energy, van der Waals (vdW) energy, non-polar solvation energy and entropy change are obtained by multivariate linear fitting to the corresponding experimental values. A comparison between the traditional MM/PBSA method and this method shows that the correlation coefficient is improved from 0.46 to 0.72 and the slope of the regression line increases from 0.10 to 1.00. More importantly, the mean absolute error (MAE) is significantly reduced from 22.52 to 1.59 kcal mol
−1
. Furthermore, the numerical stability of this method is validated on a test set with a similar correlation coefficient, slope and MAE to those of the training set. Based on the above advantages, the Δ
G
PBSA_IE
method can be a powerful tool for a reliable and accurate estimation of binding free energy and plays a significant role in a detailed energetic investigation of protein-ligand interaction.
Modifying the energy term and considering the entropic contribution by IE method significantly improve the accuracy of predicted binding free energy in MM/PBSA method.
pKa is an important property in the lead optimization process since the charge state of a molecule in physiologic pH plays a critical role in its biological activity, solubility, membrane ...permeability, metabolism, and toxicity. Accurate and fast estimation of small molecule pKa is vital during the drug discovery process. We present MolGpKa, a web server for pKa prediction using a graph-convolutional neural network model. The model works by learning pKa related chemical patterns automatically and building reliable predictors with learned features. ACD/pKa data for 1.6 million compounds from the ChEMBL database was used for model training. We found that the performance of the model is better than machine learning models built with human-engineered fingerprints. Detailed analysis shows that the substitution effect on pKa is well learned by the model. MolGpKa is a handy tool for the rapid estimation of pKa during the ligand design process. The MolGpKa server is freely available to researchers and can be accessed at https://xundrug.cn/molgpka.
Reactive molecular dynamics (MD) simulation makes it possible to study the reaction mechanism of complex reaction systems at the atomic level. However, the analysis of MD trajectories which contain ...thousands of species and reaction pathways has become a major obstacle to the application of reactive MD simulation in large-scale systems. Here, we report the development and application of the Reaction Network Generator (ReacNetGenerator) method. It can automatically extract the reaction network from the reaction trajectory without any predefined reaction coordinates and elementary reaction steps. Molecular species can be automatically identified from the cartesian coordinates of atoms and the hidden Markov model is used to filter the trajectory noises which makes the analysis process easier and more accurate. The ReacNetGenerator has been successfully used to analyze the reactive MD trajectories of the combustion of methane and 4-component surrogate fuel for rocket propellant 3 (RP-3), and it has great advantages in terms of efficiency and accuracy compared to traditional manual analysis.
The ReacNetGenerator program can automatically extract reaction information from the reactive MD trajectory and construct reaction networks.
Efficient and reliable calculation of protein–ligand binding free energy is a grand challenge in computational biology and is of critical importance in drug design and many other molecular ...recognition problems. The main challenge lies in the calculation of entropic contribution to protein–ligand binding or interaction systems. In this report, we present a new interaction entropy method which is theoretically rigorous, computationally efficient, and numerically reliable for calculating entropic contribution to free energy in protein–ligand binding and other interaction processes. Drastically different from the widely employed but extremely expensive normal mode method for calculating entropy change in protein–ligand binding, the new method calculates the entropic component (interaction entropy or −TΔS) of the binding free energy directly from molecular dynamics simulation without any extra computational cost. Extensive study of over a dozen randomly selected protein–ligand binding systems demonstrated that this interaction entropy method is both computationally efficient and numerically reliable and is vastly superior to the standard normal mode approach. This interaction entropy paradigm introduces a novel and intuitive conceptual understanding of the entropic effect in protein–ligand binding and other general interaction systems as well as a practical method for highly efficient calculation of this effect.
Entropy effects play an important role in drug-target interactions, but the entropic contribution to ligand-binding affinity is often neglected by end-point binding free energy calculation methods, ...such as MM/GBSA and MM/PBSA, due to the expensive computational cost of normal mode analysis (NMA). Here, we systematically investigated entropy effects on the prediction power of MM/GBSA and MM/PBSA using >1500 protein-ligand systems and six representative AMBER force fields. Two computationally efficient methods, including NMA based on truncated structures and the interaction entropy approach, were used to estimate the entropic contributions to ligand-target binding free energies. In terms of the overall accuracy, we found that, for the minimized structures, in most cases the inclusion of the conformational entropies predicted by truncated NMA (enthalpynmode_min_9Å) compromises the overall accuracy of MM/GBSA and MM/PBSA compared with the enthalpies calculated based on the minimized structures (enthalpymin). However, for the MD trajectories, the binding free energies can be improved by the inclusion of the conformation entropies predicted by either truncated-NMA for a relatively high dielectric constant (εin = 4) or the interaction entropy method for εin = 1-4. In terms of reproducing the absolute binding free energies, the binding free energies estimated by including the truncated-NMA entropies based on the MD trajectories (ΔGnmode_md_9Å) give the lowest average absolute deviations against the experimental data among all the tested strategies for both MM/GBSA and MM/PBSA. Although the inclusion of the truncated NMA based on the MD trajectories (ΔGnmode_md_9Å) for a relatively high dielectric constant gave the overall best result and the lowest average absolute deviations against the experimental data (for the ff03 force field), it needs too much computational time. Alternatively, considering that the interaction entropy method does not incur any additional computational cost and can give comparable (at high dielectric constant, εin = 4) or even better (at low dielectric constant, εin = 1-2) results than the truncated-NMA entropy (ΔGnmode_md_9Å), the interaction entropy approach is recommended to estimate the entropic component for MM/GBSA and MM/PBSA based on MD trajectories, especially for a diverse dataset. Furthermore, we compared the predictions of MM/GBSA with six different AMBER force fields. The results show that the ff03 force field (ff03 for proteins and gaff with AM1-BCC charges for ligands) performs the best, but the predictions given by the tested force fields are comparable, implying that the MM/GBSA predictions are not very sensitive to force fields.
In recent years, machine-learning-based scoring functions have significantly improved the scoring power. However, many of these methods do not perform well in distinguishing the native structure from ...docked decoy poses due to the lack of decoy structural information in their training data. Here, we developed a machine-learning model, named DeepBSP, that can directly predict the root mean square deviation (rmsd) of a ligand docking pose with reference to its native binding pose. Unlike the binding affinity, the rmsd between the docking poses with reference to their native structures can be straightforwardly determined. By training on a generated data set with 11,925 native complexes and more than 165,000 docked poses, our model shows excellent docking power on our test set and also on the CASF-2016 docking decoy set compared to other major scoring functions. Thus, by combining molecular dockings that generate many poses with the application of DeepBSP, one can more accurately predict the best binding pose that is closest to the native complex structure. This DeepBSP model shall be very useful in picking out poses close to their natives from many poses generated from a dock application.