Protein‐specific force field derived from the fragment molecular orbital method can improve protein–ligand binding interactions

Accurate computational estimate of the protein–ligand binding affinity is of central importance in rational drug design. To improve accuracy of the molecular mechanics (MM) force field (FF) for protein–ligand simulations, we use a protein‐specific FF derived by the fragment molecular orbital (FMO) method and by the restrained electrostatic potential (RESP) method. Applying this FMO‐RESP method to two proteins, dodecin, and lysozyme, we found that protein‐specific partial charges tend to differ more significantly from the standard AMBER charges for isolated charged atoms. We did not see the dependence of partial charges on the secondary structure. Computing the binding affinities of dodecin with five ligands by MM PBSA protocol with the FMO‐RESP charge set as well as with the standard AMBER charges, we found that the former gives better correlation with experimental affinities than the latter. While, for lysozyme with five ligands, both charge sets gave similar and relatively accurate estimates of binding affinities. © 2013 Wiley Periodicals, Inc.     

Three‐body expansion of the fragment molecular orbital method combined with density‐functional tight‐binding

The three‐body fragment molecular orbital (FMO3) method is formulated for density‐functional tight‐binding (DFTB). The energy, analytic gradient, and Hessian are derived in the gas phase, and the energy and analytic gradient are also derived for polarizable continuum model. The accuracy of FMO3‐DFTB is evaluated for five proteins, sodium cation in explicit solvent, and three isomers of polyalanine. It is shown that FMO3‐DFTB is considerably more accurate than FMO2‐DFTB. Molecular dynamics simulations for sodium cation in water are performed for 100 ps, yielding radial distribution functions and coordination numbers. © 2017 Wiley Periodicals, Inc.     

A complete set of self‐consistent charge density‐functional tight‐binding parametrization of zinc chalcogenides (ZnX; X=O,S, Se,and Te)

We have developed a complete set of self‐consistent charge density‐functional tight‐binding parameters for Zn? X (X = Zn, O, S, Se, Te, Cd, H, C, and N). The transferability of the derived parameters has been tested against Pseudo Potential‐Perdew, Burke and Ernzerhof (PP‐PBE) calculations and experimental values (whenever available) for corresponding bulk systems (e.g., hexagonal close packing, zinc‐blende, and wurtzite(wz)), various kinds of nanostructures (such as nanowires, surfaces, and nanoclusters), and also some small molecular systems. Our results show that the derived parameters reproduce the structural and energetic properties of the above‐mentioned systems very well. With the derived parameter set, one can study zinc‐chalcogenide nanostructures of relatively large size which was otherwise prohibited by other methods. The Zn‐Cd parametrization developed in this article will help in studying large semiconductor hetero‐nanostructures of Zn and Cd chalcogenides such as ZnX/CdX core/shell nanoparticles, nanotubes, nanowires, and nanoalloys. © 2012 Wiley Periodicals, Inc.     

VoteDock: Consensus docking method for prediction of protein–ligand interactions

Molecular recognition plays a fundamental role in all biological processes, and that is why great efforts have been made to understand and predict protein–ligand interactions. Finding a molecule that can potentially bind to a target protein is particularly essential in drug discovery and still remains an expensive and time‐consuming task. In silico, tools are frequently used to screen molecular libraries to identify new lead compounds, and if protein structure is known, various protein–ligand docking programs can be used. The aim of docking procedure is to predict correct poses of ligand in the binding site of the protein as well as to score them according to the strength of interaction in a reasonable time frame. The purpose of our studies was to present the novel consensus approach to predict both protein–ligand complex structure and its corresponding binding affinity. Our method used as the input the results from seven docking programs (Surflex, LigandFit, Glide, GOLD, FlexX, eHiTS, and AutoDock) that are widely used for docking of ligands. We evaluated it on the extensive benchmark dataset of 1300 protein–ligands pairs from refined PDBbind database for which the structural and affinity data was available. We compared independently its ability of proper scoring and posing to the previously proposed methods. In most cases, our method is able to dock properly approximately 20% of pairs more than docking methods on average, and over 10% of pairs more than the best single program. The RMSD value of the predicted complex conformation versus its native one is reduced by a factor of 0.5 Å. Finally, we were able to increase the Pearson correlation of the predicted binding affinity in comparison with the experimental value up to 0.5. © 2010 Wiley Periodicals, Inc. J Comput Chem 32: 568–581, 2011     

Rapid QM/MM approach for biomolecular systems under periodic boundary conditions: Combination of the density‐functional tight‐binding theory and particle mesh Ewald method

A quantum mechanical/molecular mechanical (QM/MM) approach based on the density‐functional tight‐binding (DFTB) theory is a useful tool for analyzing chemical reaction systems in detail. In this study, an efficient QM/MM method is developed by the combination of the DFTB/MM and particle mesh Ewald (PME) methods. Because the Fock matrix, which is required in the DFTB calculation, is analytically obtained by the PME method, the Coulomb energy is accurately and rapidly computed. For assessing the performance of this method, DFTB/MM calculations and molecular dynamics simulation are conducted for a system consisting of two amyloid‐β(1‐16) peptides and a zinc ion in explicit water under periodic boundary conditions. As compared with that of the conventional Ewald summation method, the computational cost of the Coulomb energy by utilizing the present approach is drastically reduced, i.e., 166.5 times faster. Furthermore, the deviation of the electronic energy is less than . © 2016 Wiley Periodicals, Inc.     

Parallel implementation of efficient charge–charge interaction evaluation scheme in periodic divide‐and‐conquer density‐functional tight‐binding calculations

A low‐computational‐cost algorithm and its parallel implementation for periodic divide‐and‐conquer density‐functional tight‐binding (DC‐DFTB) calculations are presented. The developed algorithm enables rapid computation of the interaction between atomic partial charges, which is the bottleneck for applications to large systems, by means of multipole‐ and interpolation‐based approaches for long‐ and short‐range contributions. The numerical errors of energy and forces with respect to the conventional Ewald‐based technique can be under the control of the multipole expansion order, level of unit cell replication, and interpolation grid size. The parallel performance of four different evaluation schemes combining previous approaches and the proposed one are assessed using test calculations of a cubic water box on the K computer. The largest benchmark system consisted of 3,295,500 atoms. DC‐DFTB energy and forces for this system were obtained in only a few minutes when the proposed algorithm was activated and parallelized over 16,000 nodes in the K computer. The high performance using a single node workstation was also confirmed. In addition to liquid water systems, the feasibility of the present method was examined by testing solid systems such as diamond form of carbon, face‐centered cubic form of copper, and rock salt form of sodium chloride. © 2017 Wiley Periodicals, Inc.     

Two‐dimensional replica‐exchange method for predicting protein–ligand binding structures

We have developed a two‐dimensional replica‐exchange method for the prediction of protein–ligand binding structures. The first dimension is the umbrella sampling along the reaction coordinate, which is the distance between a protein binding pocket and a ligand. The second dimension is the solute tempering, in which the interaction between a ligand and a protein and water is weakened. The second dimension is introduced to make a ligand follow the umbrella potential more easily and enhance the binding events, which should improve the sampling efficiency. As test cases, we applied our method to two protein‐ligand complex systems (MDM2 and HSP 90‐alpha). Starting from the configuration in which the protein and the ligand are far away from each other in each system, our method predicted the ligand binding structures in excellent agreement with the experimental data from Protein Data Bank much faster with the improved sampling efficiency than the replica‐exchange umbrella sampling method that we have previously developed. © 2013 Wiley Periodicals, Inc.     

Coarse‐grained molecular dynamics simulations of protein–ligand binding

Coarse‐grained molecular dynamics (CGMD) simulations with the MARTINI force field were performed to reproduce the protein–ligand binding processes. We chose two protein–ligand systems, the levansucrase–sugar (glucose or sucrose), and LinB–1,2‐dichloroethane systems, as target systems that differ in terms of the size and shape of the ligand‐binding pocket and the physicochemical properties of the pocket and the ligand. Spatial distributions of the Coarse‐grained (CG) ligand molecules revealed potential ligand‐binding sites on the protein surfaces other than the real ligand‐binding sites. The ligands bound most strongly to the real ligand‐binding sites. The binding and unbinding rate constants obtained from the CGMD simulation of the levansucrase–sucrose system were approximately 10 times greater than the experimental values; this is mainly due to faster diffusion of the CG ligand in the CG water model. We could obtain dissociation constants close to the experimental values for both systems. Analysis of the ligand fluxes demonstrated that the CG ligand molecules entered the ligand‐binding pockets through specific pathways. The ligands tended to move through grooves on the protein surface. Thus, the CGMD simulations produced reasonable results for the two different systems overall and are useful for studying the protein–ligand binding processes. © 2014 Wiley Periodicals, Inc.     

Three pillars for achieving quantum mechanical molecular dynamics simulations of huge systems: Divide‐and‐conquer,density‐functional tight‐binding,and massively parallel computation

The linear‐scaling divide‐and‐conquer (DC) quantum chemical methodology is applied to the density‐functional tight‐binding (DFTB) theory to develop a massively parallel program that achieves on‐the‐fly molecular reaction dynamics simulations of huge systems from scratch. The functions to perform large scale geometry optimization and molecular dynamics with DC‐DFTB potential energy surface are implemented to the program called DC‐DFTB‐K. A novel interpolation‐based algorithm is developed for parallelizing the determination of the Fermi level in the DC method. The performance of the DC‐DFTB‐K program is assessed using a laboratory computer and the K computer. Numerical tests show the high efficiency of the DC‐DFTB‐K program, a single‐point energy gradient calculation of a one‐million‐atom system is completed within 60 s using 7290 nodes of the K computer. © 2016 Wiley Periodicals, Inc.     

Assessment and acceleration of binding energy calculations for protein–ligand complexes by the fragment molecular orbital method

In the field of drug discovery, it is important to accurately predict the binding affinities between target proteins and drug applicant molecules. Many of the computational methods available for evaluating binding affinities have adopted molecular mechanics‐based force fields, although they cannot fully describe protein–ligand interactions. A noteworthy computational method in development involves large‐scale electronic structure calculations. Fragment molecular orbital (FMO) method, which is one of such large‐scale calculation techniques, is applied in this study for calculating the binding energies between proteins and ligands. By testing the effects of specific FMO calculation conditions (including fragmentation size, basis sets, electron correlation, exchange‐correlation functionals, and solvation effects) on the binding energies of the FK506‐binding protein and 10 ligand complex molecule, we have found that the standard FMO calculation condition, FMO2‐MP2/6‐31G(d), is suitable for evaluating the protein–ligand interactions. The correlation coefficient between the binding energies calculated with this FMO calculation condition and experimental values is determined to be R = 0.77. Based on these results, we also propose a practical scheme for predicting binding affinities by combining the FMO method with the quantitative structure–activity relationship (QSAR) model. The results of this combined method can be directly compared with experimental binding affinities. The FMO and QSAR combined scheme shows a higher correlation with experimental data (R = 0.91). Furthermore, we propose an acceleration scheme for the binding energy calculations using a multilayer FMO method focusing on the protein–ligand interaction distance. Our acceleration scheme, which uses FMO2‐HF/STO‐3G:MP2/6‐31G(d) at R_int = 7.0 Å, reduces computational costs, while maintaining accuracy in the evaluation of binding energy. © 2015 Wiley Periodicals, Inc.     

Possible mechanism of BN fullerene formation from a boron cluster: Density‐functional tight‐binding molecular dynamics simulations

We simulate the formation of a BN fullerene from an amorphous B cluster at 2000 K by quantum mechanical molecular dynamics based on the density‐functional tight‐binding method. We run 30 trajectories 200 ps in length, where N atoms are supplied around the target cluster, which is initially an amorphous B₃₆ cluster. Most of the incident N atoms are promptly incorporated into the target cluster to form B‐N‐B bridges or NB₃ pyramidal local substructures. BN fullerene formation is initiated by alternating BN ring condensation. Spontaneous atomic rearrangement and N₂ dissociation lead to the construction of an sp² single‐shelled structure, during which the BN cluster undergoes a transition from a liquid‐like to a solid‐like state. Continual atomic rearrangement and sporadic N₂ dissociation decrease the number of defective rings in the BN cluster and increase the number of six‐membered rings, forming a more regular shell structure. The number of four‐membered rings tends to remain constant, and contributes to more ordered isolated‐tetragon‐rule ring placement. © 2016 Wiley Periodicals, Inc.     

SPOT‐Ligand: Fast and effective structure‐based virtual screening by binding homology search according to ligand and receptor similarity

Structure‐based virtual screening usually involves docking of a library of chemical compounds onto the functional pocket of the target receptor so as to discover novel classes of ligands. However, the overall success rate remains low and screening a large library is computationally intensive. An alternative to this “ab initio” approach is virtual screening by binding homology search. In this approach, potential ligands are predicted based on similar interaction pairs (similarity in receptors and ligands). SPOT‐Ligand is an approach that integrates ligand similarity by Tanimoto coefficient and receptor similarity by protein structure alignment program SPalign. The method was found to yield a consistent performance in DUD and DUD‐E docking benchmarks even if model structures were employed. It improves over docking methods (DOCK6 and AUTODOCK Vina) and has a performance comparable to or better than other binding‐homology methods (FINDsite and PoLi) with higher computational efficiency. The server is available at http://sparks-lab.org . © 2016 Wiley Periodicals, Inc.     

Role of weak C—H…O and strong N—H…O intermolecular interactions on the high‐symmetry molecular packing of trans‐cyclohexane‐1,4‐dicarboxamide

An unpredicted fourfold screw N—H…O hydrogen bond C(4) motif in a primary dicarboxamide (trans‐cyclohexane‐1,4‐dicarboxamide, C₈H₁₄N₂O₂) was investigated by single‐crystal X‐ray diffraction and IR and Raman spectroscopies. Electron‐density topology and intermolecular energy analyses determined from ab initio calculations were employed to examine the influence of weak C—H…O hydrogen‐bond interactions on the peculiar arrangement of molecules in the tetragonal P4₃2₁2 space group. In addition, the way in which the co‐operative effects of those weak bonds might modify their relative influence on molecular packing was estimated from cluster calculations. Based on the results, a structural model is proposed which helps to rationalize the unusual fourfold screw molecular arrangement.     

Molecular structure, IR spectra of 2-mercaptobenzothiazole and 2-mercaptobenzoxazole by density functional theory and ab initio Hartree–Fock calculations

Quantum chemistry calculations have been performed using Gaussian03 program to compute optimized geometry, harmonic vibrational frequency along with intensities in IR and Raman spectra and atomic charges at RHF/6-31+G*, B3LYP/6-31+G* and B3LYP/6-31++G* levels for 2-mercaptobenzothiazole (MBT, C₇H₅NS₂) and 2-mercaptobenzoxazole (MBO, C₇H₅NOS) in the ground state. The scaled harmonic vibrational frequencies have been compared with experimental FT-IR and FT-Raman spectra. The results show that the scaled theoretical vibrational frequencies is very good agreement with the experimental values. A detailed interpretation of the infrared and Raman spectra of 2-mercaptobenzothiazole and 2-mercaptobenzoxazole was reported. Comparison of calculated spectra with the experimental spectra provides important information about the ability of the computational method to describe the vibrational modes.     

Roles of K151 and D180 in L‐2‐haloacid dehalogenase from Pseudomonas sp. YL: Analysis by molecular dynamics and ab initio fragment molecular orbital calculations

L ‐2‐haloacid dehalogenase (L ‐DEX) catalyzes the hydrolytic dehalogenation of L ‐2‐haloalkanoic acids to produce the corresponding D ‐2‐hydroxyalkanoic acids. This enzyme is expected to be applicable to the bioremediation of environments contaminated with halogenated organic compounds. We analyzed the reaction mechanism of L ‐DEX from Pseudomonas sp. YL (L ‐DEX YL) by using molecular modeling. The complexes of wild‐type L ‐DEX YL and its K151A and D180A mutants with its typical substrate, L ‐2‐chloropropionate, were constructed by docking simulation. Subsequently, molecular dynamics (MD) and ab initio fragment molecular orbital (FMO) calculations of the complexes were performed. The ab initio FMO method was applied at the MP2/6‐31G level to estimate interfragment interaction energies. K151 and D180, which are experimentally shown to be important for enzyme activity, interact particularly strongly with L ‐2‐chloropropionate, catalytic water, nucleophile (D10), and with each other. Our calculations suggest that K151 stabilizes substrate orientation and balances the charge around the active site, while D180 stabilizes the rotation of the nucleophile D10, fixes catalytic water around D10, and prevents K151 from approaching D10. Further, D180 may activate catalytic water on its own or with K151, S175, and N177. These roles are consistent with the previous results. Thus, MD and ab initio FMO calculations are powerful tools for the elucidation of the mechanism of enzymatic reaction at the molecular level and can be applied to other catalytically important residues. The results obtained here will play an important role in elucidating the reaction mechanism and rational design of L ‐DEX YL with improved enzymatic activity or substrate specificity. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2009     

A Bayesian statistical approach of improving knowledge‐based scoring functions for protein–ligand interactions

Knowledge‐based scoring functions are widely used for assessing putative complexes in protein–ligand and protein–protein docking and for structure prediction. Even with large training sets, knowledge‐based scoring functions face the inevitable problem of sparse data. Here, we have developed a novel approach for handling the sparse data problem that is based on estimating the inaccuracies in knowledge‐based scoring functions. This inaccuracy estimation is used to automatically weight the knowledge‐based scoring function with an alternative, force‐field‐based potential (FFP) that does not rely on training data and can, therefore, provide an improved approximation of the interactions between rare chemical groups. The current version of STScore, a protein–ligand scoring function using our method, achieves a binding mode prediction success rate of 91% on the set of 100 complexes by Wang et al., and a binding affinity correlation of 0.514 with the experimentally determined affinities in PDBbind. The method presented here may be used with other FFPs and other knowledge‐based scoring functions and can also be applied to protein–protein docking and protein structure prediction. © 2014 Wiley Periodicals, Inc.